%I #45 Feb 23 2019 11:19:22
%S 1,2,3,7,11,39,103,367,1235,4439,15935,58509,215251,799697,2983217,
%T 11187567,42109451,159082753,602809327,2290684251,8726308317,
%U 33318661277,127479700199,488672302909,1876500180291,7217308815887,27799998949873,107228568948547
%N Number of n-bead necklace structures using a maximum of four different colored beads.
%C Turning over the necklace is not allowed. Colors may be permuted without changing the necklace structure.
%D M. R. Nester (1999). Mathematical investigations of some plant interaction designs. PhD Thesis. University of Queensland, Brisbane, Australia. [See A056391 for pdf file of Chap. 2]
%H E. N. Gilbert and J. Riordan, <a href="http://projecteuclid.org/euclid.ijm/1255631587">Symmetry types of periodic sequences</a>, Illinois J. Math., 5 (1961), 657-665.
%H Marko Riedel, <a href="https://math.stackexchange.com/questions/3121744/">Necklaces with swappable colors by Power Group Enumeration</a>
%H Marko Riedel, <a href="/A056292/a056292.maple.txt">Maple code for any necklace size, any number of swappable colors, by Power Group Enumeration.</a>
%H N. J. A. Sloane, <a href="/A000013/a000013.txt">Maple code for this and related sequences</a>
%F Use de Bruijn's generalization of Polya's enumeration theorem as discussed in reference.
%F From _Robert A. Russell_, May 29 2018: (Start)
%F a(n) = (1/n) * Sum_{d|n} phi(d) * ([d==0 mod 12] * (4*S2(n/d+3, 4) - 24*S2(n/d+2, 4) + 44*S2(n/d+1, 4) - 24*S2(n/d, 4)) + [d==6 mod 12] * (3*S2(n/d+3, 4) - 18*S2(n/d+2, 4) + 33*S2(n/d+1, 4) - 18*S2(n/d, 4)) + [d==4 mod 12 | d==8 mod 12] * (3*S2(n/d+3, 4) - 19*S2(n/d+2, 4) + 38*S2(n/d+1, 4) - 24*S2(n/d, 4)) + [d==3 mod 12 | d=9 mod 12] * (2*S2(n/d+3, 4) - 13*S2(n/d+2, 4) + 26*S2(n/d+1, 4) - 15*S2(n/d, 4)) + [d==2 mod 12 | d=10 mod 12] * (2*S2(n/d+3, 4) - 13*S2(n/d+2, 4) + 27*S2[n/d+1,4) - 18*S2(n/d, 4)) + [d mod 12 in {1,5,7,11}] * (S2(n/d+3, 4) - 8*S2(n/d+2, 4) + 20*S2(n/d+1, 4) - 15*S2(n/d, 4))), where S2(n, k) is the Stirling subset number, A008277.
%F G.f.: 1 - Sum_{d>0} (phi(d) / d) * ([d==0 mod 12] * log(1-4x^d) + [d==6 mod 12] * 3*log(1-4x^d) / 4 + [d==4 mod 12 | d==8 mod 12] * (2*log(1-4x^d) + log(1-x^d)) / 3 + [d==3 mod 12 | d=9 mod 12] * (3*log(1-4x^d) + 2*log(1-2x^d)) / 8 + [d==2 mod 12 | d=10 mod 12] * (5*log(1-4x^d) + 4*log(1-x^d)) / 12 + [d mod 12 in {1,5,7,11}] * (log(1-4x^d) + 6*log(1-2x^d) + 8*log(1-x^d)) / 24).
%F (End)
%t Adn[d_, n_] := Module[{ c, t1, t2}, t2 = 0; For[c = 1, c <= d, c++, If[Mod[d, c] == 0 , t2 = t2 + (x^c/c)*(E^(c*z) - 1)]]; t1 = E^t2; t1 = Series[t1, {z, 0, n+1}]; Coefficient[t1, z, n]*n!]; Pn[n_] := Module[{ d, e, t1}, t1 = 0; For[d = 1, d <= n, d++, If[Mod[n, d] == 0, t1 = t1 + EulerPhi[d]*Adn[d, n/d]/n]]; t1/(1 - x)]; Pnq[n_, q_] := Module[{t1}, t1 = Series[Pn[n], {x, 0, q+1}] ; Coefficient[t1, x, q]]; a[n_] := Pnq[n, 4]; Table[Print[an = a[n]]; an, {n, 1, 25}] (* _Jean-François Alcover_, Oct 04 2013, after _N. J. A. Sloane_'s Maple code *)
%t (* This program uses Gilbert and Riordan's recurrence formula, which they recommend for calculations: *)
%t Adn[d_, n_] := Adn[d, n] = If[1==n, DivisorSum[d, x^# &],
%t Expand[Adn[d, 1] Adn[d, n-1] + D[Adn[d, n-1], x] x]];
%t Table[SeriesCoefficient[DivisorSum[n, EulerPhi[#] Adn[#, n/#] &]
%t /(n (1 - x)), {x, 0, 4}], {n, 1, 40}] (* _Robert A. Russell_, Feb 24 2018 *)
%t From _Robert A. Russell_, May 29 2018: (Start)
%t Table[(1/n) DivisorSum[n, EulerPhi[#] Which[Divisible[#,12], 4 StirlingS2[n/#+3,4] - 24 StirlingS2[n/#+2,4] + 44 StirlingS2[n/#+1,4] - 24 StirlingS2[n/#,4], Divisible[#,6], 3 StirlingS2[n/#+3,4] - 18 StirlingS2[n/#+2,4] + 33 StirlingS2[n/#+1,4] - 18 StirlingS2[n/#,4], Divisible[#,4], 3 StirlingS2[n/#+3,4] - 19 StirlingS2[n/#+2,4] + 38 StirlingS2[n/#+1,4] - 24 StirlingS2[n/#,4], Divisible[#,3], 2 StirlingS2[n/#+3,4] - 13 StirlingS2[n/#+2,4] + 26 StirlingS2[n/#+1,4] - 15 StirlingS2[n/#,4], Divisible[#,2], 2 StirlingS2[n/#+3,4] - 13 StirlingS2[n/#+2,4] + 27 StirlingS2[n/#+1,4] - 18 StirlingS2[n/#,4], True, StirlingS2[n/#+3,4] - 8 StirlingS2[n/#+2,4] + 20 StirlingS2[n/#+1,4] - 15 StirlingS2[n/#,4]] &],{n, 1, 40}]
%t mx = 40; Drop[CoefficientList[Series[1 - Sum[(EulerPhi[d] / d) Which[
%t Divisible[d, 12], Log[1 - 4x^d], Divisible[d, 6],
%t 3 Log[1 - 4x^d] / 4, Divisible[d, 4] ,
%t (2 Log[1 - 4x^d] + Log[1 - x^d]) / 3, Divisible[d, 3],
%t (3 Log[1 - 4x^d] + 2 Log[1 - 2x^d]) / 8,
%t Divisible[d, 2], (5 Log[1 - 4x^d] + 4 Log[1 - x^d]) / 12,
%t True, (Log[1 - 4x^d] + 6 Log[1 - 2x^d] + 8 Log[1 - x^d]) / 24], {d, 1, mx}], {x, 0, mx}], x], 1]
%t (End)
%Y Cf. A000013, A001868.
%K nonn
%O 1,2
%A _Marks R. Nester_