OFFSET
0,2
LINKS
G. C. Greubel, Table of n, a(n) for n = 0..1000
FORMULA
G.f. A(x) = (1+G003169(x))/(1-x-x*G003169(x)), where G003169(x) is the g.f. of A003169.
a(n) ~ (sqrt(3056686 + 12607266/sqrt(17)) * ((71 + 17*sqrt(17))/16)^n) / (10201 * sqrt(Pi) * n^(3/2)). - Vaclav Kotesovec, Jan 31 2023
MATHEMATICA
f[n_]:= f[n]= If[n<2, 1, If[n==2, 3, ((324*n^2 -708*n +360)*f[n-1] -(371*n^2 -1831*n +2250)*f[n-2] + (20*n^2 -130*n +210)*f[n-3])/(16*n*(2*n-1)) ]]; (* f = A003169 *)
A[n_, k_]:= A[n, k]= If[n==0, f[k], If[k==0, 1, Sum[f[k-j]*A[n-1, j], {j, 0, k}]]]; (* A = 100324 *)
a[n_]:= a[n]= Sum[A[n-k, k], {k, 0, n}]; (* a = A100325 *)
Table[a[n], {n, 0, 40}] (* G. C. Greubel, Jan 31 2023 *)
PROG
(PARI) {a(n)=local(A=1+x+x*O(x^n)); if(n==0, 1, for(i=1, n, A=1+x*A/(2-A)^2); sum(k=0, n, polcoeff(A^(n-k+1), k)))}
(SageMath)
@CachedFunction
def f(n): # f = A003169
if (n<2): return 1
elif (n==2): return 3
else: return ((324*n^2-708*n+360)*f(n-1) - (371*n^2-1831*n+2250)*f(n-2) + (20*n^2-130*n+210)*f(n-3))/(16*n*(2*n-1))
@CachedFunction
def A(n, k): # A = 100324
if (n==0): return f(k)
elif (k==0): return 1
else: return sum( f(k-j)*A(n-1, j) for j in range(k+1) )
def T(n, k): return A(n-k, k)
def A100325(n): return sum( T(n, k) for k in range(n+1) )
[A100325(n) for n in range(41)] # G. C. Greubel, Jan 31 2023
CROSSREFS
KEYWORD
nonn
AUTHOR
Paul D. Hanna, Nov 17 2004
STATUS
approved