7
$\begingroup$

In the course of my research years ago, I came across three integer sequences related to Ramanujan's pi formulas but for level $10$. Their recurrence relations may be important. (Just like the level $6$ version was important for $\zeta(3)$.) Given the binomial numbers,

$$\text{Binomial}[n,m]=\binom n m$$

the first two sequences are,

\begin{align} A_k = \sum_{j=0}^k \binom {k+j} {k-j} \sum_{m=0}^j \binom {j}{m}^4 &= 1,3,25,267,3249,42795,\dots\\ B_k = \sum_{j=0}^k (-4)^{k-j}\binom {k+j} {k-j} \sum_{m=0}^j \binom {j}{m}^4 &= 1, -2, 10, -68, 514, -4100,\dots \end{align}

which start with $k=0$, hence $A_0 = B_0 = 1.$ The third one is slightly more complicated,

\begin{align}C_n &= \sum_{j=0}^n (-5)^{n-j}\binom {n+j} {n-j} \binom{2j}j \sum_{k=0}^j \binom j k \binom{2k}{k}^{-1} \sum_{m=0}^k \binom{k}{m}^4\\ &= 1, -1, 1, -1, 1, 23, -263, 1343, -2303,\dots \end{align}

which start with $n=0$, hence $C_0 = 1.$

The recurrence of the first one is already known and discussed in this ancient MO post. Since this involves level $10$, perhaps it is not surprising it is a $5$-term recurrence of form,

$$P_1 A_{k}+P_2 A_{k+1}+P_3 A_{k+2}+P_4 A_{k+3}+P_5 A_{k+4} =0$$

where the $P_i = P_i(k)$ are deg-$4$ polynomials in $k$.


Questions

  1. May I know the recurrence relation for $B_k$ and $C_n$, or at least a simple code that will run in an old Mathematica?
  2. Given any infinite integer sequence, is there a command in Wolfram Alpha to find the recurrence relation?
$\endgroup$
12
  • $\begingroup$ Is this for the Mathematica? I mean the software Mathematica. Or is it perhaps for the Maths.S.E? $\endgroup$
    – bmf
    Commented May 11, 2023 at 9:47
  • $\begingroup$ @bmf It is really for this group. I have a very old version of Mathematica and I need to know the command, if any. Or if WolframAlpha can do it so if i have future integer sequences, then i’ll know what to do. $\endgroup$ Commented May 11, 2023 at 9:54
  • $\begingroup$ Perhaps you should stress that you want code for this, because due to the absence of code and the explicit request you can understand the confusion $\endgroup$
    – bmf
    Commented May 11, 2023 at 10:01
  • 2
    $\begingroup$ @TitoPiezasIII Maybe, for finding the recurrence formula formally without any guesses from the first few values, the function you might be looking for is DifferenceRootReduce, for example DifferenceRootReduce[(n + 1)/(n - 1), n] works in wolframalpha. $\endgroup$ Commented May 11, 2023 at 10:27
  • 1
    $\begingroup$ For finding the formula as a guess from the first few values, you might be looking for FindSequenceFunction for example, FindSequenceFunction[{1, 1, 2, 3, 5, 8, 13}, n] works in wolfram alpha $\endgroup$ Commented May 11, 2023 at 10:27

1 Answer 1

10
+50
$\begingroup$

The question is asking to find linear recurrences with polynomial coefficients for number sequences. As far as I know, Mathematica still does not have this functionality builtin. That is the context of my question 265844 "Finding polynomial relations between power series and sequences".

For this important special case, I have written a simple function that can find linear recurrences with polynomial coefficients. No doubt the code could be improved and made more robust.

findrecur::usage = "findrecur[nterm, ncoef, a, an] tries to find a linear recurrence"<>
" with number of terms (nterm) and polynomial coeffcients up to degree (ncoef-1)"<>
" given a list of numbers (a) whose name is (an).";

findrecur[nterm_Integer, ncoef_Integer, a_List, an_:"a"] := Module[
{nt=nterm, nc=ncoef, nv=Length[a], c, equList, soln, polList, test, gcd},

Print["Sequence ", an, ": number of terms, coefficients, values = ", {nt, nc, nv}];

(* Create list of equations to solve using nv sequence values list a[[]] *)
equList = Table[ 0 == Sum[a[[n+i]] * n^e* c[e,i], {e,0,nc-1}, {i,0,nt-1}], {n,nv-nt}];

(* Solve the linear equations for the polynomial coefficients *)
soln = Solve[equList,Flatten[Table[c[e,i], {e,0,nc-1}, {i,0,nt-1}],1]][[1]] //Quiet;

(* Create the list of polynomials with common factors removed *)
polList = Table[ Sum[ n^e*c[e,i], {e,0,nc-1}], {i,0,nt-1}] /.soln /.{c[0,0]->1};
polList = polList /.{c[_,_]->0} //Factor;

(* Verify the recurrence using the sequence values to be sure *)
If[And@@Table[ 0==Sum[ polList[[i+1]]*a[[n+i]], {i,0,nt-1}], {n,nv-nt}],
    gcd = PolynomialGCD@@polList; If[gcd=!=0, polList/=gcd;
 
(* Return the linear recurrence to be equated to zero *)
    Sum[ polList[[i+1]]*Subscript[an, n+i], {i, 0, nt-1}], "Fail"]]];

This function can be used for the given sequence $A_k$ as follows:

(* Define the sequence function *)
Ak[k_Integer] := Sum[ Binomial[k+j,k-j]*Sum[Binomial[j,m]^4, {m,0,j}], {j,0,k}];

(* Output the first few values as a check *)
Print["A(0..)", " = ",Table[ Ak[k], {k,0,5}]];

(* Call the function to find a polynomial recurrence *)
0 == findrecur[5, 5, Table[Ak[k], {k,1,30}], "A"] //TeXForm

(* A(0..) = {3,25,267,3249,42795} *)
(* Sequence A: number of terms, coefficients, values = {5,5,30} *)

$$ 0=(n+3) (n+1)^3 \text{A}_n\\ +\left(-20 n^4-152 n^3-420 n^2-508 n-229\right) \text{A}_{n+1}\\+\left(38 n^4+380 n^3+1416 n^2+2330 n+1431\right) \text{A}_{n+2}\\+\left(-20 n^4-248 n^3-1140 n^2-2292 n-1689\right) \text{A}_{n+3}\\ +(n+2)(n+4)^3 \text{A}_{n+4} $$

This returns a $5$ term recurrence with degree $4$ coefficients. An alternative recurrence for this sequence is returned using:

(* Call the function to find a polynomial recurrence *)
0 == findrecur[6, 4, Table[Ak[k], {k,1,29}], "A"] //TeXForm

(* Sequence A: number of terms, coefficients, values = {6,4,29} *)

$$ 0=(n+1)^3 \text{A}_n\\-3 \left(7 n^3+33 n^2+53 n+29\right) \text{A}_{n+1}\\+2 \left(29 n^3+219 n^2+549 n+460\right) \text{A}_{n+2}\\-2 \left(29 n^3+303 n^2+1053 n+1214\right) \text{A}_{n+3}\\+3 \left(7 n^3+93 n^2+413 n+613\right) \text{A}_{n+4}\\-(n+5)^3 \text{A}_{n+5} $$

N.B. The sequence values used begin with $A_1$ due to the $0^0$ feature of Mathematica.

N.B. I edited the TeXForm output a bit to improve readability.

You can use similar code for the $B_k$ and $C_k$ sequences by changing the sequence definition. The results for both of them:

(* Define the sequence function *)
Bk[k_Integer] := Sum[(-4)^(k-j)*Binomial[k+j,k-j]*
   Sum[Binomial[j, m]^4, {m,0,j}],{j,0,k}];

(* Output the first few values as a check *)
Print["B(0..)", " = ",Table[ Bk[k], {k,0,5}]];

(* Call the function to find a polynomial recurrence *)
0 == findrecur[5, 5, Table[Bk[k], {k,1,29}], "B"] //TeXForm

(* B(0..) = {1,-2,10,-68,514,-4100} *)
(* Sequence B: number of terms, coefficients, values = {5,5,29} *)

$$ 0=256 (n+3) (n+1)^3 \text{B}_n\\ +32 \left(10 n^4+79 n^3+230 n^2+296 n+143\right) \text{B}_{n+1}\\ +4\left(32 n^4+320 n^3+1199 n^2+1995 n+1244\right) \text{B}_{n+2}\\+2 \left(10 n^4+121 n^3+545 n^2+1079 n+788\right) \text{B}_{n+3}\\+(n+2) (n+4)^3 \text{B}_{n+4} $$

(* Call the function to find a polynomial recurrence *)
0 == findrecur[6, 4, Table[Bk[k], {k,1,29}], "B"] //TeXForm

(* Sequence B: number of terms, coefficients, values = {6,4,29} *)

$$ 0=1024(n+1)^3 \text{B}_n\\ +128 \left(12 n^3+63 n^2+113 n+69\right) \text{B}_{n+1}\\+16 \left(52 n^3+402 n^2+1047 n+920\right) \text{B}_{n+2}\\+4 \left(52 n^3+534 n^2+1839 n+2122\right) \text{B}_{n+3}\\+2 \left(12 n^3+153 n^2+653 n+933\right) \text{B}_{n+4}\\+(n+5)^3 \text{B}_{n+5} $$

(* Define the sequence function *)
Ck[n_Integer] := Sum[ (-5)^(n-j)*Binomial[n+j,n-j]*Binomial[2*j,j]*
  Sum[ Binomial[j,k]/Binomial[2*k,k]*
     Sum[ Binomial[k, m]^4, {m, 0, k}],{k, 0, j}],{j, 0, n}];

(* Output the first few values as a check *)
Print["C(1..)", " = ",Table[ Ck[k], {k,0,5}]];

(* Call the function to find a polynomial recurrence *)
0 == findrecur[7, 4, Table[Ck[k], {k,1,34}], "C"] //TeXForm

(* C(1..) = {1,-1,1,-1,1,23} *)
(* Sequence C: number of terms, coefficients, values = {7,4,34} *)

$$ 0=15625 (n+1)^3 \text{C}_n\\ +625\left(18 n^3+99 n^2+185 n+117\right) \text{C}_{n+1}\\ +25\left(183 n^3+1464 n^2+3966 n+3638\right) \text{C}_{n+2}\\ +2 (2 n+7) \left(279 n^2+1953 n+3523\right) \text{C}_{n+3}\\ +\left(183 n^3+2379 n^2+10371 n+15157\right) \text{C}_{n+4}\\+\left(18 n^3+279 n^2+1445 n+2501\right) \text{C}_{n+5}\\+(n+6)^3 \text{C}_{n+6} $$


I think it is worth improving the previous versions of my code. Here is the latest iteration:

findseqrecur::usage =
"findseqrecur[nt, nc, a, ao, an, n, no] finds linear recurrences with"<>
" number of terms (nt) and polynomial coeffcients up to degree (nc-1)"<>
" given a list of numbers (a ={an[ao], an[ao+1], ...) indexed by (n)"<>
" of the form sum_{i..nt} p_i(n)an_{n+no+i}.";

findseqrecur[nt_Integer, nc_Integer, a_List, ao_:1, an_:"a", n_:n, no_:0] :=
Module[{ nv=Length[a], io=1, c, equs, sols, vars, polys, recurs},

Print["Sequence ", an, ": number of terms, coefficients, values = ", {nt, nc, nv}];

(* Create list of equations to solve using nv sequence values list a[[]] *)
equs = Table[ 0 == Sum[a[[k+i]]*(k+ao)^e*c[e,i], {e,0,nc-1}, {i,1,nt}], {k,0,nv-nt}];

(* Solve the linear equations for the polynomial coefficients *)
sols = Solve[equs, Flatten[Table[c[e, i], {e,0,nc-1}, {i,1,nt}],1]][[1]] //Quiet;

(* Find a non-zero coefficients in solution *)
vars = Sort[Variables[sols[[All, 2]]]];

polys = Table[ Table[ Sum[ (n+ao)^e * c[e, i], {e, 0, nc-1}], {i, 1, nt}]
      /. sols /. {cvar -> 1} /. {c[_, _] -> 0}, {cvar, vars}];
      
recurs = Table[ Which[
    0 != Normal@Sum[ a[[n+i]] * poly[[i]], {i, 1, nt}, {n, 0, nv-nt}], Nothing,
    0 === (gcd = PolynomialGCD@@poly), Nothing,
    True, Sum[ Factor[poly[[i]]/gcd/.{n -> n+no+1-ao}] *Subscript[an, n+no+i],
       {i, 1, nt}] ], {poly, polys}];
        
DeleteDuplicates[ Flatten@recurs, (#1 == #2 || Simplify[#1 == -#2])&]];

This version returns a list of recursions instead of only a single one. For example:

(* Define the sequence function *)
Ak[k_Integer] := Sum[ Binomial[k+j,k-j]*Sum[Binomial[j,m]^4, {m,0,j}], {j,0,k}];

(* Output the first few values as a check *)
Print["A(0..) = ",Table[ Ak[k], {k,0,5}]];

(* Call the function to find a polynomial recurrence *)
result = (0 == #&) /@ findseqrecur[6, 4, Array[Ak, 28, 1], 1, "A", k, -4];

(* Verify the result recursion equations *)
And @@ Flatten @ Table[ Simplify[result/.{Subscript["A",i_]->Ak[i]}], {k,4}]

(* Output the result equations in TeX form *)
result //TeXForm

(* A(0..) = {1,3,25,267,3249,42795} *)
(* True *)

(* Sequence A: number of terms, coefficients, values = {6,4,28} *)

$$ \left\{0=+(k-2)^3 \text{A}_{k-3} \\ -3 \left(7 k^3-30 k^2+44 k-22\right) \text{A}_{k-2} \\ +2 \left(29 k^3-42 k^2+18 k+1\right) \text{A}_{k-1} \\ -2 \left(29 k^3+42 k^2+18 k-1\right) \text{A}_k \\ +3 \left(7 k^3+30 k^2+44 k+22\right) \text{A}_{k+1} \\ -(k+2)^3 \text{A}_{k+2}\right\} $$

Notice the skew-symmetry. The code could be further improved I think.

N.B. If the number of sequence values In parameter a is not sufficient to determine the coefficients uniquely, then the function may return large coefficients. You should call the function again with more sequence values and compare the returned coefficients to confirm that they are valid. Also, verify the results as in the previous example. As a very simple exercise try the code

Table[ findseqrecur[2, 2, Array[Factorial, nv]], {nv, 4}]

which shows the effect of using bigger values of nv. This also shows a slight glitch in my code which I may be able to fix later.

$\endgroup$
7
  • $\begingroup$ Thanks! I finally found the closed-form of the last three level 10 sequences. (I actually found them late 2019, then Covid-19 interrupted.) $\endgroup$ Commented May 12, 2023 at 3:24
  • $\begingroup$ Beautiful, it works! (However, I had to use the labels "ord" and "deg" because of the commands "Order" and "Degree". And my version didn't have the command "Quiet", but I figured the code will work without it.) $\endgroup$ Commented May 12, 2023 at 4:13
  • $\begingroup$ You may be interested in this new MO post. It uses these three sequences (plus a fourth one) to tie together four eta quotients. $\endgroup$ Commented May 12, 2023 at 9:47
  • $\begingroup$ @TitoPiezasIII About "order" and "degree". Yes, the symbols "Order" and "Degree" already have built-in meaning in Mathematica, but the Wolfram Language is case-sensitive so using "order" and "degree" should not cause any problem. I will look at you MO post. $\endgroup$
    – Somos
    Commented May 12, 2023 at 10:19
  • $\begingroup$ @TitoPiezasIII Thanks for that comment! I have more comments in the code and I slightly modified the code. In the case you mention, I suggest to try order = 3; degree = 4; values = 11; as the minimum. No doubt, the code should be more robust. $\endgroup$
    – Somos
    Commented May 16, 2023 at 14:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.