ore_algebra.ore_operator_1_1

Univariate operators over univariate rings

Special classes for operators living in algebras with one generator with base rings that also have one generator.

Classes

UnivariateDifferenceOperatorOverUnivariateRing(…) Element of an Ore algebra K(x)[F], where F is the forward difference operator F f(x) = f(x+1) - f(x)
UnivariateDifferentialOperatorOverUnivariateRing(…) Element of an Ore algebra K(x)[D], where D acts as derivation d/dx on K(x).
UnivariateEulerDifferentialOperatorOverUnivariateRing(…) Element of an Ore algebra K(x)[T], where T is the Euler differential operator T = x*d/dx
UnivariateOreOperatorOverUnivariateRing(…) Element of an Ore algebra with a single generator and a commutative rational function field as base ring.
UnivariateQDifferentialOperatorOverUnivariateRing(…) Element of an Ore algebra K(x)[J], where J is the Jackson q-differentiation J f(x) = (f(q*x) - f(x))/(q*(x-1))
UnivariateQRecurrenceOperatorOverUnivariateRing(…) Element of an Ore algebra K(x)[S], where S is the shift x->q*x for some q in K.
UnivariateRecurrenceOperatorOverUnivariateRing(…) Element of an Ore algebra K(x)[S], where S is the shift x->x+1.
class ore_algebra.ore_operator_1_1.UnivariateDifferenceOperatorOverUnivariateRing(parent, *data, **kwargs)

Element of an Ore algebra K(x)[F], where F is the forward difference operator F f(x) = f(x+1) - f(x)

indicial_polynomial(*args, **kwargs)

Computes the indicial polynomial of self at (a root of) p.

The indicial polynomial is a polynomial in the given variable var with coefficients in the fraction field of the base ring’s base ring.

The precise meaning of this polynomial may depend on the parent of self. A minimum requirement is that if self has a rational solution whose denominator contains sigma.factorial(p, e) but neither sigma(p, -1)*sigma.factorial(p, e) nor sigma.factorial(p, e + 1), then -e is a root of this polynomial.

Applied to p=1/x, the maximum integer root of the output should serve as a degree bound for the polynomial solutions of self.

This method is a stub. Depending on the particular subclass, restrictions on p may apply.

spread(p=0)

Returns the spread of this operator.

This is the set of integers i such that sigma(self[0], i) and sigma(self[r], -r) have a nontrivial common factor, where sigma is the shift of the parent’s algebra and r is the order of self.

If the optional argument p is given, the method is applied to gcd(self[0], p) instead of self[0].

The output set contains infty if the constant coefficient of self is zero.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']; A.<Sx> = OreAlgebra(R, 'Sx');
sage: ((x+5)*Sx - x).spread()
[4]
sage: ((x+5)*Sx - x).lclm((x+19)*Sx - x).spread()
[3, 4, 17, 18]
symmetric_product(other, solver=None)

Returns the symmetric product of self and other.

The symmetric product of two operators A and B is a minimal order operator C such that for all “functions” f and g with A.f=B.g=0 we have C.(fg)=0.

The method requires that a product rule is associated to the Ore algebra where self and other live. (See docstring of OreAlgebra for information about product rules.)

If no solver is specified, the the Ore algebra’s solver is used.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: A.<Dx> = OreAlgebra(R, 'Dx')
sage: (Dx - 1).symmetric_product(x*Dx - 1)
x*Dx - x - 1
sage: (x*Dx - 1).symmetric_product(Dx - 1)
x*Dx - x - 1
sage: ((x+1)*Dx^2 + (x-1)*Dx + 8).symmetric_product((x-1)*Dx^2 + (2*x+3)*Dx + (8*x+5))
(29*x^8 - 4*x^7 - 55*x^6 - 34*x^5 - 23*x^4 + 80*x^3 + 95*x^2 - 42*x - 46)*Dx^4 + (174*x^8 + 150*x^7 + 48*x^6 - 294*x^5 - 864*x^4 - 646*x^3 + 232*x^2 + 790*x + 410)*Dx^3 + (783*x^8 + 1661*x^7 - 181*x^6 - 1783*x^5 - 3161*x^4 - 3713*x^3 + 213*x^2 + 107*x - 1126)*Dx^2 + (1566*x^8 + 5091*x^7 + 2394*x^6 + 2911*x^5 - 10586*x^4 - 23587*x^3 - 18334*x^2 - 2047*x + 5152)*Dx + 2552*x^8 + 3795*x^7 + 8341*x^6 + 295*x^5 - 6394*x^4 - 24831*x^3 - 35327*x^2 - 23667*x - 13708

sage: A.<Sx> = OreAlgebra(R, 'Sx')
sage: (Sx - 2).symmetric_product(x*Sx - (x+1))
x*Sx - 2*x - 2
sage: (x*Sx - (x+1)).symmetric_product(Sx - 2)
x*Sx - 2*x - 2
sage: ((x+1)*Sx^2 + (x-1)*Sx + 8).symmetric_product((x-1)*Sx^2 + (2*x+3)*Sx + (8*x+5))
(-8*x^8 - 13*x^7 + 300*x^6 + 1640*x^5 + 3698*x^4 + 4373*x^3 + 2730*x^2 + 720*x)*Sx^4 + (16*x^8 + 34*x^7 - 483*x^6 - 1947*x^5 - 2299*x^4 - 2055*x^3 - 4994*x^2 - 4592*x)*Sx^3 + (-64*x^8 + 816*x^7 + 1855*x^6 - 21135*x^5 - 76919*x^4 - 35377*x^3 + 179208*x^2 + 283136*x + 125440)*Sx^2 + (1024*x^7 + 1792*x^6 - 39792*x^5 - 250472*x^4 - 578320*x^3 - 446424*x^2 + 206528*x + 326144)*Sx - 32768*x^6 - 61440*x^5 + 956928*x^4 + 4897984*x^3 + 9390784*x^2 + 7923200*x + 2329600
to_D(alg)

Returns a differential operator which annihilates every power series (about the origin) whose coefficient sequence is annihilated by self. The output operator may not be minimal.

INPUT:

  • alg – the Ore algebra in which the output should be expressed. The algebra must satisfy alg.base_ring().base_ring() == self.base_ring().base_ring() and alg.is_D() is not False. Instead of an algebra object, also a string can be passed as argument. This amounts to specifying an Ore algebra over self.base_ring() with the standard derivation with respect to self.base_ring().gen().

EXAMPLES:

sage: from ore_algebra import *
sage: Rn.<n> = ZZ['n']; Rx.<x> = ZZ['x']
sage: A.<Fn> = OreAlgebra(Rn, 'Fn')
sage: B.<Dx> = OreAlgebra(Rx, 'Dx')
sage: Fn.to_D(B)
(-x + 1)*Dx - 1
sage: ((n+1)*Fn - 1).to_D(B)
(-x^2 + x)*Dx^2 + (-4*x + 1)*Dx - 2
sage: (x*Dx-1).to_F(A).to_D(B)
x*Dx - 1
to_S(alg)

Returns the differential operator corresponding to self

INPUT:

  • alg – the Ore algebra in which the output should be expressed. The algebra must satisfy alg.base_ring().base_ring() == self.base_ring().base_ring() and alg.is_S() is not False. Instead of an algebra object, also a string can be passed as argument. This amounts to specifying an Ore algebra over self.base_ring() with a standard shift with respect to self.base_ring().gen().

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: A.<Fx> = OreAlgebra(R, 'Fx')
sage: (Fx^4).to_S(OreAlgebra(R, 'Sx'))
Sx^4 - 4*Sx^3 + 6*Sx^2 - 4*Sx + 1
sage: (Fx^4).to_S('Sx')
Sx^4 - 4*Sx^3 + 6*Sx^2 - 4*Sx + 1
to_T(alg)

Returns a differential operator, expressed in terms of the Euler derivation, which annihilates every power series (about the origin) whose coefficient sequence is annihilated by self. The output operator may not be minimal.

INPUT:

  • alg – the Ore algebra in which the output should be expressed. The algebra must satisfy alg.base_ring().base_ring() == self.base_ring().base_ring() and alg.is_T() is not False. Instead of an algebra object, also a string can be passed as argument. This amounts to specifying an Ore algebra over self.base_ring() with the Euler derivation with respect to self.base_ring().gen().

EXAMPLES:

sage: from ore_algebra import *
sage: Rn.<n> = ZZ['n']; Rx.<x> = ZZ['x']
sage: A.<Fn> = OreAlgebra(Rn, 'Fn')
sage: B.<Tx> = OreAlgebra(Rx, 'Tx')
sage: Fn.to_T(B)
(-x + 1)*Tx - x
sage: ((n+1)*Fn - 1).to_T(B)
(-x + 1)*Tx^2 - 3*x*Tx - 2*x
sage: (x*Tx-1).to_F(A).to_T(B)
x*Tx^2 + (x - 1)*Tx
to_list(*args, **kwargs)

Computes the terms of some sequence annihilated by self.

INPUT:

  • init – a vector (or list or tuple) of initial values. The components must be elements of self.base_ring().base_ring().fraction_field(). If the length is more than self.order(), we do not check whether the given terms are consistent with self.
  • n – desired number of terms.
  • start (optional) – index of the sequence term which is represented by the first entry of init. Defaults to zero.
  • append (optional) – if True, the computed terms are appended to init list. Otherwise (default), a new list is created.
  • padd (optional) – if True, the vector of initial values is implicitely prolonged to the left (!) by zeros if it is too short. Otherwise (default), the method raises a ValueError if init is too short.

OUTPUT:

A list of n terms whose k th component carries the sequence term with index start+k. Terms whose calculation causes an error are represented by None.

EXAMPLES:

sage: from ore_algebra import *
sage: R = ZZ['x']['n']; x = R('x'); n = R('n')
sage: A.<Sn> = OreAlgebra(R, 'Sn')
sage: L = ((n+2)*Sn^2 - x*(2*n+3)*Sn + (n+1))
sage: L.to_list([1, x], 5)
[1, x, (3*x^2 - 1)/2, (5*x^3 - 3*x)/2, (35*x^4 - 30*x^2 + 3)/8]
sage: polys = L.to_list([1], 5, padd=True)
sage: polys
[1, x, (3*x^2 - 1)/2, (5*x^3 - 3*x)/2, (35*x^4 - 30*x^2 + 3)/8]
sage: L.to_list([polys[3], polys[4]], 8, start=3)
[(5*x^3 - 3*x)/2,
 (35*x^4 - 30*x^2 + 3)/8,
 (63*x^5 - 70*x^3 + 15*x)/8,
 (231*x^6 - 315*x^4 + 105*x^2 - 5)/16,
 (429*x^7 - 693*x^5 + 315*x^3 - 35*x)/16,
 (6435*x^8 - 12012*x^6 + 6930*x^4 - 1260*x^2 + 35)/128,
 (12155*x^9 - 25740*x^7 + 18018*x^5 - 4620*x^3 + 315*x)/128,
 (46189*x^10 - 109395*x^8 + 90090*x^6 - 30030*x^4 + 3465*x^2 - 63)/256]
sage: ((n-5)*Sn - 1).to_list([1], 10)
[1, 1/-5, 1/20, 1/-60, 1/120, -1/120, None, None, None, None]
class ore_algebra.ore_operator_1_1.UnivariateDifferentialOperatorOverUnivariateRing(parent, *data, **kwargs)

Element of an Ore algebra K(x)[D], where D acts as derivation d/dx on K(x).

annihilator_of_composition(a, solver=None)

Returns an operator L which annihilates all the functions f(a(x)) where f runs through the functions annihilated by self. The output operator is not necessarily of smallest possible order.

INPUT:

  • a – either an element of the base ring of the parent of self, or an element of an algebraic extension of this ring.
  • solver (optional) – a callable object which applied to a matrix with polynomial entries returns its kernel.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: K.<y> = R.fraction_field()['y']
sage: K.<y> = R.fraction_field().extension(y^3 - x^2*(x+1))
sage: A.<Dx> = OreAlgebra(R, 'Dx')
sage: (x*Dx-1).annihilator_of_composition(y) # ann for x^(2/3)*(x+1)^(1/3)
(3*x^2 + 3*x)*Dx - 3*x - 2
sage: (x*Dx-1).annihilator_of_composition(y + 2*x) # ann for 2*x + x^(2/3)*(x+1)^(1/3)
(3*x^3 + 3*x^2)*Dx^2 - 2*x*Dx + 2
sage: (Dx - 1).annihilator_of_composition(y) # ann for exp(x^(2/3)*(x+1)^(1/3))
(-243*x^6 - 810*x^5 - 999*x^4 - 540*x^3 - 108*x^2)*Dx^3 + (-162*x^3 - 270*x^2 - 108*x)*Dx^2 + (162*x^2 + 180*x + 12)*Dx + 243*x^6 + 810*x^5 + 1080*x^4 + 720*x^3 + 240*x^2 + 32*x
annihilator_of_integral()

Returns an operator L which annihilates all the indefinite integrals int f where f runs through the functions annihilated by self. The output operator is not necessarily of smallest possible order.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: A.<Dx> = OreAlgebra(R, 'Dx')
sage: ((x-1)*Dx - 2*x).annihilator_of_integral()
(x - 1)*Dx^2 - 2*x*Dx
sage: _.annihilator_of_associate(Dx)
(x - 1)*Dx - 2*x
finite_singularities()

Returns a list of all the finite singularities of this operator.

OUTPUT:

For each finite singularity of the operator, the output list contains a pair (p, u) where

  • p is an irreducible polynomial, representing the finite singularity rootof(p)+ZZ
  • u is a list of pairs (v, dim, bound), where v is an integer that appears as valuation growth among the solutions of the operator, and bound is a polynomial (or rational function) such that all the solutions of valuation growth v can be written f/bound*Gamma(x-rootof(p))^v where f has minimal valuation almost everywhere. dim is a bound for the number of distinct hypergeometric solutions that may have this local behaviour at rootof(p)+ZZ.

This is a generic implementation for the case of shift and q-shift recurrences. Subclasses for other kinds of operators may need to override this method.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']; A.<Sx> = OreAlgebra(R)
sage: (x^2*(x+1)*Sx + 3*(x+1/2)).finite_singularities()
[(x + 1/2, [[1, 1, 1]]), (x, [[-3, 1, x^4 - 3*x^3 + 3*x^2 - x]])]

sage: C.<q> = ZZ[]; R.<x> = C['x']; A.<Qx> = OreAlgebra(R)
sage: ((q^2*x-1)*Qx-(x-1)).finite_singularities()
[(-x + 1, [[0, 1, q*x^2 + (-q - 1)*x + 1]])]
generalized_series_solutions(n=5, base_extend=True, ramification=True, exp=True)

Returns the generalized series solutions of this operator.

These are solutions of the form

exp(int_0^x frac{p(t^{-1/s})}t dt)*q(x^{1/s},log(x))

where

  • s is a positive integer (the object’s “ramification”)
  • p is in K[x] (the object’s “exponential part”)
  • q is in K[[x]][y] with xnmid q unless q is zero (the object’s “tail”)
  • K is some algebraic extension of the base ring’s base ring.

An operator of order r has exactly r linearly independent solutions of this form. This method computes them all, unless the flags specified in the arguments rule out some of them.

At present, the method only works for operators where the base ring’s base ring is either QQ or a number field (i.e., no finite fields, no formal parameters).

INPUT:

  • n (default: 5) – minimum number of terms in the series expansions to be computed in addition to those needed to separate all solutions from each other.
  • base_extend (default: True) – whether or not the coefficients of the solutions may belong to an algebraic extension of the base ring’s base ring.
  • ramification (default: True) – whether or not the exponential parts of the solutions may involve fractional exponents.
  • exp (default: True) – set this to False if you only want solutions that have no exponential part (viz deg(p)leq0). If set to a positive rational number alpha, the method returns all those solutions whose exponential part involves only terms x^{-i/r} with i/r<alpha.

OUTPUT:

  • a list of ContinuousGeneralizedSeries objects forming a fundamental system for this operator.

Note

  • Different solutions may require different algebraic extensions. Thus in the list returned by this method, the coefficient fields of different series typically do not coincide.
  • If a solution involves an algebraic extension of the coefficient field, then all its conjugates are solutions, too. But only one representative is listed in the output.

ALGORITHM:

  • Ince, Ordinary Differential Equations, Chapters 16 and 17
  • Kauers/Paule, The Concrete Tetrahedron, Section 7.3

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = QQ['x']; A.<Dx> = OreAlgebra(R, 'Dx')
sage: L = (6+6*x-3*x^2) - (10*x-3*x^2-3*x^3)*Dx + (4*x^2-6*x^3+2*x^4)*Dx^2
sage: L.generalized_series_solutions()
[x^3*(1 + 3/2*x + 7/4*x^2 + 15/8*x^3 + 31/16*x^4 + O(x^5)), x^(1/2)*(1 + 3/2*x + 7/4*x^2 + 15/8*x^3 + 31/16*x^4 + O(x^5))]
sage: list(map(L, _))
[0, 0]

sage: L = (1-24*x+96*x^2) + (15*x-117*x^2+306*x^3)*Dx + (9*x^2-54*x^3)*Dx^2
sage: L.generalized_series_solutions(3)
[x^(-1/3)*(1 + x + 8/3*x^2 + O(x^3)), x^(-1/3)*((1 + x + 8/3*x^2 + O(x^3))*log(x) + x - 59/12*x^2 + O(x^3))]
sage: list(map(L, _))
[0, 0]

sage: L = 216*(1+x+x^3) + x^3*(36-48*x^2+41*x^4)*Dx - x^7*(6+6*x-x^2+4*x^3)*Dx^2
sage: L.generalized_series_solutions(3)
[exp(3*x^(-2))*x^(-2)*(1 + 91/12*x^2 + O(x^3)), exp(-2*x^(-3) + x^(-1))*x^2*(1 + 41/3*x + 2849/36*x^2 + O(x^3))]
sage: list(map(L, _))
[0, 0]

sage: L = 9 - 49*x - 2*x^2 + 6*x^2*(7 + 5*x)*Dx + 36*(-1 + x)*x^3*Dx^2
sage: L.generalized_series_solutions()
[exp(x^(-1/2))*x^(4/3)*(1 + x^(2/2) + x^(4/2)), exp(-x^(-1/2))*x^(4/3)*(1 + x^(2/2) + x^(4/2))]
sage: L.generalized_series_solutions(ramification=False)
[]

sage: L = 2*x^3*Dx^2 + 3*x^2*Dx-1
sage: L.generalized_series_solutions()
[exp(a_0*x^(-1/2))]
sage: _[0].base_ring()
Number Field in a_0 with defining polynomial x^2 - 2
indicial_polynomial(p, var='alpha')

Compute the indicial polynomial of this operator at (a root of) p.

If x is the generator of the base ring, the input may be either irreducible polynomial in x or the rational function 1/x.

The output is a univariate polynomial in the given variable var with coefficients in the base ring’s base ring. It has the following property: for every nonzero series solution of self in rising powers of p, i.e. p_0 p^alpha + p_1 p^{alpha+1} + …, the minimal exponent alpha is a root of the indicial polynomial. The converse may not hold.

When p has degree one (but not in general), the degree of the indicial polynomial is equal to the order of the operator if and only if the root of p is an ordinary or regular singular point of the operator.

INPUT:

  • p – an irreducible polynomial in the base ring of the operator algebra, or 1/x.
  • var (optional) – the variable name to use for the indicial polynomial.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']; A.<Dx> = OreAlgebra(R, 'Dx');
sage: L = (x*Dx-5).lclm((x^2+1)*Dx - 7*x).lclm(Dx - 1)
sage: L.indicial_polynomial(x).factor()
5 * 2^2 * (alpha - 5) * (alpha - 1) * alpha
sage: L.indicial_polynomial(1/x).factor()
(-1) * 2 * (alpha - 7) * (alpha - 5)
sage: L.indicial_polynomial(x^2+1).factor()
5 * 7 * (alpha - 1) * alpha * (2*alpha - 7)

The indicial polynomial at p is not always the same as the indicial polynomial at a root of p:

sage: from ore_algebra.examples import cbt
sage: dop = cbt.dop[4]; dop
(-z^3 + 6*z^2 - 5*z + 1)*Dz^5 + (2*z^3 - 18*z^2 + 40*z - 15)*Dz^4 +
(-z^3 + 16*z^2 - 54*z + 41)*Dz^3 + (-4*z^2 + 22*z - 24)*Dz^2 +
(-2*z + 3)*Dz
sage: lc = dop.leading_coefficient()
sage: dop.indicial_polynomial(lc)
alpha^4 - 6*alpha^3 + 11*alpha^2 - 6*alpha
sage: K.<s> = QQ.extension(lc)
sage: z = dop.base_ring().gen()
sage: dop.change_ring(K['z']).indicial_polynomial(z-s)
7*alpha^5 + (-3*s - 50)*alpha^4 + (18*s + 125)*alpha^3 +
(-33*s - 130)*alpha^2 + (18*s + 48)*alpha

TESTS:

sage: A(x^3 - 2).indicial_polynomial(x^2 + 1)
1

sage: P.<x> = QQ[]; Q.<y> = Frac(P)[]; Dops.<Dy> = OreAlgebra(Q)
sage: dop = ((x+1)*(y*Dy)^3-x)*((y*Dy)^2+2*x*y*Dy+1)
sage: dop.indicial_polynomial(y).factor()
(x + 1) * (alpha^2 + 2*x*alpha + 1) * (alpha^3 - x/(x + 1))
sage: dop = ((((3*x - 5)/(-2*x^2 - 3/4*x - 2/41))*y^2
....:     + ((-1/8*x^2 + 903/2*x)/(-x^2 - 1))*y
....:     + (-1/59*x^2 - 1/6*x - 5)/(-8*x^2 + 1/2))*Dy^2
....:     + (((x^2 - 1/2*x + 2)/(x + 1/3))*y^2
....:     + ((2*x^2 + 19*x + 1/2)/(-5*x^2 + 21/4*x - 1/2))*y
....:     + (1/5*x^2 - 26*x - 3)/(-x^2 - 1/3*x + 1/3))*Dy
....:     + ((3*x^2 + 2/5*x + 1/2)/(-139*x^2 + 2))*y^2
....:     + ((1/2*x^2 + 1/20*x + 1)/(4/3*x^2 + 1/6*x + 4))*y
....:     + (3/5*x - 3)/(-1/2*x^2))
sage: dop.indicial_polynomial(y)
(1/472*x^2 + 1/48*x + 5/8)*alpha^2 + (-1/472*x^2 - 1/48*x - 5/8)*alpha
sage: dop.indicial_polynomial(dop.leading_coefficient())
alpha

sage: Pol.<u> = QQ[]
sage: Dop.<Du> = OreAlgebra(Pol)
sage: dop = ((-96040000*u^18 + 64038100*u^17 - 256116467*u^16 +
....: 224114567*u^15 - 32034567*u^14 + 128040267*u^13 +
....: 448194834*u^12 - 352189134*u^11 + 352189134*u^10 -
....: 448194834*u^9 - 128040267*u^8 + 32034567*u^7 - 224114567*u^6 +
....: 256116467*u^5 - 64038100*u^4 + 96040000*u^3)*Du^3 +
....: (240100000*u^17 + 96010600*u^16 - 288129799*u^15 +
....: 1008488600*u^14 - 2641222503*u^13 + 2593354404*u^12 -
....: 2977470306*u^11 + 1776857604*u^10 + 720290202*u^9 -
....: 1632885804*u^8 + 2977475205*u^7 - 2737326204*u^6 +
....: 1680832301*u^5 - 1056479200*u^4 + 288124900*u^3 -
....: 48020000*u^2)*Du^2 + (-480200000*u^16 - 672221200*u^15 +
....: 4033758398*u^14 - 5186718602*u^13 + 13062047620*u^12 -
....: 10757577620*u^11 + 11813792216*u^10 - 7971790408*u^9 +
....: 2977494796*u^8 - 2593079996*u^7 - 384081598*u^6 -
....: 2304950206*u^5 - 191923200*u^4 - 1344540400*u^3 - 96049800*u^2
....: + 96040000*u)*Du + 480200000*u^15 + 1152421200*u^14 -
....: 8931857198*u^13 + 6916036404*u^12 - 18344443640*u^11 +
....: 7588296828*u^10 - 16615302196*u^9 + 673240380*u^8 -
....: 14694120024*u^7 + 3650421620*u^6 - 8356068006*u^5 +
....: 4802156800*u^4 - 1248500400*u^3 - 96059600*u^2 + 96049800*u -
....: 96040000)
sage: dop.indicial_polynomial(70*u^2 + 69*u + 70)
alpha^3 - 3*alpha^2 + 2*alpha
local_basis_expansions(point, order=None, ring=None)

Generalized series expansions of the local basis.

INPUT:

  • point - Point where the local basis is to be computed

  • order (optional) - Number of terms to compute, starting from each “leftmost” valuation of a group of solutions with valuations differing by integers. (Thus, the absolute truncation order will be the same for all solutions in such a group, with some solutions having more actual coefficients computed that others.)

    The default is to choose the truncation order in such a way that the structure of the basis is apparent, and in particular that logarithmic terms appear if logarithms are involved at all in that basis. The corresponding order may be very large in some cases.

  • ring (optional) - Ring into which to coerce the coefficients of the expansion

EXAMPLES:

sage: from ore_algebra import *
sage: Dops, x, Dx = DifferentialOperators(QQ, 'x')

sage: (Dx - 1).local_basis_expansions(0)
[1 + x + 1/2*x^2 + 1/6*x^3]

sage: from ore_algebra.examples import ssw
sage: ssw.dop[1,0,0].local_basis_expansions(0)
[t^(-4) + 24*t^(-2)*log(t) - 48*log(t) - 96*t^2*log(t) - 88*t^2,
 t^(-2),
 1 + 2*t^2]

sage: dop = (x^2*(x^2-34*x+1)*Dx^3 + 3*x*(2*x^2-51*x+1)*Dx^2
....:     + (7*x^2-112*x+1)*Dx + (x-5))
sage: dop.local_basis_expansions(0, order=3)
[1/2*log(x)^2 + 5/2*x*log(x)^2 + 12*x*log(x) + 73/2*x^2*log(x)^2
+ 210*x^2*log(x) + 72*x^2,
log(x) + 5*x*log(x) + 12*x + 73*x^2*log(x) + 210*x^2,
1 + 5*x + 73*x^2]

sage: roots = dop.leading_coefficient().roots(AA)
sage: basis = dop.local_basis_expansions(roots[1][0], order=3)
sage: basis
[1 - (-239/12*a+169/6)*(x - 0.02943725152285942?)^2,
 (x - 0.02943725152285942?)^(1/2) - (-203/32*a+9)*(x - 0.02943725152285942?)^(3/2) + (-24031/160*a+1087523/5120)*(x - 0.02943725152285942?)^(5/2),
 (x - 0.02943725152285942?) - (-55/6*a+13)*(x - 0.02943725152285942?)^2]
sage: basis[0].base_ring()
Number Field in a with defining polynomial y^2 - 2
sage: RR(basis[0].base_ring().gen())
-1.41421356237309
sage: basis[0][-1]
(239/12*a - 169/6, (x - 0.02943725152285942?)^2)

sage: dop.local_basis_expansions(roots[1][0], order=3, ring=QQbar)
[1 - 56.33308678393081?*(x - 0.02943725152285942?)^2,
 (x - 0.02943725152285942?)^(1/2) - 17.97141728630432?*(x - 0.02943725152285942?)^(3/2) + 424.8128741711741?*(x - 0.02943725152285942?)^(5/2),
 (x - 0.02943725152285942?) - 25.96362432175337?*(x - 0.02943725152285942?)^2]

TESTS:

sage: (4*x^2*Dx^2 + (-x^2+8*x-11)).local_basis_expansions(0, 2)
[x^(-1.232050807568878?) + (-4/11*a+2/11)*x^(-0.2320508075688773?),
x^2.232050807568878? - (-4/11*a-2/11)*x^3.232050807568878?]

sage: ((27*x^2+4*x)*Dx^2 + (54*x+6)*Dx + 6).local_basis_expansions(0, 2)
[x^(-1/2) + 3/8*x^(1/2), 1 - x]

sage: dop = (Dx^3 + ((24*x^2 - 4*x - 12)/(8*x^3 - 8*x))*Dx^2 +
....:   ((32*x^2 + 32*x - 16)/(32*x^4 + 32*x^3 - 32*x^2 - 32*x))*Dx)
sage: dop.local_basis_expansions(0, 3)
[1, x^(1/2) - 1/6*x^(3/2) + 3/40*x^(5/2), x - 1/6*x^2]

Thanks to Armin Straub for this example:

sage: dop = ((81*x^4 + 14*x^3 + x^2)*Dx^3
....:       + (486*x^3 + 63*x^2 + 3*x)*Dx^2
....:       + (567*x^2 + 48*x + 1)*Dx + 81*x + 3)
sage: dop.local_basis_expansions(QQbar((4*sqrt(2)*I-7)/81), 2)
[1,
 (x + 0.0864197530864198? - 0.06983770678385654?*I)^(1/2) + (365/96*a^3+365/96*a+13/3)*(x + 0.0864197530864198? - 0.06983770678385654?*I)^(3/2),
 (x + 0.0864197530864198? - 0.06983770678385654?*I)]

and to Emre Sertöz for this one:

sage: ode = (Dx^2 + (2*x - 7/4)/(x^2 - 7/4*x + 3/4)*Dx
....:       + 3/16/(x^2 - 7/4*x + 3/4))
sage: ode.local_basis_expansions(1, 3)[1]
1 - 3/4*(x - 1) + 105/64*(x - 1)^2
local_basis_monomials(point)

Return the leading logarithmic monomials of a local basis of solutions.

INPUT:

point - a regular point of this operator

OUTPUT:

A list of expressions of the form (x-point)^λ*log(x-point)^k/k! where λ is a root of the indicial polynomial of the operator at point, and k is a nonnegative integer less than the multiplicity of that root.

If point is an ordinary point, the output is [1, x, x^2, ...]. More generally, a solution of the operator is characterized by the coefficients in its logarithmic power series expansion at point of the monomials returned by this method. The basis of solutions consisting of the local solutions in which exactly one of the monomials appears (with a coefficient equal to one), ordered as in the output of this method, is used in several functions of this package to specify vectors of “generalized initial values” at regular singular points. (The order is essentially that of asymptotic dominance as x tends to point.) Note however that this basis may not coincide with the one computed by generalized_series_solutions().

EXAMPLES:

sage: from ore_algebra import DifferentialOperators
sage: Dops, x, Dx = DifferentialOperators()
sage: ((x+1)*Dx^4+Dx-x).local_basis_monomials(0)
[1, x, x^2, x^3]
sage: ((x^2 + 1)*Dx^2 + 2*x*Dx).local_basis_monomials(i)
[log(x - I), 1]
sage: (4*x^2*Dx^2 + (-x^2+8*x-11)).local_basis_monomials(0)
[x^(-1.232050807568878?), x^2.232050807568878?]
sage: (x^3*Dx^4+3*x^2*Dx^3+x*Dx^2+x*Dx+1).local_basis_monomials(0)
[1, 1/2*x*log(x)^2, x*log(x), x]

A local basis whose elements all start with pure monomials (without logarithmic part) can nevertheless involve logarithms. In particular, the leading monomials are not enough to decide if a given solution is analytic:

sage: dop = (x^2 - x)*Dx^2 + (x - 1)*Dx + 1
sage: dop.local_basis_monomials(1)
[1, x - 1]
sage: dop.annihilator_of_composition(1 + x).generalized_series_solutions(3)
[x*(1 - x + 5/6*x^2 + O(x^3)),
 (x - x^2 + O(x^3))*log(x) - 1 + 1/2*x^2 + O(x^3)]

TESTS:

sage: ((x+1/3)*Dx^4+Dx-x).local_basis_monomials(-1/3)
[1, x + 1/3, 1/9*(3*x + 1)^2, 1/27*(3*x + 1)^3]

sage: ((x^2 - 2)^3*Dx^4+Dx-x).local_basis_monomials(sqrt(2))
[1, (x - sqrt(2))^0.978..., (x - sqrt(2))^2.044...,
(x - sqrt(2))^2.977...]

sage: dop = (Dx^3 + ((24*x^2 - 4*x - 12)/(8*x^3 - 8*x))*Dx^2 +
....:   ((32*x^2 + 32*x - 16)/(32*x^4 + 32*x^3 - 32*x^2 - 32*x))*Dx)
sage: dop.local_basis_monomials(0)
[1, sqrt(x), x]
numerical_solution(ini, path, eps=1e-16, post_transform=None, **kwds)

Evaluate an analytic solution of this operator at a point of its Riemann surface.

INPUT:

  • ini (iterable) - initial values, in number equal to the order r of the operator
  • path - a path on the complex plane, specified as a list of vertices z_0, dots, z_n
  • eps (floating-point number or ball, default 1e-16) - approximate target accuracy
  • post_transform (default: identity) - differential operator to be applied to the solutions, see examples below

OUTPUT:

A real or complex ball enclosing the value at z_n of the solution y defined in the neighborhood of z_0 by the initial values ini and extended by analytic continuation along path.

When z_0 is an ordinary point, the initial values are defined as the first r coefficients of the power series expansion at z_0 of the desired solution f. In other words, ini must be equal to

\[[f(z_0), f'(z_0), f''(z_0)/2, \dots, f^{(r-1)}(z_0)/(r-1)!].\]

Generalized initial conditions at regular singular points are also supported. If z_0 is a regular point, the entries of ini are interpreted as the coefficients of the monomials (z-z_0)^n log(z-z_0)^k/k! returned by local_basis_monomials() in the logarithmic series expansion of f at z_0. This definition reduces to the previous one when z_0 is an ordinary point.

The accuracy parameter eps is used as an indication of the absolute error the code should aim for. The diameter of the result will typically be of the order of magnitude of eps, but this is not guaranteed to be the case. (It is a bug, however, if the returned ball does not contain the exact result.)

See ore_algebra.analytic for more information.

EXAMPLES:

First a very simple example:

sage: from ore_algebra import DifferentialOperators
sage: Dops, x, Dx = DifferentialOperators()
sage: (Dx - 1).numerical_solution(ini=[1], path=[0, 1], eps=1e-50)
[2.7182818284590452353602874713526624977572470936999...]

Evaluation points can be complex and can depend on symbolic constants:

sage: (Dx - 1).numerical_solution([1], [0, i + pi])
[12.5029695888765...] + [19.4722214188416...]*I

Here, we use a more complicated analytic continuation path in order to evaluate the branch of the complex arctangent function obtained by turning around its singularity at i once:

sage: dop = (x^2 + 1)*Dx^2 + 2*x*Dx
sage: dop.numerical_solution([0, 1], [0, i+1, 2*i, i-1, 0])
[3.14159265358979...] + [+/- ...]*I

In some cases, this method is also able to compute limits of solutions at regular singular points. This only works when all solutions of the differential equation tend to finite values at the evaluation point:

sage: dop = (x - 1)^2*Dx^3 + Dx + 1
sage: dop.local_basis_monomials(1)
[1,
(x - 1)^(1.500000000000000? - 0.866025403784439?*I),
(x - 1)^(1.500000000000000? + 0.866025403784439?*I)]
sage: dop.numerical_solution(ini=[1, 0, 0], path=[0, 1])
[0.6898729110219401...] + [+/- ...]*I

sage: dop = -(x+1)*(x-1)^3*Dx^2 + (x+3)*(x-1)^2*Dx - (x+3)*(x-1)
sage: dop.local_basis_monomials(1)
[x - 1, (x - 1)^2]
sage: dop.numerical_solution([1,0], [0,1])
0

sage: (Dx*x*Dx).numerical_solution(ini=[1,0],path=[1,0])
Traceback (most recent call last):
...
ValueError: solution may not have a finite limit at evaluation
point (try using numerical_transition_matrix())

The post_transform parameter can be used to compute derivatives or linear combinations of derivatives of the solution. Here, we use this feature to evaluate the tenth derivative of the Airy Ai function:

sage: ini = [1/(3^(2/3)*gamma(2/3)), -1/(3^(1/3)*gamma(1/3))]
sage: (Dx^2-x).numerical_solution(ini, [0,2], post_transform=Dx^10)
[2.34553207877...]
sage: airy_ai(10, 2.)
2.345532078777...

A similar, slightly more complicated example:

sage: (Dx^2 - x).numerical_solution(ini, [0, 2],
....:                               post_transform=1/x + x*Dx)
[-0.08871870365567...]
sage: t = SR.var('t')
sage: (airy_ai(t)/t + t*airy_ai_prime(t))(t=2.)
-0.08871870365567...

Some notable examples of incorrect input:

sage: (Dx - 1).numerical_solution([1], [])
Traceback (most recent call last):
...
ValueError: empty path

sage: ((x - 1)*Dx + 1).numerical_solution([1], [0, 2])
Traceback (most recent call last):
...
ValueError: Step 0 --> 2 passes through or too close to singular
point 1 (to compute the connection to a singular point, make it a
vertex of the path)

sage: Dops.zero().numerical_solution([], 1)
Traceback (most recent call last):
...
ValueError: operator must be nonzero

sage: (Dx - 1).numerical_solution(ini=[], path=[0, 1])
Traceback (most recent call last):
...
ValueError: incorrect initial values: []

sage: (Dx - 1).numerical_solution([1], ["a"])
Traceback (most recent call last):
...
TypeError: unexpected value for point: 'a'
numerical_transition_matrix(path, eps=1e-16, **kwds)

Compute a transition matrix along a path drawn in the complex plane.

INPUT:

  • path - a path on the complex plane, specified as a list of vertices z_0, dots, z_n
  • eps (floating-point number or ball) - target accuracy

OUTPUT:

When self is an operator of order r, this method returns an r×r matrix of real or complex balls. The returned matrix maps a vector of “initial values at z_0” (i.e., the coefficients of the decomposition of a solution in a certain canonical local basis at z_0) to “initial values at z_n” that define the same solution, extended by analytic continuation along the path path.

The “initial values” are the coefficients of the monomials returned by local_basis_monomials() in the local logarithmic power series expansions of the solution at the corresponding point. When z_i is an ordinary point, the corresponding vector of initial values is simply

\[[f(z_i), f'(z_i), f''(z_i)/2, \dots, f^{(r-1)}(z_i)/(r-1)!].\]

The accuracy parameter eps is used as an indication of the absolute error that the code should aim for. The diameter of each entry of the result will typically be of the order of magnitude of eps, but this is not guaranteed to be the case. (It is a bug, however, if the returned ball does not contain the exact result.)

See ore_algebra.analytic for more information.

EXAMPLES:

We can compute exp(1) as the only entry of the transition matrix from 0 to 1 for the differential equation y’ = y:

sage: from ore_algebra import DifferentialOperators
sage: Dops, x, Dx = DifferentialOperators()
sage: (Dx - 1).numerical_transition_matrix([0, 1])
[[2.7182818284590452 +/- 3.54e-17]]

Now consider a second-order operator that annihilates arctan(x) and the constants. A basis of solutions is formed of the constant 1, of the form 1 + O(x^2) as x to 0, and the arctangent function, of the form x + O(x^2). Accordingly, the entries of the transition matrix from the origin to 1 + i are the values of these two functions and their first derivatives:

sage: dop = (x^2 + 1)*Dx^2 + 2*x*Dx
sage: dop.numerical_transition_matrix([0, 1+i], 1e-10)
[ [1.00...] + [+/- ...]*I  [1.017221967...] + [0.4023594781...]*I]
[ [+/- ...] + [+/- ...]*I  [0.200000000...] + [-0.400000000...]*I]

By making loops around singular points, we can compute local monodromy matrices:

sage: dop.numerical_transition_matrix([0, i + 1, 2*i, i - 1, 0])
[ [1.00...] + [+/- ...]*I  [3.141592653589793...] + [+/-...]*I]
[ [+/- ...] + [+/- ...]*I  [1.000000000000000...] + [+/-...]*I]

Then we compute a connection matrix to the singularity itself:

sage: dop.numerical_transition_matrix([0, i], 1e-10)
[            ...           [+/-...] + [-0.50000000...]*I]
[ ...1.000000...  [0.7853981634...] + [0.346573590...]*I]

Note that a path that crosses the branch cut of the complex logarithm yields a different result:

sage: dop.numerical_transition_matrix([0, i - 1, i], 1e-10)
[     [+/-...] + [+/-...]*I         [+/-...] + [-0.5000000000...]*I]
[ [1.00000...] + [+/-...]*I [-2.356194490...] + [0.3465735902...]*I]

In general, if the operator has rational coefficients, its singular points are algebraic numbers. In connection problems such as the above, they need to be specified exactly. Here is a way to do it:

sage: dop = (x^2 - 2)*Dx^2 + x + 1
sage: dop.numerical_transition_matrix([0, 1, QQbar(sqrt(2))], 1e-10)
[         [2.49388146...] + [+/-...]*I          [2.40894178...] + [+/-...]*I]
[[-0.203541775...] + [6.68738570...]*I  [0.204372067...] + [6.45961849...]*I]

The operator itself may be defined over a number field (with a complex embedding):

sage: K.<zeta7> = CyclotomicField(7)
sage: (Dx - zeta7).numerical_transition_matrix([0, 1])
[[1.32375209616333...] + [1.31434281345999...]*I]

Some notable examples of incorrect input:

sage: (Dx - 1).numerical_transition_matrix([])
Traceback (most recent call last):
...
ValueError: empty path

sage: ((x - 1)*Dx + 1).numerical_transition_matrix([0, 2])
Traceback (most recent call last):
...
ValueError: Step 0 --> 2 passes through or too close to singular
point 1 (to compute the connection to a singular point, make it a
vertex of the path)

sage: Dops.zero().numerical_transition_matrix([0, 1])
Traceback (most recent call last):
...
ValueError: operator must be nonzero
power_series_solutions(n=5)

Computes the first few terms of the power series solutions of this operator.

The method raises an error if Sage does not know how to factor univariate polynomials over the base ring’s base ring.

The base ring has to have characteristic zero.

INPUT:

  • n – minimum number of terms to be computed

OUTPUT:

A list of power series of the form x^alpha + … with pairwise distinct exponents alpha and coefficients in the base ring’s base ring’s fraction field. All expansions are computed up to order k where k is obtained by adding the maximal alpha to the maximum of n and the order of self.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: A.<Dx> = OreAlgebra(R, 'Dx')
sage: ((1-x)*Dx - 1).power_series_solutions(10) # geometric series
[1 + x + x^2 + x^3 + x^4 + x^5 + x^6 + x^7 + x^8 + O(x^9)]
sage: (Dx - 1).power_series_solutions(5) # exp(x)
[1 + x + 1/2*x^2 + 1/6*x^3 + O(x^4)]
sage: (Dx^2 - Dx + x).power_series_solutions(5) # a 2nd order equation
[x + 1/2*x^2 + 1/6*x^3 - 1/24*x^4 + O(x^5), 1 - 1/6*x^3 - 1/24*x^4 + O(x^5)]
sage: (2*x*Dx - 1).power_series_solutions(5) # sqrt(x) is not a power series
[]
spread(p=0)

Returns the spread of this operator.

This is the set of integers i such that sigma(self[0], i) and sigma(self[r], -r) have a nontrivial common factor, where sigma is the shift of the parent’s algebra and r is the order of self.

If the optional argument p is given, the method is applied to gcd(self[0], p) instead of self[0].

The output set contains infty if the constant coefficient of self is zero.

This method is a stub and may not be implemented for every algebra.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']; A.<Sx> = OreAlgebra(R, 'Sx');
sage: ((x+5)*Sx - x).spread()
[4]
sage: ((x+5)*Sx - x).lclm((x+19)*Sx - x).spread()
[3, 4, 17, 18]
to_F(alg)

Returns a difference operator annihilating the coefficient sequence of every power series (about the origin) annihilated by self.

INPUT:

  • alg – the Ore algebra in which the output should be expressed. The algebra must satisfy alg.base_ring().base_ring() == self.base_ring().base_ring() and alg.is_F() is not False. Instead of an algebra object, also a string can be passed as argument. This amounts to specifying an Ore algebra over self.base_ring() with the forward difference with respect to self.base_ring().gen().

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: A.<Dx> = OreAlgebra(R, 'Dx')
sage: R2.<n> = ZZ['n']
sage: A2.<Sn> = OreAlgebra(R2, 'Fn')
sage: (Dx - 1).to_F(A2)
(n + 1)*Fn + n
sage: ((1+x)*Dx^2 + Dx).to_F(A2)
(n^2 + n)*Fn + 2*n^2 + n
sage: ((x^3+x^2-x)*Dx + (x^2+1)).to_F(A2)
(-n - 1)*Fn^2 + (-n - 1)*Fn + n + 1
to_S(alg)

Returns a recurrence operator annihilating the coefficient sequence of every power series (about the origin) annihilated by self.

INPUT:

  • alg – the Ore algebra in which the output should be expressed. The algebra must satisfy alg.base_ring().base_ring() == self.base_ring().base_ring() and alg.is_S() is not False. Instead of an algebra object, also a string can be passed as argument. This amounts to specifying an Ore algebra over self.base_ring() with the standard shift with respect to self.base_ring().gen().

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: A.<Dx> = OreAlgebra(R, 'Dx')
sage: R2.<n> = ZZ['n']
sage: A2.<Sn> = OreAlgebra(R2, 'Sn')
sage: (Dx - 1).to_S(A2)
(n + 1)*Sn - 1
sage: ((1+x)*Dx^2 + Dx).to_S(A2)
(n^2 + n)*Sn + n^2
sage: ((x^3+x^2-x)*Dx + (x^2+1)).to_S(A2)
(-n - 1)*Sn^2 + (n + 1)*Sn + n + 1
sage: ((x+1)*Dx^3 + Dx^2).to_S(A2)
(n^3 - n)*Sn + n^3 - 2*n^2 + n
to_T(alg)

Rewrites self in terms of the eulerian derivation x*d/dx.

If the base ring of the target algebra is not a field, the operator returned by the method may not correspond exactly to self, but only to a suitable left-multiple by a term x^k.

INPUT:

  • alg – the Ore algebra in which the output should be expressed. The algebra must satisfy alg.base_ring().base_ring() == self.base_ring().base_ring() and alg.is_T() is not False. Instead of an algebra object, also a string can be passed as argument. This amounts to specifying an Ore algebra over self.base_ring() with an euler derivation with respect to self.base_ring().gen().

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: R2.<y> = ZZ['y']
sage: A.<Dx> = OreAlgebra(R, 'Dx')
sage: (Dx^4).to_T(OreAlgebra(R2, 'Ty'))
Ty^4 - 6*Ty^3 + 11*Ty^2 - 6*Ty
sage: (Dx^4).to_T('Tx').to_D(A)
x^4*Dx^4
sage: _.to_T('Tx')
Tx^4 - 6*Tx^3 + 11*Tx^2 - 6*Tx
class ore_algebra.ore_operator_1_1.UnivariateEulerDifferentialOperatorOverUnivariateRing(parent, *data, **kwargs)

Element of an Ore algebra K(x)[T], where T is the Euler differential operator T = x*d/dx

power_series_solutions(*args, **kwargs)

Computes the first few terms of the power series solutions of this operator.

The method raises an error if Sage does not know how to factor univariate polynomials over the base ring’s base ring.

The base ring has to have characteristic zero.

INPUT:

  • n – minimum number of terms to be computed

OUTPUT:

A list of power series of the form x^alpha + … with pairwise distinct exponents alpha and coefficients in the base ring’s base ring’s fraction field. All expansions are computed up to order k where k is obtained by adding the maximal alpha to the maximum of n and the order of self.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: A.<Dx> = OreAlgebra(R, 'Dx')
sage: ((1-x)*Dx - 1).power_series_solutions(10) # geometric series
[1 + x + x^2 + x^3 + x^4 + x^5 + x^6 + x^7 + x^8 + O(x^9)]
sage: (Dx - 1).power_series_solutions(5) # exp(x)
[1 + x + 1/2*x^2 + 1/6*x^3 + O(x^4)]
sage: (Dx^2 - Dx + x).power_series_solutions(5) # a 2nd order equation
[x + 1/2*x^2 + 1/6*x^3 - 1/24*x^4 + O(x^5), 1 - 1/6*x^3 - 1/24*x^4 + O(x^5)]
sage: (2*x*Dx - 1).power_series_solutions(5) # sqrt(x) is not a power series
[]
spread(p=0)

Returns the spread of this operator.

This is the set of integers i such that sigma(self[0], i) and sigma(self[r], -r) have a nontrivial common factor, where sigma is the shift of the parent’s algebra and r is the order of self.

If the optional argument p is given, the method is applied to gcd(self[0], p) instead of self[0].

The output set contains infty if the constant coefficient of self is zero.

This method is a stub and may not be implemented for every algebra.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']; A.<Sx> = OreAlgebra(R, 'Sx');
sage: ((x+5)*Sx - x).spread()
[4]
sage: ((x+5)*Sx - x).lclm((x+19)*Sx - x).spread()
[3, 4, 17, 18]
symmetric_product(other, solver=None)

Returns the symmetric product of self and other.

The symmetric product of two operators A and B is a minimal order operator C such that for all “functions” f and g with A.f=B.g=0 we have C.(fg)=0.

The method requires that a product rule is associated to the Ore algebra where self and other live. (See docstring of OreAlgebra for information about product rules.)

If no solver is specified, the the Ore algebra’s solver is used.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: A.<Dx> = OreAlgebra(R, 'Dx')
sage: (Dx - 1).symmetric_product(x*Dx - 1)
x*Dx - x - 1
sage: (x*Dx - 1).symmetric_product(Dx - 1)
x*Dx - x - 1
sage: ((x+1)*Dx^2 + (x-1)*Dx + 8).symmetric_product((x-1)*Dx^2 + (2*x+3)*Dx + (8*x+5))
(29*x^8 - 4*x^7 - 55*x^6 - 34*x^5 - 23*x^4 + 80*x^3 + 95*x^2 - 42*x - 46)*Dx^4 + (174*x^8 + 150*x^7 + 48*x^6 - 294*x^5 - 864*x^4 - 646*x^3 + 232*x^2 + 790*x + 410)*Dx^3 + (783*x^8 + 1661*x^7 - 181*x^6 - 1783*x^5 - 3161*x^4 - 3713*x^3 + 213*x^2 + 107*x - 1126)*Dx^2 + (1566*x^8 + 5091*x^7 + 2394*x^6 + 2911*x^5 - 10586*x^4 - 23587*x^3 - 18334*x^2 - 2047*x + 5152)*Dx + 2552*x^8 + 3795*x^7 + 8341*x^6 + 295*x^5 - 6394*x^4 - 24831*x^3 - 35327*x^2 - 23667*x - 13708

sage: A.<Sx> = OreAlgebra(R, 'Sx')
sage: (Sx - 2).symmetric_product(x*Sx - (x+1))
x*Sx - 2*x - 2
sage: (x*Sx - (x+1)).symmetric_product(Sx - 2)
x*Sx - 2*x - 2
sage: ((x+1)*Sx^2 + (x-1)*Sx + 8).symmetric_product((x-1)*Sx^2 + (2*x+3)*Sx + (8*x+5))
(-8*x^8 - 13*x^7 + 300*x^6 + 1640*x^5 + 3698*x^4 + 4373*x^3 + 2730*x^2 + 720*x)*Sx^4 + (16*x^8 + 34*x^7 - 483*x^6 - 1947*x^5 - 2299*x^4 - 2055*x^3 - 4994*x^2 - 4592*x)*Sx^3 + (-64*x^8 + 816*x^7 + 1855*x^6 - 21135*x^5 - 76919*x^4 - 35377*x^3 + 179208*x^2 + 283136*x + 125440)*Sx^2 + (1024*x^7 + 1792*x^6 - 39792*x^5 - 250472*x^4 - 578320*x^3 - 446424*x^2 + 206528*x + 326144)*Sx - 32768*x^6 - 61440*x^5 + 956928*x^4 + 4897984*x^3 + 9390784*x^2 + 7923200*x + 2329600
to_D(alg)

Returns the differential operator corresponding to self

INPUT:

  • alg – the Ore algebra in which the output should be expressed. The algebra must satisfy alg.base_ring().base_ring() == self.base_ring().base_ring() and alg.is_D() is not False. Instead of an algebra object, also a string can be passed as argument. This amounts to specifying an Ore algebra over self.base_ring() with the standard derivation with respect to self.base_ring().gen().

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: A.<Tx> = OreAlgebra(R, 'Tx')
sage: (Tx^4).to_D(OreAlgebra(R, 'Dx'))
x^4*Dx^4 + 6*x^3*Dx^3 + 7*x^2*Dx^2 + x*Dx
sage: (Tx^4).to_D('Dx').to_T(A)
Tx^4
to_F(alg)

Returns a difference operator annihilating the coefficient sequence of every power series (about the origin) annihilated by self.

INPUT:

  • alg – the Ore algebra in which the output should be expressed. The algebra must satisfy alg.base_ring().base_ring() == self.base_ring().base_ring() and alg.is_F() is not False. Instead of an algebra object, also a string can be passed as argument. This amounts to specifying an Ore algebra over self.base_ring() with the forward difference with respect to self.base_ring().gen().

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: A.<Tx> = OreAlgebra(R, 'Tx')
sage: R2.<n> = ZZ['n']
sage: A2.<Fn> = OreAlgebra(R2, 'Fn')
sage: (Tx - 1).to_F(A2)
n - 1
sage: ((1+x)*Tx^2 + Tx).to_F(A2)
(n^2 + 3*n + 2)*Fn + 2*n^2 + 3*n + 2
sage: ((x^3+x^2-x)*Tx + (x^2+1)).to_F(A2)
Fn^3 + (-n + 1)*Fn^2 + (-n + 1)*Fn + n + 1
to_S(alg)

Returns a recurrence operator annihilating the coefficient sequence of every power series (at the origin) annihilated by self.

INPUT:

  • alg – the Ore algebra in which the output should be expressed. The algebra must satisfy alg.base_ring().base_ring() == self.base_ring().base_ring() and alg.is_S() is not False. Instead of an algebra object, also a string can be passed as argument. This amounts to specifying an Ore algebra over self.base_ring() with the standard shift with respect to self.base_ring().gen().

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: A.<Tx> = OreAlgebra(R, 'Tx')
sage: R2.<n> = ZZ['n']
sage: A2.<Sn> = OreAlgebra(R2, 'Sn')
sage: (Tx - 1).to_S(A2)
n - 1
sage: ((1+x)*Tx^2 + Tx).to_S(A2)
(n^2 + 3*n + 2)*Sn + n^2
sage: ((x^3+x^2-x)*Tx + (x^2+1)).to_S(A2)
Sn^3 + (-n - 2)*Sn^2 + (n + 2)*Sn + n
class ore_algebra.ore_operator_1_1.UnivariateOreOperatorOverUnivariateRing(parent, *data, **kwargs)

Element of an Ore algebra with a single generator and a commutative rational function field as base ring.

TESTS:

sage: from ore_algebra import OreAlgebra
sage: R.<C> = OreAlgebra(GF(2)['x'])
sage: type(C)
<class 'ore_algebra.ore_operator_1_1.UnivariateOreOperatorOverUnivariateRing'>
sage: C.list()
[0, 1]
associate_solutions(D, p)

If self is P, this returns a list of pairs (M, m) such that D*M = p + m*P

INPUT:

  • D – a first order operator with the same parent as self. Depending on the algebra, this operator may be constrained to certain choices. For example, for differential operators, it can only be D (corresponding to integration), and for recurrence operators, it can only be S - 1 (corresponding to summation).
  • p – a nonzero base ring element

OUTPUT:

  • M – an operator of order self.order() - 1 with rational function coefficients.
  • m – a nonzero rational function.

Intended application: Express indefinite sums or integrals of holonomic functions in terms of the summand/integrand. For example, with D=S-1 and P=S^2-S-1 and p some polynomial, the output M is such that

sum_{k=0}^n p(k) F_k = const + M(F_n)

where F_k denotes the Fibonacci sequence. The rational function m does not appear in the closed form, it can be regarded as a certificate.

The method returns the empty list if and only if no nontrivial solutions exist.

This function may not be implemented for every algebra.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = QQ['x']; A.<Dx> = OreAlgebra(R, 'Dx');
sage: L = x*Dx^2 + Dx; p = 1  ## L(log(x)) == 0
sage: L.associate_solutions(Dx, p)
[(-x^2*Dx + x, -x)]
sage: (M, m) = _[0]
sage: Dx*M == p + m*L  ## this implies int(log(x)) == M(log(x)) = x*log(x) - x
True

sage: R.<x> = QQ['x']; A.<Dx> = OreAlgebra(R, 'Dx');
sage: L = x^2*Dx^2 + x*Dx + (x^2 - 1); p = 1  ## L(bessel(x)) == 0
sage: L.associate_solutions(Dx, p)
[(-Dx - 1/x, -1/x^2)]
sage: (M, m) = _[0]
sage: Dx*M == p + m*L  ## this implies int(bessel(x)) == -bessel'(x) -1/x*bessel(x)
True

sage: R.<n> = QQ['n']; A.<Sn> = OreAlgebra(R, 'Sn');
sage: L = Sn^2 - Sn - 1; p = 1  ## L(fib(n)) == 0
sage: L.associate_solutions(Sn - 1, p)
[(Sn, 1)]
sage: (M, m) = _[0]
sage: (Sn-1)*M == p + m*L  ## this implies sum(fib(n)) == fib(n+1)
True

sage: R.<n> = QQ['n']; A.<Sn> = OreAlgebra(R, 'Sn');
sage: L = Sn^3 - 2*Sn^2 - 2*Sn + 1; p = 1  ## L(fib(n)^2) == 0
sage: L.associate_solutions(Sn - 1, p)
[(1/2*Sn^2 - 1/2*Sn - 3/2, 1/2)]
sage: (M, m) = _[0]
sage: (Sn-1)*M == p + m*L  ## this implies sum(fib(n)^2) == 1/2*fib(n+2)^2 - 1/2*fib(n+1)^2 - 3/2*fib(n)^2
True
center(oBound, dBound)

Returns a Q-vector space of Ore polynomials that commute with this operator.

INPUT:

  • oBound – The maximal order of the operators in the center.
  • dBound – The maximal coefficient degree of the operators in the center.

OUTPUT:

A subspace of Q^((oBound+1)*(dBound+1)). Each entry of a vector corresponds to a coefficient of an Ore polynomial that commutes with self. To translate a vector to its corresponding Ore polynomial, call _listToOre

Note: This method only works for operators over Q[n].

degree()

Returns the maximum degree among the coefficients of self

The degree of the zero operator is -1.

If the base ring is not a polynomial ring, this causes an error.

desingularize(m=-1)

Returns a left multiple of self whose coefficients are polynomials and whose leading coefficient does not contain unnecessary factors.

INPUT:

  • m (optional) – If the order of self is r, the output operator will have order r+m. In order to ensure that all removable factors of the leading coefficient are removed in the output, m has to be chosen sufficiently large. If no m is given, a generic upper bound is determined. This feature may not be available for every class.

OUTPUT:

A left multiple of self whose coefficients are polynomials, whose order is m more than self, and whose leading coefficient has as low a degree as possible under these conditions.

The output is not unique. With low probability, the leading coefficient degree in the output may not be minimal.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<n> = ZZ['n']
sage: A.<Sn> = OreAlgebra(R, 'Sn')
sage: P = (-n^3 - 2*n^2 + 6*n + 9)*Sn^2 + (6*n^3 + 8*n^2 - 20*n - 30)*Sn - 8*n^3 - 12*n^2 + 20*n + 12
sage: Q = P.desingularize()
sage: Q.order()
3
sage: Q.leading_coefficient().degree()
1
dispersion(p=0)

Returns the dispersion of this operator.

This is the maximum nonnegative integer i such that sigma(self[0], i) and sigma(self[r], -r) have a nontrivial common factor, where sigma is the shift of the parent’s algebra and r is the order of self.

An output -1 indicates that there are no such integers i at all.

If the optional argument p is given, the method is applied to gcd(self[0], p) instead of self[0].

The output is infty if the constant coefficient of self is zero.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']; A.<Sx> = OreAlgebra(R, 'Sx');
sage: ((x+5)*Sx - x).dispersion()
4
finite_singularities()

Returns a list of all the finite singularities of this operator.

OUTPUT:

For each finite singularity of the operator, the output list contains a pair (p, u) where

  • p is an irreducible polynomial, representing the finite singularity rootof(p)+ZZ
  • u is a list of pairs (v, dim, bound), where v is an integer that appears as valuation growth among the solutions of the operator, and bound is a polynomial (or rational function) such that all the solutions of valuation growth v can be written f/bound*Gamma(x-rootof(p))^v where f has minimal valuation almost everywhere. dim is a bound for the number of distinct hypergeometric solutions that may have this local behaviour at rootof(p)+ZZ.

This is a generic implementation for the case of shift and q-shift recurrences. Subclasses for other kinds of operators may need to override this method.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']; A.<Sx> = OreAlgebra(R)
sage: (x^2*(x+1)*Sx + 3*(x+1/2)).finite_singularities()
[(x + 1/2, [[1, 1, 1]]), (x, [[-3, 1, x^4 - 3*x^3 + 3*x^2 - x]])]

sage: C.<q> = ZZ[]; R.<x> = C['x']; A.<Qx> = OreAlgebra(R)
sage: ((q^2*x-1)*Qx-(x-1)).finite_singularities()
[(-x + 1, [[0, 1, q*x^2 + (-q - 1)*x + 1]])]
indicial_polynomial(p, var='alpha')

Computes the indicial polynomial of self at (a root of) p.

The indicial polynomial is a polynomial in the given variable var with coefficients in the fraction field of the base ring’s base ring.

The precise meaning of this polynomial may depend on the parent of self. A minimum requirement is that if self has a rational solution whose denominator contains sigma.factorial(p, e) but neither sigma(p, -1)*sigma.factorial(p, e) nor sigma.factorial(p, e + 1), then -e is a root of this polynomial.

Applied to p=1/x, the maximum integer root of the output should serve as a degree bound for the polynomial solutions of self.

This method is a stub. Depending on the particular subclass, restrictions on p may apply.

left_factors(order=1, early_termination=False, infolevel=0)

Returns a list of left-hand factors of this operator.

This is a convenience method which simply returns the adjoints of the right factors of the adjoint of self. See docstring of adjoint and right_factors for further information. The method works only in algebras for which adjoint and right_factors are implemented.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ[]; A.<Sx> = OreAlgebra(R)
sage: (((x+1)*Sx + (x+5))*(2*x*Sx + 3)).left_factors()
[[(-x - 1)*Sx - x - 5]]
sage: (((x+1)*Sx + (x+5))*(2*x*Sx + 3)).right_factors()
[[x*Sx + 3/2]]
newton_polygon(p)

Computes the Newton polygon of self at (a root of) p.

INPUT:

  • p – polynomial at whose root the Newton polygon is to be determined. p must be an element of the parent’s base ring (or its fraction field). The value p=1/x represents the point at infinity.

OUTPUT:

A list of pairs (gamma, q) such that gamma is a slope in the Newton polygon and q is the associated polynomial, as elements of the base ring.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ[]; A.<Dx> = OreAlgebra(R); 
sage: L = (x^3*Dx - 1+x).lclm(x*Dx^2-1)
sage: L.newton_polygon(x)
[(1/2, x^2 - 1), (3, -x + 1)]
sage: L.newton_polygon(~x)
[(-2, -x - 1), (-1/2, x^2 - 1)]
sage: A.<Sx> = OreAlgebra(R); L = (x*Sx - 5).lclm(Sx-x^3); L.newton_polygon(~x)
[(-1, -x + 5), (3, x - 1)]

Depending on the algebra in which this operator lives, restrictions on p may apply.

polynomial_solutions(rhs=(), degree=None, solver=None)

Computes the polynomial solutions of this operator.

INPUT:

  • rhs (optional) – a list of base ring elements
  • degree (optional) – bound on the degree of interest.
  • solver (optional) – a callable for computing the right kernel of a matrix over the base ring’s base ring.

OUTPUT:

A list of tuples (p, c_0,…,c_r) such that self(p) == c_0*rhs[0] + … + c_r*rhs[r], where p is a polynomial and c_0,…,c_r are constants.

Note

  • Even if no rhs is given, the output will be a list of tuples [(p1,), (p2,),...] and not just a list of plain polynomials.
  • If no degree is given, a basis of all the polynomial solutions is returned. This feature may not be implemented for all algebras.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<n> = ZZ['n']; A.<Sn> = OreAlgebra(R, 'Sn')
sage: L = 2*Sn^2 + 3*(n-7)*Sn + 4
sage: L.polynomial_solutions((n^2+4*n-8, 4*n^2-5*n+3))
[(-70*n + 231, 242, -113)]
sage: L(-70*n + 231)
-210*n^2 + 1533*n - 2275
sage: 242*(n^2+4*n-8) - 113*(4*n^2-5*n+3)
-210*n^2 + 1533*n - 2275

sage: R.<x> = ZZ['x']; A.<Dx> = OreAlgebra(R, 'Dx')
sage: L = (x*Dx - 19).lclm( x*Dx - 4 )
sage: L.polynomial_solutions()
[(x^4,), (x^19,)]
radical()

Computes the radical of an Ore polynomial P, i.e. an operator L and an integer k such that P=L^k and k is maximal among all the integers for which such an L exists.

OUTPUT:

A tuple (L,k) such that self is equal to L^k and there is no larger integer k’ for which such an L exists.

Note: This method only works for operators over Q[x].

rational_solutions(rhs=(), denominator=None, degree=None, solver=None)

Computes the rational solutions of this operator.

INPUT:

  • rhs (optional) – a list of base ring elements
  • denominator (optional) – bound on the degree of interest.
  • degree (optional) – bound on the degree of interest.
  • solver (optional) – a callable for computing the right kernel of a matrix over the base ring’s base ring.

OUTPUT:

A list of tuples (r, c_0,…,c_r) such that self(r) == c_0*rhs[0] + … + c_r*rhs[r], where r is a rational function and c_0,…,c_r are constants.

Note

  • Even if no rhs is given, the output will be a list of tuples [(p1,), (p2,),...] and not just a list of plain rational functions.
  • If no denominator is given, a basis of all the rational solutions is returned. This feature may not be implemented for all algebras.
  • If no degree is given, a basis of all the polynomial solutions is returned. This feature may not be implemented for all algebras.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']; A.<Dx> = OreAlgebra(R, 'Dx')
sage: L = ((x+3)*Dx + 2).lclm(x*Dx + 3).symmetric_product((x+4)*Dx-2)
sage: L.rational_solutions()
[((x^2 + 8*x + 16)/x^3,), ((x^2 + 8*x + 16)/(x^2 + 6*x + 9),)]
sage: L.rational_solutions((1, x))
[((7*x^5 + 21*x^4 + 73*x^2 + 168*x + 144)/(x^5 + 6*x^4 + 9*x^3), 5184, 756),
 ((4*x^2 + 14*x + 1)/(x^2 + 6*x + 9), 2592, 378),
 ((7*x^2 + 24*x)/(x^2 + 6*x + 9), 4608, 672)]
sage: L(_[0][0]) == _[0][1] + _[0][2]*x
True

sage: (x*(x*Dx-5)).rational_solutions([1])
[(1/x, -6), (x^5, 0)]

sage: R.<n> = ZZ['n']; A.<Sn> = OreAlgebra(R, 'Sn');
sage: L = ((n+3)*Sn - n).lclm((2*n+5)*Sn - (2*n+1))
sage: L.rational_solutions()
[((-4*n^3 - 8*n^2 + 3)/(4*n^5 + 20*n^4 + 35*n^3 + 25*n^2 + 6*n),), (1/(4*n^2 + 8*n + 3),)]

sage: L = (2*n^2 - n - 2)*Sn^2 + (-n^2 - n - 1)*Sn + n^2 - 14
sage: y = (-n + 1)/(n^2 + 2*n - 2)
sage: L.rational_solutions((L(y),))
[((-n + 1)/(n^2 + 2*n - 2), 1)]          
right_factors(order=1, early_termination=False, infolevel=0)

Returns a list of right hand factors of this operator.

INPUT:

  • order (default=1) – only determine right factors of at most this order
  • early_termination (optional) – if set to True, the search for factors will be aborted as soon as one factor has been found. A list containing this single factor will be returned (or the empty list if there are no first order factors). If set to False (default), a complete list will be computed.
  • infolevel (optional) – nonnegative integer specifying the amount of progress reports that should be printed during the calculation. Defaults to 0 for no output.

OUTPUT:

A list of bases for all vector spaces of first-order operators living in the parent of self of which self is a left multiple.

Note that this implementation does not construct factors that involve algebraic extensions of the constant field.

This is a generic implementation for the case of shift and q-shift recurrences. Subclasses for other kinds of operators may need to override this method.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<n> = ZZ['n']; A.<Sn> = OreAlgebra(R, 'Sn');
sage: L = (-25*n^6 - 180*n^5 - 584*n^4 - 1136*n^3 - 1351*n^2 - 860*n - 220)*Sn^2 + (50*n^6 + 560*n^5 + 2348*n^4 + 5368*n^3 + 7012*n^2 + 4772*n + 1298)*Sn - 200*n^5 - 1540*n^4 - 5152*n^3 - 8840*n^2 - 7184*n - 1936
sage: L.right_factors()
[[(n^2 + 6/5*n + 9/25)*Sn - 2*n^2 - 32/5*n - 128/25], [(n^2 + 2*n + 1)*Sn - 4*n - 2]]
sage: ((Sn - n)*(n*Sn - 1)).right_factors()
[[n*Sn - 1]]
sage: ((Sn - n).lclm(n*Sn - 1)).right_factors()
[[n*Sn - 1], [Sn - n]]
sage: (Sn^2 - 2*Sn + 1).right_factors()
[[Sn - 1, n*Sn - n - 1]]

sage: R.<x> = QQ['x']; A.<Qx> = OreAlgebra(R, q=2) 
sage: ((2*x+3)*Qx - (8*x+3)).lclm(x*Qx-2*(x+5)).right_factors()
[[(x + 3/2)*Qx - 4*x - 3/2], [x*Qx - 2*x - 10]]
sage: (((2*x-1)*Qx-(x-1)).lclm(Qx-(x-3))).right_factors()
[[(x - 1/2)*Qx - 1/2*x + 1/2], [Qx - x + 3]]
sage: (((2*x-1)*Qx-(x-1))*(Qx-(x-3))).right_factors()
[[Qx - x + 3]]
sage: (((2*x-1)*Qx-(x-1))*(x^2*Qx-(x-3))).right_factors()
[[x^2*Qx - x + 3]]
singularities(backwards=False)

Returns the integer singularities of the operator self.

INPUT:

  • backwards (default False) – boolean value that decides whether the singularities of the leading coefficient are returned (when backwards is False) or those of the coefficient with minimal degree (regarding Sn or Dx)

OUTPUT:

  • If backwards is False, a set containing the roots of the leading coefficient of the annihilator of self shifted by its order are returned
  • If backwards is True, a set containing the roots of the coefficient with minimal degree (regarding Sn or Dx respectively) are returned; shifted by the degree of this coefficient

EXAMPLES:

sage: from ore_algebra import OreAlgebra
sage: A = OreAlgebra(QQ['n'],'Sn')
sage: a = A("(n-3)*(n+2)*Sn^3 + n^2*Sn^2 - (n-1)*(n+5)*Sn")
sage: a.singularities()
{1, 6}
sage: a.singularities(True)
{-4, 2}

return the integer singularities of the Ore Operator self, i.e. the roots of the leading coefficient shifted by the order of the operator if backwards``is false; when``backwards is true then the roots of the smallest non-zero term (concerning the degree) are returned (shifted by the degree of this term)

spread(p=0)

Returns the spread of this operator.

This is the set of integers i such that sigma(self[0], i) and sigma(self[r], -r) have a nontrivial common factor, where sigma is the shift of the parent’s algebra and r is the order of self.

If the optional argument p is given, the method is applied to gcd(self[0], p) instead of self[0].

The output set contains infty if the constant coefficient of self is zero.

This method is a stub and may not be implemented for every algebra.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']; A.<Sx> = OreAlgebra(R, 'Sx');
sage: ((x+5)*Sx - x).spread()
[4]
sage: ((x+5)*Sx - x).lclm((x+19)*Sx - x).spread()
[3, 4, 17, 18]
class ore_algebra.ore_operator_1_1.UnivariateQDifferentialOperatorOverUnivariateRing(parent, *data, **kwargs)

Element of an Ore algebra K(x)[J], where J is the Jackson q-differentiation J f(x) = (f(q*x) - f(x))/(q*(x-1))

annihilator_of_integral()

Returns an operator L which annihilates all the indefinite q-integrals int_q f where f runs through the functions annihilated by self. The output operator is not necessarily of smallest possible order.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['q'].fraction_field()['x']
sage: A.<Jx> = OreAlgebra(R, 'Jx')
sage: ((x-1)*Jx - 2*x).annihilator_of_integral()
(x - 1)*Jx^2 - 2*x*Jx
sage: _.annihilator_of_associate(Jx)
(x - 1)*Jx - 2*x
power_series_solutions(n=5)

Computes the first few terms of the power series solutions of this operator.

The method raises an error if Sage does not know how to factor univariate polynomials over the base ring’s base ring.

The base ring has to have characteristic zero.

INPUT:

  • n – minimum number of terms to be computed

OUTPUT:

A list of power series of the form x^{alpha} + … with pairwise distinct exponents alpha and coefficients in the base ring’s base ring’s fraction field. All expansions are computed up to order k where k is obtained by adding the maximal alpha to the maximum of n and the order of self.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = QQ['x']
sage: A.<Jx> = OreAlgebra(R, 'Jx', q=2)
sage: (Jx-1).lclm((1-x)*Jx-1).power_series_solutions()
[x^2 + x^3 + 3/5*x^4 + 11/35*x^5 + O(x^6), 1 + x - 2/7*x^3 - 62/315*x^4 - 146/1395*x^5 + O(x^6)]
spread(p=0)

Returns the spread of this operator.

This is the set of integers i such that sigma(self[0], i) and sigma(self[r], -r) have a nontrivial common factor, where sigma is the shift of the parent’s algebra and r is the order of self.

If the optional argument p is given, the method is applied to gcd(self[0], p) instead of self[0].

The output set contains infty if the constant coefficient of self is zero.

This method is a stub and may not be implemented for every algebra.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']; A.<Sx> = OreAlgebra(R, 'Sx');
sage: ((x+5)*Sx - x).spread()
[4]
sage: ((x+5)*Sx - x).lclm((x+19)*Sx - x).spread()
[3, 4, 17, 18]
symmetric_product(other, solver=None)

Returns the symmetric product of self and other.

The symmetric product of two operators A and B is a minimal order operator C such that for all “functions” f and g with A.f=B.g=0 we have C.(fg)=0.

The method requires that a product rule is associated to the Ore algebra where self and other live. (See docstring of OreAlgebra for information about product rules.)

If no solver is specified, the the Ore algebra’s solver is used.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: A.<Dx> = OreAlgebra(R, 'Dx')
sage: (Dx - 1).symmetric_product(x*Dx - 1)
x*Dx - x - 1
sage: (x*Dx - 1).symmetric_product(Dx - 1)
x*Dx - x - 1
sage: ((x+1)*Dx^2 + (x-1)*Dx + 8).symmetric_product((x-1)*Dx^2 + (2*x+3)*Dx + (8*x+5))
(29*x^8 - 4*x^7 - 55*x^6 - 34*x^5 - 23*x^4 + 80*x^3 + 95*x^2 - 42*x - 46)*Dx^4 + (174*x^8 + 150*x^7 + 48*x^6 - 294*x^5 - 864*x^4 - 646*x^3 + 232*x^2 + 790*x + 410)*Dx^3 + (783*x^8 + 1661*x^7 - 181*x^6 - 1783*x^5 - 3161*x^4 - 3713*x^3 + 213*x^2 + 107*x - 1126)*Dx^2 + (1566*x^8 + 5091*x^7 + 2394*x^6 + 2911*x^5 - 10586*x^4 - 23587*x^3 - 18334*x^2 - 2047*x + 5152)*Dx + 2552*x^8 + 3795*x^7 + 8341*x^6 + 295*x^5 - 6394*x^4 - 24831*x^3 - 35327*x^2 - 23667*x - 13708

sage: A.<Sx> = OreAlgebra(R, 'Sx')
sage: (Sx - 2).symmetric_product(x*Sx - (x+1))
x*Sx - 2*x - 2
sage: (x*Sx - (x+1)).symmetric_product(Sx - 2)
x*Sx - 2*x - 2
sage: ((x+1)*Sx^2 + (x-1)*Sx + 8).symmetric_product((x-1)*Sx^2 + (2*x+3)*Sx + (8*x+5))
(-8*x^8 - 13*x^7 + 300*x^6 + 1640*x^5 + 3698*x^4 + 4373*x^3 + 2730*x^2 + 720*x)*Sx^4 + (16*x^8 + 34*x^7 - 483*x^6 - 1947*x^5 - 2299*x^4 - 2055*x^3 - 4994*x^2 - 4592*x)*Sx^3 + (-64*x^8 + 816*x^7 + 1855*x^6 - 21135*x^5 - 76919*x^4 - 35377*x^3 + 179208*x^2 + 283136*x + 125440)*Sx^2 + (1024*x^7 + 1792*x^6 - 39792*x^5 - 250472*x^4 - 578320*x^3 - 446424*x^2 + 206528*x + 326144)*Sx - 32768*x^6 - 61440*x^5 + 956928*x^4 + 4897984*x^3 + 9390784*x^2 + 7923200*x + 2329600
to_Q(alg)

Returns a q-recurrence operator which annihilates the coefficient sequence of every power series (about the origin) annihilated by self. The output operator may not be minimal.

INPUT:

  • alg – the Ore algebra in which the output should be expressed. The algebra must satisfy alg.base_ring().base_ring() == self.base_ring().base_ring() and alg.is_Q() == self.parent().is_J(). Instead of an algebra object, also a string can be passed as argument. This amounts to specifying an Ore algebra over self.base_ring() with the q-shift with respect to self.base_ring().gen().

EXAMPLES:

sage: from ore_algebra import *
sage: Rn.<n> = ZZ['n']; Rx.<x> = ZZ['x']
sage: A.<Jx> = OreAlgebra(Rx, 'Jx', q=2)
sage: B.<Qn> = OreAlgebra(Rn, 'Qn', q=2)
sage: (Jx - 1).to_Q(B)
(2*n - 1)*Qn - 1
sage: ((x+1)*Jx - 1).to_Q(B)
(4*n - 1)*Qn^2 + (2*n - 2)*Qn
sage: (n*Qn-1).to_J(A).to_Q(B) % (n*Qn - 1)
0 
class ore_algebra.ore_operator_1_1.UnivariateQRecurrenceOperatorOverUnivariateRing(parent, *data, **kwargs)

Element of an Ore algebra K(x)[S], where S is the shift x->q*x for some q in K.

annihilator_of_composition(a, solver=None)

Returns an operator L which annihilates all the sequences f(a(n)) where f runs through the functions annihilated by self. The output operator is not necessarily of smallest possible order.

INPUT:

  • a – a polynomial u*x+v where x is the generator of the base ring, u and v are integers.
  • solver (optional) – a callable object which applied to a matrix with polynomial entries returns its kernel.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = QQ['x']
sage: A.<Qx> = OreAlgebra(R, 'Qx', q=3)
sage: L = (x+3)*Qx^2 - (5*x+3)*Qx + 2*x-1
sage: data = L.to_list([1,2], 11)
sage: data
[1, 2, 15/4, 115/12, 1585/48, 19435/144, 2387975/4032, 188901875/70848, 488427432475/40336128, 1461633379710215/26500836096, 14580926901721431215/57983829378048]
sage: L2 = L.annihilator_of_composition(2*x)
sage: L2.to_list([1,15/4], 5)
[1, 15/4, 1585/48, 2387975/4032, 488427432475/40336128]
sage: Lrev = L.annihilator_of_composition(10 - x)
sage: Lrev.to_list([data[10], data[9]], 11)
[14580926901721431215/57983829378048, 1461633379710215/26500836096, 488427432475/40336128, 188901875/70848, 2387975/4032, 19435/144, 1585/48, 115/12, 15/4, 2, 1]
annihilator_of_sum()

Returns an operator L which annihilates all the indefinite sums sum_{k=0}^n a_k where a_n runs through the sequences annihilated by self. The output operator is not necessarily of smallest possible order.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['q'].fraction_field()['x']
sage: A.<Qx> = OreAlgebra(R, 'Qx')
sage: ((x+1)*Qx - x).annihilator_of_sum()
(q*x + 1)*Qx^2 + (-2*q*x - 1)*Qx + q*x
spread(p=0)

Returns the spread of this operator.

This is the set of integers i such that sigma(self[0], i) and sigma(self[r], -r) have a nontrivial common factor, where sigma is the shift of the parent’s algebra and r is the order of self.

If the optional argument p is given, the method is applied to gcd(self[0], p) instead of self[0].

The output set contains infty if the constant coefficient of self is zero.

This method is a stub and may not be implemented for every algebra.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']; A.<Sx> = OreAlgebra(R, 'Sx');
sage: ((x+5)*Sx - x).spread()
[4]
sage: ((x+5)*Sx - x).lclm((x+19)*Sx - x).spread()
[3, 4, 17, 18]
to_J(alg)

Returns a q-differential operator which annihilates every power series (about the origin) whose coefficient sequence is annihilated by self. The output operator may not be minimal.

INPUT:

  • alg – the Ore algebra in which the output should be expressed. The algebra must satisfy alg.base_ring().base_ring() == self.base_ring().base_ring() and alg.is_J() is not False. Instead of an algebra object, also a string can be passed as argument. This amounts to specifying an Ore algebra over self.base_ring() with the q-derivation with respect to self.base_ring().gen().

EXAMPLES:

sage: from ore_algebra import *
sage: Rn.<n> = ZZ['n']; Rx.<x> = ZZ['x']
sage: A.<Qn> = OreAlgebra(Rn, 'Qn', q=2)
sage: B.<Jx> = OreAlgebra(Rx, 'Jx', q=2)
sage: (Qn - 1).to_J(B)
(-2*x + 1)*Jx - 1
sage: ((n+1)*Qn - 1).to_J(B)
2*x*Jx^2 + (-4*x + 4)*Jx - 2
sage: (x*Jx-1).to_Q(A).to_J(B) % (x*Jx - 1)
0
to_list(init, n, start=0, append=False, padd=False)

Computes the terms of some sequence annihilated by self.

INPUT:

  • init – a vector (or list or tuple) of initial values. The components must be elements of self.base_ring().base_ring().fraction_field(). If the length is more than self.order(), we do not check whether the given terms are consistent with self.
  • n – desired number of terms.
  • start (optional) – index of the sequence term which is represented by the first entry of init. Defaults to zero.
  • append (optional) – if True, the computed terms are appended to init list. Otherwise (default), a new list is created.
  • padd (optional) – if True, the vector of initial values is implicitely prolonged to the left (!) by zeros if it is too short. Otherwise (default), the method raises a ValueError if init is too short.

OUTPUT:

A list of n terms whose k th component carries the sequence term with index start+k. Terms whose calculation causes an error are represented by None.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = QQ['x']; A.<Qx> = OreAlgebra(R, 'Qx', q=3)
sage: (Qx^2-x*Qx + 1).to_list([1,1], 10)
[1, 1, 0, -1, -9, -242, -19593, -4760857, -3470645160, -7590296204063]
sage: (Qx^2-x*Qx + 1)(_)
[0, 0, 0, 0, 0, 0, 0, 0]
class ore_algebra.ore_operator_1_1.UnivariateRecurrenceOperatorOverUnivariateRing(parent, *data, **kwargs)

Element of an Ore algebra K(x)[S], where S is the shift x->x+1.

annihilator_of_composition(a, solver=None)

Returns an operator L which annihilates all the sequences f(floor(a(n))) where f runs through the functions annihilated by self. The output operator is not necessarily of smallest possible order.

INPUT:

  • a – a polynomial u*x+v where x is the generator of the base ring, u and v are integers or rational numbers. If they are rational, the base ring of the parent of self must contain QQ.
  • solver (optional) – a callable object which applied to a matrix with polynomial entries returns its kernel.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = QQ['x']
sage: A.<Sx> = OreAlgebra(R, 'Sx')
sage: ((2+x)*Sx^2-(2*x+3)*Sx+(x+1)).annihilator_of_composition(2*x+5) 
(16*x^3 + 188*x^2 + 730*x + 936)*Sx^2 + (-32*x^3 - 360*x^2 - 1340*x - 1650)*Sx + 16*x^3 + 172*x^2 + 610*x + 714
sage: ((2+x)*Sx^2-(2*x+3)*Sx+(x+1)).annihilator_of_composition(1/2*x)
(x^2 + 11*x + 30)*Sx^6 + (-3*x^2 - 25*x - 54)*Sx^4 + (3*x^2 + 17*x + 26)*Sx^2 - x^2 - 3*x - 2
sage: ((2+x)*Sx^2-(2*x+3)*Sx+(x+1)).annihilator_of_composition(100-x)
(-x + 99)*Sx^2 + (2*x - 199)*Sx - x + 100
annihilator_of_interlacing(*other)

Returns an operator L which annihilates any sequence which can be obtained by interlacing sequences annihilated by self and the operators given in the arguments.

More precisely, if self and the operators given in the arguments are denoted L_1,L_2,dots,L_m, and if f_1(n),dots,f_m(n) are some sequences such that L_i annihilates f_i(n), then the output operator L annihilates sequence f_1(0),f_2(0),dots,f_m(0),f_1(1),f_2(1),dots,f_m(1),dots, the interlacing sequence of f_1(n),dots,f_m(n).

The output operator is not necessarily of smallest possible order.

The other operators must be coercible to the parent of self.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = QQ['x']
sage: A.<Sx> = OreAlgebra(R, 'Sx')
sage: (x*Sx - (x+1)).annihilator_of_interlacing(Sx - (x+1), Sx + 1)
(x^3 + 17/2*x^2 + 5/2*x - 87/2)*Sx^9 + (-1/3*x^4 - 11/2*x^3 - 53/2*x^2 - 241/6*x + 14)*Sx^6 + (7/2*x^2 + 67/2*x + 205/2)*Sx^3 + 1/3*x^4 + 13/2*x^3 + 77/2*x^2 + 457/6*x + 45
annihilator_of_sum()

Returns an operator L which annihilates all the indefinite sums sum_{k=0}^n a_k where a_n runs through the sequences annihilated by self. The output operator is not necessarily of smallest possible order.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: A.<Sx> = OreAlgebra(R, 'Sx')
sage: ((x+1)*Sx - x).annihilator_of_sum() # constructs L such that L(H_n) == 0
(x + 2)*Sx^2 + (-2*x - 3)*Sx + x + 1
forward_matrix_bsplit(n, start=0)

Uses division-free binary splitting to compute a product of n consecutive companion matrices of self.

If self annihilates some sequence c of order r, this allows rapidly computing c_n, ldots, c_{n+r-1} (or just c_n) without generating all the intermediate values.

INPUT:

  • n – desired number of terms to move forward
  • start (optional) – starting index. Defaults to zero.

OUTPUT:

A pair (M, Q) where M is an r by r matrix and Q is a scalar, such that M / Q is the product of the companion matrix at n consecutive indices.

We have Q [c_{s+n}, ldots, c_{s+r-1+n}]^T = M [c_s, c_{s+1}, ldots, c_{s+r-1}]^T, where s is the initial position given by start.

EXAMPLES:

sage: from ore_algebra import *
sage: R = ZZ
sage: Rx.<x> = R[]
sage: Rxk.<k> = Rx[]
sage: Rxks = OreAlgebra(Rxk, 'Sk')
sage: ann = Rxks([1+k, -3*x - 2*k*x, 2+k])
sage: initial = Matrix([[1], [x]])
sage: M, Q = ann.forward_matrix_bsplit(5)
sage: (M * initial).change_ring(QQ['x']) / Q
[               63/8*x^5 - 35/4*x^3 + 15/8*x]
[231/16*x^6 - 315/16*x^4 + 105/16*x^2 - 5/16]

sage: Matrix([[legendre_P(5, x)], [legendre_P(6, x)]])
[               63/8*x^5 - 35/4*x^3 + 15/8*x]
[231/16*x^6 - 315/16*x^4 + 105/16*x^2 - 5/16]


sage: Sk = Rxks.gen()
sage: (Sk^2 - 1).forward_matrix_param_rectangular(1, 10)
(
[1 0]
[0 1], 1
)

TODO: this should detect if the base coefficient ring is QQ (etc.) and then switch to ZZ (etc.) internally.

forward_matrix_param_rectangular(value, n, start=0, m=None)

Assuming the coefficients of self are in R[x][k], computes the nth forward matrix with the parameter x evaluated at value, using rectangular splitting with a step size of m.

TESTS:

sage: from sage.all import Matrix, randrange
sage: from ore_algebra import *
sage: R = ZZ
sage: Rx = R['x']; x = Rx.gen()
sage: Rxk = Rx['k']; k = Rxk.gen()
sage: Rxks = OreAlgebra(Rxk, 'Sk')
sage: V = QQ
sage: Vks = OreAlgebra(V['k'], 'Sk')
sage: for i in range(1000): # long time (2.5 s)
....:     A = Rxks.random_element(randrange(1,4))
....:     r = A.order()
....:     v = V.random_element()
....:     initial = [V.random_element() for i in range(r)]
....:     start = randrange(0,5)
....:     n = randrange(0,30)
....:     m = randrange(0,10)
....:     B = Vks(list(A.polynomial()(x=v)))
....:     M, Q = A.forward_matrix_param_rectangular(v, n, m=m, start=start)
....:     if Q != 0:
....:         V1 = M * Matrix(initial).transpose() / Q
....:         values = B.to_list(initial, n + r, start)
....:         V2 = Matrix(values[-r:]).transpose()
....:         if V1 != V2:
....:             raise ValueError
generalized_series_solutions(n=5, dominant_only=False, real_only=False, infolevel=0)

Returns the generalized series solutions of this operator.

These are solutions of the form

(x/e)^{x u/v}rho^xexpbigl(c_1 x^{1/m} +…+ c_{v-1} x^{1-1/m}bigr)x^alpha p(x^{-1/m},log(x))

where

  • e is Euler’s constant (2.71…)
  • v is a positive integer
  • u is an integer; the term (x/e)^(v/u) is called the “superexponential part” of the solution
  • rho is an element of an algebraic extension of the coefficient field K (the algebra’s base ring’s base ring); the term rho^x is called the “exponential part” of the solution
  • c_1,…,c_{v-1} are elements of K(rho); the term exp(…) is called the “subexponential part” of the solution
  • m is a positive integer multiple of v, it is called the object’s “ramification”
  • alpha is an element of some algebraic extension of K(rho); the term n^alpha is called the “polynomial part” of the solution (even if alpha is not an integer)
  • p is an element of K(rho)(alpha)[[x]][y]. It is called the “expansion part” of the solution.

An operator of order r has exactly r linearly independent solutions of this form. This method computes them all, unless the flags specified in the arguments rule out some of them.

Generalized series solutions are asymptotic expansions of sequences annihilated by the operator.

At present, the method only works for operators where K is some field which supports coercion to QQbar.

INPUT:

  • n (default: 5) – minimum number of terms in the expansions parts to be computed.
  • dominant_only (default: False) – if set to True, only compute solution(s) with maximal growth.
  • real_only (default: False) – if set to True, only compute solution(s) where rho,c_1,…,c_{v-1},alpha are real.
  • infolevel (default: 0) – if set to a positive integer, the methods prints some messages about the progress of the computation.

OUTPUT:

  • a list of DiscreteGeneralizedSeries objects forming a fundamental system for this operator.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<n> = QQ['n']; A.<Sn> = OreAlgebra(R, 'Sn')
sage: (Sn - (n+1)).generalized_series_solutions()
[(n/e)^n*n^(1/2)*(1 + 1/12*n^(-1) + 1/288*n^(-2) - 139/51840*n^(-3) - 571/2488320*n^(-4) + O(n^(-5)))]
sage: list(map(Sn - (n+1), _))
[0]

sage: L = ((n+1)*Sn - n).annihilator_of_sum().symmetric_power(2)
sage: L.generalized_series_solutions()
[1 + O(n^(-5)),
 (1 + O(n^(-5)))*log(n) + 1/2*n^(-1) - 1/12*n^(-2) + 1/120*n^(-4) + O(n^(-5)),
 (1 + O(n^(-5)))*log(n)^2 + (n^(-1) - 1/6*n^(-2) + 1/60*n^(-4) + O(n^(-5)))*log(n) + 1/4*n^(-2) - 1/12*n^(-3) + 1/144*n^(-4) + O(n^(-5))]
sage: list(map(L, _))
[0, 0, 0]

sage: L = n^2*(1-2*Sn+Sn^2) + (n+1)*(1+Sn+Sn^2)
sage: L.generalized_series_solutions() # long time
[exp(3.464101615137755?*I*n^(1/2))*n^(1/4)*(1 - 2.056810333988042?*I*n^(-1/2) - 1107/512*n^(-2/2) + (0.?e-19 + 1.489453749877895?*I)*n^(-3/2) + 2960239/2621440*n^(-4/2) + (0.?e-19 - 0.926161373412572?*I)*n^(-5/2) - 16615014713/46976204800*n^(-6/2) + (0.?e-20 + 0.03266142931818572?*I)*n^(-7/2) + 16652086533741/96207267430400*n^(-8/2) + (0.?e-20 - 0.1615093987591473?*I)*n^(-9/2) + O(n^(-10/2))), exp(-3.464101615137755?*I*n^(1/2))*n^(1/4)*(1 + 2.056810333988042?*I*n^(-1/2) - 1107/512*n^(-2/2) + (0.?e-19 - 1.489453749877895?*I)*n^(-3/2) + 2960239/2621440*n^(-4/2) + (0.?e-19 + 0.926161373412572?*I)*n^(-5/2) - 16615014713/46976204800*n^(-6/2) + (0.?e-20 - 0.03266142931818572?*I)*n^(-7/2) + 16652086533741/96207267430400*n^(-8/2) + (0.?e-20 + 0.1615093987591473?*I)*n^(-9/2) + O(n^(-10/2)))]

sage: L = guess([(-3)^k*(k+1)/(2*k+4) - 2^k*k^3/(k+3) for k in range(500)], A)
sage: L.generalized_series_solutions()
[2^n*n^2*(1 - 3*n^(-1) + 9*n^(-2) - 27*n^(-3) + 81*n^(-4) + O(n^(-5))), (-3)^n*(1 - n^(-1) + 2*n^(-2) - 4*n^(-3) + 8*n^(-4) + O(n^(-5)))]
sage: L.generalized_series_solutions(dominant_only=True)
[(-3)^n*(1 - n^(-1) + 2*n^(-2) - 4*n^(-3) + 8*n^(-4) + O(n^(-5)))]

TESTS:

sage: rop = (-8 -12*Sn + (n^2+5*n+6)*Sn^3)
sage: rop
(n^2 + 5*n + 6)*Sn^3 - 12*Sn - 8
sage: rop.generalized_series_solutions(1) # long time
[(n/e)^(-2/3*n)*2^n*exp(3*n^(1/3))*n^(-2/3)*(1 + 3/2*n^(-1/3) + 9/8*n^(-2/3) + O(n^(-3/3))),
(n/e)^(-2/3*n)*(-1.000000000000000? + 1.732050807568878?*I)^n*exp((-1.500000000000000? + 2.598076211353316?*I)*n^(1/3))*n^(-2/3)*(1 + (-0.750000000000000? - 1.299038105676658?*I)*n^(-1/3) + (-0.562500000000000? + 0.974278579257494?*I)*n^(-2/3) + O(n^(-3/3))),
(n/e)^(-2/3*n)*(-1.000000000000000? - 1.732050807568878?*I)^n*exp((-1.500000000000000? - 2.598076211353316?*I)*n^(1/3))*n^(-2/3)*(1 + (-0.750000000000000? + 1.299038105676658?*I)*n^(-1/3) + (-0.562500000000000? - 0.974278579257494?*I)*n^(-2/3) + O(n^(-3/3)))]
spread(p=0)

Returns the spread of this operator.

This is the set of integers i such that sigma(self[0], i) and sigma(self[r], -r) have a nontrivial common factor, where sigma is the shift of the parent’s algebra and r is the order of self.

If the optional argument p is given, the method is applied to gcd(self[0], p) instead of self[0].

The output set contains infty if the constant coefficient of self is zero.

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']; A.<Sx> = OreAlgebra(R, 'Sx');
sage: ((x+5)*Sx - x).spread()
[4]
sage: ((x+5)*Sx - x).lclm((x+19)*Sx - x).spread()
[3, 4, 17, 18]
to_D(alg)

Returns a differential operator which annihilates every power series whose coefficient sequence is annihilated by self. The output operator may not be minimal.

INPUT:

  • alg – the Ore algebra in which the output should be expressed. The algebra must satisfy alg.base_ring().base_ring() == self.base_ring().base_ring() and alg.is_D() is not False. Instead of an algebra object, also a string can be passed as argument. This amounts to specifying an Ore algebra over self.base_ring() with the standard derivation with respect to self.base_ring().gen().

EXAMPLES:

sage: from ore_algebra import *
sage: Rn.<n> = ZZ['n']; Rx.<x> = ZZ['x']
sage: A.<Sn> = OreAlgebra(Rn, 'Sn')
sage: B.<Dx> = OreAlgebra(Rx, 'Dx')
sage: (Sn - 1).to_D(B)
(-x + 1)*Dx - 1
sage: ((n+1)*Sn - 1).to_D(B)
x*Dx^2 + (-x + 1)*Dx - 1
sage: (x*Dx-1).to_S(A).to_D(B)
x*Dx - 1
to_F(alg)

Returns the difference operator corresponding to self

INPUT:

  • alg – the Ore algebra in which the output should be expressed. The algebra must satisfy alg.base_ring().base_ring() == self.base_ring().base_ring() and alg.is_F() is not False. Instead of an algebra object, also a string can be passed as argument. This amounts to specifying an Ore algebra over self.base_ring() with the forward difference with respect to self.base_ring().gen().

EXAMPLES:

sage: from ore_algebra import *
sage: R.<x> = ZZ['x']
sage: A.<Sx> = OreAlgebra(R, 'Sx')
sage: (Sx^4).to_F(OreAlgebra(R, 'Fx'))
Fx^4 + 4*Fx^3 + 6*Fx^2 + 4*Fx + 1
sage: (Sx^4).to_F('Fx').to_S(A)
Sx^4
to_T(alg)

Returns a differential operator, expressed in terms of the Euler derivation, which annihilates every power series (about the origin) whose coefficient sequence is annihilated by self. The output operator may not be minimal.

INPUT:

  • alg – the Ore algebra in which the output should be expressed. The algebra must satisfy alg.base_ring().base_ring() == self.base_ring().base_ring() and alg.is_T() is not False. Instead of an algebra object, also a string can be passed as argument. This amounts to specifying an Ore algebra over self.base_ring() with the Euler derivation with respect to self.base_ring().gen().

EXAMPLES:

sage: from ore_algebra import *
sage: Rn.<n> = ZZ['n']; Rx.<x> = ZZ['x']
sage: A.<Sn> = OreAlgebra(Rn, 'Sn')
sage: B.<Tx> = OreAlgebra(Rx, 'Tx')
sage: (Sn - 1).to_T(B)
(-x + 1)*Tx - x
sage: ((n+1)*Sn - 1).to_T(B)
Tx^2 - x*Tx - x
sage: (x*Tx-1).to_S(A).to_T(B)
x*Tx^2 + (x - 1)*Tx
to_list(init, n, start=0, append=False, padd=False)

Computes the terms of some sequence annihilated by self.

INPUT:

  • init – a vector (or list or tuple) of initial values. The components must be elements of self.base_ring().base_ring().fraction_field(). If the length is more than self.order(), we do not check whether the given terms are consistent with self.
  • n – desired number of terms.
  • start (optional) – index of the sequence term which is represented by the first entry of init. Defaults to zero.
  • append (optional) – if True, the computed terms are appended to init list. Otherwise (default), a new list is created.
  • padd (optional) – if True, the vector of initial values is implicitely prolonged to the left (!) by zeros if it is too short. Otherwise (default), the method raises a ValueError if init is too short.

OUTPUT:

A list of n terms whose k th component carries the sequence term with index start+k. Terms whose calculation causes an error are represented by None.

EXAMPLES:

sage: from ore_algebra import *
sage: R = ZZ['x']['n']; x = R('x'); n = R('n')
sage: A.<Sn> = OreAlgebra(R, 'Sn')
sage: L = ((n+2)*Sn^2 - x*(2*n+3)*Sn + (n+1))
sage: L.to_list([1, x], 5)
[1, x, (3*x^2 - 1)/2, (5*x^3 - 3*x)/2, (35*x^4 - 30*x^2 + 3)/8]
sage: polys = L.to_list([1], 5, padd=True)
sage: polys
[1, x, (3*x^2 - 1)/2, (5*x^3 - 3*x)/2, (35*x^4 - 30*x^2 + 3)/8]
sage: L.to_list([polys[3], polys[4]], 8, start=3)
[(5*x^3 - 3*x)/2,
 (35*x^4 - 30*x^2 + 3)/8,
 (63*x^5 - 70*x^3 + 15*x)/8,
 (231*x^6 - 315*x^4 + 105*x^2 - 5)/16,
 (429*x^7 - 693*x^5 + 315*x^3 - 35*x)/16,
 (6435*x^8 - 12012*x^6 + 6930*x^4 - 1260*x^2 + 35)/128,
 (12155*x^9 - 25740*x^7 + 18018*x^5 - 4620*x^3 + 315*x)/128,
 (46189*x^10 - 109395*x^8 + 90090*x^6 - 30030*x^4 + 3465*x^2 - 63)/256]
sage: ((n-5)*Sn - 1).to_list([1], 10)
[1, 1/-5, 1/20, 1/-60, 1/120, -1/120, None, None, None, None]