95 def weak_normalizer(a, d, DE, z=None):  argument
116     dn, ds = splitfactor(d, DE)
120     g = gcd(dn, dn.diff(DE.t))
124     a1, b = gcdex_diophantine(d.quo(d1).as_poly(DE.t), d1.as_poly(DE.t),
125         a.as_poly(DE.t))
126     r = (a - Poly(z, DE.t)*derivation(d1, DE)).as_poly(DE.t).resultant(
127         d1.as_poly(DE.t))
131         return (Poly(1, DE.t), (a, d))
135     q = reduce(mul, [gcd(a - Poly(n, DE.t)*derivation(d1, DE), d1) for n in N],
136         Poly(1, DE.t))
138     dq = derivation(q, DE)
146 def normal_denom(fa, fd, ga, gd, DE):  argument
162     dn, ds = splitfactor(fd, DE)
163     en, es = splitfactor(gd, DE)
166     h = en.gcd(en.diff(DE.t)).quo(p.gcd(p.diff(DE.t)))
176     ba = a*fa - dn*derivation(h, DE)*fd
183 def special_denom(a, ba, bd, ca, cd, DE, case='auto'):  argument
208         case = DE.case
211         p = Poly(DE.t, DE.t)
213         p = Poly(DE.t**2 + 1, DE.t)
217         return (a, B, C, Poly(1, DE.t))
222     nb = order_at(ba, p, DE.t) - order_at(bd, p, DE.t)
223     nc = order_at(ca, p, DE.t) - order_at(cd, p, DE.t)
230             dcoeff = DE.d.quo(Poly(DE.t, DE.t))
231             with DecrementLevel(DE):  # We are guaranteed to not have problems,
233                 alphaa, alphad = frac_in(-ba.eval(0)/bd.eval(0)/a.eval(0), DE.t)
234                 etaa, etad = frac_in(dcoeff, DE.t)
235                 A = parametric_log_deriv(alphaa, alphad, etaa, etad, DE)
242             dcoeff = DE.d.quo(Poly(DE.t**2+1, DE.t))
243             with DecrementLevel(DE):  # We are guaranteed to not have problems,
245 …          alphaa, alphad = frac_in(im(-ba.eval(sqrt(-1))/bd.eval(sqrt(-1))/a.eval(sqrt(-1))), DE.t)
246 …            betaa, betad = frac_in(re(-ba.eval(sqrt(-1))/bd.eval(sqrt(-1))/a.eval(sqrt(-1))), DE.t)
247                 etaa, etad = frac_in(dcoeff, DE.t)
249                 if recognize_log_derivative(Poly(2, DE.t)*betaa, betad, DE):
250 …A = parametric_log_deriv(alphaa*Poly(sqrt(-1), DE.t)*betad+alphad*betaa, alphad*betad, etaa, etad,…
260     B = ba*pN.quo(bd) + Poly(n, DE.t)*a*derivation(p, DE).quo(p)*pN
268 def bound_degree(a, b, cQ, DE, case='auto', parametric=False):  argument
291         case = DE.case
293     da = a.degree(DE.t)
294     db = b.degree(DE.t)
298         dc = max([i.degree(DE.t) for i in cQ])
300         dc = cQ.degree(DE.t)
302     alpha = cancel(-b.as_poly(DE.t).LC().as_expr()/
303         a.as_poly(DE.t).LC().as_expr())
316         etaa, etad = frac_in(DE.d, DE.T[DE.level - 1])
318         t1 = DE.t
319         with DecrementLevel(DE):
320             alphaa, alphad = frac_in(alpha, DE.t)
325                         DE)
338                 A = is_log_deriv_k_t_radical_in_field(alphaa, alphad, DE)
342                         beta = -(a*derivation(z, DE).as_poly(t1) +
344                         betaa, betad = frac_in(beta, DE.t)
347                                 [(etaa, etad)], DE)
358             etaa, etad = frac_in(DE.d.quo(Poly(DE.t, DE.t)), DE.T[DE.level - 1])
359             with DecrementLevel(DE):
360                 alphaa, alphad = frac_in(alpha, DE.t)
361                 A = parametric_log_deriv(alphaa, alphad, etaa, etad, DE)
370         delta = DE.d.degree(DE.t)
371         lam = DE.d.LC()
384 def spde(a, b, c, n, DE):  argument
401     zero = Poly(0, DE.t)
403     alpha = Poly(1, DE.t)
404     beta = Poly(0, DE.t)
418         if a.degree(DE.t) == 0:
424         b += derivation(a, DE)
425         c = z - derivation(r, DE)
426         n -= a.degree(DE.t)
431 def no_cancel_b_large(b, c, n, DE):  argument
445     q = Poly(0, DE.t)
448         m = c.degree(DE.t) - b.degree(DE.t)
452         p = Poly(c.as_poly(DE.t).LC()/b.as_poly(DE.t).LC()*DE.t**m, DE.t,
456         c = c - derivation(p, DE) - b*p
461 def no_cancel_b_small(b, c, n, DE):  argument
477     q = Poly(0, DE.t)
483             m = c.degree(DE.t) - DE.d.degree(DE.t) + 1
489             p = Poly(c.as_poly(DE.t).LC()/(m*DE.d.as_poly(DE.t).LC())*DE.t**m,
490                 DE.t, expand=False)
492             if b.degree(DE.t) != c.degree(DE.t):
494             if b.degree(DE.t) == 0:
495                 return (q, b.as_poly(DE.T[DE.level - 1]),
496                     c.as_poly(DE.T[DE.level - 1]))
497             p = Poly(c.as_poly(DE.t).LC()/b.as_poly(DE.t).LC(), DE.t,
502         c = c - derivation(p, DE) - b*p
508 def no_cancel_equal(b, c, n, DE):  argument
524     q = Poly(0, DE.t)
525     lc = cancel(-b.as_poly(DE.t).LC()/DE.d.as_poly(DE.t).LC())
532         m = max(M, c.degree(DE.t) - DE.d.degree(DE.t) + 1)
537         u = cancel(m*DE.d.as_poly(DE.t).LC() + b.as_poly(DE.t).LC())
541             p = Poly(c.as_poly(DE.t).LC()/u*DE.t**m, DE.t, expand=False)
543             if c.degree(DE.t) != DE.d.degree(DE.t) - 1:
546                 p = c.as_poly(DE.t).LC()/b.as_poly(DE.t).LC()
550         c = c - derivation(p, DE) - b*p
555 def cancel_primitive(b, c, n, DE):  argument
570     with DecrementLevel(DE):
571         ba, bd = frac_in(b, DE.t)
572         A = is_log_deriv_k_t_radical_in_field(ba, bd, DE)
586     if n < c.degree(DE.t):
589     q = Poly(0, DE.t)
591         m = c.degree(DE.t)
594         with DecrementLevel(DE):
595             a2a, a2d = frac_in(c.LC(), DE.t)
596             sa, sd = rischDE(ba, bd, a2a, a2d, DE)
597         stm = Poly(sa.as_expr()/sd.as_expr()*DE.t**m, DE.t, expand=False)
600         c -= b*stm + derivation(stm, DE)
605 def cancel_exp(b, c, n, DE):  argument
620     eta = DE.d.quo(Poly(DE.t, DE.t)).as_expr()
622     with DecrementLevel(DE):
623         etaa, etad = frac_in(eta, DE.t)
624         ba, bd = frac_in(b, DE.t)
625         A = parametric_log_deriv(ba, bd, etaa, etad, DE)
640     if n < c.degree(DE.t):
643     q = Poly(0, DE.t)
645         m = c.degree(DE.t)
650         with DecrementLevel(DE):
652             a1a, a1d = frac_in(a1, DE.t)
653             a1a = a1a*etad + etaa*a1d*Poly(m, DE.t)
656             a2a, a2d = frac_in(c.LC(), DE.t)
658             sa, sd = rischDE(a1a, a1d, a2a, a2d, DE)
659         stm = Poly(sa.as_expr()/sd.as_expr()*DE.t**m, DE.t, expand=False)
662         c -= b*stm + derivation(stm, DE)  # deg(c) becomes smaller
666 def solve_poly_rde(b, cQ, n, DE, parametric=False):  argument
679     if not b.is_zero and (DE.case == 'base' or
680             b.degree(DE.t) > max(0, DE.d.degree(DE.t) - 1)):
683             return prde_no_cancel_b_large(b, cQ, n, DE)
684         return no_cancel_b_large(b, cQ, n, DE)
686     elif (b.is_zero or b.degree(DE.t) < DE.d.degree(DE.t) - 1) and \
687             (DE.case == 'base' or DE.d.degree(DE.t) >= 2):
690             return prde_no_cancel_b_small(b, cQ, n, DE)
692         R = no_cancel_b_small(b, cQ, n, DE)
699             with DecrementLevel(DE):
700                 b0, c0 = b0.as_poly(DE.t), c0.as_poly(DE.t)
705                 y = solve_poly_rde(b0, c0, n, DE).as_poly(DE.t)
708     elif DE.d.degree(DE.t) >= 2 and b.degree(DE.t) == DE.d.degree(DE.t) - 1 and \
709             n > -b.as_poly(DE.t).LC()/DE.d.as_poly(DE.t).LC():
713         if not b.as_poly(DE.t).LC().is_number:
720         R = no_cancel_equal(b, cQ, n, DE)
727             y = solve_poly_rde(b, C, m, DE)
736             if DE.case == 'exp':
740                 return cancel_exp(b, cQ, n, DE)
742             elif DE.case == 'primitive':
746                 return cancel_primitive(b, cQ, n, DE)
750                     "cases are not yet implemented (%s)." % DE.case)
759 def rischDE(fa, fd, ga, gd, DE):  argument
774     _, (fa, fd) = weak_normalizer(fa, fd, DE)
775     a, (ba, bd), (ca, cd), hn = normal_denom(fa, fd, ga, gd, DE)
776     A, B, C, hs = special_denom(a, ba, bd, ca, cd, DE)
784         n = bound_degree(A, B, C, DE)
792     B, C, m, alpha, beta = spde(A, B, C, n, DE)
796         y = solve_poly_rde(B, C, m, DE)