1
\$\begingroup\$

I recently played around with prime factorization using a prime wheel. Much of the heavy lifting is done by the PrimeWheel class. I think the implementation is nice to anyone who is interested in wheels but a little confused on how to use them. Note the PrimeWheel class could also be used for a fast but simple IsPrime primality test by modifying GetPrimeDivisors to return a false the instant the divisor is found.

PrimeWheel Class

public class PrimeWheel
{
    public PrimeWheel()
    {
        Product = GetProduct();
    }

    // You may play around by adding 13 and maybe 17 but performance may actually degrade
    // because it will take longer for GetSieveFlags to run.  For 32 bit integers,
    // I find 11 to be in the sweet spot as the last prime in the wheel.
    // If last base prime = 11, then product = 2310 and possible remainders = 339.
    // If last base prime = 13, then product = 30030 and possible remainders = 3243.
    // Performance example for first 100K numbers: 
    //      If last is 11, it took 1.3 seconds (Release x64) but 
    //      if last is 13, that number jumped to 7.9 seconds.
    public static IList<int> BasePrimes => new[] { 2, 3, 5, 7, 11 };
    public int Product { get; } = 1;

    private int GetProduct()
    {
        var product = 1;
        foreach (var prime in BasePrimes)
        {
            product *= prime;
        }
        return product;
    }

    private bool[] SieveFlags = null;

    private bool[] GetSieveFlags()
    {
        var flags = new bool[Product + 1];
        flags[0] = false;
        flags[1] = false;
        for (var i = 2; i < flags.Length; i++)
        {
            // Assume flags are true and turn them off below.
            flags[i] = true;
        }

        // If last base prime = 11, then product = 2310 and sqroot = 48.
        var sqroot = (int)Math.Sqrt(Product);

        // Output primes up to the square root of wheel.Product
        for (var i = 2; i <= sqroot; i++)
        {
            // If this bit has already been turned off, then its associated number is composite. 
            if (!flags[i]) { continue; }
            // Any multiples of number are composite and their respective bits should be turned off.
            for (var j = i * i; j < flags.Length; j += i)
            {
                flags[j] = false;
            }
        }

        return flags;
    }

    public IEnumerable<int> Primes
    {
        get
        {
            if (SieveFlags == null)
            {
                SieveFlags = GetSieveFlags();
            }
            for (var i = 2; i < SieveFlags.Length; i++)
            {
                if (SieveFlags[i])
                {
                    yield return i;
                }
            }
        }
    }

    public IEnumerable<int> PrimesLargerThanBase
    {
        get
        {
            if (SieveFlags == null)
            {
                SieveFlags = GetSieveFlags();
            }
            for (var i = BasePrimes.Last() + 1; i < SieveFlags.Length; i++)
            {
                if (SieveFlags[i])
                {
                    yield return i;
                }
            }
        }
    }

    public IEnumerable<int> PossibleRemainders
    {
        get
        {
            // Builds a list of possible remainders when rotating thru the wheel.
            // This eliminates any remainders ABSOLUTELY known to be composite.
            yield return 1;  // 1 is always a possible remainder 
            foreach (var prime in PrimesLargerThanBase)
            {
                yield return prime;
            }
        }
    }
}

A static GetPrimeDivisors method relies upon the PrimeWheel. I originally used the ValueTuple package to return IEnumerable<(int Prime, int Power)>. As the method is very fast and returns at most 9 items (2*3*5*7*11*17*19*23*29), I changed the output type to be IDictionary<int, int> where the Key`` is the prime divisor and theValue` is the power that the divisor is raised to.

I could have combined the first and second passes below, since creating the underlying sieve is quite fast. But I chose to break it up to defer creating the sieve until it is absolutely needed.

GetPrimeDivisors method

public static class WheelFactorization
{
    public static IDictionary<int, int> GetPrimeDivisors(int number)
    {
        if (number < 2)
        {
            throw new ArgumentException("Input value must be > 1.", nameof(number));
        }

        var divisors = new Dictionary<int, int>();

        // Local Action delegate may reduce 'number'.
        Action<int> checkDivisor = (prime) =>
        {
            var power = 0;
            while (number % prime == 0)
            {
                ++power;
                number /= prime;
            }
            if (power > 0)
            {
                divisors.Add(prime, power);
            }
        };


        // FIRST PASS:  Find divisors within PrimeWheel.BasePrimes.
        // This is a very fast check before the need to construct a wheel instance. 
        // This pass may be considered a "quarter turn" of the first rotation of the wheel.
        foreach (var prime in PrimeWheel.BasePrimes)
        {
            checkDivisor(prime);
            if (number == 1)
            {
                return divisors;
            }
        }

        // Much of the heavy lifting is done inside the Wheel class.
        var wheel = new PrimeWheel();

        // SECOND PASS:  Find other prime divisors <= wheel.Product.
        // This relies upon building a Sieve of Eratosthenes internal to the wheel.
        // This pass completes the first full rotation of the wheel.
        foreach (var prime in wheel.PrimesLargerThanBase)
        {
            checkDivisor(prime);
            if (number == 1)
            {
                return divisors;
            }
        }

        // THIRD PASS:  Find prime divisors > wheel.Product.
        // If the last base prime is 11, the Product is 2310 and there will be 309 possible remainders
        // to be checked.  Of the 309 remainders, 308 are primes ranging from 13 to 2309.
        // This means this Third Pass will not be reached unless number is composed of at least 1 prime
        // greater than 2309 (including number itself if number is a prime).
        // Here the wheel will continue to rotate until it can't turn anymore.
        var sqroot = (int)Math.Sqrt(number);
        var remainders = wheel.PossibleRemainders.ToList();
        for (var factor = wheel.Product; factor <= sqroot; factor += wheel.Product)
        {
            foreach (var remainder in remainders)
            {
                var candidate = factor + remainder;
                checkDivisor(candidate);
                if (number == 1)
                {
                    return divisors;
                }
            }
        }

        // If this point is reached, number is a prime and should output itself.
        divisors.Add(number, 1);
        return divisors;
    }
}

Simple Tests

Here is the code for a brief test I put together:

TriggerJitter();  

TestAssorted();

const int rangeSize = 250000;
var sizeK = rangeSize / 1000.0;
TestRange($"Test First {sizeK}K - Show Aggregate Time", 2, rangeSize);
TestRange($"Test Last {sizeK}K - Show Aggregate Time", int.MaxValue - rangeSize, rangeSize);

Methods called by tests

private static void TriggerJitter()
{
    // For all the associated methods to be fully processed by Jitter,
    // we need to pass a number that will force Wheel.PossibleRemainders to execute.
    var number = int.MaxValue;
    var elapsed = new double[2];
    for (var i = 0; i < elapsed.Length; i++)
    {
        var sw = Stopwatch.StartNew();
        var divisors = WheelFactorization.GetPrimeDivisors(number);
        sw.Stop();
        elapsed[i] = sw.Elapsed.TotalMilliseconds;
    }
    Console.WriteLine();
    Console.WriteLine($"Jitter triggered for {number}");
    Console.WriteLine($"    First  Run: {elapsed[0]} milliseconds");
    Console.WriteLine($"    Second Run: {elapsed[1]}");
    Console.WriteLine($"    Difference: {(elapsed[0] - elapsed[1])}");
    Console.WriteLine();
}

private static void TestAssorted()
{
    Console.WriteLine("Test Assorted - Show Individual Times for assorted numbers");

    var numbers = new int[]
    {
        2,
        3,
        Magic.ProductFirst5Primes - 1,
        Magic.ProductFirst5Primes,
        Magic.ProductFirst5Primes + 1,
        (int)Math.Pow(2, 22),
        Magic.ProductSecond5Primes,
        Magic.Largest15BitPrime,
        Magic.Largest15BitPrime * Magic.Largest15BitPrime,
        Magic.Largest31BitPrime
    };

    foreach (var number in numbers)
    {
        var sw = Stopwatch.StartNew();
        var divisors = WheelFactorization.GetPrimeDivisors(number);
        sw.Stop();
        var indent = "    ";
        Console.WriteLine($"{indent}Prime Divisors for {number}");
        indent += indent;
        foreach (var factor in divisors)
        {
            Console.WriteLine($"{indent}{factor.Key} raised to {factor.Value}");
        }
        Console.WriteLine($"{indent}Milliseconds: {sw.Elapsed.TotalMilliseconds}");
    }
    Console.WriteLine();
}

private static void TestRange(string title, int first, int rangeSize)
{
    // For this test we don't display the divisors returned but do show the total elapsed milliseconds.
    Console.WriteLine(title);
    var last = first + rangeSize;
    var sw = Stopwatch.StartNew();
    for (var number = first; number <= last; number++)
    {
        // Break needed if last == int.MaxValue because number++ wraps to negatives.
        if (number < 2) { break; }
        var divisors = WheelFactorization.GetPrimeDivisors(number);
    }
    sw.Stop();
    Console.WriteLine($"    Milliseconds: {sw.Elapsed.TotalMilliseconds}");
    Console.WriteLine($"    Average: {(sw.Elapsed.TotalMilliseconds) / rangeSize}");
    Console.WriteLine();
}

Magic number class

public static class Magic
{
    // For more about Magic Numbers: https://msdn.microsoft.com/en-us/library/ee621251.aspx
    public const int Largest15BitPrime = 32749;
    public const int Largest31BitPrime = int.MaxValue;

    public const int ProductFirst5Primes = 2 * 3 * 5 * 7 * 11;
    public const int ProductSecond5Primes = 13 * 17 * 19 * 23 * 29;
}

Sample Console Output

Jitter triggered for 2147483647
    First  Run: 6.9883 milliseconds
    Second Run: 0.0877
    Difference: 6.9006

Test Assorted - Show Individual Times for assorted numbers
    Prime Divisors for 2
        2 raised to 1
        Milliseconds: 0.0012
    Prime Divisors for 3
        3 raised to 1
        Milliseconds: 0.0009
    Prime Divisors for 2309
        2309 raised to 1
        Milliseconds: 0.0153
    Prime Divisors for 2310
        2 raised to 1
        3 raised to 1
        5 raised to 1
        7 raised to 1
        11 raised to 1
        Milliseconds: 0.0063
    Prime Divisors for 2311
        2311 raised to 1
        Milliseconds: 0.0327
    Prime Divisors for 4194304
        2 raised to 22
        Milliseconds: 0.0009
    Prime Divisors for 2800733
        13 raised to 1
        17 raised to 1
        19 raised to 1
        23 raised to 1
        29 raised to 1
        Milliseconds: 0.0042
    Prime Divisors for 32749
        32749 raised to 1
        Milliseconds: 0.0294
    Prime Divisors for 1072497001
        32749 raised to 2
        Milliseconds: 0.0489
    Prime Divisors for 2147483647
        2147483647 raised to 1
        Milliseconds: 0.0895

Test First 250K - Show Aggregate Time
    Milliseconds: 3807.9156
    Average: 0.0152316624

Test Last 250K - Show Aggregate Time
    Milliseconds: 7898.3017
    Average: 0.0315932068

Concerns/Not Concerns

I have no concerns about the sieve. I have many different versions - some that use BitArrays, track odd numbers only, or have an offset. I am quite familiar with various things I could try. For a tiny sieve of 2311 Booleans, I expressly did not want to transform a number into an index or vice versa. Thus, I will be quite resistant to altering the sieve just to be different.

An early version of the PrimeWheel constructor accepted a seed integer. This made sure the Product was never greater than the input seed, which meant different instances of the wheel could vary with different last known primes. The version shown above does not vary between instances in that Product is always 2130. This means it could that PrimeWheel could be a static class, but the code feels cleaner without PrimeWheel being a static class. But I am open to other opinions on this.

Given both its quickness and a max of 9 items, I don’t see the need for GetPrimeDivisors to be an IEnumerable. Yet I’m not totally committed to dictionary either. If I changed to IEnumerable or IList, I would probably switch back to a tuple because I find IEnumerable<KeyValuePair<int, int>> to be ugly.

Some purists may argue that a property getter should not generate data but only echo data already stored in private fields. In general,I tend to agree with that philosophy but make room for exceptions. Since the entire GetPrimeDivisors runs in 1/10th of a 1 millisecond, and the PrimeWheel.Primes property is a tiny fraction of that, I am willing to make an exception that GetSieveFlags can be called in a property getter.

\$\endgroup\$
  • \$\begingroup\$ Could you add the class definition for the GetPrimeDivisors method? I'd be easier to copy/paste it ;-) \$\endgroup\$ – t3chb0t Jun 9 '18 at 22:50
  • \$\begingroup\$ @t3chb0t Added as requested, but not much to it as you can see. \$\endgroup\$ – Rick Davin Jun 9 '18 at 23:49
  • 1
    \$\begingroup\$ For the downvoter(s), is there anything about the code or question that can be improved? \$\endgroup\$ – Rick Davin Jun 11 '18 at 11:12
  • \$\begingroup\$ I'd like to know that too. I find this is a very interesting question and I was the first one to give you +1 for it ;-) \$\endgroup\$ – t3chb0t Jun 11 '18 at 11:14

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.