OEIS/Holonomic Project

This project tries to define as many existing OEIS sequences as possible by holonomic recurrences.

Whereas the usual linear recurrence has constant coefficients, a holonomic recurrence has coefficients that are polynomials in n (the sequence index). Holonomic recurrences can be derived from linear differential equations, and are applicable to a wide class of series - see Wikipedia: Holonomic function.

Example

A000957 Fine's sequence (or Fine numbers): 
  number of relations of valence >= 1 on an n-set; 
  also number of ordered rooted trees with n edges having root of even degree.
0, 1, 0, 1, 2, 6, 18, 57, 186, 622, 2120, 7338, 25724, 91144, 325878, ...
G.f.: (1-sqrt(1-4*x))/(3-sqrt(1-4*x))
D-finite with recurrence: 2*n*a(n) +(12-7*n)*a(n-1) +2*(3-2*n)*a(n-2)=0. - R. J. Mathar

Types of recurrences

There is a hierarchy of sequence types which are included in this project, with increasing complexity and computational effort:

• constant
• periodic
  • finite (followed by infinitely many zeroes), polynomials
  • continued fractions of sqrt(n)
• recurrent with constant coefficients (cf. OEIS Index entries for linear recurrences with constant coefficients)
  • defined by a rational generating function (fraction of two polynomials in x), among them sequences for Coexter groups and (lattice) coordination sequences
• holonomic recurrent (D-finite, hypergeometric functions, algebraic equations)

Initial terms

For all types the recurrence frequently starts only after some initial terms which do not participate in the recurrence equation.

Exponential generating function

The project also handles e.g.f.s, in which case the terms are multiplied by n! during the generation.of the sequence.

Embedding in jOEIS

jOEIS is a project which aims to implement all OEIS sequences in Java. Currently (Feb. 2021) there are more than 104,000 sequences implemented. They are compiled in a single jar file, and any sequence can be called repeatedly in a uniform way to produce the next terms.

In jOEIS there is a series of basic Java classes which generate sequences of the holonomic type. They can all be found under https://github.com/archmageirvine/joeis/tree/master/src/irvine/oeis:

BlockMultAddSequence
ContinuedFractionOfSqrtSequence
CoordinationSequence
CoxeterSequence
FiniteSequence
GeneratingFunctionSequence
HolonomicRecurrence
LatticeCoordinationSequence
LinearRecurrence
PaddingSequence
PeriodicSequence
PrependSequence

Currently the following numbers of sequences are covered by this project:

    36 block
  2911 conti
   946 coord
  2251 coxet
  4601 finit
  8520 gener
  5902 holon
   340 latti
 20035 linea
   336 paddi
  1202 perio
   700 prepe
------------
 47780 total

Program for holonomic recurrences

Recurrences are attractive because it is rather simple to implement them in any programming language. A CAS may have special functions for them (e.g. LinearRecurrence and RecurrenceTable in Mathematica). The algorithm is straightforward, at it needs only fixed space even in the holonomic case. In jOEIS, HolonomicRecurrence.java uses a ring buffer of size k for the recurrence variables a(n), a(n-1), ... a(n-k+1), and a second array where the powers of the index n are maintained. The length of this array corresponds with the degree of the equation (the highest power of n). The program is not too big, and could be rewritten in other languages, for example Python.

Parameterization

The program uses these parameters:

• recurrence equation
• initial terms
• index of the first term to be returned (OEIS offset1)

Equation

Given that the recurrence equation is written in the form of an annihilator:

p[0]*1 + p[1]*a[n-k+1] + p[2]*a[n-k+2] + ... + p[k-1]*a[n-k+k-1] + p[k]*a[n] = 0

then the corresponding parameter is represented by the array

[p[0], p[1], p[2], ... p[k]]

where each polynomial p[i] is represented by an array [c0, c1, c2 ... cm] that defines the exponents of n. The initial example equation

2*n*a(n) +(12-7*n)*a(n-1) +2*(3-2*n)*a(n-2)=0

is represented by:

[[0],[6,-4],[12,-7],[0,2]]

Constant term

The first polynomial allows for a constant term, but normally the homogenous form is used with p[0] = [0].

Index shifting

Usually the last recurrence variable is a(n), but it is also possible to shift the recurrence index, for example to have

 2*a(n) - 2*a(n+1) - a(n+2) = 0

Such a shifting is specified by an additional parameter dist which is 0 by default. The example would get dist=2. In the future, all equations will be recomputed such that the constant term as well as the shift become always zero.

Initial terms

For an order (number of variables - 1) k a recurrence normally needs k initial terms to be computable. With the current implementation it is possible

• that there are more initial terms, in which case the recurrence only applies to the last k terms,
• that there are fewer initial terms, in which case a suitable number of zeroes is implied before the first term.

Black startable sequences

As a special possibility a recurrence can be started "out of nothing", that is the initial terms can be omitted altogether. The implementation then implies a leading sequence term one, and a suitable number of zeroes before. We call this property black startable (an expression used for power stations), and there is a surprising number of sequences that have it.

Database table

The main result of the project is a database table with suitable parameter sets for the OEIS sequences that are covered. The table is create with the following MySQL statements:

--  Table for holonomic recurrences in the OEIS
DROP    TABLE  IF EXISTS holref;
CREATE  TABLE            holref
   ( aseqno   VARCHAR(10) NOT NULL  -- A-number
   , callcode VARCHAR(32)           -- "holos", generate calls for HolonomicRecurrence
   , offset   INT                   -- offset
   , matrix   VARCHAR(8192)         -- coefficients of polynomials in n
   , init     VARCHAR(4096)         -- initial terms for the recurrence
   , shift    INT                   -- optional shift (0)
   , gftype   INT                   -- 0 for o.g.f, 1 for e.g.f
   , keyword  VARCHAR(128)          -- list of keywords for source sequence type etc.
   , PRIMARY KEY(aseqno, callcode)
   );
COMMIT;

The tab-separated data loaded into the database table can be found under holref.txt.

Parameters of hypergeometric series

There is a number of OEIS sequences that can be generated by a recurrence

f(n) * a(n) = g(n) * a(n-1)

For example, formulae like

A217364 a(n) = 2^n*binomial(5*n, n)/(4*n+1).
A274136 a(n) = (n+1)*(2*n+2)!/(n+2)
A275655 a(n) = binomial(6*n,3*n)*binomial(2*n,n).

are of this type. In these cases, the parameters for HolonomicRecurrence can be computed by the following Maple procedure:

with(StringTools):
holo := proc(as,f); local b, num, den, fd, str, filename;
  b := simplify(simplify(-f(n)/f(n-1), GAMMA));
  printf("***** b = %a\n", b);
  num := seq(coeff(numer(b),n,j),j=0..degree(numer(b)));
  den := seq(coeff(denom(b),n,j),j=0..degree(denom(b)));
  str := cat("[[0],", DeleteSpace(convert([num],string)), ",", DeleteSpace(convert([den],string)), "]");
  printf("%a\tholos\t0\t\%s\t[%a]\t0\t\n", as, str, f(0));
end;

We obtain, for example:

holo(A217364, proc(n) 2^n*binomial(5*n, n)/(4*n+1) end);
***** b = (-3125*n^4+6250*n^3-4375*n^2+1250*n-120)/(128*n^4-64*n^3-8*n^2+4*n)
A217364 holos   0       [[0],[-120,1250,-4375,6250,-3125],[0,4,-8,-64,128]]     [1]    0

Now we follow the Maple file above by an empty line and a pattern specification line:

holo($(PARM0),proc(n) $(PARM3) end proc);

Then it forms a pattern file for the Perl program callmaple.pl which merges the pattern file with a list of recurrence expressions in "seq4" format, and the output is in the same format with all recurrence expressions replaced by the parameters suitable for HolonomicRecurrence.