the 10001st prime number

Question

Problem description:

By listing the first six prime numbers: 2, 3, 5, 7, 11, and 13, we can see that the 6th prime is 13.

What is the 10001st prime number?

Prime number:

A prime number is a whole number greater than 1 whose only factors are 1 and itself. A factor is a whole number that can be divided evenly into another number.

Example:

The first 25 prime numbers (all the prime numbers less than 100) are:

2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97

The best solution is with Sieve Of Eratosthenes

Sieve Of Eratosthenes:

In mathematics, the sieve of Eratosthenes is an ancient algorithm for finding all prime numbers up to any given limit.

It does so by iteratively marking as composite (i.e., not prime) the multiples of each prime, starting with the first prime number, 2.

The multiples of a given prime are generated as a sequence of numbers starting from that prime, with constant difference between them that is equal to that prime.

This is the sieve's key distinction from using trial division to sequentially test each candidate number for divisibility by each prime.

Once all the multiples of each discovered prime have been marked as composites, the remaining unmarked numbers are primes.

The sieve of Eratosthenes is one of the most efficient ways to find all primes smaller than n when n is smaller than 10 million or so.

Sieve Of Eratosthenes(Step by Step):

1- Create a list of consecutive integers from 2 through n: (2, 3, 4, ..., n).

2- Initially, let p equal 2, the smallest prime number.

3- Enumerate the multiples of p by counting in increments of p from 2p to n, and mark them in the list (these will be 2p, 3p, 4p, ...; the p itself should not be marked).

4- Find the smallest number in the list greater than p that is not marked. If there was no such number, stop. Otherwise, let p now equal this new number (which is the next prime), and repeat from step 3.

5- When the algorithm terminates, the numbers remaining not marked in the list are all the primes below n.

My Solution

This is my solution for problem 7 of Project Euler using Python:

 def sieve_of_eratosthenes(n): 
  '''method for finding all primes up to 
     a given natural n.
  '''
    is_prime = [True]*n
    is_prime[0] = False
    is_prime[1] = False
    
    for i in range(2,int(math.sqrt(n)+1)):  
        index = i*2
        while index < n:
            is_prime[index] = False
            index = index+i
    prime = []
    for i in range(n):
        if is_prime[i] == True:
            prime.append(i)
    return prime

How could my code be improved?

Answer 1

Initial improvements

Don't use int(math.sqrt(n)). Python's integers can be very large, and math.sqrt(n) needs to first convert the integer to a floating point value before the square-root is calculated, which can lead to inaccurate results when n exceeds \$2^{53}\$. While you are likely nowhere near that limit, the function math.isqrt(n) (introduced in Python 3.8) avoids the conversion to floating point, and is actually shorter to type.

Your while loop ...

        index = i*2
        while index < n:
            is_prime[index] = False
            index = index+i

... can be replaced with a for loop using Python's range object:

        for index in range(i * 2, n, i):
            is_prime[index] = False

This should be faster than your while loop, due to the looping condition and increment all handled internally in the range() object.

The final for loop used to collect primes could be replaced with list comprehension using the enumerate function:

    primes = [i for i, flag in enumerate(is_prime) if flag]

Functional improvements

2 is the only even prime number. It is worth while handling it separately, and only testing odd numbers in the sieve.

Since only odd primes are being considered in the sieve, you can double your "crossing out" increment to 2p.

Also, crossing out multiples 2p, 3p, 4p, ... is inefficient, since all odd multiples below \$p^2\$ will have been crossed out during prior steps.

As pointed out by vnp, you're missing "step 4". You don't need to cross out multiples of composite numbers, since they've already been crossed out.

    is_prime = [False, True] * (n // 2)
    is_prime[1] = False
    is_prime[2] = True

    for i in range(3, math.isqrt(n) + 1, 2):
        if is_prime[i]:
            for index in range(i * i, n, 2 * i):
                is_prime[index] = False

Tool changes

https://pypi.org/project/bitarray/

There is a module called bitarray which can be installed. Instead of using a list to store the is_prime flags, the bitarray can store the flags packed into individual bits of one long buffer. This is a phenomenal memory optimization.

Moreover, using bitarray further simplifies the crossing-out of multiples to a slice assignment operation:

import bitarray

...

    is_prime = bitarray.bitarray(n)
    is_prime.setall(False)
    is_prime[2] = True
    is_prime[3::2] = True

    ...

    for i in range(3, math.isqrt(n) + 1, 2):
        if is_prime[i]:
            is_prime[i*i::2*i] = False

    ...

Answer 2

The code does not implement Step 4. Before going into elimination loop, you have to test is_prime[index]. It it is False, don't bother.

Step 3 is suboptimal. Every multiple of \$p\$ less than \$p^2\$ has a factor less than \$p\$, and therefore was already marked as non-prime. It is safe to start elimination directly from \$p^2\$.

This observation also allows to remove a call to math.sqrt, which looks very out of place in this problem.

Answer 3

Your code is clear in intent and structure. Some suggestions on performance:

• only keep an index of odd numbers. This halfs memory footprint and looping.
• start applying the sieve (is_prime[index] = False) from index * index. All other multiples of index are already False (for example 7 * index is already sieved by 7.)
• the loop at the end producing the primes seems a bit clunky. 4 lines for returning a filter on an index/range. Could be something like: primes = [index for index in range(n) if is_prime[index]]

Finally: you could use type hinting… That seems to be the current standard in Python.