Sieve of Eratosthenes: making it quicker

Question

I was thinking about doing problem 23 of Project Euler. It includes the difficulty where you have to factorize primes. I had already done that in problems I solved earlier, but it was only necessary to have a prime list until a relative small number.

So I had already implemented the algorithm when I took a look at this problem. I stored numbers in a list, flagged the ones that weren't prime numbers and removed them from the list, al of this in Python. It worked well for small lists, but when I went big, the memory usage was just too much.

I still wanted to use the algorithm of Eratosthenes, but I still needed a solution to my memory issue. So I decided, I don't know if it is the best solution, to write the list to a file and to constantly rewrite this file.

I succeeded in making this algorithm, but it's not really fast. To find primes up to the number 20 000 it took my program 405 seconds, or 6'45".

Is there any way I can improve this implementation? No file rewriting, cashing data in the RAM, etc.?

SieveCash.py (Python 3)

import tempfile
import os
import sys

def createSieve(outputFile, maxp, logging=0):
    """Creates a sieve and stores one prime per line in file outputFile."""
    # Write a list of all numbers to check to a file
    myNumFile = createBeginNumberList(maxp, outputFile)

    # Do the looping
    rfh = open(myNumFile, "r")
    try:
        totalRead = 0
        for i in range(maxp+1):
            # log if necessary
            if logging > 0 and i % logging == 0:
                print("Checking number", i, "(" + str(int(i/(maxp+1)*100)) + "%)")

            # do thingy
            bytesRead, myNum = readNumber(rfh)
            if bytesRead == 0:
                break

            totalRead += bytesRead

            if myNum != 0:
                rfh.close()
                rfh = None

                rewriteFileWithoutMultiples(myNumFile, myNum, totalRead)

                rfh = open(myNumFile, "r")
                rfh.seek(totalRead)

    finally:
        if rfh != None:
            rfh.close()

    # Remove zeros
    removeZerosInFile(outputFile)

    return myNumFile


def createBeginNumberList(maxn, fileName):
    """@return value: the file path"""

    tfh = open(fileName, "w")
    try:
        writeLineToFile(tfh, "0")
        writeLineToFile(tfh, "0")
        for i in range(2, maxn+1):
            writeLineToFile(tfh, str(i))
    finally:
        tfh.close()

    return fileName

def rewriteFileWithoutMultiples(inputFilename, multiple, startPos):
    rfh = open(inputFilename, "r")
    try:
        writeFile = tempfile.mkstemp()
        os.close(writeFile[0])
        wfh = open(writeFile[1], "w")
    except:
        rfh.close()
        raise

    # Rewrite file
    try:
        readBytes = 0
        while readBytes < startPos:
            b, aNum = readNumber(rfh)
            readBytes += b
            writeLineToFile(wfh, str(aNum))

        while True:
            b, aNum = readNumber(rfh)
            if aNum == None:
                break

            if aNum % multiple == 0:
                writeLineToFile(wfh, "0")
            else:
                writeLineToFile(wfh, str(aNum))
    finally:
        rfh.close()
        wfh.close()

    # Copy
    copyFile(writeFile[1], inputFilename)

def copyFile(source, dest):
    rfh = open(source, "r")
    try:
        wfh = open(dest, "w")
    except:
        rfh.close()
        raise

    try:
        aLine = rfh.readline()
        while len(aLine) != 0:
            wfh.write(aLine)
            aLine = rfh.readline()
    finally:
        rfh.close()
        wfh.close()

def removeZerosInFile(inputFile):
    rfh = open(inputFile, "r")
    try:
        writeFile = tempfile.mkstemp()
        os.close(writeFile[0])
        wfh = open(writeFile[1], "w")
    except:
        rfh.close()
        raise

    # Rewrite file
    try:
        while True:
            b, aNum = readNumber(rfh)
            if aNum == None:
                break

            if aNum != 0:
                writeLineToFile(wfh, str(aNum))

    finally:
        rfh.close()
        wfh.close()

    # Copy
    copyFile(writeFile[1], inputFile)



def writeLineToFile(fh, l):
    try:
        fh.write(l)
        fh.write('\n')
    except:
        fh.close()
        raise

def readNumber(fh):
    mynum = fh.readline()
    numb = len(mynum)
    mynum = mynum.strip('\n')

    if len(mynum) == 0:
        return (0, None)

    return (numb, int(mynum))

if __name__ == "__main__":
    if len(sys.argv) < 2:
        print("Usage:", sys.argv[0], "MAX_PRIME [LOG_INTERVAL]")
        sys.exit(1)

    maxp = int(sys.argv[1])
    logging = int(maxp / 100)
    try:
        logging = int(sys.argv[2])
    except:
        pass

    print("Creating seive until number", maxp, "with logging interval", logging)
    createSieve(os.path.expanduser("~/Desktop/Sieve-" + str(maxp) + ".txt"),
                maxp, logging)

Program Output (On an iMac 2.16 GHz Intel Core 2 Duo, 1GB RAM)

imac-van-ief2:23 ief2$ SECONDS=0; python3 SieveCash.py 20000 2000; echo $SECONDS
Creating seive until number 20000 with logging interval 2000
Checking number 0 (0%)
Checking number 2000 (9%)
Checking number 4000 (19%)
Checking number 6000 (29%)
Checking number 8000 (39%)
Checking number 10000 (49%)
Checking number 12000 (59%)
Checking number 14000 (69%)
Checking number 16000 (79%)
Checking number 18000 (89%)
Checking number 20000 (99%)
405

Answer 1

Writing to disk is going to be really slow. Don't do that.

If you are running out of memory, try using the standard array module or the numpy library. Both of those provide efficient arrays for holding large numbers of values. Your only running out of memory because python's lists don't store their values in a memory efficient manner.

If you are running out of memory for some reason even taking that into account, post that code here and we'll be glad to point out how to improve the memory usage.

Answer 2

If space is the problem, a bit-array may help. Straightforward and unoptimised (I don't know Python so well),

import sys
import array

def sieve(n):
    sieveBits = (n-1) // 2
    sieveInts = (sieveBits+31) // 32
    sieveBound = int(n**0.5) // 2
    arr = array.array('I')
    arr.extend((0,)*sieveInts)
    for i in xrange(sieveBound):
        if (arr[i >> 5] & (1 << (i&31))) == 0:
            for j in xrange(2*i*(i+3)+3,sieveBits,2*i+3):
                arr[j >> 5] |= 1 << (j&31)
    return arr


def primes(n):
    arr = sieve(n)
    primes = [2] + [2*i+3 for i in xrange((n-1)//2) if arr[i >> 5] & (1 << (i & 31)) == 0]
    return primes

if __name__ == "__main__":
    if len(sys.argv) > 1:
        up = int(sys.argv[1])
    else:
        up = 100000
    print len(primes(up))

It's dog-slow (4.6 seconds to sieve up to 10,000,000), but the sieving doesn't take much memory.

Answer 3

Here is my version of the sieve. Ideone finds the 2262 primes less than 20000 in 0.02 seconds.

def sieve(n):
    m = (n-1) // 2
    b = [True]*m
    i,p,ps = 0,3,[2]
    while p*p < n:
        if b[i]:
            ps.append(p)
            j = 2*i*i + 6*i + 3
            while j < m:
                b[j] = False
                j = j + 2*i + 3
        i+=1; p+=2
    while i < m:
        if b[i]:
            ps.append(p)
        i+=1; p+=2
    return ps

Answer 4

One minor improvement is to use xrange instead of range. In your case, you don't need range, and xrange has better memory performance, which might make a little difference.