login
a(1) = 1; a(n+1) = Sum_{d|n, gcd(d, n/d) = 1} a(n/d) * a(d).
2

%I #8 Mar 09 2020 11:32:39

%S 1,1,2,4,8,16,36,72,144,288,592,1184,2384,4768,9608,19248,38496,76992,

%T 154272,308544,617152,1234448,2470080,4940160,9880608,19761216,

%U 39527200,79054400,158109088,316218176,632456976,1264913952,2529827904,5059658176,10119393344,20238787264

%N a(1) = 1; a(n+1) = Sum_{d|n, gcd(d, n/d) = 1} a(n/d) * a(d).

%H Robert Israel, <a href="/A333051/b333051.txt">Table of n, a(n) for n = 1..3320</a>

%p a[1]:= 1:

%p for n from 1 to 40 do

%p P:= ifactors(n)[2];

%p k:= nops(P);

%p t:= 0;

%p for S in combinat:-powerset(k) do

%p d:= mul(P[i][1]^P[i][2],i=S);

%p t:= t + a[d]*a[n/d]

%p od;

%p a[n+1]:= t

%p od:

%p seq(a[i],i=1..41); # _Robert Israel_, Mar 09 2020

%t a[1] = 1; a[n_] := a[n] = Sum[If[GCD[(n - 1)/d, d] == 1, a[(n - 1)/d] a[d], 0], {d, Divisors[n - 1]}]; Table[a[n], {n, 1, 36}]

%Y Cf. A038044, A122698.

%K nonn

%O 1,3

%A _Ilya Gutkovskiy_, Mar 06 2020