login
Coefficients of the eigenfunction of a sequence transformation.
4

%I #12 Jan 18 2022 20:06:19

%S 1,3,6,45,126,750,2796,19389,75894,449562,2027796,12211794,57895596,

%T 332787324,1677545304,9766642077,50378641830,286825948194,

%U 1529968671492,8729259097158,47374697101572,269062276076868,1484430536591592

%N Coefficients of the eigenfunction of a sequence transformation.

%C G.f. A(x) satisfies A(x^2) = (A(x/2)-1)/x - A(x/2)^2/2.

%C B(x) := 1/(2*x) - x*A(x^2) satisfies B(x)^2 + 1 = B(2*x^2).

%C Define f(n, c) := x - Sum_{k>=0} a(k)/(2*x)^(2*k+1) where x = c^(2^n). Then A003095(n+1) = A004019(n) + 1 = f(n, 1.502836801...). Also, A062013(n) = f(n, 1.78050350...). - _Michael Somos_, Jun 07 2021

%e G.f. = A(x) = 1 + 3*x + 6*x^2 + 45*x^3 + 126*x^4 + 750*x^5 + 2796*x^6 + ...

%e B(x) = 1/(2*x) - x - 3*x^3 - 6*x^5 - 45*x^7 - 126*x^9 - 750*x^11 - ... - _Michael Somos_, Jul 11 2019

%t a[ n_] := If[n < 0, 0, Module[{A = 1 + O[x], m = 2}, While[m < n + 2, m *= 2; A = (Normal[ 1/x - Sqrt[ 1/x^2 - 2/x - 2*(Normal[A] /. x -> x^2) + O[x]^(m - 2)]] /. x -> 2*x) + O[x]^(m - 1) //PowerExpand]; SeriesCoefficient[A, n]]]; (* _Michael Somos_, Jun 07 2021 *)

%o (PARI) {a(n) = my(A, m); if( n < 0, 0, m=2; A = 1 + O(x); while( m < n+2, m*=2; A = subst(1/x - sqrt(2*(subst((1/2)/x - A, x, x^2) - 1/x)), x, 2*x)); polcoeff(A, n))};

%Y Cf. A003095, A004019, A062013.

%K nonn

%O 0,2

%A _Michael Somos_, Oct 04 2003