Prime Number Determination Using Wheel Factorization

• Download source code - 26.47 KB
• Download demo - 12.42 KB

Introduction

For a positive integer n, we have to determine if n is prime, where n is small (i.e., it can be assigned to a ulong type in C#). Thus, n will be constrained to be less than 2^64-1. The algorithm presented here will make some notable improvements over the brute force method by using Wheel Factorization.

Background

Prime numbers are interesting numbers to the math community. And, the search for prime numbers seems to be the hobby (or profession) of many mathematicians. Prime numbers are also important in the RSA encryption algorithm. There is even a prize money for finding very large primes.

The Algorithm

Our goal is to determine if a positive integer is prime or not. The brute force way would look something like this:

// Warning: Slow code ahead, do not use.
// The brute force way of determining if a number is prime.
public bool IsPrime(ulong primeSuspect)
{
    if (primeSuspect < 2) return false;

    if (primeSuspect == 2) return true;

    for (ulong divisor = 2; divisor < primeSuspect; divisor++)
    {
        // If no remainder after dividing by the divisor
        if (primeSuspect % divisor == 0) 
        {
           return false;  // It is not prime
        }  
    }
    // If we did not find a divisor, it is prime
    return true; 
}

The above code is inefficient because it often checks a number as a possible divisor that has already been eliminated because it is a multiple of a smaller number. No need to check 4 as a divisor if 2 was already checked. And, no need to check 9 as a divisor if 3 has already been checked, and so on. Another problem is that it checks all numbers up to our prime candidate as a divisor. We will fix this first.

So, our first improvement to the algorithm will be to only search for factors for our prime candidate up to the square root of our prime candidate. If there is a factor to our candidate greater than its square root, then there must also be a factor less than its square root such that when these two factors are multiplied, it gives us our original candidate.

Our next improvement will be to use a method called Wheel Factorization. If we already know all the primes less than or equal to the square root of our candidate, this would be optimal. However, searching for all these lesser primes comes at the price of more computation time, and negates much of our gains. So, I'll use Wheel Factorization to speed up the search.

In Wheel Factorization, you start with the first few primes. In this example, I will use 2, 3, and 5, the first three primes, to make it simple. (In the downloaded code, I use more than the first three primes, which will give us some more improvement.) This gives us a Wheel Factorization of 30, the product of the first three primes (2*3*5). You then make a list of the integers from 1 to 30, and eliminate all the numbers in the list that are multiples of 2, 3, or 5. This gives us this list: {1, 7, 11, 13, 17, 19, 23, 29} of sieved numbers. These sieved numbers give us a pattern of numbers that repeat and are not multiples of 2, 3, or 5. Thus, if you add 30, or 60, or 90, etc. to each of these numbers in the list, none are divisible by 2, 3, or 5. I will make a small modification to this sieved list of numbers to make the loop simpler. I will remove the 1, and add it to 30 at the tail of the list. So now, our list of numbers is {7, 11, 13, 17, 19, 23, 29, 31}. This is so I can do a pass = 0 and don't have to divide by 1.

So, here is the heart of the program (simplified for this article by using just the first three primes to create the sieve):

private static ulong[] aSieve30 = new ulong[]
         {7, 11, 13, 17, 19, 23, 29, 31};

// Find the first divisor greater than 1 of our candidatePrime.
public static ulong FirstDivisor(ulong candidatePrime)
{
    if (candidatePrime == 0)
        throw new ArgumentException ("Zero is an invalid parameter!");

    // A List of the first three primes
    List<ulong> firstPrimes = 
           new List<ulong>(new ulong[] {2, 3, 5}); 

    WheelFactor = 30;  // The product of the primes in firstPrimes

    if (candidatePrime == 1)
    {  
        // 1 is not considered a prime or a composite number.
        return 0; // So return any number other than 1.
    }
    foreach (ulong prime in firstPrimes)
    {
        if (candidatePrime % prime == 0) return prime;
    }

    // No need to search beyond the square root for divisors
    ulong theSqrt = (ulong)Math.Sqrt((double)candidatePrime); 

    for (ulong pass = 0; pass < theSqrt; pass += WheelFactor)
    {
        foreach (ulong sieve in aSieve30)
        {
            if (candidatePrime % (pass + sieve) == 0)
            {
                  return pass + sieve;
            }
        }
    }
    // If we got this far our number is a prime
    return candidatePrime;
}

public static bool IsPrime(ulong primeSuspect)
{
   if (primeSuspect == 0) return false;
   return (FirstDivisor(primeSuspect) == primeSuspect);
}

As you can see from the for loop above, it is incremented by our WheelFactor, which is 30 in this case. And, the inner loop checks all 8 primes in our sieved list. Therefore, just 8 of 30 numbers are checked, which is more than a 73% improvement over the brute force way. Unfortunately, increasing the number of primes in our first primes list has a diminishing improvement on our search. The downloaded code uses the first 8 primes, which gives about an 83% improvement over the brute force way.

The figure illustrates a wheel factorization of 30. After performing trial divisions of 2, 3, and 5, then you only have to do trial divisions for those spokes of the wheel that are white. The spokes of the wheel that are red have been eliminated for consideration as possible divisors.

Conclusion

In Wheel Factorization, you get some good performance improvement in determining if a number is prime by skipping all the multiples of the first few primes.

Reference

• The University of Tennessee at Martin's Primes Pages

History

• 19 Nov 2008 - Original submission
• 31 Jan 2009 - The BuildSieve() method was improved