OEIS/Generating Functions and Recurrences

From tehowiki
Revision as of 17:22, 25 May 2019 by imported>Gfis
Jump to navigation Jump to search

Linear Recurrences

Ordinary generating function o.g.f.

An o.g.f. consists of a fraction of two univariate polynomials (the numerator and the denominator polynomial). This rational function is expanded into a power series by polynomial division, beginning with the lowest power of the variable (we always use x), and the correspondant coefficients of the powers x are the terms of the integer sequences defined by the o.g.f.

Example for the polynomial division:

 1 / (1-x-x^2) = 1 + x + 2x^2 + 3x^3 + 5x^5 ... (Fibonnacci, A000045)
-1 + x + x^2    <---- denominator * (-1)
------------    <---- add
 0 + x +  x^2  
    -x +  x^2 + x^3
    --------------
     0 + 2x^2 + x^3
        -2x^2 + 2x^3 + 2x^4
       -------------------
          0   + 3x^2 + 2x^4
               -3x^3 + 3x^4 + 3x^5
               -------------------
                 0   + 5x^4 + 3x^5 ...

Corresponding Mathematica functions:

 CoefficientList[Series[numerator_polynomial / denominator_polynomial, {x, 0, number of terms}], x]
 CoefficientList[Series[-(x/(-1 + x + x^2)), {x, 0, 20}], x]
 FindGeneratingFunction[{terms}]

Linear Recurrence with constant coefficients

There is a one-to-one relation between an o.g.f. and a linear recurrence. The denominator's coefficients (without the constant 1) are the signature of the linear recurrence. The number of elements in the signature is the order of the linear recurrence.

In our example the linear recurrence is the following:

 a(0) = 1, a(1) = 1, a(n) = a(n - 1) + a(n - 2) for n > 1. 

The number of initial terms must be greater than or equal to the order.

Corresponding Mathematica functions:

 LinearRecurrence[{signature}, {initial terms}, number of terms]
 LinearRecurrence[{1, 1}, {1, 1}, 8] for the Fibonnacci sequence
 FindLinearRecurrence[{terms}]

Linear recurrences can be implemented very easily by additions and multiplications only. An array of length order may be used as a rolling buffer where all indexes are taken mod order.

Initial terms

OEIS sequences often start out with some initial terms before the recurrence really starts. For example:

A049221 Number of horizontally convex n-ominoes in which the top row has exactly 1 square, 
        which is not above the rightmost square in the second row.
Offset 1
1, 0, 1, 4, 14, 46, 148, 474, 1518, 4864, 15590, 49974, 160196
a(n) = 5a(n-1) - 7a(n-2) + 4a(n-3) for n >= 6
linear recurrence with constant coefficients, signature (5,-7,4)

The condition n>=6 and the order 3 show that the recurrence only starts at a(3) = 1, and that a(1) = 1 and a(2) = 0 were prefixed to the sequence. The second parameter of the Mathematica call

LinearRecurrence[{5, -7, 4}, {1, 0, 1, 4, 14}, 40]

is convenient because such initial terms may simply be written before the terms belonging to the recurrence.

How do we see such initial terms in the o.g.f.? The numerator and the denominator polynomials should have the same degree. If the numerator has a higher degree, we can transform the fraction into mixed form:

G.f.: x*(1-x)^2*(1-3x+x^2)/(1-5x+7x^2-4x^3) 
  = (x-5x^2+6x^3-5x^4+x^5)/(1-5x+7x^2-4x^3)
  = 0  + 1x + 0x^2 + x^3(1-x+x^2)/(1-5x+7x^2-4x^3)
    a(0) a(1) a(2)   a(3 ...)

Now we can directly recognize the additional initial terms.

Offset

The power series of generated by an o.g.f. usually has the general form

c0*x^0 + c1*x^1 + c2*x^2 + c3*x^3 ...

The integer sequence defined by this o.g.f. then simply is

c0, c1, c2, c3, ...

All OEIS sequences are assigned an offset that relates them to the context for which they are defined. The most natural offset for a sequence defined by an o.g.f. is 0, that is the sequence starts with a(0) = c0.



Usually the expansion of an o.g.f. starts

Computation of the o.g.f. from the linear recurrence

Example (from youtube.com):

a(n) = 2*a(n-1) + 4*a(n-2); a(0) = 1, a(1) = 3; [1, 3, 10, 32, ...]
a(n) - 2*a(n-1) - 4*a(n-2) = 0
 A       = 1 + 3x + 10x^2 + 32x^3 + ...
-(2x*A   =     2x +  6x^2 + 20x^3 + ...)
-(4x^2*A =           4x^2 + 12x^3 + ...)
----------------------------------------
(1 - 2x - 4x^2)*A = 1+x
A = (1 + x) / (1 - 2x - 4x^2) 

Example (from youtube.com):

a(n) = a(n-1) + 2*a(n-2) + 3; a(0) = 2, a(1) = 2; [2, 2, 9, 16, 37, 72, 149, 296, ...]
a(n) - a(n-1) - 2*a(n-2) = 3
 A       = 2 + 2x + 9x^2 + 16x^3 + 37x^4 + ...
-(x*A    =     2x + 2x^2 +  9x^3 + 16x^4 + ...)
-(2x^2*A =          4x^2 +  4x^3 + 18x^4 + ...)
----------------------------------------
(1 - x - 2x^2)*A = 2 + 3x^2 +3x^3 + 3x^4 + ... = 2 + 3x^2/(1-x) = (2 - 2x + 3x^2)/(1-x)
A = (2 - 2x + 3x^2) / ((1 - x - 2x^2)*(1 - x))

Example:

A173652: 0, 0, 0,  1,   4,   8,  33,  68, 245, 549, 1815 ...
                                            x^3  x^4  x^5  x^6  x^7  x^8  x^9  x^10
LinearRecurrence[{0,    9,     3,    -17,    -8,     6,     7,    -1,    1}, {0, 0, 0, 1, 4, 8, 33, 68, 245}, 31]
denominator: (1 + 0x^1 - 9x^2 - 3x^3 + 17x^4 + 8x^5 - 6x^6 - 7x^7 + x^8 - x^9)
numerator: x^3

Coxeter Groups

The main entry is:

A169452 
Number of reduced words of length n in Coxeter group on 7 generators S_i 
with relations (S_i)^2 = (S_i S_j)^33 = I.
1, 7, 42, 252, 1512, 9072, 54432, 326592, 1959552
signature(5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,-15).
G.f.:(t^33+2*t^32+2*t^31+2*t^30
+2*t^29+2*t^28+2*t^27+2*t^26+2*t^25+2*t^24+2*t^23+2*t^22+2*t^21+2*t^20
+2*t^19+2*t^18+2*t^17+2*t^16+2*t^15+2*t^14+2*t^13+2*t^12+2*t^11+2*t^10
+2*t^9+2*t^8+2*t^7+2*t^6+2*t^5+2*t^4+2*t^3+2*t^2+2*t+1) / 
(15*t^33-5*t^32-5*t^31-5*t^30
-5*t^29-5*t^28-5*t^27-5*t^26-5*t^25-5*t^24-5*t^23-5*t^22-5*t^21-5*t^20
-5*t^19-5*t^18-5*t^17-5*t^16-5*t^15-5*t^14-5*t^13-5*t^12-5*t^11-5*t^10
-5*t^9-5*t^8-5*t^7-5*t^6-5*t^5-5*t^4-5*t^3-5*t^2-5*t+1)
  • Harvey P. Dale gave this Mathematica program for the g.f. on Aug 16 2014:
coxG[{pwr_, c1_, c2_, trms_:20}] := Module[{
num=Total[  2t^Range[pwr-1]]+   t^pwr+1, 
den=Total[c2*t^Range[pwr-1]]+c1*t^pwr+1}, 
CoefficientList[Series[num/den, {t, 0, trms}], t]]; 
coxG[{33, 15, -5, 30}] (* for A169452 *)
  • with the parameters:
pwr  = largest exponent in the g.f. 
    and of "(S_i S_j)" in the name,
c1   = first coefficient in the denominator of the g.f., 
     = triangular(-c2) = binomial(-c2 + 1, 2) 
c2   = - second coefficient in the denominator of the g.f., 
     = 2 - (numbers of generators in the name);
trms = number of terms desired (with a default number of 20)