%I #6 Oct 10 2016 04:31:59
%S 1,4,33,424,7505,170496,4744873,156529024,5974216641,258970009600,
%T 12566664261041,674795685758976,39720422453156497,2543022838953017344,
%U 175923061842374645625,13076498369827187163136,1039320236257785348449537,87954586779787961844105216,7895887532418683295505005121,749448035808323155802521600000,74989090946223628553344278643281,7888932153987131087072869161631744
%N E.g.f.: A(x) = x*exp(A(x) - A(x)^2) + A(x)^2.
%F E.g.f. A(x) satisfies:
%F (1) exp(A(x) - A(x)^2) = LambertW(-x)/(-x).
%F (2) A(x) = -LambertW(-x) + A(x)^2.
%F (3) A(x) = C( -LambertW(-x) ), where C(x) = x + C(x)^2 is a g.f. of the Catalan numbers.
%F a(n) = Sum_{k=1..n} A000108(k-1) * n^(n-k) * k! * binomial(n-1,k-1), where A000108 is the Catalan numbers.
%F a(n) ~ 2^(2*n - 1/2) * n^(n-1) / (sqrt(3) * exp(3*n/4)). - _Vaclav Kotesovec_, Oct 10 2016
%e E.g.f.: A(x) = x + 4*x^2/2! + 33*x^3/3! + 424*x^4/4! + 7505*x^5/5! + 170496*x^6/6! + 4744873*x^7/7! + 156529024*x^8/8! + 5974216641*x^9/9! + 258970009600*x^10/10! +...
%e such that
%e A(x) - A(x)^2 = x + 2*x^2/2! + 9*x^3/3! + 64*x^4/4! + 625*x^5/5! + 7776*x^6/6! + 117649*x^7/7! +...+ n^(n-1)*x^n/n! +...
%e which equals -LambertW(-x).
%e RELATED SERIES.
%e A(x)^2 = 2*x^2/2! + 24*x^3/3! + 360*x^4/4! + 6880*x^5/5! + 162720*x^6/6! + 4627224*x^7/7! + 154431872*x^8/8! + 5931169920*x^9/9! + 257970009600*x^10/10! +...
%e exp(A(x)) = 1 + x + 5*x^2/2! + 46*x^3/3! + 629*x^4/4! + 11556*x^5/5! + 268537*x^6/6! + 7578040*x^7/7! + 252168009*x^8/8! + 9677553040*x^9/9! + 421010089901*x^10/10! +...
%e exp(A(x)^2) = 1 + 2*x^2/2! + 24*x^3/3! + 372*x^4/4! + 7360*x^5/5! + 179400*x^6/6! + 5228664*x^7/7! + 177953552*x^8/8! + 6940738368*x^9/9! + 305570622240*x^10/10! +...
%t Rest[CoefficientList[Series[(1 - Sqrt[1 + 4*LambertW[-x]])/2, {x, 0, 20}], x] * Range[0, 20]!] (* _Vaclav Kotesovec_, Oct 10 2016 *)
%o (PARI) {a(n) = sum(k=1,n, n^(n-k) * (2*k-2)!/(k-1)!^2 * (n-1)!/(n-k)! )}
%o for(n=1,25,print1(a(n),", "))
%o (PARI) {a(n) = my(A=x); for(i=0,n, A = x*exp(A - A^2 +x*O(x^n)) + A^2 ); n!*polcoeff(A,n)}
%o for(n=1,25,print1(a(n),", "))
%K nonn
%O 1,2
%A _Paul D. Hanna_, Oct 09 2016