%I #36 May 17 2018 21:02:30
%S 1,5,87,2971,163121,12962661,1395857215,194634226067,33990369362241,
%T 7247035915622821,1848636684656077991,555005864462114884875,
%U 193458213840943964983537,77399534126148191747554181,35196002960227350045891984591,18037244723394790042393195636291
%N a(n) = 1/exp(1) * Sum_{k>=0} (k+n)!^2 / k!^3.
%C From _Peter Luschny_, Mar 27 2011: (Start)
%C Let B_{n}(x) = sum_{j>=0}(exp(j!/(j-n)!*x-1)/j!) then a(n) = 3! [x^3] taylor(B_{n}(x)), where [x^3] denotes the coefficient of x^3 in the Taylor series for B_{n}(x).
%C a(n) is column 3 of the square array representation of A090210. (End)
%H G. C. Greubel, <a href="/A069948/b069948.txt">Table of n, a(n) for n = 0..240</a>
%H K. A. Penson, P. Blasiak, A. Horzela, G. H. E. Duchamp and A. I. Solomon, <a href="http://arxiv.org/abs/0904.0369">Laguerre-type derivatives: Dobinski relations and combinatorial identities</a>, J. Math. Phys. vol. 50, 083512 (2009)
%F Integral representation as n-th moment of a positive function on a positive halfaxis (solution of the Stieltjes moment problem), in Maple notation: a(n)=int(x^n*2*BesselK(0,2*sqrt(x))*hypergeom([],[1,1],x)/exp(1), x=0..infinity), n=0,1... Special values of the hypergeometric function of type 2F2: a(n)=exp(-1)*GAMMA(n+1)^2*hypergeom([n+1, n+1], [1, 1], 1). - _Karol A. Penson_ and G. H. E. Duchamp (gduchamp2(AT)free.fr), Jan 09 2007
%F Recurrence: (8*n-7)*a(n) = (24*n^3 + 3*n^2 - 26*n + 4)*a(n-1) - (n-1)^2*(24*n^3 - 85*n^2 + 66*n + 13)*a(n-2) + (n-1)^2*(8*n+1)*(n-2)^4*a(n-3). - _Vaclav Kotesovec_, Jul 30 2013
%F a(n) ~ n^(2*n+1/3)*exp(n^(1/3) + 3*n^(2/3) - 2*n - 2/3)/sqrt(3) * (1 + 41/(54*n^(1/3)) + 13769/(29160*n^(2/3))). - _Vaclav Kotesovec_, Jul 30 2013
%p A069948 := proc(n) exp(-x)*n!^2*hypergeom([n+1,n+1],[1,1],x);
%p round(evalf(subs(x=1,%),99)) end:
%p seq(A069948(n),n=0..13); # _Peter Luschny_, Mar 30 2011
%p # second Maple program:
%p a:= n-> sum((k+n)!^2/k!^3, k=0..infinity)/exp(1):
%p seq(a(n), n=0..15); # _Alois P. Heinz_, May 17 2018
%t f[n_] := f[n] = Sum[(k + n)!^3/((k + n)!*(k!^3)*E), {k, 0, Infinity}]; Table[ f[n], {n, 0, 13}] (* or *)
%t Table[n!^2*HypergeometricPFQ[{n + 1, n + 1}, {1, 1}, 1]/Exp[1], {n, 0, 13}] (* _Robert G. Wilson v_, Jan 11 2007 *)
%o (PARI) {default(realprecision, 200)}; for(n=0,30, print1(round(exp(-1)*(n!)^2*sum(k=0,500, binomial(n+k, k)^2/k!)), ", ")) \\ _G. C. Greubel_, May 17 2018
%Y Cf. A000110, A020556, A069223, A090210.
%K nonn
%O 0,2
%A _Robert G. Wilson v_, May 02 2002
%E More terms from _Robert G. Wilson v_, Jan 11 2007