I would like to recursively replicate each stage of the Euler product proof. This does it:

ap[r1_, r2_, e_] :=  With[{r = Range[r1, r2]}, Pick[r, PrimeOmega[r], e]]
ap2 = ap[1, 10^3, 1];
β = 10^3; 
a = Table[1/n, {n, β}]


b = Reverse@Complement[a, 1/ap2[[1]] a]



c1 = Reverse@Sort[Join[Complement[b, 1/ap2[[2]] b], -Complement[1/ap2[[2]] b, b]]];
c2 = DeleteCases[Table[If[Denominator[c1[[α]]] > β, 0, c1[[α]]], {α, Length@c1}], 0]



d1 = Reverse@Sort[Join[Complement[c2, 1/ap2[[3]] c2], -Complement[1/ap2[[3]] c2, c2]]];
d2 = DeleteCases[Table[If[Denominator[d1[[α]]] > β, 0, d1[[α]]], {α, Length@d1}], 0];


I can do it by naming each stage, but am finding difficulty coding it using Fold.

Note: The "minus Complement" is not necessary for the primes, but for other Euler products, it is necessary.


Hi Martin: Here is a partial solution using FoldList, along with a helper function fz defined as

fz[z_, p_Integer] := 

I use fz to multiply the current expansion z by 1-p^-s, for prime p, then mess with the expression to get a list of single integers to the power -s. The initial value for z is Table[1/k^s,{k,1,n}], where n is a maximum. In the question, you used n=1000.

Then use EulerProductRecursion[n] to get a table of successive iterations.

EulerProductRecursion[n_Integer] :=
   Module[{z = Table[1/k^s, {k, 1, n}], 
           p = Prime[Range[PrimePi[Sqrt[n]]]]},
      Total[FoldList[fz[#1, #2] &, z, p], {2}]]

The prime cutoff using PrimePi[Sqrt[n]] is not fundamental to the solution, and can be set to simply PrimePi[n], but this latter setting is boring for p>Sqrt[n].

  • $\begingroup$ thanks for your answer - will try to have a play around with it :) $\endgroup$ – martin Dec 30 '14 at 4:25

I actually went for this in the end:

recursionNumber = 10; range = 100;
f[list_, a_] :=  list - Take[Riffle[Table[0, {x, plist[[a]] range}], list/plist[[a]], 
plist[[a]]], range]
plist = Prime[Range[range]];
list = DeleteCases[Fold[f, Table[1/n, {n, range}], Range[recursionNumber]], 0];
Plus @@ Join[{1}, Rest@Table[DisplayForm[FractionBox[Numerator[list[[n]]], 
Denominator[list[[n]]^s]]], {n, Length@list}]]

as it was more what I was after (should really have been a bit clearer in question), but hats off to KennyColnago for his excellent answer.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.