Coverage for local/lib/python2.7/site-packages/sage/rings/cfinite_sequence.py : 73%

# -*- coding: utf-8 -*- 

r""" 

C-Finite Sequences 

C-finite infinite sequences satisfy homogenous linear recurrences with constant coefficients: 

.. MATH:: 

a_{n+d} = c_0a_n + c_1a_{n+1} + \cdots + c_{d-1}a_{n+d-1}, \quad d>0. 

CFiniteSequences are completely defined by their ordinary generating function (o.g.f., which 

is always a :mod:`fraction <sage.rings.fraction_field_element>` of 

:mod:`polynomials <sage.rings.polynomial.polynomial_element>` over `\mathbb{Z}` or `\mathbb{Q}` ). 

EXAMPLES:: 

sage: fibo = CFiniteSequence(x/(1-x-x^2)) # the Fibonacci sequence 

sage: fibo 

C-finite sequence, generated by x/(-x^2 - x + 1) 

sage: fibo.parent() 

The ring of C-Finite sequences in x over Rational Field 

sage: fibo.parent().category() 

Category of commutative rings 

sage: C.<x> = CFiniteSequences(QQ); 

sage: fibo.parent() == C 

True 

sage: C 

The ring of C-Finite sequences in x over Rational Field 

sage: C(x/(1-x-x^2)) 

C-finite sequence, generated by x/(-x^2 - x + 1) 

sage: C(x/(1-x-x^2)) == fibo 

True 

sage: var('y') 

y 

sage: CFiniteSequence(y/(1-y-y^2)) 

C-finite sequence, generated by y/(-y^2 - y + 1) 

sage: CFiniteSequence(y/(1-y-y^2)) == fibo 

False 

Finite subsets of the sequence are accessible via python slices:: 

sage: fibo[137] #the 137th term of the Fibonacci sequence 

19134702400093278081449423917 

sage: fibo[137] == fibonacci(137) 

True 

sage: fibo[0:12] 

[0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89] 

sage: fibo[14:4:-2] 

[377, 144, 55, 21, 8] 

They can be created also from the coefficients and start values of a recurrence:: 

sage: r = C.from_recurrence([1,1],[0,1]) 

sage: r == fibo 

True 

Given enough values, the o.g.f. of a C-finite sequence 

can be guessed:: 

sage: r = C.guess([0,1,1,2,3,5,8]) 

sage: r == fibo 

True 

.. SEEALSO:: 

:func:`fibonacci`, :class:`BinaryRecurrenceSequence` 

AUTHORS: 

- Ralf Stephan (2014): initial version 

REFERENCES: 

.. [GK82] Daniel H. Greene and Donald E. (1982), "2.1.1 Constant 

coefficients - A) Homogeneous equations", Mathematics for the Analysis 

of Algorithms (2nd ed.), Birkhauser, p. 17. 

.. [KP11] Manuel Kauers and Peter Paule. The Concrete Tetrahedron. 

Springer-Verlag, 2011. 

.. [SZ94] Bruno Salvy and Paul Zimmermann. - Gfun: a Maple package for 

the manipulation of generating and holonomic functions in one variable. 

- Acm transactions on mathematical software, 20.2:163-177, 1994. 

.. [Z11] Doron Zeilberger. "The C-finite ansatz." The Ramanujan Journal 

(2011): 1-10. 

""" 

#***************************************************************************** 

# Copyright (C) 2014 Ralf Stephan <gtrwst9@gmail.com> 

# 

# This program is free software: you can redistribute it and/or modify 

# it under the terms of the GNU General Public License as published by 

# the Free Software Foundation, either version 2 of the License, or 

# (at your option) any later version. 

# http://www.gnu.org/licenses/ 

#***************************************************************************** 

from six.moves import range 

from six import add_metaclass 

from sage.categories.fields import Fields 

from sage.misc.inherit_comparison import InheritComparisonClasscallMetaclass 

from sage.rings.ring import CommutativeRing 

from sage.rings.integer import Integer 

from sage.rings.integer_ring import ZZ 

from sage.rings.rational_field import QQ 

from sage.arith.all import gcd 

from sage.rings.polynomial.polynomial_ring_constructor import PolynomialRing 

from sage.rings.polynomial.polynomial_ring import PolynomialRing_general 

from sage.rings.laurent_series_ring import LaurentSeriesRing 

from sage.rings.power_series_ring import PowerSeriesRing 

from sage.rings.fraction_field import FractionField 

from sage.structure.element import FieldElement 

from sage.structure.unique_representation import UniqueRepresentation 

from sage.interfaces.gp import Gp 

from sage.misc.all import sage_eval 

_gp = None 

def CFiniteSequences(base_ring, names = None, category = None): 

r""" 

Return the ring of C-Finite sequences. 

The ring is defined over a base ring (`\mathbb{Z}` or `\mathbb{Q}` ) 

and each element is represented by its ordinary generating function (ogf) 

which is a rational function over the base ring. 

INPUT: 

- ``base_ring`` -- the base ring to construct the fraction field 

representing the C-Finite sequences 

- ``names`` -- (optional) the list of variables. 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: C 

The ring of C-Finite sequences in x over Rational Field 

sage: C.an_element() 

C-finite sequence, generated by (-x + 2)/(-x^2 - x + 1) 

sage: C.category() 

Category of commutative rings 

sage: C.one() 

Finite sequence [1], offset = 0 

sage: C.zero() 

Constant infinite sequence 0. 

sage: C(x) 

Finite sequence [1], offset = 1 

sage: C(1/x) 

Finite sequence [1], offset = -1 

sage: C((-x + 2)/(-x^2 - x + 1)) 

C-finite sequence, generated by (-x + 2)/(-x^2 - x + 1) 

TESTS:: 

sage: TestSuite(C).run() 

""" 

if isinstance(base_ring, PolynomialRing_general): 

polynomial_ring = base_ring 

base_ring = polynomial_ring.base_ring() 

if names is None: 

names = ['x'] 

elif len(names)>1: 

raise NotImplementedError("Multidimensional o.g.f. not implemented.") 

if category is None: 

category = Fields() 

if not(base_ring in (QQ, ZZ)): 

raise ValueError("O.g.f. base not rational.") 

polynomial_ring = PolynomialRing(base_ring, names) 

return CFiniteSequences_generic(polynomial_ring, category) 

@add_metaclass(InheritComparisonClasscallMetaclass) 

class CFiniteSequence(FieldElement): 

r""" 

Create a C-finite sequence given its ordinary generating function. 

INPUT: 

- ``ogf`` -- a rational function, the ordinary generating function 

(can be a an element from the symbolic ring, fraction field or polynomial 

ring) 

OUTPUT: 

- A CFiniteSequence object 

EXAMPLES:: 

sage: CFiniteSequence((2-x)/(1-x-x^2)) # the Lucas sequence 

C-finite sequence, generated by (-x + 2)/(-x^2 - x + 1) 

sage: CFiniteSequence(x/(1-x)^3) # triangular numbers 

C-finite sequence, generated by x/(-x^3 + 3*x^2 - 3*x + 1) 

Polynomials are interpreted as finite sequences, or recurrences of degree 0:: 

sage: CFiniteSequence(x^2-4*x^5) 

Finite sequence [1, 0, 0, -4], offset = 2 

sage: CFiniteSequence(1) 

Finite sequence [1], offset = 0 

This implementation allows any polynomial fraction as o.g.f. by interpreting 

any power of `x` dividing the o.g.f. numerator or denominator as a right or left shift 

of the sequence offset:: 

sage: CFiniteSequence(x^2+3/x) 

Finite sequence [3, 0, 0, 1], offset = -1 

sage: CFiniteSequence(1/x+4/x^3) 

Finite sequence [4, 0, 1], offset = -3 

sage: P = LaurentPolynomialRing(QQ.fraction_field(), 'X') 

sage: X=P.gen() 

sage: CFiniteSequence(1/(1-X)) 

C-finite sequence, generated by 1/(-X + 1) 

The o.g.f. is always normalized to get a denominator constant coefficient of `+1`:: 

sage: CFiniteSequence(1/(x-2)) 

C-finite sequence, generated by -1/2/(-1/2*x + 1) 

The given ``ogf`` is used to create an appropriate parent: it can 

be a symbolic expression, a polynomial , or a fraction field element 

as long as it can be coerced into a proper fraction field over the 

rationals:: 

sage: var('x') 

x 

sage: f1 = CFiniteSequence((2-x)/(1-x-x^2)) 

sage: P.<x> = QQ[] 

sage: f2 = CFiniteSequence((2-x)/(1-x-x^2)) 

sage: f1 == f2 

True 

sage: f1.parent() 

The ring of C-Finite sequences in x over Rational Field 

sage: f1.ogf().parent() 

Fraction Field of Univariate Polynomial Ring in x over Rational Field 

sage: CFiniteSequence(log(x)) 

Traceback (most recent call last): 

... 

TypeError: unable to convert log(x) to a rational 

TESTS:: 

sage: P.<x> = QQ[] 

sage: CFiniteSequence(0.1/(1-x)) 

C-finite sequence, generated by 1/10/(-x + 1) 

sage: CFiniteSequence(pi/(1-x)) 

Traceback (most recent call last): 

... 

TypeError: unable to convert -pi to a rational 

sage: P.<x,y> = QQ[] 

sage: CFiniteSequence(x*y) 

Traceback (most recent call last): 

... 

NotImplementedError: Multidimensional o.g.f. not implemented. 

""" 

@staticmethod 

def __classcall_private__(cls, ogf): 

r""" 

Ensures that elements created by :class:`CFiniteSequence` have the same 

parent than the ones created by the parent itself and follow the category 

framework (they should be instance of :class:`CFiniteSequences` automatic 

element class). 

This method is called before the ``__init__`` method, it checks the 

o.g.f to create the appropriate parent. 

INPUT: 

- ``ogf`` - a rational function 

TESTS:: 

sage: f1 = CFiniteSequence((2-x)/(1-x-x^2)) 

sage: f1 

C-finite sequence, generated by (-x + 2)/(-x^2 - x + 1) 

sage: C.<x> = CFiniteSequences(QQ); 

sage: f2 = CFiniteSequence((2-x)/(1-x-x^2)) 

sage: f2 

C-finite sequence, generated by (-x + 2)/(-x^2 - x + 1) 

sage: f3 = C((2-x)/(1-x-x^2)) 

sage: f3 

C-finite sequence, generated by (-x + 2)/(-x^2 - x + 1) 

sage: f1 == f2 and f2 == f3 

True 

sage: f1.parent() == f2.parent() and f2.parent() == f3.parent() 

True 

sage: type(f1) 

<class 'sage.rings.cfinite_sequence.CFiniteSequences_generic_with_category.element_class'> 

sage: type(f1) == type(f2) and type(f2) == type(f3) 

True 

sage: CFiniteSequence(log(x)) 

Traceback (most recent call last): 

... 

TypeError: unable to convert log(x) to a rational 

sage: CFiniteSequence(pi) 

Traceback (most recent call last): 

... 

TypeError: unable to convert pi to a rational 

sage: var('y') 

y 

sage: f4 = CFiniteSequence((2-y)/(1-y-y^2)) 

sage: f4 

C-finite sequence, generated by (-y + 2)/(-y^2 - y + 1) 

sage: f4 == f1 

False 

sage: f4.parent() == f1.parent() 

False 

sage: f4.parent() 

The ring of C-Finite sequences in y over Rational Field 

""" 

br = ogf.base_ring() 

if not(br in (QQ, ZZ)): 

br = QQ # if the base ring of the o.g.f is not QQ, we force it to QQ and see if the o.g.f converts nicely 

# trying to figure out the ogf variables 

variables = [] 

if not ogf in br: 

if hasattr(ogf, 'variables'): 

variables = ogf.variables() 

elif hasattr(ogf.parent(), 'gens'): 

variables = ogf.parent().gens() 

# for some reason, fraction field elements don't have the variables 

# method, but symbolic elements don't have the gens method so we check both 

if len(variables)==0: 

parent = CFiniteSequences(QQ) # if we cannot find variables, we create the default parent (with x) 

else: 

parent = CFiniteSequences(QQ, variables) 

return parent(ogf) # if ogf cannot be converted to a fraction field, this will break and raise the proper error 

def __init__(self, parent, ogf): 

r""" 

Initialize the C-Finite sequence. 

The ``__init__`` method can only be called by the :class:`CFiniteSequences` 

class. By Default, a class call reaches the ``__classcall_private__`` 

which first creates a proper parent and then call the ``__init__``. 

INPUT: 

- ``ogf`` -- the ordinary generating function, a fraction of polynomials over the rationals 

- ``parent`` -- the parent of the C-Finite sequence, an occurrence of :class:`CFiniteSequences` 

OUTPUT: 

- A CFiniteSequence object 

TESTS:: 

sage: C.<x> = CFiniteSequences(QQ); 

sage: C((2-x)/(1-x-x^2)) # indirect doctest 

C-finite sequence, generated by (-x + 2)/(-x^2 - x + 1) 

""" 

br = parent.base_ring() 

ogf = parent.fraction_field()(ogf) 

P = parent.polynomial_ring() 

num = ogf.numerator() 

den = ogf.denominator() 

FieldElement.__init__(self, parent) 

if den == 1: 

self._c = [] 

self._off = num.valuation() 

self._deg = 0 

if ogf == 0: 

self._a = [0] 

else: 

self._a = P((num / (P.gen()) ** self._off)).list() 

else: 

# Transform the ogf numerator and denominator to canonical form 

# to get the correct offset, degree, and recurrence coeffs and 

# start values. 

self._off = 0 

self._deg = 0 

x = P.gen() 

if num.constant_coefficient() == 0: 

self._off = num.valuation() 

num = P(num / x ** self._off) 

elif den.constant_coefficient() == 0: 

self._off = -den.valuation() 

den = P(den * x ** self._off) 

f = den.constant_coefficient() 

num = P(num / f) 

den = P(den / f) 

f = gcd(num, den) 

num = P(num / f) 

den = P(den / f) 

self._deg = den.degree() 

self._c = [-den.list()[i] for i in range(1, self._deg + 1)] 

if self._off >= 0: 

num = x ** self._off * num 

else: 

den = x ** (-self._off) * den 

# determine start values (may be different from _get_item_ values) 

alen = max(self._deg, num.degree() + 1) 

R = LaurentSeriesRing(br, parent.variable_name(), default_prec=alen) 

rem = num % den 

if den != 1: 

self._a = R(num / den).list() 

self._aa = R(rem / den).list()[:self._deg] # needed for _get_item_ 

else: 

self._a = num.list() 

if len(self._a) < alen: 

self._a.extend([0] * (alen - len(self._a))) 

ogf = num / den 

self._ogf = ogf 

def _repr_(self): 

""" 

Return textual definition of sequence. 

TESTS:: 

sage: CFiniteSequence(1/x^5) 

Finite sequence [1], offset = -5 

sage: CFiniteSequence(x^3) 

Finite sequence [1], offset = 3 

""" 

if self._deg == 0: 

if self.ogf() == 0: 

return 'Constant infinite sequence 0.' 

else: 

return 'Finite sequence ' + str(self._a) + ', offset = ' + str(self._off) 

else: 

return 'C-finite sequence, generated by ' + str(self.ogf()) 

def __hash__(self): 

r""" 

Hash value for C finite sequence. 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: hash(C((2-x)/(1-x-x^2))) # random 

42 

""" 

return hash(self.parent()) ^ hash(self._ogf) 

def _add_(self, other): 

""" 

Addition of C-finite sequences. 

TESTS:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: r = C(1/(1-2*x)) 

sage: r[0:5] # a(n) = 2^n 

[1, 2, 4, 8, 16] 

sage: s = C.from_recurrence([1],[1]) 

sage: (r + s)[0:5] # a(n) = 2^n + 1 

[2, 3, 5, 9, 17] 

sage: r + 0 == r 

True 

sage: (r + x^2)[0:5] 

[1, 2, 5, 8, 16] 

sage: (r + 3/x)[-1] 

3 

sage: r = CFiniteSequence(x) 

sage: r + 0 == r 

True 

sage: CFiniteSequence(0) + CFiniteSequence(0) 

Constant infinite sequence 0. 

""" 

return CFiniteSequence(self.ogf() + other.numerator() / other.denominator()) 

def _sub_(self, other): 

""" 

Subtraction of C-finite sequences. 

TESTS:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: r = C(1/(1-2*x)) 

sage: r[0:5] # a(n) = 2^n 

[1, 2, 4, 8, 16] 

sage: s = C.from_recurrence([1],[1]) 

sage: (r - s)[0:5] # a(n) = 2^n + 1 

[0, 1, 3, 7, 15] 

""" 

return CFiniteSequence(self.ogf() - other.numerator() / other.denominator()) 

def _mul_(self, other): 

""" 

Multiplication of C-finite sequences. 

TESTS:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: r = C.guess([1,2,3,4,5,6]) 

sage: (r*r)[0:6] # self-convolution 

[1, 4, 10, 20, 35, 56] 

sage: r = C(x) 

sage: r*1 == r 

True 

sage: r*-1 

Finite sequence [-1], offset = 1 

sage: C(0) * C(1) 

Constant infinite sequence 0. 

""" 

return CFiniteSequence(self.ogf() * other.numerator() / other.denominator()) 

def _div_(self, other): 

""" 

Division of C-finite sequences. 

TESTS:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: r = C.guess([1,2,3,4,5,6]) 

sage: (r/2)[0:6] 

[1/2, 1, 3/2, 2, 5/2, 3] 

sage: s = C(x) 

sage: s/(s*-1 + 1) 

C-finite sequence, generated by x/(-x + 1) 

""" 

return CFiniteSequence(self.ogf() / (other.numerator() / other.denominator())) 

def coefficients(self): 

""" 

Return the coefficients of the recurrence representation of the 

C-finite sequence. 

OUTPUT: 

- A list of values 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: lucas = C((2-x)/(1-x-x^2)) # the Lucas sequence 

sage: lucas.coefficients() 

[1, 1] 

""" 

return self._c 

def __eq__(self, other): 

""" 

Compare two CFiniteSequences. 

EXAMPLES:: 

sage: f = CFiniteSequence((2-x)/(1-x-x^2)) 

sage: f2 = CFiniteSequence((2-x)/(1-x-x^2)) 

sage: f == f2 

True 

sage: f == (2-x)/(1-x-x^2) 

False 

sage: (2-x)/(1-x-x^2) == f 

False 

sage: C.<x> = CFiniteSequences(QQ) 

sage: r = C.from_recurrence([1,1],[2,1]) 

sage: s = C.from_recurrence([-1],[1]) 

sage: r == s 

False 

sage: r = C.from_recurrence([-1],[1]) 

sage: s = C(1/(1+x)) 

sage: r == s 

True 

""" 

if not isinstance(other, CFiniteSequence): 

return False 

return self.ogf() == other.ogf() 

def __ne__(self, other): 

""" 

Compare two CFiniteSequences. 

EXAMPLES:: 

sage: f = CFiniteSequence((2-x)/(1-x-x^2)) 

sage: f2 = CFiniteSequence((2-x)/(1-x-x^2)) 

sage: f != f2 

False 

sage: f != (2-x)/(1-x-x^2) 

True 

sage: (2-x)/(1-x-x^2) != f 

True 

sage: C.<x> = CFiniteSequences(QQ) 

sage: r = C.from_recurrence([1,1],[2,1]) 

sage: s = C.from_recurrence([-1],[1]) 

sage: r != s 

True 

sage: r = C.from_recurrence([-1],[1]) 

sage: s = C(1/(1+x)) 

sage: r != s 

False 

""" 

return not self.__eq__(other) 

def __getitem__(self, key): 

r""" 

Return a slice of the sequence. 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: r = C.from_recurrence([3,3],[2,1]) 

sage: r[2] 

9 

sage: r[101] 

16158686318788579168659644539538474790082623100896663971001 

sage: r = C(1/(1-x)) 

sage: r[5] 

1 

sage: r = C(x) 

sage: r[0] 

0 

sage: r[1] 

1 

sage: r = C(0) 

sage: r[66] 

0 

sage: lucas = C.from_recurrence([1,1],[2,1]) 

sage: lucas[5:10] 

[11, 18, 29, 47, 76] 

sage: r = C((2-x)/x/(1-x-x*x)) 

sage: r[0:4] 

[1, 3, 4, 7] 

sage: r = C(1-2*x^2) 

sage: r[0:4] 

[1, 0, -2, 0] 

sage: r[-1:4] # not tested, python will not allow this! 

[0, 1, 0 -2, 0] 

sage: r = C((-2*x^3 + x^2 + 1)/(-2*x + 1)) 

sage: r[0:5] # handle ogf > 1 

[1, 2, 5, 8, 16] 

sage: r[-2] 

0 

sage: r = C((-2*x^3 + x^2 - x + 1)/(2*x^2 - 3*x + 1)) 

sage: r[0:5] 

[1, 2, 5, 9, 17] 

sage: s=C((1-x)/(-x^2 - x + 1)) 

sage: s[0:5] 

[1, 0, 1, 1, 2] 

sage: s=C((1+x^20+x^40)/(1-x^12)/(1-x^30)) 

sage: s[0:20] 

[1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0] 

sage: s=C(1/((1-x^2)*(1-x^6)*(1-x^8)*(1-x^12))) 

sage: s[999998] 

289362268629630 

""" 

if isinstance(key, slice): 

m = max(key.start, key.stop) 

return [self[ii] for ii in range(*key.indices(m + 1))] 

elif isinstance(key, (int, Integer)): 

from sage.matrix.constructor import Matrix 

d = self._deg 

if (self._off <= key and key < self._off + len(self._a)): 

return self._a[key - self._off] 

elif d == 0: 

return 0 

(quo, rem) = self.numerator().quo_rem(self.denominator()) 

wp = quo[key - self._off] 

if key < self._off: 

return wp 

A = Matrix(QQ, 1, d, self._c) 

B = Matrix.identity(QQ, d - 1) 

C = Matrix(QQ, d - 1, 1, 0) 

if quo == 0: 

V = Matrix(QQ, d, 1, self._a[:d][::-1]) 

else: 

V = Matrix(QQ, d, 1, self._aa[:d][::-1]) 

M = Matrix.block([[A], [B, C]], subdivide=False) 

return wp + list(M ** (key - self._off) * V)[d - 1][0] 

else: 

raise TypeError("invalid argument type") 

def ogf(self): 

""" 

Return the ordinary generating function associated with 

the CFiniteSequence. 

This is always a fraction of polynomials in the base ring. 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: r = C.from_recurrence([2],[1]) 

sage: r.ogf() 

1/(-2*x + 1) 

sage: C(0).ogf() 

0 

""" 

return self._ogf 

def numerator(self): 

r""" 

Return the numerator of the o.g.f of ``self``. 

EXAMPLES:: 

sage: f = CFiniteSequence((2-x)/(1-x-x^2)); f 

C-finite sequence, generated by (-x + 2)/(-x^2 - x + 1) 

sage: f.numerator() 

-x + 2 

""" 

return self.ogf().numerator() 

def denominator(self): 

r""" 

Return the numerator of the o.g.f of ``self``. 

EXAMPLES:: 

sage: f = CFiniteSequence((2-x)/(1-x-x^2)); f 

C-finite sequence, generated by (-x + 2)/(-x^2 - x + 1) 

sage: f.denominator() 

-x^2 - x + 1 

""" 

return self.ogf().denominator() 

def recurrence_repr(self): 

""" 

Return a string with the recurrence representation of 

the C-finite sequence. 

OUTPUT: 

- A string 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: C((2-x)/(1-x-x^2)).recurrence_repr() 

'Homogenous linear recurrence with constant coefficients of degree 2: a(n+2) = a(n+1) + a(n), starting a(0...) = [2, 1]' 

sage: C(x/(1-x)^3).recurrence_repr() 

'Homogenous linear recurrence with constant coefficients of degree 3: a(n+3) = 3*a(n+2) - 3*a(n+1) + a(n), starting a(1...) = [1, 3, 6]' 

sage: C(1).recurrence_repr() 

'Finite sequence [1], offset 0' 

sage: r = C((-2*x^3 + x^2 - x + 1)/(2*x^2 - 3*x + 1)) 

sage: r.recurrence_repr() 

'Homogenous linear recurrence with constant coefficients of degree 2: a(n+2) = 3*a(n+1) - 2*a(n), starting a(0...) = [1, 2, 5, 9]' 

sage: r = CFiniteSequence(x^3/(1-x-x^2)) 

sage: r.recurrence_repr() 

'Homogenous linear recurrence with constant coefficients of degree 2: a(n+2) = a(n+1) + a(n), starting a(3...) = [1, 1, 2, 3]' 

""" 

if self._deg == 0: 

return 'Finite sequence %s, offset %d' % (str(self._a), self._off) 

else: 

if self._c[0] == 1: 

cstr = 'a(n+%d) = a(n+%d)' % (self._deg, self._deg - 1) 

elif self._c[0] == -1: 

cstr = 'a(n+%d) = -a(n+%d)' % (self._deg, self._deg - 1) 

else: 

cstr = 'a(n+%d) = %s*a(n+%d)' % (self._deg, str(self._c[0]), self._deg - 1) 

for i in range(1, self._deg): 

j = self._deg - i - 1 

if self._c[i] < 0: 

if self._c[i] == -1: 

cstr = cstr + ' - a(n+%d)' % (j,) 

else: 

cstr = cstr + ' - %d*a(n+%d)' % (-(self._c[i]), j) 

elif self._c[i] > 0: 

if self._c[i] == 1: 

cstr = cstr + ' + a(n+%d)' % (j,) 

else: 

cstr = cstr + ' + %d*a(n+%d)' % (self._c[i], j) 

cstr = cstr.replace('+0', '') 

astr = ', starting a(%s...) = [' % str(self._off) 

maxwexp = self.numerator().quo_rem(self.denominator())[0].degree() + 1 

for i in range(maxwexp + self._deg): 

astr = astr + str(self[self._off + i]) + ', ' 

astr = astr[:-2] + ']' 

return 'Homogenous linear recurrence with constant coefficients of degree ' + str(self._deg) + ': ' + cstr + astr 

def series(self, n): 

""" 

Return the Laurent power series associated with the 

CFiniteSequence, with precision `n`. 

INPUT: 

- `n` -- a nonnegative integer 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: r = C.from_recurrence([-1,2],[0,1]) 

sage: s = r.series(4); s 

x + 2*x^2 + 3*x^3 + 4*x^4 + O(x^5) 

sage: type(s) 

<type 'sage.rings.laurent_series_ring_element.LaurentSeries'> 

""" 

R = LaurentSeriesRing(QQ, 'x', default_prec=n) 

return R(self.ogf()) 

def closed_form(self, n = 'n'): 

""" 

Return a symbolic expression in ``n``, which equals the n-th term of 

the sequence. 

It is a well-known property of C-finite sequences ``a_n`` that they 

have a closed form of the type: 

.. MATH:: 

a_n = \sum_{i=1}^d c_i(n) \cdot r_i^n, 

where ``r_i`` are the roots of the characteristic equation and 

``c_i(n)`` is a polynomial (whose degree equals the multiplicity of 

``r_i`` minus one). This is a natural generalization of Binet's 

formula for Fibonacci numbers. See, for instance, [KP, Theorem 4.1]. 

Note that if the o.g.f. has a polynomial part, that is, if the 

numerator degree is not strictly less than the denominator degree, 

then this closed form holds only when ``n`` exceeds the degree of that 

polynomial part. In that case, the returned expression will differ 

from the sequence for small ``n``. 

EXAMPLES:: 

sage: CFiniteSequence(1/(1-x)).closed_form() 

1 

sage: CFiniteSequence(x^2/(1-x)).closed_form() 

1 

sage: CFiniteSequence(1/(1-x^2)).closed_form() 

1/2*(-1)^n + 1/2 

sage: CFiniteSequence(1/(1+x^3)).closed_form() 

1/3*(-1)^n + 1/3*(1/2*I*sqrt(3) + 1/2)^n + 1/3*(-1/2*I*sqrt(3) + 1/2)^n 

sage: CFiniteSequence(1/(1-x)/(1-2*x)/(1-3*x)).closed_form() 

9/2*3^n - 4*2^n + 1/2 

Binet's formula for the Fibonacci numbers:: 

sage: CFiniteSequence(x/(1-x-x^2)).closed_form() 

sqrt(1/5)*(1/2*sqrt(5) + 1/2)^n - sqrt(1/5)*(-1/2*sqrt(5) + 1/2)^n 

sage: [_.subs(n=k).full_simplify() for k in range(6)] 

[0, 1, 1, 2, 3, 5] 

sage: CFiniteSequence((4*x+3)/(1-2*x-5*x^2)).closed_form() 

1/2*(sqrt(6) + 1)^n*(7*sqrt(1/6) + 3) - 1/2*(-sqrt(6) + 1)^n*(7*sqrt(1/6) - 3) 

Examples with multiple roots:: 

sage: CFiniteSequence(x*(x^2+4*x+1)/(1-x)^5).closed_form() 

1/4*n^4 + 1/2*n^3 + 1/4*n^2 

sage: CFiniteSequence((1+2*x-x^2)/(1-x)^4/(1+x)^2).closed_form() 

1/12*n^3 - 1/8*(-1)^n*(n + 1) + 3/4*n^2 + 43/24*n + 9/8 

sage: CFiniteSequence(1/(1-x)^3/(1-2*x)^4).closed_form() 

4/3*(n^3 - 3*n^2 + 20*n - 36)*2^n + 1/2*n^2 + 19/2*n + 49 

sage: CFiniteSequence((x/(1-x-x^2))^2).closed_form() 

1/5*(n - sqrt(1/5))*(1/2*sqrt(5) + 1/2)^n + 1/5*(n + sqrt(1/5))*(-1/2*sqrt(5) + 1/2)^n 

""" 

from sage.arith.all import binomial 

from sage.rings.qqbar import QQbar 

from sage.symbolic.ring import SR 

n = SR(n) 

expr = SR(0) 

R = FractionField(PolynomialRing(QQbar, self.parent().variable_name())) 

ogf = R(self.ogf()) 

__, parts = ogf.partial_fraction_decomposition(decompose_powers=False) 

for part in parts: 

denom = part.denominator().factor() 

denom_base, denom_exp = denom[0] 

# denominator is of the form (x+b)^{m+1} 

m = denom_exp - 1 

b = denom_base.constant_coefficient() 

# check that the partial fraction decomposition was indeed done correctly 

# (that is, there is only one factor, of degree 1, and monic) 

assert len(denom) == 1 and len(denom_base.list()) == 2 and denom_base.list()[1] == 1 and denom.unit() == 1 

r = SR((-1/b).radical_expression()) 

c = SR(0) 

for (k, a) in enumerate(part.numerator().list()): 

a = QQbar(a) 

c += binomial(n+m-k,m) * SR(((-1)**k*a*b**(k-m-1)).radical_expression()) 

expr += c.expand() * r**n 

return expr 

class CFiniteSequences_generic(CommutativeRing, UniqueRepresentation): 

r""" 

The class representing the ring of C-Finite Sequences 

TESTS:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: from sage.rings.cfinite_sequence import CFiniteSequences_generic 

sage: isinstance(C,CFiniteSequences_generic) 

True 

sage: type(C) 

<class 'sage.rings.cfinite_sequence.CFiniteSequences_generic_with_category'> 

sage: C 

The ring of C-Finite sequences in x over Rational Field 

""" 

Element = CFiniteSequence 

def __init__(self, polynomial_ring, category): 

r""" 

Create the ring of CFiniteSequences over ``base_ring`` 

INPUT: 

- ``base_ring`` -- the base ring for the o.g.f (either ``QQ`` or ``ZZ``) 

- ``names`` -- an iterable of variables (shuould contain only one variable) 

- ``category`` -- the category of the ring (default: ``Fields()``) 

TESTS:: 

sage: C.<y> = CFiniteSequences(QQ); C 

The ring of C-Finite sequences in y over Rational Field 

sage: C.<x> = CFiniteSequences(QQ); C 

The ring of C-Finite sequences in x over Rational Field 

sage: C.<x> = CFiniteSequences(ZZ); C 

The ring of C-Finite sequences in x over Integer Ring 

sage: C.<x,y> = CFiniteSequences(ZZ) 

Traceback (most recent call last): 

... 

NotImplementedError: Multidimensional o.g.f. not implemented. 

sage: C.<x> = CFiniteSequences(CC) 

Traceback (most recent call last): 

... 

ValueError: O.g.f. base not rational. 

""" 

base_ring = polynomial_ring.base_ring() 

self._polynomial_ring = polynomial_ring 

self._fraction_field = FractionField(self._polynomial_ring) 

CommutativeRing.__init__(self,base_ring, self._polynomial_ring.gens(), category) 

def _repr_(self): 

r""" 

Return the string representation of ``self`` 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: C 

The ring of C-Finite sequences in x over Rational Field 

""" 

return "The ring of C-Finite sequences in {} over {}".format(self.gen(), self.base_ring()) 

def _element_constructor_(self, ogf): 

r""" 

Construct a C-Finite Sequence 

INPUT: 

- ``ogf`` -- the ordinary generating function, a fraction of polynomials over the rationals 

TESTS:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: C((2-x)/(1-x-x^2)) 

C-finite sequence, generated by (-x + 2)/(-x^2 - x + 1) 

sage: C(x/(1-x)^3) 

C-finite sequence, generated by x/(-x^3 + 3*x^2 - 3*x + 1) 

sage: C(x^2-4*x^5) 

Finite sequence [1, 0, 0, -4], offset = 2 

sage: C(x^2+3/x) 

Finite sequence [3, 0, 0, 1], offset = -1 

sage: C(1/x + 4/x^3) 

Finite sequence [4, 0, 1], offset = -3 

sage: P = LaurentPolynomialRing(QQ.fraction_field(), 'X') 

sage: X = P.gen() 

sage: C(1/(1-X)) 

C-finite sequence, generated by 1/(-x + 1) 

sage: C = CFiniteSequences(QQ) 

sage: C(x) 

Finite sequence [1], offset = 1 

""" 

ogf = self.fraction_field()(ogf) 

return self.element_class(self, ogf) 

def ngens(self): 

r""" 

Return the number of generators of ``self`` 

EXAMPLES:: 

sage: from sage.rings.cfinite_sequence import CFiniteSequences 

sage: C.<x> = CFiniteSequences(QQ); 

sage: C.ngens() 

1 

""" 

return 1 

def gen(self,i=0): 

r""" 

Return the i-th generator of ``self``. 

INPUT: 

- ``i`` -- an integer (default:0) 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ); 

sage: C.gen() 

x 

sage: x == C.gen() 

True 

TESTS:: 

sage: C.gen(2) 

Traceback (most recent call last): 

... 

ValueError: The ring of C-Finite sequences in x over Rational Field has only one generator (i=0) 

""" 

if i!= 0: 

raise ValueError("{} has only one generator (i=0)".format(self)) 

return self.polynomial_ring().gen() 

def an_element(self): 

r""" 

Return an element of C-Finite Sequences. 

OUTPUT: 

The Lucas sequence. 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ); 

sage: C.an_element() 

C-finite sequence, generated by (-x + 2)/(-x^2 - x + 1) 

""" 

x = self.gen() 

return self((2-x)/(1-x-x**2)) 

def __contains__(self, x): 

""" 

Return True if x is an element of ``CFiniteSequences`` or 

canonically coerces to this ring. 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ); 

sage: x in C 

True 

sage: 1/x in C 

True 

sage: 5 in C 

True 

sage: pi in C 

False 

sage: Cy.<y> = CFiniteSequences(QQ); 

sage: y in C 

False 

sage: y in Cy 

True 

""" 

if x.parent() == self: 

return True 

try: 

self._coerce_(x) 

except TypeError: 

return False 

return True 

def fraction_field(self): 

r""" 

Return the fraction field used to represent the elements of ``self``. 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ); 

sage: C.fraction_field() 

Fraction Field of Univariate Polynomial Ring in x over Rational Field 

""" 

return self._fraction_field 

def polynomial_ring(self): 

r""" 

Return the polynomial ring used to represent the elements of ``self``. 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ); 

sage: C.polynomial_ring() 

Univariate Polynomial Ring in x over Rational Field 

""" 

return self._polynomial_ring 

def _coerce_map_from_(self, S): 

""" 

A coercion from `S` exists, if `S` coerces into ``self``'s fraction 

field. 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ); 

sage: C.has_coerce_map_from(C.fraction_field()) 

True 

sage: C.has_coerce_map_from(QQ) 

True 

sage: C.has_coerce_map_from(QQ[x]) 

True 

sage: C.has_coerce_map_from(ZZ) 

True 

""" 

if self.fraction_field().has_coerce_map_from(S): 

return True 

def from_recurrence(self, coefficients, values): 

""" 

Create a C-finite sequence given the coefficients $c$ and 

starting values $a$ of a homogenous linear recurrence. 

.. MATH:: 

a_{n+d} = c_0a_n + c_1a_{n+1} + \cdots + c_{d-1}a_{n+d-1}, \quad d\ge0. 

INPUT: 

- ``coefficients`` -- a list of rationals 

- ``values`` -- start values, a list of rationals 

OUTPUT: 

- A CFiniteSequence object 

EXAMPLES:: 

sage: C.<x> = CFiniteSequences(QQ) 

sage: C.from_recurrence([1,1],[0,1]) # Fibonacci numbers 

C-finite sequence, generated by x/(-x^2 - x + 1) 

sage: C.from_recurrence([-1,2],[0,1]) # natural numbers 

C-finite sequence, generated by x/(x^2 - 2*x + 1) 

sage: r = C.from_recurrence([-1],[1]) 

sage: s = C.from_recurrence([-1],[1,-1]) 

sage: r == s 

True 

sage: r = C(x^3/(1-x-x^2)) 

sage: s = C.from_recurrence([1,1],[0,0,0,1,1]) 

sage: r == s 

True 

sage: C.from_recurrence(1,1) 

Traceback (most recent call last): 

... 

ValueError: Wrong type for recurrence coefficient list. 

""" 

if not isinstance(coefficients, list): 

raise ValueError("Wrong type for recurrence coefficient list.") 

if not isinstance(values, list): 

raise ValueError("Wrong type for recurrence start value list.") 

deg = len(coefficients) 

co = coefficients[::-1] 

co.extend([0] * (len(values) - deg)) 

R = self.polynomial_ring() 

x = R.gen() 

den = -1 + sum([x ** (n + 1) * co[n] for n in range(deg)]) 

num = -values[0] + sum([x ** n * (-values[n]