|
%I #15 Mar 29 2023 17:05:36
%S 1,1,2,12,96,1000,13500,221718,4301808,97747200,2555001360,
%T 75526842600,2503943418240,92407030642056,3759862792921872,
%U 167429488088545200,8120958429706093440,426777425467443381120,24161214872571127574400,1467122583066982481802816
%N Expansion of e.g.f.: Sum_{n>=0} log(1 + x/(1-x)^n)^n / n!.
%H G. C. Greubel, <a href="/A219119/b219119.txt">Table of n, a(n) for n = 0..300</a>
%F E.g.f.: Sum_{n>=0} binomial(1/(1-x)^n, n) * x^n.
%F E.g.f.: Sum_{n>=0} x^n/n! * Product_{k=0..n-1} (1/(1-x)^n - k).
%F E.g.f.: Sum_{n>=0} x^n/n! * Sum_{k=0..n} Stirling1(n,k) / (1-x)^(n*k).
%e E.g.f.: A(x) = 1 + x + 2*x^2/2! + 12*x^3/3! + 96*x^4/4! + 1000*x^5/5! + ...
%e where the g.f. satisfies the identities:
%e A(x) = 1 + log(1+x/(1-x)) + log(1+x/(1-x)^2)^2/2! + log(1+x/(1-x)^3)^3/3! + log(1+x/(1-x)^4)^4/4! + log(1+x/(1-5*x)^5)^5/5! + ...
%e A(x) = 1 + binomial(1/(1-x),1)*x + binomial(1/(1-x)^2,2)*x^2 + binomial(1/(1-x)^3,3)*x^3 + binomial(1/(1-x)^4,4)*x^4 + binomial(1/(1-x)^5,5)*x^5 + ...
%e A(x) = 1 + x/(1-x) + x^2/(1-x)^4*(1-(1-x)^2)/2! + x^3/(1-x)^9*(1-(1-x)^3)*(1-2*(1-x)^3)/3! + x^4/(1-x)^16*(1-(1-x)^4)*(1-2*(1-x)^4)*(1-3*(1-x)^4)/4! + ...
%t m:= 50;
%t f[x_, m_]:= Sum[Product[(1/(1-x)^n -j), {j,0,n-1}]*x^n/n!, {n,0,m+1}];
%t CoefficientList[Series[f[x,m], {x,0,m}], x]*Range[0,m]!] (* _G. C. Greubel_, Feb 02 2023 *)
%o (PARI) {a(n)=n!*polcoeff(sum(m=0,n,log(1+x/(1-x+x*O(x^n))^m)^m/m!),n)}
%o (PARI) {a(n)=n!*polcoeff(sum(m=0,n,binomial(1/(1-x+x*O(x^n))^m,m)*x^m), n)}
%o for(n=0,30,print1(a(n),", "))
%o (PARI) {a(n)=n!*polcoeff(sum(m=0,n, x^m/m!*prod(k=0,m-1, (1/(1-x)^m-k+x*O(x^n)))),n)}
%o (PARI) {Stirling1(n, k)=n!*polcoeff(binomial(x, n), k)}
%o {a(n)=local(A=1+x); A=sum(m=0, n, sum(k=0, m, Stirling1(m, k)/(1-x+x*O(x^n))^(m*k))*x^m/m!); n!*polcoeff(A, n)}
%o for(n=0, 30, print1(a(n), ", "))
%o (Magma)
%o m:=50;
%o f:= func< x | 1 + (&+[(&*[1/(1-x)^n -j: j in [0..n-1]])*x^n/Factorial(n) : n in [1..m+2]]) >;
%o R<x>:=PowerSeriesRing(Integers(), m);
%o Coefficients(R!(Laplace( f(x) ))); // _G. C. Greubel_, Feb 02 2023
%o (SageMath)
%o m=50
%o def f(x): return sum(product(1/(1-x)^n - k for k in range(n))*x^n/factorial(n) for n in range(m+2))
%o def A219119_list(prec):
%o P.<x> = PowerSeriesRing(QQ, prec)
%o return P( f(x) ).egf_to_ogf().list()
%o A219119_list(m) # _G. C. Greubel_, Feb 02 2023
%Y Cf. A216839.
%K sign
%O 0,3
%A _Paul D. Hanna_, Nov 13 2012
| 