|
%I #13 Nov 02 2019 03:09:40
%S 1,1,4,38,532,9329,190312,4340296,108043128,2890318936,82209697588,
%T 2467155342740,77676395612884,2554497746708964,87449858261161216,
%U 3107829518797739032,114399270654847628768,4353537522757357068296
%N G.f. satisfies: A(x/A(x)^4) = 1 + x; thus A(x) = 1 + series_reversion(x/A(x)^4).
%H Robert Israel, <a href="/A120974/b120974.txt">Table of n, a(n) for n = 0..250</a>
%F G.f. satisfies: A(x) = 1 + x*B(x)^4 = 1 + (1 + x*C(x)^4 )^4 where B(x) and C(x) satisfy: C(x) = B(x)*B(A(x)-1), B(x) = A(A(x)-1), B(A(x)-1) = A(B(x)-1), B(x/A(x)^4) = A(x), B(x) = A(x*B(x)^4) and B(x) is g.f. of A120975.
%p A:= x -> 1:
%p for m from 1 to 30 do
%p Ap:= unapply(A(x)+c*x^m,x);
%p S:= series(Ap(x/Ap(x)^4)-1-x, x, m+1);
%p cs:= solve(convert(S,polynom),c);
%p A:= subs(c=cs, eval(Ap));
%p od:
%p seq(coeff(A(x),x,m),m=0..30);# _Robert Israel_, Oct 25 2019
%t nmax = 17; sol = {a[0] -> 1};
%t Do[A[x_] = Sum[a[k] x^k, {k, 0, n}] /. sol; eq = CoefficientList[ A[x/A[x]^4] - 1 - x + O[x]^(n + 1), x] == 0 /. sol; sol = sol ~Join~ Solve[eq][[1]], {n, 1, nmax}];
%t sol /. Rule -> Set;
%t a /@ Range[0, nmax] (* _Jean-François Alcover_, Nov 02 2019 *)
%o (PARI) {a(n)=local(A=[1,1]);for(i=2,n,A=concat(A,0); A[ #A]=-Vec(subst(Ser(A),x,x/Ser(A)^4))[ #A]);A[n+1]}
%Y Cf. A120975; variants: A120970, A120972, A120976, A030266, A067145, A107096.
%K nonn
%O 0,3
%A _Paul D. Hanna_, Jul 20 2006
|
