nthprime

Finds either the n'th prime, or compute the number of primes that are less than some given value.

I have absolutely no reason for having written this code, except that it was interesting to write. In an idle moment, I wondered how one might efficiently compute the n'th prime from the prime sequence, given the index N.

Yes, one might use the primes function, which is reasonably fast for small enough sets of primes. But primes by itself does not solve this problem. Primes returns the entire list of primes less than or equal to some value. So you might call primes, only to ind that you had generated too few primes to get the specific prime number that you wanted. Worse, suppose you wanted to find P(1e8)? It is very inefficient to generate the entire list of 100,000,000 primes, only needing the last element of that list.

A closely related question is to ask how many primes are less than a given value. We can get this by numel(primes(K)), but if the number K is very large, perhaps as large as 2^32, the call to primes will take far too long to execute.

The nthprime function solves both of these problems efficiently. For example, what is P(12345678)?

nthprime(123456)
ans =
1632899

See that it is prime.

isprime(1632899)
ans =
1

Also, we can verify that it is the 123456'th prime in the complete sequence of prime numbers.

numel(primes(1632899))
ans =
123456

Find the primes numbered [1 10 100 1000 10000 100000 1000000 10000000 100000000]

p = nthprime(10.^(0:8))
p =
2 29 541 7919 104729 1299709 15485863 179424673 2038074743

nthprime is fairly efficient. For example, there are 7603553 primes less than 2^27, but computing the entire list of those primes can be a large task. If for some reason, you only wanted to know the last prime on that list, the cost of calling nthprime is far less than the time to compute the entire list.

tic,p = primes(2^27);toc
Elapsed time is 4.517565 seconds.

tic,plast = nthprime(7603553);toc
Elapsed time is 0.003868 seconds.

Finally, nthprime can work with primes as high as 2^32, whereas the primes function runs out of steam far below that number. There are roughly 2e8 primes below 2^32. To be exact, we can use nthprime to tell us how many there are:

nthprime(2^32,1)
ans =
203280221

Thus the second argument allows us to specify the operation of nthprime. The call nthprime(K,1) is the same as numel(primes(K)). But it is far more efficient for large values of K.

Again, I have absolutely no target for this code, no reason to have written it, except for working out how I might solve the problem itself. The code was written using a relatively small database of primes to localize the problem, then it applies a simple prime sieve on top of that.