%I #25 Jan 14 2020 14:09:34
%S 1,5,23,120,724,5012,39332,345832,3371976,36135792,422379792,
%T 5349561984,72996193152,1067779243008,16670798231040,276718772067840,
%U 4866610479828480,90401487246167040,1768784607499944960,36360467544043008000,783508616506603008000
%N Sum (number of cycles)^2 over all n! permutations of [1..n].
%C Sum (number of cycles) over all n! permutations of [1..n] gives A000254.
%C a(n) equals -1 times the coefficient of x of the characteristic polynomial of the n X n matrix whose (i,j)-entry is equal to i+1 if i=j, and is equal to 1 otherwise. [_John M. Campbell_, May 24 2011]
%H N. J. A. Sloane and Alois P. Heinz, <a href="/A151881/b151881.txt">Table of n, a(n) for n = 1..450</a> (first 30 terms from N. J. A. Sloane)
%F a(n) = (-1)^(n+1)*(Stirling1(n+1,2)-2*Stirling1(n+1,3)). - _Vladeta Jovovic_, Jul 22 2009
%F a(n) = n!*(H(n)+H(n)^2-H2(n)), where H(n)=sum(k=1..n,1/k), H2(n)=sum(k=1..n,1/k^2). - _Vladimir Reshetnikov_, Oct 20 2015
%F E.g.f.: (log(1-x)^2-log(1-x))/(1-x). - _Vladimir Reshetnikov_, Oct 20 2015
%F a(n) = Sum_{k=0..n} |Stirling1(n,k)|*k^2. - _Peter Luschny_, Jan 14 2020
%p with(combinat): with(numtheory):
%p M:=30;
%p for n from 1 to M do
%p p:=partition(n); s:=0:
%p for k from 1 to nops(p) do
%p # get next partition of n
%p # convert partition to list of sizes of parts
%p q:=convert(p[k], multiset);
%p for i from 1 to n do a(i):=0: od:
%p for i from 1 to nops(q) do a(q[i][1]):=q[i][2]: od:
%p # get number of parts:
%p nump := add(a(i), i=1..n);
%p # get multiplicity:
%p c:=1: for i from 1 to n do c:=c*a(i)!*i^a(i): od:
%p prop:=nump^2;
%p s:=s + (n!/c)*prop;
%p od;
%p lprint(n, s);
%p A[n]:=s;
%p od:
%p [seq(A[n], n=1..M)];
%p # Alternatively after Reshetnikov:
%p a := n -> n!*(add(1/k,k=1..n)^2 + add(1/k-1/k^2,k=1..n)):
%p seq(a(n), n=1..19); # _Peter Luschny_, Oct 21 2015
%t Table[-Coefficient[CharacteristicPolynomial[ Array[KroneckerDelta[#1, #2]((((#1+1)))-1)+1&,{n,n}],x],x], {n,1,10}] (* _John M. Campbell_, May 24 2011 *)
%t Table[n! (HarmonicNumber[n] + HarmonicNumber[n]^2 - HarmonicNumber[n, 2]), {n, 1, 20}] (* _Vladimir Reshetnikov_, Oct 20 2015 *)
%Y Cf. A000254, A151882.
%K nonn
%O 1,2
%A _N. J. A. Sloane_, Jul 22 2009