def test_order_at():
a = Poly(t**4, t)
b = Poly((t**2 + 1)**3*t, t)
c = Poly((t**2 + 1)**6*t, t)
d = Poly((t**2 + 1)**10*t**10, t)
e = Poly((t**2 + 1)**100*t**37, t)
p1 = Poly(t, t)
p2 = Poly(1 + t**2, t)
assert order_at(a, p1, t) == 4
assert order_at(b, p1, t) == 1
assert order_at(c, p1, t) == 1
assert order_at(d, p1, t) == 10
assert order_at(e, p1, t) == 37
assert order_at(a, p2, t) == 0
assert order_at(b, p2, t) == 3
assert order_at(c, p2, t) == 6
assert order_at(d, p1, t) == 10
assert order_at(e, p2, t) == 100
assert order_at(Poly(0, t), Poly(t, t), t) == oo
assert order_at_oo(Poly(t**2 - 1, t), Poly(t + 1), t) == \
order_at_oo(Poly(t - 1, t), Poly(1, t), t) == -1
assert order_at_oo(Poly(0, t), Poly(1, t), t) == oo
python类Poly()的实例源码
def test_special_denom():
# TODO: add more tests here
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t, t)]})
assert special_denom(Poly(1, t), Poly(t**2, t), Poly(1, t), Poly(t**2 - 1, t),
Poly(t, t), DE) == \
(Poly(1, t), Poly(t**2 - 1, t), Poly(t**2 - 1, t), Poly(t, t))
# assert special_denom(Poly(1, t), Poly(2*x, t), Poly((1 + 2*x)*t, t), DE) == 1
# issue 3940
# Note, this isn't a very good test, because the denominator is just 1,
# but at least it tests the exp cancellation case
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(-2*x*t0, t0),
Poly(I*k*t1, t1)]})
DE.decrement_level()
assert special_denom(Poly(1, t0), Poly(I*k, t0), Poly(1, t0), Poly(t0, t0),
Poly(1, t0), DE) == \
(Poly(1, t0), Poly(I*k, t0), Poly(t0, t0), Poly(1, t0))
def test_bound_degree():
# Base
DE = DifferentialExtension(extension={'D': [Poly(1, x)]})
assert bound_degree(Poly(1, x), Poly(-2*x, x), Poly(1, x), DE) == 0
# Primitive (see above test_bound_degree_fail)
# TODO: Add test for when the degree bound becomes larger after limited_integrate
# TODO: Add test for db == da - 1 case
# Exp
# TODO: Add tests
# TODO: Add test for when the degree becomes larger after parametric_log_deriv()
# Nonlinear
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t**2 + 1, t)]})
assert bound_degree(Poly(t, t), Poly((t - 1)*(t**2 + 1), t), Poly(1, t), DE) == 0
def test_spde():
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t**2 + 1, t)]})
raises(NonElementaryIntegralException, lambda: spde(Poly(t, t), Poly((t - 1)*(t**2 + 1), t), Poly(1, t), 0, DE))
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t, t)]})
assert spde(Poly(t**2 + x*t*2 + x**2, t), Poly(t**2/x**2 + (2/x - 1)*t, t),
Poly(t**2/x**2 + (2/x - 1)*t, t), 0, DE) == \
(Poly(0, t), Poly(0, t), 0, Poly(0, t), Poly(1, t))
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t0/x**2, t0), Poly(1/x, t)]})
assert spde(Poly(t**2, t), Poly(-t**2/x**2 - 1/x, t),
Poly((2*x - 1)*t**4 + (t0 + x)/x*t**3 - (t0 + 4*x**2)/(2*x)*t**2 + x*t, t), 3, DE) == \
(Poly(0, t), Poly(0, t), 0, Poly(0, t),
Poly(t0*t**2/2 + x**2*t**2 - x**2*t, t))
DE = DifferentialExtension(extension={'D': [Poly(1, x)]})
assert spde(Poly(x**2 + x + 1, x), Poly(-2*x - 1, x), Poly(x**5/2 +
3*x**4/4 + x**3 - x**2 + 1, x), 4, DE) == \
(Poly(0, x), Poly(x/2 - S(1)/4, x), 2, Poly(x**2 + x + 1, x), Poly(5*x/4, x))
assert spde(Poly(x**2 + x + 1, x), Poly(-2*x - 1, x), Poly(x**5/2 +
3*x**4/4 + x**3 - x**2 + 1, x), n, DE) == \
(Poly(0, x), Poly(x/2 - S(1)/4, x), -2 + n, Poly(x**2 + x + 1, x), Poly(5*x/4, x))
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(1, t)]})
raises(NonElementaryIntegralException, lambda: spde(Poly((t - 1)*(t**2 + 1)**2, t), Poly((t - 1)*(t**2 + 1), t), Poly(1, t), 0, DE))
def test_solve_poly_rde_no_cancel():
# deg(b) large
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(1 + t**2, t)]})
assert solve_poly_rde(Poly(t**2 + 1, t), Poly(t**3 + (x + 1)*t**2 + t + x + 2, t),
oo, DE) == Poly(t + x, t)
# deg(b) small
DE = DifferentialExtension(extension={'D': [Poly(1, x)]})
assert solve_poly_rde(Poly(0, x), Poly(x/2 - S(1)/4, x), oo, DE) == \
Poly(x**2/4 - x/4, x)
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t**2 + 1, t)]})
assert solve_poly_rde(Poly(2, t), Poly(t**2 + 2*t + 3, t), 1, DE) == \
Poly(t + 1, t, x)
# deg(b) == deg(D) - 1
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t**2 + 1, t)]})
assert no_cancel_equal(Poly(1 - t, t),
Poly(t**3 + t**2 - 2*x*t - 2*x, t), oo, DE) == \
(Poly(t**2, t), 1, Poly((-2 - 2*x)*t - 2*x, t))
def test_is_deriv_k():
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(1/x, t1), Poly(1/(x + 1), t2)],
'L_K': [1, 2], 'E_K': [], 'L_args': [x, x + 1], 'E_args': []})
assert is_deriv_k(Poly(2*x**2 + 2*x, t2), Poly(1, t2), DE) == \
([(t1, 1), (t2, 1)], t1 + t2, 2)
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(1/x, t1), Poly(t2, t2)],
'L_K': [1], 'E_K': [2], 'L_args': [x], 'E_args': [x]})
assert is_deriv_k(Poly(x**2*t2**3, t2), Poly(1, t2), DE) == \
([(x, 3), (t1, 2)], 2*t1 + 3*x, 1)
# TODO: Add more tests, including ones with exponentials
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(2/x, t1)],
'L_K': [1], 'E_K': [], 'L_args': [x**2], 'E_args': []})
assert is_deriv_k(Poly(x, t1), Poly(1, t1), DE) == \
([(t1, S(1)/2)], t1/2, 1)
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(2/(1 + x), t0)],
'L_K': [1], 'E_K': [], 'L_args': [x**2 + 2*x + 1], 'E_args': []})
assert is_deriv_k(Poly(1 + x, t0), Poly(1, t0), DE) == \
([(t0, S(1)/2)], t0/2, 1)
def test_legendre_poly():
raises(ValueError, lambda: legendre_poly(-1, x))
assert legendre_poly(1, x, polys=True) == Poly(x)
assert legendre_poly(0, x) == 1
assert legendre_poly(1, x) == x
assert legendre_poly(2, x) == Q(3, 2)*x**2 - Q(1, 2)
assert legendre_poly(3, x) == Q(5, 2)*x**3 - Q(3, 2)*x
assert legendre_poly(4, x) == Q(35, 8)*x**4 - Q(30, 8)*x**2 + Q(3, 8)
assert legendre_poly(5, x) == Q(63, 8)*x**5 - Q(70, 8)*x**3 + Q(15, 8)*x
assert legendre_poly(6, x) == Q(
231, 16)*x**6 - Q(315, 16)*x**4 + Q(105, 16)*x**2 - Q(5, 16)
assert legendre_poly(1).dummy_eq(x)
assert legendre_poly(1, polys=True) == Poly(x)
def test_laguerre_poly():
raises(ValueError, lambda: laguerre_poly(-1, x))
assert laguerre_poly(1, x, polys=True) == Poly(-x + 1)
assert laguerre_poly(0, x) == 1
assert laguerre_poly(1, x) == -x + 1
assert laguerre_poly(2, x) == Q(1, 2)*x**2 - Q(4, 2)*x + 1
assert laguerre_poly(3, x) == -Q(1, 6)*x**3 + Q(9, 6)*x**2 - Q(18, 6)*x + 1
assert laguerre_poly(4, x) == Q(
1, 24)*x**4 - Q(16, 24)*x**3 + Q(72, 24)*x**2 - Q(96, 24)*x + 1
assert laguerre_poly(5, x) == -Q(1, 120)*x**5 + Q(25, 120)*x**4 - Q(
200, 120)*x**3 + Q(600, 120)*x**2 - Q(600, 120)*x + 1
assert laguerre_poly(6, x) == Q(1, 720)*x**6 - Q(36, 720)*x**5 + Q(450, 720)*x**4 - Q(2400, 720)*x**3 + Q(5400, 720)*x**2 - Q(4320, 720)*x + 1
assert laguerre_poly(0, x, a) == 1
assert laguerre_poly(1, x, a) == -x + a + 1
assert laguerre_poly(2, x, a) == x**2/2 + (-a - 2)*x + a**2/2 + 3*a/2 + 1
assert laguerre_poly(3, x, a) == -x**3/6 + (a/2 + Q(
3)/2)*x**2 + (-a**2/2 - 5*a/2 - 3)*x + a**3/6 + a**2 + 11*a/6 + 1
assert laguerre_poly(1).dummy_eq(-x + 1)
assert laguerre_poly(1, polys=True) == Poly(-x + 1)
def leading_term(expr, *args):
"""
Returns the leading term in a sympy Polynomial.
expr: sympy.Expr
any sympy expression
args: list
list of input symbols for univariate/multivariate polynomials
"""
expr = sympify(str(expr))
P = Poly(expr, *args)
return LT(P)
# ...
# ...
def test_order_at():
a = Poly(t**4, t)
b = Poly((t**2 + 1)**3*t, t)
c = Poly((t**2 + 1)**6*t, t)
d = Poly((t**2 + 1)**10*t**10, t)
e = Poly((t**2 + 1)**100*t**37, t)
p1 = Poly(t, t)
p2 = Poly(1 + t**2, t)
assert order_at(a, p1, t) == 4
assert order_at(b, p1, t) == 1
assert order_at(c, p1, t) == 1
assert order_at(d, p1, t) == 10
assert order_at(e, p1, t) == 37
assert order_at(a, p2, t) == 0
assert order_at(b, p2, t) == 3
assert order_at(c, p2, t) == 6
assert order_at(d, p1, t) == 10
assert order_at(e, p2, t) == 100
assert order_at(Poly(0, t), Poly(t, t), t) == oo
assert order_at_oo(Poly(t**2 - 1, t), Poly(t + 1), t) == \
order_at_oo(Poly(t - 1, t), Poly(1, t), t) == -1
assert order_at_oo(Poly(0, t), Poly(1, t), t) == oo
def test_special_denom():
# TODO: add more tests here
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t, t)]})
assert special_denom(Poly(1, t), Poly(t**2, t), Poly(1, t), Poly(t**2 - 1, t),
Poly(t, t), DE) == \
(Poly(1, t), Poly(t**2 - 1, t), Poly(t**2 - 1, t), Poly(t, t))
# assert special_denom(Poly(1, t), Poly(2*x, t), Poly((1 + 2*x)*t, t), DE) == 1
# issue 3940
# Note, this isn't a very good test, because the denominator is just 1,
# but at least it tests the exp cancellation case
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(-2*x*t0, t0),
Poly(I*k*t1, t1)]})
DE.decrement_level()
assert special_denom(Poly(1, t0), Poly(I*k, t0), Poly(1, t0), Poly(t0, t0),
Poly(1, t0), DE) == \
(Poly(1, t0), Poly(I*k, t0), Poly(t0, t0), Poly(1, t0))
def test_bound_degree():
# Base
DE = DifferentialExtension(extension={'D': [Poly(1, x)]})
assert bound_degree(Poly(1, x), Poly(-2*x, x), Poly(1, x), DE) == 0
# Primitive (see above test_bound_degree_fail)
# TODO: Add test for when the degree bound becomes larger after limited_integrate
# TODO: Add test for db == da - 1 case
# Exp
# TODO: Add tests
# TODO: Add test for when the degree becomes larger after parametric_log_deriv()
# Nonlinear
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t**2 + 1, t)]})
assert bound_degree(Poly(t, t), Poly((t - 1)*(t**2 + 1), t), Poly(1, t), DE) == 0
def test_spde():
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t**2 + 1, t)]})
raises(NonElementaryIntegralException, lambda: spde(Poly(t, t), Poly((t - 1)*(t**2 + 1), t), Poly(1, t), 0, DE))
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t, t)]})
assert spde(Poly(t**2 + x*t*2 + x**2, t), Poly(t**2/x**2 + (2/x - 1)*t, t),
Poly(t**2/x**2 + (2/x - 1)*t, t), 0, DE) == \
(Poly(0, t), Poly(0, t), 0, Poly(0, t), Poly(1, t))
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t0/x**2, t0), Poly(1/x, t)]})
assert spde(Poly(t**2, t), Poly(-t**2/x**2 - 1/x, t),
Poly((2*x - 1)*t**4 + (t0 + x)/x*t**3 - (t0 + 4*x**2)/(2*x)*t**2 + x*t, t), 3, DE) == \
(Poly(0, t), Poly(0, t), 0, Poly(0, t),
Poly(t0*t**2/2 + x**2*t**2 - x**2*t, t))
DE = DifferentialExtension(extension={'D': [Poly(1, x)]})
assert spde(Poly(x**2 + x + 1, x), Poly(-2*x - 1, x), Poly(x**5/2 +
3*x**4/4 + x**3 - x**2 + 1, x), 4, DE) == \
(Poly(0, x), Poly(x/2 - S(1)/4, x), 2, Poly(x**2 + x + 1, x), Poly(5*x/4, x))
assert spde(Poly(x**2 + x + 1, x), Poly(-2*x - 1, x), Poly(x**5/2 +
3*x**4/4 + x**3 - x**2 + 1, x), n, DE) == \
(Poly(0, x), Poly(x/2 - S(1)/4, x), -2 + n, Poly(x**2 + x + 1, x), Poly(5*x/4, x))
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(1, t)]})
raises(NonElementaryIntegralException, lambda: spde(Poly((t - 1)*(t**2 + 1)**2, t), Poly((t - 1)*(t**2 + 1), t), Poly(1, t), 0, DE))
DE = DifferentialExtension(extension={'D': [Poly(1, x)]})
assert spde(Poly(x**2 - x, x), Poly(1, x), Poly(9*x**4 - 10*x**3 + 2*x**2, x), 4, DE) == (Poly(0, x), Poly(0, x), 0, Poly(0, x), Poly(3*x**3 - 2*x**2, x))
assert spde(Poly(x**2 - x, x), Poly(x**2 - 5*x + 3, x), Poly(x**7 - x**6 - 2*x**4 + 3*x**3 - x**2, x), 5, DE) == \
(Poly(1, x), Poly(x + 1, x), 1, Poly(x**4 - x**3, x), Poly(x**3 - x**2, x))
def test_solve_poly_rde_cancel():
# exp
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t, t)]})
assert cancel_exp(Poly(2*x, t), Poly(2*x, t), 0, DE) == \
Poly(1, t)
assert cancel_exp(Poly(2*x, t), Poly((1 + 2*x)*t, t), 1, DE) == \
Poly(t, t)
# TODO: Add more exp tests, including tests that require is_deriv_in_field()
# primitive
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(1/x, t)]})
# If the DecrementLevel context manager is working correctly, this shouldn't
# cause any problems with the further tests.
raises(NonElementaryIntegralException, lambda: cancel_primitive(Poly(1, t), Poly(t, t), oo, DE))
assert cancel_primitive(Poly(1, t), Poly(t + 1/x, t), 2, DE) == \
Poly(t, t)
assert cancel_primitive(Poly(4*x, t), Poly(4*x*t**2 + 2*t/x, t), 3, DE) == \
Poly(t**2, t)
# TODO: Add more primitive tests, including tests that require is_deriv_in_field()
def test_is_log_deriv_k_t_radical_in_field():
# NOTE: any potential constant factor in the second element of the result
# doesn't matter, because it cancels in Da/a.
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(1/x, t)]})
assert is_log_deriv_k_t_radical_in_field(Poly(5*t + 1, t), Poly(2*t*x, t), DE) == \
(2, t*x**5)
assert is_log_deriv_k_t_radical_in_field(Poly(2 + 3*t, t), Poly(5*x*t, t), DE) == \
(5, x**3*t**2)
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(-t/x**2, t)]})
assert is_log_deriv_k_t_radical_in_field(Poly(-(1 + 2*t), t),
Poly(2*x**2 + 2*x**2*t, t), DE) == \
(2, t + t**2)
assert is_log_deriv_k_t_radical_in_field(Poly(-1, t), Poly(x**2, t), DE) == \
(1, t)
assert is_log_deriv_k_t_radical_in_field(Poly(1, t), Poly(2*x**2, t), DE) == \
(2, 1/t)
def test_legendre_poly():
raises(ValueError, lambda: legendre_poly(-1, x))
assert legendre_poly(1, x, polys=True) == Poly(x)
assert legendre_poly(0, x) == 1
assert legendre_poly(1, x) == x
assert legendre_poly(2, x) == Q(3, 2)*x**2 - Q(1, 2)
assert legendre_poly(3, x) == Q(5, 2)*x**3 - Q(3, 2)*x
assert legendre_poly(4, x) == Q(35, 8)*x**4 - Q(30, 8)*x**2 + Q(3, 8)
assert legendre_poly(5, x) == Q(63, 8)*x**5 - Q(70, 8)*x**3 + Q(15, 8)*x
assert legendre_poly(6, x) == Q(
231, 16)*x**6 - Q(315, 16)*x**4 + Q(105, 16)*x**2 - Q(5, 16)
assert legendre_poly(1).dummy_eq(x)
assert legendre_poly(1, polys=True) == Poly(x)
def _split_reaction_monom(reaction, species):
"""Split a reaction into separate reactions, one
for each monomial in rate (compare to split_reactions).
:rtype: list of Reactions.
"""
ratenumer, ratedenom = reaction.rate.cancel().as_numer_denom()
ratenumer = ratenumer.expand()
species = map(sp.Symbol, species)
ratendict = sp.Poly(ratenumer, *species).as_dict()
if len(ratendict) > 1:
reactions = []
i = 0
for degrees in ratendict:
i = i + 1
ratenpart = sp.Mul(*[species[r]**degrees[r] for r in range(len(species))]) * ratendict[degrees]
reactions.append(Reaction(reaction.reactionid + "_" + str(i), \
reaction.reactant, \
reaction.product, \
ratenpart / ratedenom))
return reactions
return [reaction]
def Lx(a,n):
x=symbols('x')
return Matrix(n, 1, lambda i,j: Poly((reduce(mul, ((x-a[k] if k!=i else 1) for k in range(0,n)), 1)).expand(basic=True), x))
def main():
setup_matplotlib_rc()
expr = get_CRAM_from_cache(degree, prec)
c = 0.6*degree
# Get the translated approximation on [-1, 1]. This is similar logic from CRAM_exp().
n, d = map(Poly, fraction(expr))
inv = -c*(t + 1)/(t - 1)
p, q = map(lambda i: Poly(i, t), fraction(inv))
n, d = n.transform(p, q), d.transform(p, q)
rat_func = n/d.TC()/(d/d.TC())
rat_func = rat_func.evalf(prec)
plt.clf()
fig, ax = plt.subplots()
fig.set_size_inches(4, 4)
plot_in_terminal(expr - exp(-t), (0, 100), prec=prec, points=points,
axes=ax)
# Hide grid lines
ax.grid(False)
# Hide axes ticks
ax.set_xticks([])
ax.set_yticks([])
plt.axis('off')
plt.tight_layout()
plt.savefig('logo.png', transparent=True, bbox_inches='tight', pad_inches=0)
def test_gosper_normal():
assert gosper_normal(4*n + 5, 2*(4*n + 1)*(2*n + 3), n) == \
(Poly(S(1)/4, n), Poly(n + S(3)/2), Poly(n + S(1)/4))
def test_solve_poly_system():
assert solve_poly_system([x - 1], x) == [(S.One,)]
assert solve_poly_system([y - x, y - x - 1], x, y) is None
assert solve_poly_system([y - x**2, y + x**2], x, y) == [(S.Zero, S.Zero)]
assert solve_poly_system([2*x - 3, 3*y/2 - 2*x, z - 5*y], x, y, z) == \
[(Rational(3, 2), Integer(2), Integer(10))]
assert solve_poly_system([x*y - 2*y, 2*y**2 - x**2], x, y) == \
[(0, 0), (2, -sqrt(2)), (2, sqrt(2))]
assert solve_poly_system([y - x**2, y + x**2 + 1], x, y) == \
[(-I*sqrt(S.Half), -S.Half), (I*sqrt(S.Half), -S.Half)]
f_1 = x**2 + y + z - 1
f_2 = x + y**2 + z - 1
f_3 = x + y + z**2 - 1
a, b = sqrt(2) - 1, -sqrt(2) - 1
assert solve_poly_system([f_1, f_2, f_3], x, y, z) == \
[(0, 0, 1), (0, 1, 0), (1, 0, 0), (a, a, a), (b, b, b)]
solution = [(1, -1), (1, 1)]
assert solve_poly_system([Poly(x**2 - y**2), Poly(x - 1)]) == solution
assert solve_poly_system([x**2 - y**2, x - 1], x, y) == solution
assert solve_poly_system([x**2 - y**2, x - 1]) == solution
assert solve_poly_system(
[x + x*y - 3, y + x*y - 4], x, y) == [(-3, -2), (1, 2)]
raises(NotImplementedError, lambda: solve_poly_system([x**3 - y**3], x, y))
raises(PolynomialError, lambda: solve_poly_system([1/x], x))
def test_weak_normalizer():
a = Poly((1 + x)*t**5 + 4*t**4 + (-1 - 3*x)*t**3 - 4*t**2 + (-2 + 2*x)*t, t)
d = Poly(t**4 - 3*t**2 + 2, t)
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t, t)]})
r = weak_normalizer(a, d, DE, z)
assert r == (Poly(t**5 - t**4 - 4*t**3 + 4*t**2 + 4*t - 4, t),
(Poly((1 + x)*t**2 + x*t, t), Poly(t + 1, t)))
assert weak_normalizer(r[1][0], r[1][1], DE) == (Poly(1, t), r[1])
r = weak_normalizer(Poly(1 + t**2), Poly(t**2 - 1, t), DE, z)
assert r == (Poly(t**4 - 2*t**2 + 1, t), (Poly(-3*t**2 + 1, t), Poly(t**2 - 1, t)))
assert weak_normalizer(r[1][0], r[1][1], DE, z) == (Poly(1, t), r[1])
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(1 + t**2)]})
r = weak_normalizer(Poly(1 + t**2), Poly(t, t), DE, z)
assert r == (Poly(t, t), (Poly(0, t), Poly(1, t)))
assert weak_normalizer(r[1][0], r[1][1], DE, z) == (Poly(1, t), r[1])
def test_bound_degree_fail():
# Primitive
DE = DifferentialExtension(extension={'D': [Poly(1, x),
Poly(t0/x**2, t0), Poly(1/x, t)]})
assert bound_degree(Poly(t**2, t), Poly(-(1/x**2*t**2 + 1/x), t),
Poly((2*x - 1)*t**4 + (t0 + x)/x*t**3 - (t0 + 4*x**2)/2*x*t**2 + x*t,
t), DE) == 3
def test_rischDE():
# TODO: Add more tests for rischDE, including ones from the text
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t, t)]})
DE.decrement_level()
assert rischDE(Poly(-2*x, x), Poly(1, x), Poly(1 - 2*x - 2*x**2, x),
Poly(1, x), DE) == \
(Poly(x + 1, x), Poly(1, x))
def test_prde_normal_denom():
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(1 + t**2, t)]})
fa = Poly(1, t)
fd = Poly(x, t)
G = [(Poly(t, t), Poly(1 + t**2, t)), (Poly(1, t), Poly(x + x*t**2, t))]
assert prde_normal_denom(fa, fd, G, DE) == \
(Poly(x, t), (Poly(1, t), Poly(1, t)), [(Poly(x*t, t),
Poly(t**2 + 1, t)), (Poly(1, t), Poly(t**2 + 1, t))], Poly(1, t))
G = [(Poly(t, t), Poly(t**2 + 2*t + 1, t)), (Poly(x*t, t),
Poly(t**2 + 2*t + 1, t)), (Poly(x*t**2, t), Poly(t**2 + 2*t + 1, t))]
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t, t)]})
assert prde_normal_denom(Poly(x, t), Poly(1, t), G, DE) == \
(Poly(t + 1, t), (Poly((-1 + x)*t + x, t), Poly(1, t)), [(Poly(t, t),
Poly(1, t)), (Poly(x*t, t), Poly(1, t)), (Poly(x*t**2, t),
Poly(1, t))], Poly(t + 1, t))
def test_prde_special_denom():
a = Poly(t + 1, t)
ba = Poly(t**2, t)
bd = Poly(1, t)
G = [(Poly(t, t), Poly(1, t)), (Poly(t**2, t), Poly(1, t)), (Poly(t**3, t), Poly(1, t))]
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t, t)]})
assert prde_special_denom(a, ba, bd, G, DE) == \
(Poly(t + 1, t), Poly(t**2, t), [(Poly(t, t), Poly(1, t)),
(Poly(t**2, t), Poly(1, t)), (Poly(t**3, t), Poly(1, t))], Poly(1, t))
G = [(Poly(t, t), Poly(1, t)), (Poly(1, t), Poly(t, t))]
assert prde_special_denom(Poly(1, t), Poly(t**2, t), Poly(1, t), G, DE) == \
(Poly(1, t), Poly(t**2 - 1, t), [(Poly(t**2, t), Poly(1, t)),
(Poly(1, t), Poly(1, t))], Poly(t, t))
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(-2*x*t0, t0)]})
DE.decrement_level()
G = [(Poly(t, t), Poly(t**2, t)), (Poly(2*t, t), Poly(t, t))]
assert prde_special_denom(Poly(5*x*t + 1, t), Poly(t**2 + 2*x**3*t, t), Poly(t**3 + 2, t), G, DE) == \
(Poly(5*x*t + 1, t), Poly(0, t), [(Poly(t, t), Poly(t**2, t)),
(Poly(2*t, t), Poly(t, t))], Poly(1, x))
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly((t**2 + 1)*2*x, t)]})
G = [(Poly(t + x, t), Poly(t*x, t)), (Poly(2*t, t), Poly(x**2, x))]
assert prde_special_denom(Poly(5*x*t + 1, t), Poly(t**2 + 2*x**3*t, t), Poly(t**3, t), G, DE) == \
(Poly(5*x*t + 1, t), Poly(0, t), [(Poly(t + x, t), Poly(x*t, t)),
(Poly(2*t, t, x), Poly(x**2, t, x))], Poly(1, t))
assert prde_special_denom(Poly(t + 1, t), Poly(t**2, t), Poly(t**3, t), G, DE) == \
(Poly(t + 1, t), Poly(0, t), [(Poly(t + x, t), Poly(x*t, t)), (Poly(2*t, t, x),
Poly(x**2, t, x))], Poly(1, t))
def test_prde_linear_constraints():
DE = DifferentialExtension(extension={'D': [Poly(1, x)]})
G = [(Poly(2*x**3 + 3*x + 1, x), Poly(x**2 - 1, x)), (Poly(1, x), Poly(x - 1, x)),
(Poly(1, x), Poly(x + 1, x))]
assert prde_linear_constraints(Poly(1, x), Poly(0, x), G, DE) == \
((Poly(2*x, x), Poly(0, x), Poly(0, x)), Matrix([[1, 1, -1], [5, 1, 1]]))
G = [(Poly(t, t), Poly(1, t)), (Poly(t**2, t), Poly(1, t)), (Poly(t**3, t), Poly(1, t))]
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t, t)]})
assert prde_linear_constraints(Poly(t + 1, t), Poly(t**2, t), G, DE) == \
((Poly(t, t), Poly(t**2, t), Poly(t**3, t)), Matrix())
G = [(Poly(2*x, t), Poly(t, t)), (Poly(-x, t), Poly(t, t))]
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(1/x, t)]})
prde_linear_constraints(Poly(1, t), Poly(0, t), G, DE) == \
((Poly(0, t), Poly(0, t)), Matrix([[2*x, -x]]))
def test_constant_system():
A = Matrix([[-(x + 3)/(x - 1), (x + 1)/(x - 1), 1],
[-x - 3, x + 1, x - 1],
[2*(x + 3)/(x - 1), 0, 0]])
u = Matrix([(x + 1)/(x - 1), x + 1, 0])
DE = DifferentialExtension(extension={'D': [Poly(1, x)]})
assert constant_system(A, u, DE) == \
(Matrix([[1, 0, 0],
[0, 1, 0],
[0, 0, 0],
[0, 0, 1]]), Matrix([0, 1, 0, 0]))
def test_prde_no_cancel():
# b large
DE = DifferentialExtension(extension={'D': [Poly(1, x)]})
assert prde_no_cancel_b_large(Poly(1, x), [Poly(x**2, x), Poly(1, x)], 2, DE) == \
([Poly(x**2 - 2*x + 2, x), Poly(1, x)], Matrix([[1, 0, -1, 0],
[0, 1, 0, -1]]))
assert prde_no_cancel_b_large(Poly(1, x), [Poly(x**3, x), Poly(1, x)], 3, DE) == \
([Poly(x**3 - 3*x**2 + 6*x - 6, x), Poly(1, x)], Matrix([[1, 0, -1, 0],
[0, 1, 0, -1]]))
# b small
# XXX: Is there a better example of a monomial with D.degree() > 2?
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(t**3 + 1, t)]})
# My original q was t**4 + t + 1, but this solution implies q == t**4
# (c1 = 4), with some of the ci for the original q equal to 0.
G = [Poly(t**6, t), Poly(x*t**5, t), Poly(t**3, t), Poly(x*t**2, t), Poly(1 + x, t)]
assert prde_no_cancel_b_small(Poly(x*t, t), G, 4, DE) == \
([Poly(t**4/4 - x/12*t**3 + x**2/24*t**2 + (-S(11)/12 - x**3/24)*t + x/24, t),
Poly(x/3*t**3 - x**2/6*t**2 + (-S(1)/3 + x**3/6)*t - x/6, t), Poly(t, t),
Poly(0, t), Poly(0, t)], Matrix([[1, 0, -1, 0, 0, 0, 0, 0, 0, 0],
[0, 1, -S(1)/4, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 1, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 1, 0, 0, 0, 0, 0],
[1, 0, 0, 0, 0, -1, 0, 0, 0, 0],
[0, 1, 0, 0, 0, 0, -1, 0, 0, 0],
[0, 0, 1, 0, 0, 0, 0, -1, 0, 0],
[0, 0, 0, 1, 0, 0, 0, 0, -1, 0],
[0, 0, 0, 0, 1, 0, 0, 0, 0, -1]]))
# TODO: Add test for deg(b) <= 0 with b small
def test_limited_integrate_reduce():
DE = DifferentialExtension(extension={'D': [Poly(1, x), Poly(1/x, t)]})
assert limited_integrate_reduce(Poly(x, t), Poly(t**2, t), [(Poly(x, t),
Poly(t, t))], DE) == \
(Poly(t, t), Poly(-1/x, t), Poly(t, t), 1, (Poly(x, t), Poly(1, t)),
[(Poly(-x*t, t), Poly(1, t))])
