Hallo,
Die tupledDivisors[] laufen besser, wenn man die rekursive Formulierung
vereinfacht und den Gedanken implementiert, dass Divisortupeln, denen man
selbst wg. der Graduierung der Tupel eine 1 vorangestellt hat, von
branch[] nicht zu behandeln sind:
In[1]:= Remove[tupledDivisors, branch, leaf]
leaf[o_Integer /; PrimeQ[o]] := leaf[o] = {{1, o}}
leaf[o_Integer /; ! PrimeQ[o]] := leaf[o] =(* Skip {1,o} terms *)
Rest[Take[Transpose[{#, Reverse[#]} &[#]], Ceiling[Length[#]/2]] &[
Divisors[o]]]
branch[l_List] :=
Block[ (* handle distinct members of l only *) {z =
Flatten[Thread /@
Transpose[{Most[FoldList[Plus, 1, Length /@ #]],
leaf /@ First /@ #} &[Split[l]]], 1]},
Sort /@
MapThread[
Flatten[ReplacePart[#1, Rule @@ #2]] &, {Table[l, {Length[z]}], z}]
]
tupledDivisors[A_Integer?Positive, P_Integer?Positive] :=
Nest[((Join[{1}, #] & /@ #) \[Union]
Select[Flatten[branch /@ Select[#, First[#] > 1 &], 1],
First[#] > 1 &]) &, {{P}}, A - 1]
In[6]:= Length[tupledDivisors[12, 8 5 7 18 36]] // Timing
Out[6]= {0.858, 6756}
Zeichen der immer noch vorhandenen Mehrfacherzeugung von Tupeln ist das
Union[]: Einem effizienten Algorithmus sollte ein Join[] reichen. Und die
Tests noch:
In[29]:= (* This is http://oeis.org/A000041 *)
With[{N = 29,(* jede Primzahl funktioniert hier *) x = 17},
Table[Length[tupledDivisors[N, x^A]], {A, N}] ==
Table[Length[IntegerPartitions[A]], {A, N}]
]
Out[29]= True
In[30]:= (* This is http://oeis.org/A000110 *)
With[{N = 10},
Table[Length[
tupledDivisors[12, Times @@ Table[Prime[j], {j, A}]]], {A, N}] ==
BellB[Range[N]]
] // Timing
Out[30]= {20.046, True}
In[31]:= (* This is http://oeis.org/A001055, Test von Peter *)
Clear[pPSoll]
pPSoll = {1, 1, 1, 2, 1, 2, 1, 3, 2, 2, 1, 4, 1, 2, 2, 5, 1, 4, 1, 4,
2, 2, 1, 7, 2, 2, 3, 4, 1, 5, 1, 7, 2, 2, 2, 9, 1, 2, 2, 7, 1, 5,
1, 4, 4, 2, 1, 12, 2, 4, 2, 4, 1, 7, 2, 7, 2, 2, 1, 11, 1, 2, 4,
11, 2, 5, 1, 4, 2, 5, 1, 16, 1, 2, 4, 4, 2, 5, 1, 12, 5, 2, 1, 11,
2, 2, 2, 7, 1, 11, 2, 4, 2, 2, 2, 19, 1, 4, 4, 9};
In[33]:= (* Test auf A001055 *)
Table[tupledDivisors[Total[FactorInteger[n][[All, 2]]], n] //
Length, {n, 100}] == pPSoll
Out[33]= True
Gruss
Udo.