Is there a way to find the approximate value of the nth prime?
Asked Answered
S

7

17

Is there a function which will return the approximate value of the n th prime? I think this would be something like an approximate inverse prime counting function. For example, if I gave this function 25 it would return a number around 100, or if I gave this function 1000 it would return a number around 8000. I don't care if the number returned is prime or not, but I do want it to be fast (so no generating the first n prime numbers to return the n th.)

I would like this so that I can generate the first n prime numbers using a sieve (Eratosthenes or Atkin). Therefore, the approximation for n th the would ideally never underestimate the value of the actual n th prime.

(Update: see my answer for a good method of finding the upper bound of the n th prime number.)

Succory answered 25/6, 2009 at 8:6 Comment(3)
en.wikipedia.org/wiki/…Poodle
just do your sieve segmented. and, definitely Eratosthenes.Chyle
The absolute best work on this subject today is Christian Axler's dissertation.Swelling
O
13

Tighter bounds:

static const unsigned short primes_small[] = {0,2,3,5,7,11};

static unsigned long nth_prime_upper(unsigned long n) {
  double fn = (double) n;
  double flogn, flog2n, upper;
  if (n < 6)  return primes_small[n];
  flogn  = log(n);
  flog2n = log(flogn);

  if      (n >= 688383)    /* Dusart 2010 page 2 */
    upper = fn * (flogn + flog2n - 1.0 + ((flog2n-2.00)/flogn));
  else if (n >= 178974)    /* Dusart 2010 page 7 */
    upper = fn * (flogn + flog2n - 1.0 + ((flog2n-1.95)/flogn));
  else if (n >=  39017)    /* Dusart 1999 page 14 */
    upper = fn * (flogn + flog2n - 0.9484);
  else                    /* Modified from Robin 1983 for 6-39016 _only_ */
    upper = fn * ( flogn  +  0.6000 * flog2n );

  if (upper >= (double) ULONG_MAX) {
     /* Adjust this as needed for your type and exception method */
    if (n <= 425656284035217743UL) return 18446744073709551557UL;
    fprintf(stderr, "nth_prime_upper overflow\n"; exit(-1);
  }

  return (unsigned long) ceil(upper);
}

These should not ever be less than the actual nth_prime, should work for any 64-bit input, and be an order of magnitude or more closer than the formula from Robin given earlier (or Wimblik's complicated range-limited formula). For my use I have a slightly larger small primes table so can tighten up the last estimate a bit more. Technically from the formulas we could use floor() instead of ceil() but I worry about precision.

Edit: Another option for improving this a bit is implementing good prime count bounds (e.g. Axler 2014) and doing a binary search on them. My code for this method takes ~10x longer than the above (still running in under a millisecond), but can reduce the error percentage by an order of magnitude.

If you want an estimate for the nth prime, you can do:

  • Cipolla 1902 (see Dusart 1999 page 12 or this paper. I find three terms (m=2) plus a third order correction factor to be useful, but with more terms it oscillates too much. The formula shown in the Wikipedia link is this formula (with m=2). Using the two-term inverse li or inverse Riemann R below will give better results.
  • Calculate the Dusart 2010 upper and lower bounds and average the results. Not too bad, though I suspect using a weighted average will work better as the bounds aren't equally tight.
  • li^{-1}(n) Since li(n) is a decent approximation to the prime count, the inverse is a decent nth_prime approximation. This, and all the rest, can be done fairly quickly as a binary search on the function.
  • li^{-1}(n) + li^{-1}(sqrt(n))/4 Closer, since this is getting closer to R(n)
  • R^{-1} The inverse Riemann R function is the closest average approximation I know that's reasonable.

Lastly, if you have a very fast prime count method such as one of the LMO implementations (there are three open source implementations now), you can write a fast precise nth_prime method. Computing the 10^10th prime can be done in a few milliseconds, and the 10^13th in a couple seconds (on a modern fast machine). The approximations are extremely fast at all sizes and work for far larger numbers, but everyone has a different idea of what "large" means.

Ottoman answered 22/8, 2014 at 6:12 Comment(6)
This is the best answer but I think there is little tiny bug in it. The statement flogn = log(flogn); should probably be flog2n = log(flogn);Beatty
@GregS, thank you! Fixed. I also added a paragraph about using the inverse prime count bound.Ottoman
@Ottoman This is one best works I've ever read on primes I think primes are going to obsess me. is there any chance of using prime number theorem at the expense of precision?Calkins
@DarshanPatil It is not even close to the estimates by Christian Axler, which is the best work on the subject today.Swelling
@Swelling Yes, I agree Thanks for letting me know!🤝 ✌️Calkins
I tried the solution here, and the estimate seems ok after 1000 primes, with margin < 1%. But for the first 1000 primes the estimate here is bad, with the error margin way above 1%.Swelling
S
11

Thanks for all of those answers. I suspected there was something fairly simple like that, but I couldn't find it at the time. I've done a bit more research too.

As I want it for a sieve to generate the first n prime numbers, I want the approximation to be greater than or equal to the nth prime. (Therefore, I want the upper bound of the nth prime number.)

Wikipedia gives the following upper bound for n >= 6

p_n <= n log n + n log log n   (1)

where p_n is the nth prime, and log is the natural logarithm. This is a good start, but it can overestimate by a not inconsiderable amount. This article in The College Mathematics Journal gives a tighter upper bound for n >= 7022

p_n <= n log n + n (log log n - 0.9385)   (2)

This is a much tighter bound as the following table shows

n         p_n         approx 1    error%  approx 2    error%
1         2                            
10        29          31          6.90 
100       541         613         13.31
1000      7919        8840        11.63
10000     104729      114306      9.14    104921      0.18
100000    1299709     1395639     7.38    1301789     0.16
1000000   15485863    16441302    6.17    15502802    0.11
10000000  179424673   188980382   5.33    179595382   0.10

I implemented my nth prime approximation function to use the second approximation for n >= 7022, the first approximation for 6 <= n < 7022, and an array lookup for the first 5 prime numbers.

(Although the first method isn't a very tight bound, especially for the range where I use it, I am not concerned as I want this for a sieve, and a sieve of smaller numbers is computationally cheap.)

Succory answered 1/7, 2009 at 13:0 Comment(6)
Unfortunately, this formula appears to be broken. p_8597 = 88789, but the formula gives 88759.3, which is an underestimate. It seems like it may be okay for n >= 8602.Bragi
How peculiar. That formula comes straight from this paper, which in turn is based on this French paper...Succory
From checking the primes up to 2e9 I got: p_8621 = 89009, p ∈ [88746.2, 89014.4], delta ∈ [-0.9718, -0.9407] By choosing 8621 you get the somewhat better approximation using the constant -0.9407. The next improvement is at p_15957. Yes, the paper must be wrong.Lunn
The bound actually holds from n > 39017 according to Dussart: "The kth prime is greater than k(lnk+lnlnk−1) for k≥2.", Math. Comp., 68(225), 411--415.Nalepka
@DavidJohnstone : The article link seems broken. Could please provide an updated link?Staphyloplasty
@Gaurav, here is a currently working link to the paper, An Upper Bound on the Nth Prime: maa.org/sites/default/files/jaroma03200545640.pdfDoomsday
S
4

Prime number theorem gives a number of primes below a threshold value, so it could be used to give an approximate value for the nth prime.

Sneaking answered 25/6, 2009 at 8:12 Comment(0)
M
4

As a rough estimate, you can use n*ln(n) as an approximation for the nth prime number. There is a much more complex, but more accurate method, details of which you can find on Wikipedia here.

Mide answered 25/6, 2009 at 9:16 Comment(0)
B
1

To complement Dana J's Upper bound this formula should give you a good lower bound.

P(n) = (((2 Log(3, n + 2))/(Log(2.5, 2) + Log(3, 3)) + (2 Log(3, n - 2))/(Log(3, 2) + Log(3, 3)))/2) n;
Billow answered 21/12, 2017 at 6:0 Comment(0)
P
0

An efficient implementation is probably not possible with a sieve. Think what would happen if you want to have the first 10.000 prime numbers. You probably would have to make a sieve over a huge bigger amount of numbers.

Your own implentation in this question and my answer are good ways to implement this without knowing the appr. value of a prime

Pye answered 25/6, 2009 at 10:16 Comment(0)
R
0

My Best Prime(n) Estimate

1/2*(8-8.7*n-n^2+
1/2*(2*abs(log(n)/log(3)+log(log(n)/log(2))/log(2))+
abs((log(log(3))-log(log(n))+2*n*log(log(n)/log(2))+
sqrt(((8*log(3)*log(n))/log(2)-log(log(2))+
log(log(n)))*log(log(n)/log(2))))/log(log(n)/log(2))))*(-1+
abs(log(n)/log(3)+log(log(n)/log(2))/log(2))+abs(-(1/2)+n+
sqrt(((8*log(3)*log(n))/log(2)-log(log(2))+
log(log(n)))*log(log(n)/log(2)))/(2*log(log(n)/log(2))))))

Here's my most recent more experimental formula. btw. The ten trillionth prime is 323,780,508,946,331 this formula works quite well at that scale not sure if it continues to get closer than n*ln(n)+n*(ln(ln(n))-0.9385).

1/2*(3-(8+ln(2.3))*n-n^2+1/2*(-1+
abs(-(1/2)+n+sqrt(ln(ln(n)/ln(2))*(-ln(ln(2))+ln(ln(n))+
(8*ln(3)*ln((n*ln(8*n))/ln(n)))/ln(2)))/(2*ln(ln((n*ln(8*n))/
ln(n))/ln(2))))+abs(ln(n)/ln(3)+ln(ln((n*ln(8*n))/ln(n))/ln(2))/
ln(2)))*(2*abs(ln((n*ln(8*n))/ln(n))/ln(3)+ln(ln((n*ln(8*n))/ln(n))/
ln(2))/ln(2))+abs(1/ln(ln(n)/ln(2))*(ln(ln(3))-ln(ln(n))+2*n*ln(ln(n)/
ln(2))+sqrt(((8*ln(3)*ln(n))/ln(2)-ln(ln(2))+ln(ln((n*ln(8*n))/ln(n))))*
ln(ln((n*ln(8*n))/ln(n))/ln(2)))))))
Rowboat answered 28/2, 2012 at 18:49 Comment(1)
Neither of these work. I give it 1000, getting a negative number.Swelling

© 2022 - 2024 — McMap. All rights reserved.