Restore Pow2 benchmark data reuse

This commit is contained in:
Marek Lemańczyk committed 2025-10-18 12:22:08 +02:00
1 parent b2ed0d5832
commit c6ad53fffb
3 files changed
+260 -322

No files matched your search

@@ -4,6 +4,7 @@ using BenchmarkDotNet.Attributes;
using BenchmarkDotNet.Jobs;
using System.Numerics;
using PerfectNumbers.Core;
using PerfectNumbers.Core.Gpu;
namespace EvenPerfectBitScanner.Benchmarks;
@@ -20,16 +21,18 @@ namespace EvenPerfectBitScanner.Benchmarks;
/// </remarks>
[MemoryDiagnoser]
[SimpleJob(RuntimeMoniker.Net80, launchCount: 1, warmupCount: 1, iterationCount: 5)]
public class Pow2MontgomeryModBenchmarks
public class Pow2ModBenchmarks
{
private const int SampleCount = 256;
private const ulong OneShiftedLeft63 = 1UL << 63;
private const long OneShiftedLeft60 = 1L << 60;
private const ulong Pow2WindowFallbackThreshold = 32UL;
private MontgomeryDivisorData[] _smallDivisors = Array.Empty<MontgomeryDivisorData>();
private MontgomeryDivisorData[] _largeDivisors = Array.Empty<MontgomeryDivisorData>();
private MontgomeryDivisorData[] _veryLargeDivisors = Array.Empty<MontgomeryDivisorData>();
private UInt128[] _wideModuli = Array.Empty<UInt128>();
private UInt128[] _wideExponents = Array.Empty<UInt128>();
private ulong[] _smallExponents = Array.Empty<ulong>();
private ulong[] _largeExponents = Array.Empty<ulong>();
private ulong[] _veryLargeExponents = Array.Empty<ulong>();
@@ -58,6 +61,7 @@ public class Pow2MontgomeryModBenchmarks
_largeDivisors = cache.LargeDivisors;
_veryLargeDivisors = cache.VeryLargeDivisors;
_wideModuli = cache.WideModuli;
_wideExponents = cache.WideExponents;
_smallExponents = cache.SmallExponents;
_largeExponents = cache.LargeExponents;
_veryLargeExponents = cache.VeryLargeExponents;
@@ -70,8 +74,8 @@ public class Pow2MontgomeryModBenchmarks
/// <summary>
/// Baseline right-to-left Montgomery ladder; measured 36.36 μs on the small sample set and 117.32 μs on the large one.
/// </summary>
[Benchmark(Baseline = true)]
public ulong BaselineCpu()
[Benchmark]
public ulong MontgomeryWindowedCpu()
{
GetData(out var exponents, out var divisors, out _);
ulong checksum = 0UL;
@@ -88,7 +92,7 @@ public class Pow2MontgomeryModBenchmarks
/// Baseline right-to-left Montgomery ladder; measured ??.?? μs on the small sample set and ??.?? μs on the large one.
/// </summary>
[Benchmark]
public ulong BaselineGpu()
public ulong MontgomeryWindowedGpu()
{
GetData(out var exponents, out var divisors, out _);
ulong checksum = 0UL;
@@ -105,7 +109,7 @@ public class Pow2MontgomeryModBenchmarks
/// Uses a known divisor cycle to fold the exponent first, dropping runtimes to 26.04 μs on small inputs (0.72× baseline) and 9.53 μs on large inputs (0.08× baseline).
/// </summary>
[Benchmark]
public ulong BaselineWithKnownCycleCpu()
public ulong MontgomeryWindowedWithCycleCpu()
{
GetData(out var exponents, out var divisors, out var cycles);
ulong checksum = 0UL;
@@ -122,7 +126,7 @@ public class Pow2MontgomeryModBenchmarks
/// Uses a known divisor cycle to fold the exponent first, dropping runtimes to 26.04 μs on small inputs (0.72× baseline) and 9.53 μs on large inputs (0.08× baseline).
/// </summary>
[Benchmark]
public ulong BaselineWithKnownCycleGpu()
public ulong MontgomeryWindowedWithCycleGpu()
{
GetData(out var exponents, out var divisors, out var cycles);
ulong checksum = 0UL;
@@ -271,8 +275,8 @@ public class Pow2MontgomeryModBenchmarks
/// <summary>
/// Classic binary ladder without Montgomery reduction; runs 24.75 μs on small inputs (0.68× baseline) but collapses to 295.29 μs on large ones (2.52×).
/// </summary>
[Benchmark]
public ulong BinaryMod()
[Benchmark(Baseline = true)]
public ulong StandardPow2Mod()
{
GetData(out var exponents, out var divisors, out _);
ulong checksum = 0UL;
@@ -319,6 +323,38 @@ public class Pow2MontgomeryModBenchmarks
return checksum;
}
/// <summary>
/// Current UInt128-native windowed Montgomery ladder that avoids GpuUInt128 on CPU paths.
/// </summary>
[Benchmark]
public UInt128 MontgomeryWindowedUInt128()
{
UInt128 checksum = UInt128.Zero;
for (int i = 0; i < SampleCount; i++)
{
checksum ^= _wideExponents[i].Pow2MontgomeryModWindowed(_wideModuli[i]);
}
return checksum;
}
/// <summary>
/// Previous CPU shim that reused the GPU-oriented GpuUInt128 implementation.
/// </summary>
[Benchmark]
public UInt128 MontgomeryWindowedGpuShimUInt128()
{
UInt128 checksum = UInt128.Zero;
for (int i = 0; i < SampleCount; i++)
{
checksum ^= Pow2MontgomeryModWindowedGpuShim(_wideExponents[i], _wideModuli[i]);
}
return checksum;
}
private void GetData(out ulong[] exponents, out MontgomeryDivisorData[] divisors, out ulong[] cycles)
{
switch (Scale)
@@ -350,6 +386,7 @@ public class Pow2MontgomeryModBenchmarks
MontgomeryDivisorData[] largeDivisors,
MontgomeryDivisorData[] veryLargeDivisors,
UInt128[] wideModuli,
UInt128[] wideExponents,
ulong[] smallExponents,
ulong[] largeExponents,
ulong[] veryLargeExponents,
@@ -362,6 +399,7 @@ public class Pow2MontgomeryModBenchmarks
LargeDivisors = largeDivisors;
VeryLargeDivisors = veryLargeDivisors;
WideModuli = wideModuli;
WideExponents = wideExponents;
SmallExponents = smallExponents;
LargeExponents = largeExponents;
VeryLargeExponents = veryLargeExponents;
@@ -375,6 +413,7 @@ public class Pow2MontgomeryModBenchmarks
public MontgomeryDivisorData[] LargeDivisors { get; }
public MontgomeryDivisorData[] VeryLargeDivisors { get; }
public UInt128[] WideModuli { get; }
public UInt128[] WideExponents { get; }
public ulong[] SmallExponents { get; }
public ulong[] LargeExponents { get; }
public ulong[] VeryLargeExponents { get; }
@@ -394,11 +433,13 @@ public class Pow2MontgomeryModBenchmarks
{
Random random = new(13);
Random heuristicRandom = new(113);
Random wideRandom = new(19);
var smallDivisors = new MontgomeryDivisorData[SampleCount];
var largeDivisors = new MontgomeryDivisorData[SampleCount];
var veryLargeDivisors = new MontgomeryDivisorData[SampleCount];
var wideModuli = new UInt128[SampleCount];
var wideExponents = new UInt128[SampleCount];
var smallExponents = new ulong[SampleCount];
var largeExponents = new ulong[SampleCount];
var veryLargeExponents = new ulong[SampleCount];
@@ -427,8 +468,9 @@ public class Pow2MontgomeryModBenchmarks
veryLargeExponents[i] = NextVeryLargeExponent(heuristicRandom);
veryLargeCycles[i] = CalculateCycleLengthWithHeuristics(veryLargeModulus, ref previousPrimeOrder);
UInt128 wideModulus = NextWideOddModulus(heuristicRandom);
UInt128 wideModulus = NextWideOddModulus(wideRandom);
wideModuli[i] = wideModulus;
wideExponents[i] = NextWideExponent(wideRandom);
wideCycles[i] = CalculateCycleLengthWithHeuristics(wideModulus, ref previousWidePrimeOrder);
}
@@ -437,6 +479,7 @@ public class Pow2MontgomeryModBenchmarks
largeDivisors,
veryLargeDivisors,
wideModuli,
wideExponents,
smallExponents,
largeExponents,
veryLargeExponents,
@@ -483,10 +526,37 @@ public class Pow2MontgomeryModBenchmarks
{
Span<byte> buffer = stackalloc byte[16];
random.NextBytes(buffer);
ulong low = BinaryPrimitives.ReadUInt64LittleEndian(buffer);
ulong high = BinaryPrimitives.ReadUInt64LittleEndian(buffer[8..]) | (1UL << 63);
UInt128 value = ((UInt128)high << 64) | low;
return value | UInt128.One;
ulong high = BinaryPrimitives.ReadUInt64LittleEndian(buffer[8..]);
if (high == 0UL)
{
high = 1UL;
}
UInt128 modulus = ((UInt128)high << 64) | low;
modulus |= UInt128.One;
if (modulus <= UInt128.One)
{
modulus += 2UL;
}
return modulus;
}
private static UInt128 NextWideExponent(Random random)
{
Span<byte> buffer = stackalloc byte[16];
random.NextBytes(buffer);
ulong low = BinaryPrimitives.ReadUInt64LittleEndian(buffer);
ulong high = BinaryPrimitives.ReadUInt64LittleEndian(buffer[8..]);
if (high == 0UL)
{
high = 1UL;
}
return ((UInt128)high << 64) | low;
}
private static ulong CalculateCycleLengthWithHeuristics(ulong modulus, ref ulong? previousPrimeOrder)
@@ -851,6 +921,182 @@ public class Pow2MontgomeryModBenchmarks
return (ulong)result;
}
private static UInt128 Pow2MontgomeryModWindowedGpuShim(UInt128 exponent, UInt128 modulus)
{
if (modulus == UInt128.One)
{
return UInt128.Zero;
}
if (exponent == UInt128.Zero)
{
return UInt128.One % modulus;
}
GpuUInt128 modulusGpu = new(modulus);
GpuUInt128 exponentGpu = new(exponent);
GpuUInt128 baseValue = new(2UL);
if (baseValue.CompareTo(modulusGpu) >= 0)
{
baseValue.Sub(modulusGpu);
}
if (ShouldUseSingleBit(exponentGpu))
{
GpuUInt128 singleBitResult = Pow2MontgomeryModSingleBit(exponentGpu, modulusGpu, baseValue);
return (UInt128)singleBitResult;
}
int bitLength = exponentGpu.GetBitLength();
int windowSize = GetWindowSize(bitLength);
int oddPowerCount = 1 << (windowSize - 1);
Span<GpuUInt128> oddPowers = stackalloc GpuUInt128[oddPowerCount];
InitializeOddPowers(baseValue, modulusGpu, oddPowers);
GpuUInt128 result = GpuUInt128.One;
int index = bitLength - 1;
while (index >= 0)
{
if (!IsBitSet(exponentGpu, index))
{
result.MulMod(result, modulusGpu);
index--;
continue;
}
int windowStart = index - windowSize + 1;
if (windowStart < 0)
{
windowStart = 0;
}
while (!IsBitSet(exponentGpu, windowStart))
{
windowStart++;
}
int windowBitCount = index - windowStart + 1;
for (int square = 0; square < windowBitCount; square++)
{
result.MulMod(result, modulusGpu);
}
ulong windowValue = ExtractWindowValue(exponentGpu, windowStart, windowBitCount);
int tableIndex = (int)((windowValue - 1UL) >> 1);
GpuUInt128 factor = oddPowers[tableIndex];
result.MulMod(factor, modulusGpu);
index = windowStart - 1;
}
return (UInt128)result;
}
private static bool ShouldUseSingleBit(GpuUInt128 exponent)
{
if (exponent.High == 0UL && exponent.Low <= Pow2WindowFallbackThreshold)
{
return true;
}
return false;
}
private static int GetWindowSize(int bitLength)
{
if (bitLength <= 8)
{
return Math.Max(bitLength, 1);
}
if (bitLength <= 23)
{
return 4;
}
if (bitLength <= 79)
{
return 5;
}
if (bitLength <= 239)
{
return 6;
}
if (bitLength <= 671)
{
return 7;
}
return 8;
}
private static void InitializeOddPowers(GpuUInt128 baseValue, GpuUInt128 modulus, Span<GpuUInt128> oddPowers)
{
oddPowers[0] = baseValue;
if (oddPowers.Length == 1)
{
return;
}
baseValue.MulMod(baseValue, modulus);
for (int i = 1; i < oddPowers.Length; i++)
{
oddPowers[i] = oddPowers[i - 1];
oddPowers[i].MulMod(baseValue, modulus);
}
}
private static GpuUInt128 Pow2MontgomeryModSingleBit(GpuUInt128 exponent, GpuUInt128 modulus, GpuUInt128 baseValue)
{
GpuUInt128 result = GpuUInt128.One;
while (!exponent.IsZero)
{
if ((exponent.Low & 1UL) != 0UL)
{
result.MulMod(baseValue, modulus);
}
exponent.ShiftRight(1);
if (exponent.IsZero)
{
break;
}
baseValue.MulMod(baseValue, modulus);
}
return result;
}
private static bool IsBitSet(GpuUInt128 value, int bitIndex)
{
if (bitIndex >= 64)
{
return ((value.High >> (bitIndex - 64)) & 1UL) != 0UL;
}
return ((value.Low >> bitIndex) & 1UL) != 0UL;
}
private static ulong ExtractWindowValue(GpuUInt128 exponent, int windowStart, int windowBitCount)
{
if (windowStart != 0)
{
GpuUInt128 shifted = exponent;
shifted.ShiftRight(windowStart);
ulong mask = (1UL << windowBitCount) - 1UL;
return shifted.Low & mask;
}
ulong directMask = (1UL << windowBitCount) - 1UL;
return exponent.Low & directMask;
}
private static ulong MultiplyMod(ulong x, ulong y, ulong modulus)
{
UInt128 product = (UInt128)x * y;
@@ -1,308 +0,0 @@
using System.Buffers.Binary;
using BenchmarkDotNet.Attributes;
using BenchmarkDotNet.Jobs;
using PerfectNumbers.Core;
using PerfectNumbers.Core.Gpu;
namespace EvenPerfectBitScanner.Benchmarks;
[MemoryDiagnoser]
[SimpleJob(RuntimeMoniker.Net80, launchCount: 1, warmupCount: 1, iterationCount: 5)]
public class Pow2MontgomeryModUInt128Benchmarks
{
private const int SampleCount = 256;
private const ulong Pow2WindowFallbackThreshold = 32UL;
private UInt128[] _moduli = Array.Empty<UInt128>();
private UInt128[] _exponents = Array.Empty<UInt128>();
[GlobalSetup]
public void Setup()
{
WideSampleData cache = WideSampleDataCache.Instance;
_moduli = cache.Moduli;
_exponents = cache.Exponents;
}
private sealed class WideSampleData
{
public WideSampleData(UInt128[] moduli, UInt128[] exponents)
{
Moduli = moduli;
Exponents = exponents;
}
public UInt128[] Moduli { get; }
public UInt128[] Exponents { get; }
}
private static class WideSampleDataCache
{
private static readonly Lazy<WideSampleData> Cache = new(Create);
public static WideSampleData Instance => Cache.Value;
private static WideSampleData Create()
{
Random random = new(19);
var moduli = new UInt128[SampleCount];
var exponents = new UInt128[SampleCount];
for (int i = 0; i < SampleCount; i++)
{
moduli[i] = NextWideOddModulus(random);
exponents[i] = NextWideExponent(random);
}
return new WideSampleData(moduli, exponents);
}
private static UInt128 NextWideOddModulus(Random random)
{
Span<byte> buffer = stackalloc byte[16];
random.NextBytes(buffer);
ulong low = BinaryPrimitives.ReadUInt64LittleEndian(buffer);
ulong high = BinaryPrimitives.ReadUInt64LittleEndian(buffer[8..]);
if (high == 0UL)
{
high = 1UL;
}
UInt128 modulus = ((UInt128)high << 64) | low;
modulus |= UInt128.One;
if (modulus <= UInt128.One)
{
modulus += 2UL;
}
return modulus;
}
private static UInt128 NextWideExponent(Random random)
{
Span<byte> buffer = stackalloc byte[16];
random.NextBytes(buffer);
ulong low = BinaryPrimitives.ReadUInt64LittleEndian(buffer);
ulong high = BinaryPrimitives.ReadUInt64LittleEndian(buffer[8..]);
if (high == 0UL)
{
high = 1UL;
}
return ((UInt128)high << 64) | low;
}
}
/// <summary>
/// Current UInt128-native windowed Montgomery ladder that avoids GpuUInt128 on CPU paths.
/// </summary>
[Benchmark(Baseline = true)]
public UInt128 UInt128Windowed()
{
UInt128 checksum = UInt128.Zero;
for (int i = 0; i < SampleCount; i++)
{
checksum ^= _exponents[i].Pow2MontgomeryModWindowed(_moduli[i]);
}
return checksum;
}
/// <summary>
/// Previous CPU shim that reused the GPU-oriented GpuUInt128 implementation.
/// </summary>
[Benchmark]
public UInt128 GpuStructShim()
{
UInt128 checksum = UInt128.Zero;
for (int i = 0; i < SampleCount; i++)
{
checksum ^= Pow2MontgomeryModWindowedGpuShim(_exponents[i], _moduli[i]);
}
return checksum;
}
private static UInt128 Pow2MontgomeryModWindowedGpuShim(UInt128 exponent, UInt128 modulus)
{
if (modulus == UInt128.One)
{
return UInt128.Zero;
}
if (exponent == UInt128.Zero)
{
return UInt128.One % modulus;
}
GpuUInt128 modulusGpu = new(modulus);
GpuUInt128 exponentGpu = new(exponent);
GpuUInt128 baseValue = new(2UL);
if (baseValue.CompareTo(modulusGpu) >= 0)
{
baseValue.Sub(modulusGpu);
}
if (ShouldUseSingleBit(exponentGpu))
{
GpuUInt128 singleBitResult = Pow2MontgomeryModSingleBit(exponentGpu, modulusGpu, baseValue);
return (UInt128)singleBitResult;
}
int bitLength = exponentGpu.GetBitLength();
int windowSize = GetWindowSize(bitLength);
int oddPowerCount = 1 << (windowSize - 1);
Span<GpuUInt128> oddPowers = stackalloc GpuUInt128[oddPowerCount];
InitializeOddPowers(baseValue, modulusGpu, oddPowers);
GpuUInt128 result = GpuUInt128.One;
int index = bitLength - 1;
while (index >= 0)
{
if (!IsBitSet(exponentGpu, index))
{
result.MulMod(result, modulusGpu);
index--;
continue;
}
int windowStart = index - windowSize + 1;
if (windowStart < 0)
{
windowStart = 0;
}
while (!IsBitSet(exponentGpu, windowStart))
{
windowStart++;
}
int windowBitCount = index - windowStart + 1;
for (int square = 0; square < windowBitCount; square++)
{
result.MulMod(result, modulusGpu);
}
ulong windowValue = ExtractWindowValue(exponentGpu, windowStart, windowBitCount);
int tableIndex = (int)((windowValue - 1UL) >> 1);
GpuUInt128 factor = oddPowers[tableIndex];
result.MulMod(factor, modulusGpu);
index = windowStart - 1;
}
return (UInt128)result;
}
private static bool ShouldUseSingleBit(GpuUInt128 exponent)
{
if (exponent.High == 0UL && exponent.Low <= Pow2WindowFallbackThreshold)
{
return true;
}
return false;
}
private static int GetWindowSize(int bitLength)
{
if (bitLength <= 8)
{
return Math.Max(bitLength, 1);
}
if (bitLength <= 23)
{
return 4;
}
if (bitLength <= 79)
{
return 5;
}
if (bitLength <= 239)
{
return 6;
}
if (bitLength <= 671)
{
return 7;
}
return 8;
}
private static void InitializeOddPowers(GpuUInt128 baseValue, GpuUInt128 modulus, Span<GpuUInt128> oddPowers)
{
oddPowers[0] = baseValue;
if (oddPowers.Length == 1)
{
return;
}
// Reusing baseValue to hold base^2 for the shared odd-power ladder.
baseValue.MulMod(baseValue, modulus);
for (int i = 1; i < oddPowers.Length; i++)
{
oddPowers[i] = oddPowers[i - 1];
oddPowers[i].MulMod(baseValue, modulus);
}
}
private static GpuUInt128 Pow2MontgomeryModSingleBit(GpuUInt128 exponent, GpuUInt128 modulus, GpuUInt128 baseValue)
{
GpuUInt128 result = GpuUInt128.One;
while (!exponent.IsZero)
{
if ((exponent.Low & 1UL) != 0UL)
{
result.MulMod(baseValue, modulus);
}
exponent.ShiftRight(1);
if (exponent.IsZero)
{
break;
}
// Reusing baseValue to store the squared base for the next iteration.
baseValue.MulMod(baseValue, modulus);
}
return result;
}
private static bool IsBitSet(GpuUInt128 value, int bitIndex)
{
if (bitIndex >= 64)
{
return ((value.High >> (bitIndex - 64)) & 1UL) != 0UL;
}
return ((value.Low >> bitIndex) & 1UL) != 0UL;
}
private static ulong ExtractWindowValue(GpuUInt128 exponent, int windowStart, int windowBitCount)
{
if (windowStart != 0)
{
GpuUInt128 shifted = exponent;
shifted.ShiftRight(windowStart);
ulong mask = (1UL << windowBitCount) - 1UL;
return shifted.Low & mask;
}
ulong directMask = (1UL << windowBitCount) - 1UL;
return exponent.Low & directMask;
}
}
+2 -2
View File
@@ -9,9 +9,9 @@ public static class Program
// var config = DefaultConfig.Instance.WithOptions(ConfigOptions.JoinSummary | ConfigOptions.DisableOptimizationsValidator);
// BenchmarkSwitcher.FromAssembly(typeof(Program).Assembly).Run(args, config);
Pow2MontgomeryModBenchmarks benchmarks = new()
Pow2ModBenchmarks benchmarks = new()
{
Scale = Pow2MontgomeryModBenchmarks.InputScale.Large
Scale = Pow2ModBenchmarks.InputScale.Large
};
GpuPrimeWorkLimiter.SetLimit(1);