%I #36 May 20 2021 10:56:19
%S 1,1,1,1,1,1,1,1,2,1,1,1,3,4,1,1,1,4,10,10,1,1,1,5,19,43,26,1,1,1,6,
%T 31,116,201,80,1,1,1,7,46,245,776,1038,246,1,1,1,8,64,446,2126,5620,
%U 5538,810,1,1,1,9,85,735,4751,19811,42288,30667,2704,1
%N Number A(n,k) of necklaces with n white beads and k*n black beads; square array A(n,k), n>=0, k>=0, read by antidiagonals.
%C For k>=1 is column k asymptotic to (k+1)^((k+1)*n-1/2) / (sqrt(2*Pi) * k^(k*n+1/2) * n^(3/2)). - _Vaclav Kotesovec_, Aug 22 2015
%H Alois P. Heinz, <a href="/A261494/b261494.txt">Antidiagonals n = 0..140, flattened</a>
%H F. Ruskey, <a href="http://combos.org/necklace">Necklaces, Lyndon words, De Bruijn sequences, etc.</a>
%H F. Ruskey, <a href="/A000011/a000011.pdf">Necklaces, Lyndon words, De Bruijn sequences, etc.</a> [Cached copy, with permission, pdf format only]
%H Eric Weisstein's World of Mathematics, <a href="http://mathworld.wolfram.com/Necklace.html">Necklace</a>
%H Wikipedia, <a href="https://en.wikipedia.org/wiki/Necklace_(combinatorics)">Necklace (combinatorics)</a>
%H <a href="/index/Ne#necklaces">Index entries for sequences related to necklaces</a>
%F A(n,k) = 1/((k+1)*n) * Sum_{d|n} C((k+1)*n/d,n/d) * A000010(d) for n>0, A(0,k) = 1.
%F A(n,k) = 1/((k+1)*n)*Sum_{i=1..n} C((k+1)*gcd(n,i),gcd(n,i)) = 1/((k+1)*n)*Sum_{i=1..n} C((k+1)*n/gcd(n,i),n/gcd(n,i))*phi(gcd(n,i))/phi(n/gcd(n,i)) for n >= 1, where phi = A000010. - _Richard L. Ollerton_, May 19 2021
%e A(2,2) = 3: 000011, 000101, 001001.
%e A(3,2) = 10: 000000111, 000001011, 000010011, 000100011, 001000011, 010000011, 000010101, 000100101, 001000101, 001001001.
%e Square array A(n,k) begins:
%e 1, 1, 1, 1, 1, 1, 1, ...
%e 1, 1, 1, 1, 1, 1, 1, ...
%e 1, 2, 3, 4, 5, 6, 7, ...
%e 1, 4, 10, 19, 31, 46, 64, ...
%e 1, 10, 43, 116, 245, 446, 735, ...
%e 1, 26, 201, 776, 2126, 4751, 9276, ...
%e 1, 80, 1038, 5620, 19811, 54132, 124936, ...
%p with(numtheory):
%p A:= (n, k)-> `if`(n=0, 1, add(binomial((k+1)*n/d, n/d)
%p *phi(d), d=divisors(n))/((k+1)*n)):
%p seq(seq(A(n, d-n), n=0..d), d=0..14);
%t A[n_, k_] := If[n==0, 1, DivisorSum[n, Binomial[(k+1)*n/#, n/#]*EulerPhi[#] /((k+1)*n)&]]; Table[A[n, d-n], {d, 0, 14}, {n, 0, d}] // Flatten (* _Jean-François Alcover_, Feb 19 2017, translated from Maple *)
%o (PARI) a(n,k) = if(n<1, 1, sumdiv(n, d, binomial((k + 1)*n/d, n/d) * eulerphi(d)) / ((k + 1)*n));
%o for(d=0, 14, for(n=0, d, print1(a(n, d - n),", ");); print();) \\ _Indranil Ghosh_, Mar 25 2017
%Y Columns k=0-10 give: A000012, A003239, A082936, A261497, A261498, A261499, A261500, A261501, A261502, A261503, A261504.
%Y Main diagonal gives A261495.
%Y Lower diagonal gives A261496.
%Y Cf. A000010.
%K nonn,tabl
%O 0,9
%A _Alois P. Heinz_, Aug 21 2015