3
$\begingroup$

Update

Is there any better way of generating the nth prime power?

Chip Hurst gave a great solution for a list below, so

PrimePowersUpTo[x_] := Union @@ Table[Prime[Range[PrimePi[x^(1/n)]]]^n, {n, Log2[x]}]
pp=PrimePowersUpTo[10^7];
pp[[n]]

would work, but I was wondering whether there was a more direct way?

Original question

I am trying to get a list of prime powers, but I can only seem to get one that goesx to 10^7.

xmax = 10000000;
rr = DeleteCases[Table[If[PrimePowerQ[n] == True, n, 0], {n, 1, xmax}], 0];

or

Select[Range[10], PrimePowerQ]

I tried using Artes' idea of taking log base prime as per here, but didn't get very far.

Is there any way I can go to 10^12?

$\endgroup$
  • $\begingroup$ PrimePi[10^12] is a very large number ... $\endgroup$ – Dr. belisarius Nov 26 '14 at 17:17
  • $\begingroup$ @belisarius yes, I was rather hoping I could go that high ... $\endgroup$ – martin Nov 26 '14 at 17:19
  • $\begingroup$ You can estimate reasonably PrimePi with RiemannR even if Mathematica implementation doesn't help, e.g. RiemannR[10^12] // N yields 3.76079*10^10. See also Approximation to the prime counting function. $\endgroup$ – Artes Nov 26 '14 at 18:06
  • $\begingroup$ @Artes I'm not sure I follow - how would this help me find the nth prime power? $\endgroup$ – martin Nov 26 '14 at 18:13
3
$\begingroup$

Here's my way of finding the nth prime power. I do a binary search. Unfortunately it's slow for large inputs.

lower[_] = 2;
upper[x_] := Prime[x]

PrimePower[n_] := Module[{lo = lower[n], hi = upper[n], mid, pp},
    mid = Quotient[lo + hi, 2];

    While[lo < mid < hi && (pp = PrimePowerPi[mid]) != n,
        If[pp > n,
            hi = mid,
            lo = mid
        ];
        mid = Quotient[lo + hi, 2];
    ];

    PreviousPrimePower[mid]
]

PreviousPrimePower[x_] := Max[Table[NextPrime[x^(1/n), -1]^n, {n, Log2[x]}]]

Examples:

PrimePowerPi[10^7]
665134
PrimePower[665134]
9999991
PrimePowerPi[10^12]
37607992088
PrimePower[37607992088]
999999999989

From Daniel Lichtblau:

Feel free to use this variant (divides more finely).

PrimePowerb[n_] := Module[
  {lo = lower[n], hi = upper[n], mid, ppmid, pphi, pplo}, 
  mid = Round[.95*hi];
  ppmid = PrimePowerPi[mid];
  pplo = 1;
  pphi = PrimePowerPi[hi];
  While[lo < mid < hi && ppmid != n,
   If[ppmid > n,
    hi = mid;
    pphi = ppmid;
    mid = mid + Round[(lo - mid)*(n - ppmid)/(pplo - ppmid)];
    ,
    lo = mid;
    pplo = ppmid;
    mid = mid + Round[(hi - mid)*(n - ppmid)/(pphi - ppmid)];
    ];
   If[mid == lo, mid++];
   If[mid == hi, mid--];
   ppmid = PrimePowerPi[mid];
   ];
  PreviousPrimePower[mid]]
$\endgroup$
  • $\begingroup$ thank you! this works really well :) $\endgroup$ – martin Nov 26 '14 at 18:45
  • 3
    $\begingroup$ Feel free to use this variant (divides more finely). PrimePowerb[n_] := Module[ {lo = lower[n], hi = upper[n], mid, ppmid, pphi, pplo}, mid = Round[.95*hi]; ppmid = PrimePowerPi[mid]; pplo = 1; pphi = PrimePowerPi[hi]; While[lo < mid < hi && ppmid != n, If[ppmid > n, hi = mid; pphi = ppmid; mid = mid + Round[(lo - mid)*(n - ppmid)/(pplo - ppmid)]; , lo = mid; pplo = ppmid; mid = mid + Round[(hi - mid)*(n - ppmid)/(pphi - ppmid)]; ]; If[mid == lo, mid++]; If[mid == hi, mid--]; ppmid = PrimePowerPi[mid]; ]; PreviousPrimePower[mid]] $\endgroup$ – Daniel Lichtblau Nov 26 '14 at 21:28
  • $\begingroup$ @DanielLichtblau just seen thas - thank you! :) $\endgroup$ – martin Nov 28 '14 at 18:03
3
$\begingroup$

How about

PrimePowersUpTo[x_] := Union @@ Table[Prime[Range[PrimePi[x^(1/n)]]]^n, {n, Log2[x]}]

PrimePowersUpTo[10^7] // Length
665134

Now about getting the first 10^12 prime powers. You would want to get them in batches (segmented sieve) because there are about 38 billion of them:

PrimePowerPi[x_] := Sum[PrimePi[x^(1/n)], {n, Log2[x]}]

PrimePowerPi[10^12]
37607992088
$\endgroup$
  • $\begingroup$ is there any way of generating the nth prime power? $\endgroup$ – martin Nov 26 '14 at 17:35
  • $\begingroup$ You might generate the nth prime and then subtract your way down to account for the relatively fewer powers that also near that value. Possibly you can use a tighter estimate since you have ways to approximate how many of those there are, but then you also have to count carefully so as not to include powers that are too large. $\endgroup$ – Daniel Lichtblau Nov 26 '14 at 18:04
  • $\begingroup$ @DanielLichtblau I don't suppose you could elaborate in an answer, could you? $\endgroup$ – martin Nov 26 '14 at 18:05

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.