Revert "Use standard modular math for single-use GPU exponent checks"

This commit is contained in:
Marek Lemańczyk authored and GitHub committed 2025-10-13 04:59:42 +02:00
1 parent ad6a5c6e61
commit 8b02993bac
6 files changed
+16 -71

No files matched your search

@@ -652,7 +652,7 @@ public sealed class MersenneNumberDivisorByDivisorCpuTester : IMersenneNumberDiv
private bool _disposed;
private ulong[]? _primeDeltas;
private ulong[]? _cycleRemainders;
private ulong[]? _residues;
private ulong[]? _montgomeryDeltas;
private int _bufferCapacity;
private MersenneNumberDivisorResidueGpuEvaluator? _gpuResidueEvaluator;
@@ -765,12 +765,13 @@ public sealed class MersenneNumberDivisorByDivisorCpuTester : IMersenneNumberDiv
ComputeCycleRemainders(divisorCycle, primes, deltaSpan, remainderSpan);
EnsureGpuEvaluator();
Span<ulong> residueSpan = _residues!.AsSpan(0, length);
Span<ulong> residueSpan = _montgomeryDeltas!.AsSpan(0, length);
_gpuResidueEvaluator!.ComputeResidues(primes, divisorData, residueSpan);
ulong montgomeryOne = divisorData.MontgomeryOne;
for (int i = 0; i < length; i++)
{
hits[i] = remainderSpan[i] == 0UL && residueSpan[i] == 1UL ? (byte)1 : (byte)0;
hits[i] = remainderSpan[i] == 0UL && residueSpan[i] == montgomeryOne ? (byte)1 : (byte)0;
}
}
@@ -800,7 +801,7 @@ public sealed class MersenneNumberDivisorByDivisorCpuTester : IMersenneNumberDiv
return;
}
if (requiredLength <= _bufferCapacity && _primeDeltas is not null && _cycleRemainders is not null && _residues is not null)
if (requiredLength <= _bufferCapacity && _primeDeltas is not null && _cycleRemainders is not null && _montgomeryDeltas is not null)
{
return;
}
@@ -820,7 +821,7 @@ public sealed class MersenneNumberDivisorByDivisorCpuTester : IMersenneNumberDiv
_primeDeltas = ArrayPool<ulong>.Shared.Rent(newCapacity);
_cycleRemainders = ArrayPool<ulong>.Shared.Rent(newCapacity);
_residues = ArrayPool<ulong>.Shared.Rent(newCapacity);
_montgomeryDeltas = ArrayPool<ulong>.Shared.Rent(newCapacity);
_bufferCapacity = newCapacity;
}
@@ -838,10 +839,10 @@ public sealed class MersenneNumberDivisorByDivisorCpuTester : IMersenneNumberDiv
_cycleRemainders = null;
}
if (_residues is not null)
if (_montgomeryDeltas is not null)
{
ArrayPool<ulong>.Shared.Return(_residues, clearArray: false);
_residues = null;
ArrayPool<ulong>.Shared.Return(_montgomeryDeltas, clearArray: false);
_montgomeryDeltas = null;
}
_bufferCapacity = 0;
@@ -35,7 +35,7 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
_kernelCache.GetOrAdd(accelerator, acc => acc.LoadAutoGroupedStreamKernel<Index1D, ArrayView<MontgomeryDivisorData>, ArrayView<ulong>, ArrayView<byte>>(CheckKernel));
private Action<Index1D, MontgomeryDivisorData, ArrayView<ulong>, ArrayView<ulong>> GetKernelByPrimeExponent(Accelerator accelerator) =>
_kernelExponentCache.GetOrAdd(accelerator, acc => acc.LoadAutoGroupedStreamKernel<Index1D, MontgomeryDivisorData, ArrayView<ulong>, ArrayView<ulong>>(ComputeExponentKernel));
_kernelExponentCache.GetOrAdd(accelerator, acc => acc.LoadAutoGroupedStreamKernel<Index1D, MontgomeryDivisorData, ArrayView<ulong>, ArrayView<ulong>>(ComputeMontgomeryExponentKernel));
private Action<Index1D, ArrayView<ulong>, ArrayView<byte>, byte> GetRemainderDeltaKernel(Accelerator accelerator) =>
_remainderDeltaKernelCache.GetOrAdd(accelerator, acc => acc.LoadAutoGroupedStreamKernel<Index1D, ArrayView<ulong>, ArrayView<byte>, byte>(ComputeRemainderDeltasKernel));
@@ -569,7 +569,7 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
}
ulong exponent = exponents[index];
hits[index] = exponent.Pow2ModBinaryGpu(modulus) == 1UL ? (byte)1 : (byte)0;
hits[index] = exponent.Pow2MontgomeryModWindowedGpu(divisor, keepMontgomery: false) == 1UL ? (byte)1 : (byte)0;
}
public sealed class DivisorScanSession : IMersenneNumberDivisorByDivisorTester.IDivisorScanSession
@@ -777,7 +777,7 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
mask[globalIndex] = 1;
}
private static void ComputeExponentKernel(Index1D index, MontgomeryDivisorData divisor, ArrayView<ulong> exponents, ArrayView<ulong> results)
private static void ComputeMontgomeryExponentKernel(Index1D index, MontgomeryDivisorData divisor, ArrayView<ulong> exponents, ArrayView<ulong> results)
{
ulong modulus = divisor.Modulus;
if (modulus <= 1UL || (modulus & 1UL) == 0UL)
@@ -788,7 +788,7 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
ulong exponent = exponents[index];
results[index] = exponent.Pow2ModBinaryGpu(modulus);
results[index] = exponent.Pow2MontgomeryModWindowedGpu(divisor, keepMontgomery: true);
}
public ulong DivisorLimit
@@ -109,13 +109,6 @@ internal sealed class MersenneNumberDivisorResidueGpuEvaluator : IDisposable
}
ulong exponent = exponents[idx];
ulong modulus = divisorData.Modulus;
if (modulus <= 1UL || (modulus & 1UL) == 0UL)
{
results[idx] = 0UL;
return;
}
results[idx] = exponent.Pow2ModBinaryGpu(modulus);
results[idx] = exponent.Pow2MontgomeryModWindowedGpu(divisorData, keepMontgomery: true);
}
}
@@ -811,7 +811,7 @@ internal static partial class PrimeOrderGpuHeuristics
}
ulong reduced = phi / factor;
if (reduced.Pow2ModBinaryGpu(divisor.Modulus) == 1UL)
if (reduced.Pow2MontgomeryModWindowedGpu(divisor, keepMontgomery: false) == 1UL)
{
return false;
}
@@ -1595,7 +1595,7 @@ internal static partial class PrimeOrderGpuHeuristics
private static bool Pow2EqualsOneKernel(ulong exponent, in MontgomeryDivisorData divisor)
{
return exponent.Pow2ModBinaryGpu(divisor.Modulus) == 1UL;
return exponent.Pow2MontgomeryModWindowedGpu(divisor, keepMontgomery: false) == 1UL;
}
private static ulong CalculateByDoublingKernel(ulong prime)
-31
View File
@@ -310,37 +310,6 @@ public static class ULongExtensions
return state.MulModWithNativeModulo(b, modulus);
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public static ulong Pow2ModBinaryGpu(this ulong exponent, ulong modulus)
{
if (modulus <= 1UL)
{
return 0UL;
}
ulong remainingExponent = exponent;
ulong result = 1UL % modulus;
ulong baseValue = 2UL % modulus;
while (remainingExponent != 0UL)
{
if ((remainingExponent & 1UL) != 0UL)
{
result = result.MulMod64GpuCompatible(baseValue, modulus);
}
remainingExponent >>= 1;
if (remainingExponent == 0UL)
{
break;
}
baseValue = baseValue.MulMod64GpuCompatible(baseValue, modulus);
}
return result;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public static ulong MontgomeryMultiply(this ulong a, ulong b, ulong modulus, ulong nPrime)
{
-18
View File
@@ -1,18 +0,0 @@
# GPU Modular Arithmetic in the `--mersenne=bydivisor` Path
The GPU implementation now routes single-use exponent checks through standard modular arithmetic while keeping Montgomery reduction only where repeated multiplications reuse the same modulus.
## Kernels using standard modular arithmetic
- `CheckKernel` in `MersenneNumberDivisorByDivisorGpuTester` evaluates each divisor candidate with `Pow2ModBinaryGpu` and compares the result directly against 1. 【F:PerfectNumbers.Core/Gpu/MersenneNumberDivisorByDivisorGpuTester.cs†L561-L645】【F:PerfectNumbers.Core/ULongExtensions.cs†L314-L339】
- `ComputeExponentKernel` in the same tester precomputes batched residues with the plain helper so CPU lookups receive standard-domain results. 【F:PerfectNumbers.Core/Gpu/MersenneNumberDivisorByDivisorGpuTester.cs†L780-L791】【F:PerfectNumbers.Core/ULongExtensions.cs†L314-L339】
- `ComputeResiduesKernel` inside `MersenneNumberDivisorResidueGpuEvaluator` now mirrors that logic when building residue tables for mixed CPU/GPU scans. 【F:PerfectNumbers.Core/Gpu/MersenneNumberDivisorResidueGpuEvaluator.cs†L102-L114】【F:PerfectNumbers.Core/ULongExtensions.cs†L314-L339】
- The order-search helpers in `PrimeOrderGpuHeuristics` (`TrySpecialMaxKernel`, `ExponentLoweringKernel`, and `Pow2EqualsOneKernel`) also rely on `Pow2ModBinaryGpu` for their equality checks. 【F:PerfectNumbers.Core/Gpu/PrimeOrderGpuHeuristics.cs†L798-L866】【F:PerfectNumbers.Core/Gpu/PrimeOrderGpuHeuristics.cs†L1596-L1598】
## Kernels still using Montgomery arithmetic
- `MontgomeryOddPowerKernelScan` in `GpuKernelPool` continues to build windowed odd-power tables with `MontgomeryMultiply`, which amortizes the Montgomery setup over many multiplications per modulus. 【F:PerfectNumbers.Core/Gpu/GpuKernelPool.cs†L200-L237】【F:PerfectNumbers.Core/ULongExtensions.cs†L532-L559】
## Shared helper
- `Pow2ModBinaryGpu` in `ULongExtensions` implements the GPU-friendly binary exponentiation loop that powers the standard modular path across these kernels. 【F:PerfectNumbers.Core/ULongExtensions.cs†L314-L339】