%I #27 Oct 18 2018 03:06:03
%S 1,3,9,22,51,108,221,429,810,1476,2631,4572,7802,13056,21519,34918,
%T 55935,88452,138332,213990,327852,497592,748833,1117692,1655719,
%U 2434938,3556791,5161808,7445631,10677096,15226658,21599469,30485268,42817788
%N Series expansion of (eta(q^9) / eta(q))^3 in powers of q.
%C Cubic AGM theta functions: a(q) (see A004016), b(q) (A005928), c(q) (A005882).
%H Seiichi Manyama, <a href="/A121589/b121589.txt">Table of n, a(n) for n = 1..10000</a>
%H Kevin Acres, David Broadhurst, <a href="https://arxiv.org/abs/1810.07478">Eta quotients and Rademacher sums</a>, arXiv:1810.07478 [math.NT], 2018. See Table 1 p. 10.
%F Euler transform of period 9 sequence [3, 3, 3, 3, 3, 3, 3, 3, 0, ...].
%F G.f.: x * (Product_{k>0} (1 - x^(9*k) / (1 - x^k))^3.
%F Expansion of c(q^3) / (3 * b(q)) = (c(q) / (3 * b(q^3))^3 in powers of q where b(), c() are cubic AGM functions.
%F G.f. A(x) satisfies 0 = f(A(x), A(x^2)) where f(u, v) = (u - v^2) * (u^2 - v) - 2 * u * v * ( 3 * (u + v) + 13 * u * v ).
%F G.f. A(x) satisfies 0 = f(A(x), A(x^3)) where f(u, v) = u^3 - v * (1 + 9 * v + 27 * v^2) * (1 + 9 * u + 27 * u^2).
%F G.f. A(x) satisfies 0 = f(A(x), A(x^2), A(x^3), A(x^6)) where f(u1, u2, u3, u6) = u1 * u2 - (1 + 3 * (u1 + u2)) * (u3 + u6 + 9 * u3 * u6).
%F G.f. is a period 1 Fourier series which satisfies f(-1 / (9 t)) = (1/27) g(t) where q = exp(2 Pi i t) and g() is the g.f. of A131986.
%F a(n) ~ exp(4*Pi*sqrt(n)/3) / (27 * sqrt(6) * n^(3/4)). - _Vaclav Kotesovec_, Sep 07 2015
%F a(1) = 1, a(n) = (3/(n-1))*Sum_{k=1..n-1} A116607(k)*a(n-k) for n > 1. - _Seiichi Manyama_, Apr 01 2017
%F Convolution inverse of A131986. Convolution cube of A104502. - _Michael Somos_, Nov 02 2017
%e G.f. = q + 3*q^2 + 9*q^3 + 22*q^4 + 51*q^5 + 108*q^6 + 221*q^7 + 429*q^8 + ...
%p N:= 100: # to get a(1)..a(N)
%p S:= series(q*Product(1-q^(9*k),k=1..N/9)/Product((1-q^k)^3, k=1..N),q,N+1):
%p seq(coeff(S,q,n),n=1..N); # _Robert Israel_, Nov 02 2017
%t nmax = 40; CoefficientList[Series[Product[((1-x^(9*k))/(1-x^k))^3, {k, 1, nmax}], {x, 0, nmax}], x] (* _Vaclav Kotesovec_, Sep 07 2015 *)
%t QP = QPochhammer; s = (QP[q^9]/QP[q])^3 + O[q]^40; CoefficientList[s, q] (* _Jean-François Alcover_, Nov 30 2015, adapted from PARI *)
%t a[ n_] := SeriesCoefficient[ q (QPochhammer[ q^9] / QPochhammer[ q])^3, {q, 0, n}]; (* _Michael Somos_, Nov 02 2017 *)
%o (PARI) {a(n) = my(A); if( n<1, 0, n--; A = x * O(x^n); polcoeff( (eta(x^9 + A) / eta(x + A))^3, n))};
%Y Cf. A104502, A131986.
%K nonn
%O 1,2
%A _Michael Somos_, Aug 09 2006
%E Second formula corrected by _Vaclav Kotesovec_, Sep 07 2015