%I #26 Mar 03 2024 08:14:02
%S 1,4,28,224,1948,17928,171776,1695872,17133436,176297668,1841222776,
%T 19467629120,207978652416,2241618514120,24345336854400,
%U 266168049520832,2927074607294300,32356419163487336,359330087240388828,4007079691584624576
%N G.f. satisfies: A(x) = Product_{n>=1} (1 + x^n*A(x))^2/(1 - x^n*A(x))^2.
%C Related q-series (Heine) identity:
%C 1 + Sum_{n>=1} x^n*Product_{k=0..n-1} (y+q^k)*(z+q^k)/((1-x*q^k)*(1-q^(k+1)) = Product_{n>=0} (1+x*y*q^n)*(1+x*z*q^n)/((1-x*q^n)*(1-x*y*z*q^n)); here q=x, x=x*A(x), y=z=1.
%H Vaclav Kotesovec, <a href="/A192620/b192620.txt">Table of n, a(n) for n = 0..300</a>
%F G.f. A(x) satisfies:
%F (1) A(x) = 1 + Sum_{n>=1} x^n*A(x)^n * Product_{k=1..n} (1+x^(k-1))^2/((1-x^k*A(x))*(1-x^k), due to the Heine identity.
%F (2) A(x)^(1/2) = 1 + Sum_{n>=1} x^n*A(x)^n * Product_{k=1..n} (1+x^(k-1))/(1-x^k), due to the q-binomial theorem.
%F Equals the self-convolution of A192621.
%F a(n) ~ c * d^n / n^(3/2), where d = 12.042513458183758627924432194393539477581... and c = 1.04958502741924123967536156787764354342367951743839... - _Vaclav Kotesovec_, Oct 04 2023
%F Radius of convergence r = 0.083039143238027913107320323917684421045... = 1/d and A(r) = 2.638555772492608872250287025192536127217... satisfy A(r) = 1 / Sum_{n>=1} 4*r^n/(1 - r^(2*n)*A(r)^2) and A(r) = Product_{n>=1} (1 + r^n*A(r))^2/(1 - r^n*A(r))^2. - _Paul D. Hanna_, Mar 02 2024
%e G.f.: A(x) = 1 + 4*x + 28*x^2 + 224*x^3 + 1948*x^4 + 17928*x^5 + ...
%e The g.f. A = A(x) satisfies the following relations:
%e (0) A = (1+x*A)^2/(1-x*A)^2 * (1+x^2*A)^2/(1-x^2*A)^2 * (1+x^3*A)^2/(1-x^3*A)^2 * ...
%e (1) A = 1 + 4*x*A/((1-x*A)*(1-x)) + 4*x^2*A^2*(1+x)^2/((1-x*A)*(1-x^2*A)*(1-x)*(1-x^2)) + 4*x^3*A^3*(1+x)^2*(1+x^2)^2/((1-x*A)*(1-x^2*A)*(1-x^3*A)*(1-x)*(1-x^2)*(1-x^3)) + ...
%e (2) A^(1/2) = 1 + 2*x*A/(1-x) + 2*x^2*A^2*(1+x)/((1-x)*(1-x^2)) + 2*x^3*A^3*(1+x)*(1+x^2)/((1-x)*(1-x^2)*(1-x^3)) + ...
%t nmax = 30; A[_] = 0; Do[A[x_] = Product[(1 + x^k*A[x])^2/(1 - x^k*A[x])^2, {k, 1, nmax}] + O[x]^(nmax + 1) // Normal, nmax + 1]; CoefficientList[A[x], x] (* _Vaclav Kotesovec_, Oct 04 2023 *)
%t (* Calculation of constants {d,c}: *) Chop[{1/r, Sqrt[(r*s^(3/2)*((-1 + s)*Derivative[0, 1][QPochhammer][-s, r] + Sqrt[s]*(1 + s)*Derivative[0, 1][QPochhammer][s, r]))/(2* Pi*(1 + s)*QPochhammer[s, r]* (2* s*((1 + s^2)/(-1 + s^2)^2) + (QPolyGamma[1, Log[-s]/Log[r], r] - QPolyGamma[1, Log[s]/Log[r], r])/ Log[r]^2))]} /. FindRoot[{(-1 + s)^2*(QPochhammer[-s, r]^2/((1 + s)^2*QPochhammer[s, r]^2)) == s, 1 - 4*(s/(-1 + s^2)) + (2*(QPolyGamma[0, Log[-s]/Log[r], r] - QPolyGamma[0, Log[s]/Log[r], r]))/Log[r] == 0}, {r, 1/12}, {s, 2}, WorkingPrecision -> 120]] (* _Vaclav Kotesovec_, Mar 03 2024 *)
%o (PARI) {a(n)=local(A=1+x);for(i=1,n,A=prod(k=1,n,(1+x^k*A)^2/(1-x^k*A+x*O(x^n))^2));polcoeff(A,n)}
%o (PARI) {a(n)=local(A=1+x);for(i=1,n,A=1+sum(m=1,n,x^m*A^m*prod(k=1,m,(1+x^(k-1))^2/((1-x^k*A+x*O(x^n))*(1-x^k)))));polcoeff(A,n)}
%o (PARI) {a(n)=local(A=1+x);for(i=1,n,A=(1+sum(m=1,n,x^m*A^m*prod(k=1,m,(1+x^(k-1))/(1-x^k+x*O(x^n)))))^2);polcoeff(A,n)}
%Y Cf. A190862, A192621, A192622, A192623, A192624, A192625.
%K nonn
%O 0,2
%A _Paul D. Hanna_, Jul 06 2011