I’ve created a much cleaner, better designed version to my parallel sieve. I’ve implemented most of EBrown’s micro-optimizations but also revamped the code on my own (that is not in direct response to someone’s answer).
public class Sieve32FastV2
{
private static ArgumentException BadUpperLimitException(string paramName) => new ArgumentException($"{paramName} be must greater than or equal to 2.");
public static IEnumerable<uint> Primes(int upperLimit)
{
if (upperLimit < 2)
{
throw BadUpperLimitException(nameof(upperLimit));
}
return Primes((uint)upperLimit);
}
public static IEnumerable<uint> Primes(uint upperLimit)
{
if (upperLimit < 2)
{
throw BadUpperLimitException(nameof(upperLimit));
}
var instance = new Sieve32FastV2(upperLimit);
return instance.EnumeratePrimes();
}
private Sieve32FastV2(uint upperLimit)
{
_vectors = VectorList.Create(upperLimit);
}
private readonly VectorList _vectors = null;
// Performance Tweak:
// Favor "(number & 0x01U) == 0U" over "number % 2 == 0" or "number % 2U == 0U"
// Favor "x >> 1" over "x / 2"
// Favor "x << 1" over "x * 2" or "x + x"
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static bool IsEven(uint number) => (number & 0x01U) == 0;
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static bool IsEven(int index) => (index & 0x01) == 0;
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static int Halved(int value) => value >> 1;
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static uint Halved(uint value) => value >> 1;
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static int Doubled(int value) => value << 1;
private IEnumerable<uint> EnumeratePrimes()
{
if (_vectors.UpperLimit < 2) { yield break; }
yield return 2;
if (_vectors.UpperLimit == 2) { yield break; }
// I call _vectors[0] the rootVector not just because its the very first one, but also because
// it was intentionally created so that the index of upperLimit's square root is contained within _vectors[0].
var rootVector = _vectors[0];
// Performance Tweak: Create copy of value in local variable
var rootBitIndex = _vectors.SquareRootIndex;
// The number of times a bit in all BitArray(s) are accessed:
//
// UpperLimit = int.MaxValue => 3,315,151,693 times
// UpperLimit = uint.MaxValue => 6,701,709,402 times
for (var bitIndex = 0; bitIndex <= rootBitIndex; bitIndex++)
{
if (rootVector[bitIndex])
{
var prime = rootVector.ToNumber(bitIndex);
yield return prime;
// All multiples of prime - on all vectors - are composites and should be marked as such.
MarkCompositesInParallel(prime);
}
}
// Output remaining primes found beyond rootBitIndex
// Performance Tweak: Favor "for" over "foreach"
for (var vectorIndex = 0; vectorIndex < _vectors.Count; vectorIndex++)
{
var vector = _vectors[vectorIndex];
// Performance Tweak: Favor "<= stopIndex" over "< vector.BitLength". Calc stopIndex once before loop.
var stopIndex = vector.BitLength - 1;
// Not a performance tweak: startIndex is a separate variable just for readability.
var startIndex = (vectorIndex == 0) ? rootBitIndex + 1 : 0;
for (var bitIndex = startIndex; bitIndex <= stopIndex; bitIndex++)
{
if (vector[bitIndex]) { yield return vector.ToNumber(bitIndex); }
}
}
}
private void MarkCompositesInParallel(uint prime)
{
// Performance Tweak: Favor "Parallel.For" over "Parallel.ForEach"
Parallel.For(0, _vectors.Count, vectorIndex =>
{
var vector = _vectors[vectorIndex];
// Performance Tweak: Create copy local to lambda. Cast once to bitStep before the inner loop.
var bitStep = (int)prime;
var startIndex = 0;
if (vectorIndex == 0)
{
startIndex = vector.ToIndex(prime * prime);
}
else
{
var remainder = (int)(vector.StartingNumber % prime);
// If remainder is 0, then we will use startIndex's initial value of 0.
if (remainder != 0)
{
startIndex = bitStep - remainder;
// On the full number scale, every other multiple of prime is even and should be skipped
// over for the next multiple, which is an odd number.
if (IsEven(remainder))
{
startIndex += bitStep;
}
startIndex = Halved(startIndex);
}
}
// Performance Tweak: Favor "<= stopIndex" over "< vector.BitLength". Calc stopIndex once before loop.
var stopIndex = vector.BitLength - 1;
for (int bitIndex = startIndex; bitIndex <= stopIndex; bitIndex += bitStep)
{
vector[bitIndex] = false;
}
});
}
private class VectorList : List<Vector>
{
// May not slap you in the face, but these 2 properties are readonly.
public uint UpperLimit { get; }
public int SquareRootIndex { get; }
public static VectorList Create(uint upperLimit)
{
var instance = new VectorList(upperLimit);
instance.CreateVectors();
return instance;
}
private VectorList(uint upperLimit)
{
// Any upperLimit > 2 should be odd for working with VectorList and Vector(s).
if ((upperLimit > 2) && IsEven(upperLimit))
{
upperLimit--;
}
this.UpperLimit = upperLimit;
SquareRootIndex = ToFlatIndex((uint)Math.Sqrt(upperLimit));
}
private void CreateVectors()
{
var typicalBitLength = CalcTypicalBitLength();
var typicalNumberRange = (uint)Doubled(typicalBitLength);
var count = (UpperLimit / typicalNumberRange) + 1;
for (uint i = 0, endingNumber = 1; i <= count; i++)
{
if (endingNumber >= UpperLimit) { break; }
// The first vector may have to be longer to accomodate the index of UpperLimit's square root.
var length = (i == 0) ? GetSpecialFirstLength(typicalBitLength) : typicalBitLength;
var startingNumber = endingNumber + 2;
var vector = new Vector(startingNumber, length, UpperLimit);
this.Add(vector);
endingNumber = vector.EndingNumber;
}
}
private int CalcTypicalBitLength()
{
// This is called before we have created a _vectors[0].
var length = ToFlatIndex(UpperLimit) + 1;
// Small enough values will result in 1 vector
const uint smallNumberCutoff = 10000;
if (UpperLimit < smallNumberCutoff) { return length; }
// Divide length for later parallelization over many (but not too many) vectors.
const int tinyFactor = 8;
var maxVectorCount = tinyFactor * Environment.ProcessorCount;
length = (length / maxVectorCount) + 1;
return PaddedLength(length);
}
private int GetSpecialFirstLength(int length) => (SquareRootIndex < length) ? length : PaddedLength(SquareRootIndex + 1);
private static int ToFlatIndex(uint number) => (int)Halved(number - 3);
private static int PaddedLength(int length)
{
// BitArray internally uses 32 bit int[] so align upwards to a 32 bit boundary,
// i.e. pad the end of length (in bits) to consume a full 32 bit int.
var remainder = length % 32;
return (remainder == 0) ? length : length + 32 - remainder;
}
}
// A Vector is aware of its bits, length, starting number, and ending number.
private class Vector
{
private readonly BitArray _bits = null;
public Vector(uint startNumber, int length, uint upperLimit)
{
StartingNumber = startNumber;
var endNumber = startNumber + Doubled(length - 1);
// In this constructor, endNumber is a long that could be > uint.MaxValue,
// in which case we clamp the length appropriately.
if (endNumber > upperLimit)
{
length = ToIndex(upperLimit) + 1;
}
_bits = new BitArray(length, defaultValue: true);
}
public bool this[int index] { get { return _bits[index]; } set { _bits[index] = value; } }
public int BitLength => _bits.Length;
// May not slap you in the face, but this property is readonly.
public uint StartingNumber { get; }
public uint EndingNumber => ToNumber(_bits.Length - 1);
public Func<uint, int> ToIndex => (uint number) => (int)Halved(number - StartingNumber);
public Func<int, uint> ToNumber => (int bitIndex) => (uint)Doubled(bitIndex) + StartingNumber;
}
}
Changes from Version 1
MarkCompositesInParallel has been improved for better readability (more below). The input bitIndex was never used so it was dropped from the method signature. The new method has its own bitIndex that has no relation to the original.
I ditched the names prime32 and prime31. I want to distinguish between between a uint value that is Number-related and an int value that is Index-related. So prime31 was renamed bitStep since its really for stepping over the indices of the BitArray. With that change, prime32 is now simply prime.
Anywhere there is a section of code that is a performance tweak or micro-optimization, a comment is provided.
The class level field _upperLimit was omitted. Instead it references the UpperLimit property from the private VectorList class. VectorList also has a new property named SquareRootIndex.
Where appropriate, many properties and fields are now declared as readonly, such as _vectors. For VectorList this would be UpperLimit and SquareRootIndex. For Vector this would be _bits and StartingNumber.
ToIndex and ToNumber have been moved to the Vector class. As such, I no longer have to pass a StartingNumber to each method. This reduces confusion of why did I pass 3U instead of rootVector.StartingNumber? I really liked this design though it did lead to the following…
At the risk of being the tiniest bit WET, VectorList has a one-liner ToFlatIndex function delegate that is similar to ToIndex. This is both intentional and acceptable to me. There are 2 calls to ToFlatIndex that are called before the first vector is ever created.
Additional tweaked methods were added. Instead of being function delegates, I made them methods but flagged them with [MethodImpl(MethodImplOptions.AggressiveInlining)]. Note this requires a using System.Runtime.CompilerServices;.
Performance Tweaks and Micro-Optimizations
I implemented most of EBrown’s suggestions except for replacing - 1 with + -1.
// Performance Tweak:
// Favor "(number & 0x01U) == 0U" over "number % 2 == 0" or "number % 2U == 0U"
// Favor "x >> 1" over "x / 2"
// Favor "x << 1" over "x * 2" or "x + x"
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static bool IsEven(uint number) => (number & 0x01U) == 0;
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static bool IsEven(int index) => (index & 0x01) == 0;
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static int Halved(int value) => value >> 1;
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static uint Halved(uint value) => value >> 1;
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static int Doubled(int value) => value << 1;
The respective IsEven methods are referenced only once each. Why not use them directly inline once? Since they use micro-optimizations, I like having them grouped together and commented accordingly. Plus later the code is more readable.
if (IsEven(remainder)) { }
is more readable than
if ((remainder & 0x01U) == 0U) { }
Other tweaks include:
- Having local copies of values, particular in a lambda.
- Using
forrather thanforeach. Additionallyforprovides avectorIndex. - Using
<= stopIndexrather than< bitLength.
For an upper limit of uint.MaxValue, the aggregate performance of all tweaks averaged 3.56 seconds faster over 8 timed runs. In fact the slowest tweaked version was 2.71 seconds faster than the fastest non-tweaked version.
MarkCompositesInParallel
I went through 2 revisions of this. The first was choosing better names to provide context and intent. For instance, I renamed targetNumber to be firstOddMultiple since that’s what the target number really is meant to be.
DRAFT VERSION
private void MarkCompositesInParallel(uint prime)
{
// Performance Tweak: Favor "Parallel.For" over "Parallel.ForEach"
Parallel.For(0, _vectors.Count, vectorIndex =>
{
var vector = _vectors[vectorIndex];
// initialize firstOddMultiple with a default value - but validate it below.
var firstOddMultiple = vector.StartingNumber;
if (vectorIndex == 0)
{
firstOddMultiple = prime * prime;
}
else
{
var remainder = vector.StartingNumber % prime;
if (remainder != 0U)
{
// CAUTION: Could this code block overflow past uint.MaxValue?
firstOddMultiple = vector.StartingNumber + prime - remainder;
// On the full number scale, every other multiple of prime is even and should be skipped
// over for the next multiple, which is an odd number.
if (IsEven(remainder))
{
firstOddMultiple += prime;
}
}
}
// Performance Tweak: Favor "<= stopIndex" over "< vector.BitLength". Calc stopIndex once before loop.
var stopIndex = vector.BitLength - 1;
// Performance Tweak: Create copy local to lambda. Cast once to bitStep before the inner loop.
var bitStep = (int)prime;
for (int bitIndex = vector.ToIndex(firstOddMultiple); bitIndex <= stopIndex; bitIndex += bitStep)
{
vector[bitIndex] = false;
}
});
}
Did you notice the comment // CAUTION: Could this code block overflow past uint.MaxValue?
I know that prime * prime will not overflow but I did not work out any rigorous math based on number of processors and tinyFactor to conclusively say that
firstOddMultiple = vector.StartingNumber + prime - remainder;
Will not overflow. The final version below is just as fast as the above, but the risk of overflow is mitigated. Note that remainder is changed from uint to int.
SAFER VERSION
private void MarkCompositesInParallel(uint prime)
{
// Performance Tweak: Favor "Parallel.For" over "Parallel.ForEach"
Parallel.For(0, _vectors.Count, vectorIndex =>
{
var vector = _vectors[vectorIndex];
// Performance Tweak: Create copy local to lambda. Cast once to bitStep before the inner loop.
var bitStep = (int)prime;
var startIndex = 0;
if (vectorIndex == 0)
{
startIndex = vector.ToIndex(prime * prime);
}
else
{
var remainder = (int)(vector.StartingNumber % prime);
// If remainder is 0, then we use startIndex's initial value of 0.
if (remainder != 0)
{
startIndex = bitStep - remainder;
// On the full number scale, every other multiple of prime is even and should be skipped
// over for the next multiple, which is an odd number.
if (IsEven(remainder))
{
startIndex += bitStep;
}
startIndex = Halved(startIndex);
}
}
// Performance Tweak: Favor "<= stopIndex" over "< vector.BitLength". Calc stopIndex once before loop.
var stopIndex = vector.BitLength - 1;
for (int bitIndex = startIndex; bitIndex <= stopIndex; bitIndex += bitStep)
{
vector[bitIndex] = false;
}
});
}
Misc
At one time, the Vector class had
public bool IsRoot => StartingNumber == 3U;
But as I favor for over foreach this was not needed (except when I tested a non-tweaked version that did use foreach).
I still have only one defined BadUpperLimitException but when if it gets thrown it incorporates EBrown’s suggestion of using nameof(variable).
Overall the improvements did not help with performance but did lead to what I believe is a nicer design.
Linq.Paralleletc... even withunsafecode and pointers ;-) but I couldn't make it any faster ;-( \$\endgroup\$