Clarify Pollard-Rho parity guard and seed CPU batches
This commit is contained in:
1 parent
8e3e990499
commit
7feeb08943
4 files changed
+160
-74
No files matched your search
@@ -43,16 +43,13 @@ internal sealed class MersenneCpuDivisorScanSession : IMersenneNumberDivisorByDi
|
|||||||
// Keep these remainder steppers in place so future updates continue reusing the previously computed residues.
|
// Keep these remainder steppers in place so future updates continue reusing the previously computed residues.
|
||||||
// They are critical for avoiding repeated full Montgomery exponentiation work when scanning divisors.
|
// They are critical for avoiding repeated full Montgomery exponentiation work when scanning divisors.
|
||||||
var exponentStepper = new ExponentRemainderStepper(cachedData);
|
var exponentStepper = new ExponentRemainderStepper(cachedData);
|
||||||
if (!exponentStepper.IsValidModulus)
|
|
||||||
{
|
|
||||||
throw new InvalidOperationException($"Invalid modulus {cachedData.Modulus} for divisor scans.");
|
|
||||||
}
|
|
||||||
|
|
||||||
var cycleStepper = new CycleRemainderStepper(divisorCycle);
|
var cycleStepper = new CycleRemainderStepper(divisorCycle);
|
||||||
|
|
||||||
|
bool initialUnity = exponentStepper.InitializeCpuIsUnity(primes[0]);
|
||||||
ulong remainder = cycleStepper.Initialize(primes[0]);
|
ulong remainder = cycleStepper.Initialize(primes[0]);
|
||||||
hits[0] = remainder == 0UL
|
hits[0] = remainder == 0UL
|
||||||
? (exponentStepper.ComputeNextIsUnity(primes[0]) ? (byte)1 : (byte)0)
|
? (initialUnity ? (byte)1 : (byte)0)
|
||||||
: (byte)0;
|
: (byte)0;
|
||||||
|
|
||||||
for (int i = 1; i < length; i++)
|
for (int i = 1; i < length; i++)
|
||||||
|
|||||||
@@ -1,5 +1,9 @@
|
|||||||
namespace PerfectNumbers.Core;
|
namespace PerfectNumbers.Core;
|
||||||
|
|
||||||
|
/// <summary>
|
||||||
|
/// Tracks successive Montgomery residues for a single divisor while scanning exponent candidates.
|
||||||
|
/// Callers are responsible for seeding the stepper once per divisor and then advancing exponents in ascending order.
|
||||||
|
/// </summary>
|
||||||
internal struct ExponentRemainderStepper
|
internal struct ExponentRemainderStepper
|
||||||
{
|
{
|
||||||
private readonly MontgomeryDivisorData _divisor;
|
private readonly MontgomeryDivisorData _divisor;
|
||||||
@@ -8,7 +12,6 @@ internal struct ExponentRemainderStepper
|
|||||||
private readonly ulong _montgomeryOne;
|
private readonly ulong _montgomeryOne;
|
||||||
private ulong _previousExponent;
|
private ulong _previousExponent;
|
||||||
private ulong _currentMontgomery;
|
private ulong _currentMontgomery;
|
||||||
private bool _hasState;
|
|
||||||
|
|
||||||
public ExponentRemainderStepper(in MontgomeryDivisorData divisor)
|
public ExponentRemainderStepper(in MontgomeryDivisorData divisor)
|
||||||
{
|
{
|
||||||
@@ -18,13 +21,8 @@ internal struct ExponentRemainderStepper
|
|||||||
_montgomeryOne = divisor.MontgomeryOne;
|
_montgomeryOne = divisor.MontgomeryOne;
|
||||||
_previousExponent = 0UL;
|
_previousExponent = 0UL;
|
||||||
_currentMontgomery = divisor.MontgomeryOne;
|
_currentMontgomery = divisor.MontgomeryOne;
|
||||||
_hasState = false;
|
|
||||||
}
|
}
|
||||||
|
|
||||||
public bool IsValidModulus => _modulus > 1UL && (_modulus & 1UL) != 0UL;
|
|
||||||
|
|
||||||
public bool HasState => _hasState;
|
|
||||||
|
|
||||||
public ulong PreviousExponent => _previousExponent;
|
public ulong PreviousExponent => _previousExponent;
|
||||||
|
|
||||||
public bool MatchesDivisor(in MontgomeryDivisorData divisor)
|
public bool MatchesDivisor(in MontgomeryDivisorData divisor)
|
||||||
@@ -32,29 +30,47 @@ internal struct ExponentRemainderStepper
|
|||||||
return _modulus == divisor.Modulus && _nPrime == divisor.NPrime && _montgomeryOne == divisor.MontgomeryOne;
|
return _modulus == divisor.Modulus && _nPrime == divisor.NPrime && _montgomeryOne == divisor.MontgomeryOne;
|
||||||
}
|
}
|
||||||
|
|
||||||
|
/// <summary>
|
||||||
|
/// Clears the cached Montgomery residue so the stepper can be reused for another divisor.
|
||||||
|
/// Callers must reinitialize the stepper before consuming any residues after a reset.
|
||||||
|
/// </summary>
|
||||||
public void Reset()
|
public void Reset()
|
||||||
{
|
{
|
||||||
_previousExponent = 0UL;
|
_previousExponent = 0UL;
|
||||||
_currentMontgomery = _montgomeryOne;
|
_currentMontgomery = _montgomeryOne;
|
||||||
_hasState = false;
|
|
||||||
}
|
}
|
||||||
|
|
||||||
public ulong ComputeNext(ulong exponent)
|
/// <summary>
|
||||||
|
/// Initializes the stepper for CPU calculations using the first exponent in a strictly ascending sequence.
|
||||||
|
/// The returned value is the canonical (non-Montgomery) residue corresponding to the provided exponent.
|
||||||
|
/// </summary>
|
||||||
|
/// <param name="exponent">First exponent evaluated for the current divisor. The caller must reuse the same divisor metadata.</param>
|
||||||
|
/// <returns>The canonical residue produced by <c>2^exponent mod divisor</c>.</returns>
|
||||||
|
public ulong InitializeCpu(ulong exponent)
|
||||||
{
|
{
|
||||||
if (!IsValidModulus)
|
InitializeCpuState(exponent);
|
||||||
{
|
return ReduceCurrent();
|
||||||
_hasState = false;
|
}
|
||||||
return 0UL;
|
|
||||||
}
|
|
||||||
|
|
||||||
if (!_hasState || exponent <= _previousExponent)
|
/// <summary>
|
||||||
{
|
/// Initializes the CPU stepping state and reports whether the seed exponent already yields unity.
|
||||||
_currentMontgomery = exponent.Pow2MontgomeryModWindowedCpu(_divisor, keepMontgomery: true);
|
/// </summary>
|
||||||
_previousExponent = exponent;
|
/// <param name="exponent">First exponent evaluated for the current divisor. Must match the divisor passed to the constructor.</param>
|
||||||
_hasState = true;
|
/// <returns><see langword="true"/> when the seed exponent produces a Montgomery residue equal to one.</returns>
|
||||||
return ReduceCurrent();
|
public bool InitializeCpuIsUnity(ulong exponent)
|
||||||
}
|
{
|
||||||
|
InitializeCpuState(exponent);
|
||||||
|
return _currentMontgomery == _montgomeryOne;
|
||||||
|
}
|
||||||
|
|
||||||
|
/// <summary>
|
||||||
|
/// Advances the CPU stepping state to the next exponent in a strictly ascending sequence that shares the same divisor.
|
||||||
|
/// Callers must invoke <see cref="InitializeCpu(ulong)"/> beforehand.
|
||||||
|
/// </summary>
|
||||||
|
/// <param name="exponent">Next exponent in the ascending sequence. Callers must never move backwards or reuse an older exponent.</param>
|
||||||
|
/// <returns>The canonical residue produced by <c>2^exponent mod divisor</c>.</returns>
|
||||||
|
public ulong ComputeNextCpu(ulong exponent)
|
||||||
|
{
|
||||||
ulong delta = exponent - _previousExponent;
|
ulong delta = exponent - _previousExponent;
|
||||||
// TODO: Once divisor cycle lengths are mandatory, pull the delta multiplier from the
|
// TODO: Once divisor cycle lengths are mandatory, pull the delta multiplier from the
|
||||||
// single-block divisor-cycle snapshot so we can skip the powmod entirely and reuse the
|
// single-block divisor-cycle snapshot so we can skip the powmod entirely and reuse the
|
||||||
@@ -65,23 +81,37 @@ internal struct ExponentRemainderStepper
|
|||||||
return ReduceCurrent();
|
return ReduceCurrent();
|
||||||
}
|
}
|
||||||
|
|
||||||
|
/// <summary>
|
||||||
|
/// Initializes the stepper for GPU calculations using the first exponent in a strictly ascending sequence.
|
||||||
|
/// </summary>
|
||||||
|
/// <param name="exponent">First exponent evaluated on the GPU path. The same divisor metadata must be reused for future calls.</param>
|
||||||
|
/// <returns>The canonical residue produced by <c>2^exponent mod divisor</c>.</returns>
|
||||||
|
public ulong InitializeGpu(ulong exponent)
|
||||||
|
{
|
||||||
|
InitializeGpuState(exponent);
|
||||||
|
return ReduceCurrent();
|
||||||
|
}
|
||||||
|
|
||||||
|
/// <summary>
|
||||||
|
/// Initializes the GPU stepping state and reports whether the seed exponent produces unity.
|
||||||
|
/// </summary>
|
||||||
|
/// <param name="exponent">First exponent evaluated on the GPU path. Must match the divisor supplied during construction.</param>
|
||||||
|
/// <returns><see langword="true"/> when the seed exponent produces a Montgomery residue equal to one.</returns>
|
||||||
|
public bool InitializeGpuIsUnity(ulong exponent)
|
||||||
|
{
|
||||||
|
InitializeGpuState(exponent);
|
||||||
|
return _currentMontgomery == _montgomeryOne;
|
||||||
|
}
|
||||||
|
|
||||||
|
/// <summary>
|
||||||
|
/// Advances the CPU stepping state to the next exponent and reports whether it evaluates to unity.
|
||||||
|
/// Callers must seed the stepper through <see cref="InitializeCpu(ulong)"/> or <see cref="InitializeCpuIsUnity(ulong)"/> and then
|
||||||
|
/// pass exponents in strictly ascending order.
|
||||||
|
/// </summary>
|
||||||
|
/// <param name="exponent">Next exponent in the ascending sequence.</param>
|
||||||
|
/// <returns><see langword="true"/> when the canonical residue equals one.</returns>
|
||||||
public bool ComputeNextIsUnity(ulong exponent)
|
public bool ComputeNextIsUnity(ulong exponent)
|
||||||
{
|
{
|
||||||
if (!IsValidModulus)
|
|
||||||
{
|
|
||||||
throw new InvalidOperationException($"Modulus is invalid for exponent {exponent}");
|
|
||||||
// _hasState = false;
|
|
||||||
// return false;
|
|
||||||
}
|
|
||||||
|
|
||||||
if (!_hasState || exponent <= _previousExponent)
|
|
||||||
{
|
|
||||||
_currentMontgomery = exponent.Pow2MontgomeryModWindowedCpu(_divisor, keepMontgomery: true);
|
|
||||||
_previousExponent = exponent;
|
|
||||||
_hasState = true;
|
|
||||||
return _currentMontgomery == _montgomeryOne;
|
|
||||||
}
|
|
||||||
|
|
||||||
ulong delta = exponent - _previousExponent;
|
ulong delta = exponent - _previousExponent;
|
||||||
// TODO: Reuse the divisor-cycle derived Montgomery delta once the cache exposes single-cycle
|
// TODO: Reuse the divisor-cycle derived Montgomery delta once the cache exposes single-cycle
|
||||||
// lookups so this branch also avoids recomputing powmods when the snapshot lacks the divisor.
|
// lookups so this branch also avoids recomputing powmods when the snapshot lacks the divisor.
|
||||||
@@ -91,40 +121,32 @@ internal struct ExponentRemainderStepper
|
|||||||
return _currentMontgomery == _montgomeryOne;
|
return _currentMontgomery == _montgomeryOne;
|
||||||
}
|
}
|
||||||
|
|
||||||
|
/// <summary>
|
||||||
|
/// Seeds the stepper with a precomputed Montgomery residue produced by the caller.
|
||||||
|
/// The provided exponent becomes the new baseline and all subsequent calls must continue with greater exponents.
|
||||||
|
/// </summary>
|
||||||
|
/// <param name="exponent">Exponent associated with <paramref name="montgomeryResult"/>.</param>
|
||||||
|
/// <param name="montgomeryResult">Residue in Montgomery form produced externally for the configured divisor.</param>
|
||||||
|
/// <param name="isUnity">Outputs whether the supplied residue equals Montgomery one.</param>
|
||||||
|
/// <returns><see langword="true"/> when the residue is accepted and cached.</returns>
|
||||||
public bool TryInitializeFromMontgomeryResult(ulong exponent, ulong montgomeryResult, out bool isUnity)
|
public bool TryInitializeFromMontgomeryResult(ulong exponent, ulong montgomeryResult, out bool isUnity)
|
||||||
{
|
{
|
||||||
if (!IsValidModulus)
|
|
||||||
{
|
|
||||||
_hasState = false;
|
|
||||||
isUnity = false;
|
|
||||||
return false;
|
|
||||||
}
|
|
||||||
|
|
||||||
_currentMontgomery = montgomeryResult;
|
_currentMontgomery = montgomeryResult;
|
||||||
_previousExponent = exponent;
|
_previousExponent = exponent;
|
||||||
_hasState = true;
|
|
||||||
isUnity = _currentMontgomery == _montgomeryOne;
|
isUnity = _currentMontgomery == _montgomeryOne;
|
||||||
return true;
|
return true;
|
||||||
}
|
}
|
||||||
|
|
||||||
|
/// <summary>
|
||||||
|
/// Advances the stepper using a caller-supplied Montgomery delta that represents the gap to the next exponent.
|
||||||
|
/// Callers must seed the stepper first, reuse the same divisor metadata, and continue feeding exponents in ascending order so the delta remains valid.
|
||||||
|
/// </summary>
|
||||||
|
/// <param name="exponent">Next exponent in the ascending sequence.</param>
|
||||||
|
/// <param name="montgomeryDelta">Montgomery residue delta corresponding to the exponent gap.</param>
|
||||||
|
/// <param name="isUnity">Outputs whether the resulting residue equals Montgomery one.</param>
|
||||||
|
/// <returns><see langword="true"/> when the advance succeeds.</returns>
|
||||||
public bool TryAdvanceWithMontgomeryDelta(ulong exponent, ulong montgomeryDelta, out bool isUnity)
|
public bool TryAdvanceWithMontgomeryDelta(ulong exponent, ulong montgomeryDelta, out bool isUnity)
|
||||||
{
|
{
|
||||||
if (!IsValidModulus)
|
|
||||||
{
|
|
||||||
_hasState = false;
|
|
||||||
isUnity = false;
|
|
||||||
return false;
|
|
||||||
}
|
|
||||||
|
|
||||||
if (!_hasState || exponent <= _previousExponent)
|
|
||||||
{
|
|
||||||
_currentMontgomery = exponent.Pow2MontgomeryModWindowedCpu(_divisor, keepMontgomery: true);
|
|
||||||
_previousExponent = exponent;
|
|
||||||
_hasState = true;
|
|
||||||
isUnity = _currentMontgomery == _montgomeryOne;
|
|
||||||
return true;
|
|
||||||
}
|
|
||||||
|
|
||||||
// TODO: Once the divisor-cycle cache exposes a direct Montgomery delta, multiply it here instead
|
// TODO: Once the divisor-cycle cache exposes a direct Montgomery delta, multiply it here instead
|
||||||
// of relying on the caller-provided delta so incremental scans remain in sync with the snapshot
|
// of relying on the caller-provided delta so incremental scans remain in sync with the snapshot
|
||||||
// without computing additional cycles or mutating cache state.
|
// without computing additional cycles or mutating cache state.
|
||||||
@@ -134,5 +156,17 @@ internal struct ExponentRemainderStepper
|
|||||||
return true;
|
return true;
|
||||||
}
|
}
|
||||||
|
|
||||||
|
private void InitializeCpuState(ulong exponent)
|
||||||
|
{
|
||||||
|
_currentMontgomery = exponent.Pow2MontgomeryModWindowedCpu(_divisor, keepMontgomery: true);
|
||||||
|
_previousExponent = exponent;
|
||||||
|
}
|
||||||
|
|
||||||
|
private void InitializeGpuState(ulong exponent)
|
||||||
|
{
|
||||||
|
_currentMontgomery = exponent.Pow2MontgomeryModWindowedCpu(_divisor, keepMontgomery: true);
|
||||||
|
_previousExponent = exponent;
|
||||||
|
}
|
||||||
|
|
||||||
private ulong ReduceCurrent() => _currentMontgomery.MontgomeryMultiply(1UL, _modulus, _nPrime);
|
private ulong ReduceCurrent() => _currentMontgomery.MontgomeryMultiply(1UL, _modulus, _nPrime);
|
||||||
}
|
}
|
||||||
@@ -564,7 +564,8 @@ public class MersenneDivisorCycles
|
|||||||
|
|
||||||
stepper.Reset();
|
stepper.Reset();
|
||||||
|
|
||||||
for (int j = 0; j < actual; j++)
|
evaluations[0] = stepper.InitializeCpuIsUnity(candidates[0]);
|
||||||
|
for (int j = 1; j < actual; j++)
|
||||||
{
|
{
|
||||||
evaluations[j] = stepper.ComputeNextIsUnity(candidates[j]);
|
evaluations[j] = stepper.ComputeNextIsUnity(candidates[j]);
|
||||||
}
|
}
|
||||||
@@ -582,10 +583,13 @@ public class MersenneDivisorCycles
|
|||||||
|
|
||||||
private static ulong PollardRho64(ulong n)
|
private static ulong PollardRho64(ulong n)
|
||||||
{
|
{
|
||||||
if ((n & 1UL) == 0UL)
|
// The EvenPerfectBitScanner path removes all factors of two before invoking Pollard-Rho, so this branch never executes
|
||||||
{
|
// outside synthetic tests and benchmarks. Keep the guard documented to prevent future regressions without paying the hot-
|
||||||
return 2UL;
|
// path cost.
|
||||||
}
|
// if ((n & 1UL) == 0UL)
|
||||||
|
// {
|
||||||
|
// return 2UL;
|
||||||
|
// }
|
||||||
|
|
||||||
while (true)
|
while (true)
|
||||||
{
|
{
|
||||||
|
|||||||
@@ -165,7 +165,14 @@ internal static partial class PrimeOrderCalculator
|
|||||||
|
|
||||||
ExponentRemainderStepper stepper = ThreadStaticPools.RentExponentStepper(divisorData);
|
ExponentRemainderStepper stepper = ThreadStaticPools.RentExponentStepper(divisorData);
|
||||||
int candidateCount = candidates.Length;
|
int candidateCount = candidates.Length;
|
||||||
for (int i = 0; i < candidateCount; i++)
|
|
||||||
|
if (stepper.InitializeCpuIsUnity(candidates[0]))
|
||||||
|
{
|
||||||
|
ThreadStaticPools.ReturnExponentStepper(stepper);
|
||||||
|
return false;
|
||||||
|
}
|
||||||
|
|
||||||
|
for (int i = 1; i < candidateCount; i++)
|
||||||
{
|
{
|
||||||
if (stepper.ComputeNextIsUnity(candidates[i]))
|
if (stepper.ComputeNextIsUnity(candidates[i]))
|
||||||
{
|
{
|
||||||
@@ -300,7 +307,8 @@ internal static partial class PrimeOrderCalculator
|
|||||||
|
|
||||||
stepper.Reset();
|
stepper.Reset();
|
||||||
|
|
||||||
for (int j = 0; j < actual; j++)
|
evaluations[0] = stepper.InitializeCpuIsUnity(candidates[0]);
|
||||||
|
for (int j = 1; j < actual; j++)
|
||||||
{
|
{
|
||||||
evaluations[j] = stepper.ComputeNextIsUnity(candidates[j]);
|
evaluations[j] = stepper.ComputeNextIsUnity(candidates[j]);
|
||||||
}
|
}
|
||||||
@@ -459,7 +467,12 @@ internal static partial class PrimeOrderCalculator
|
|||||||
|
|
||||||
stepper.Reset();
|
stepper.Reset();
|
||||||
|
|
||||||
for (int j = 0; j < actual; j++)
|
if (stepper.InitializeCpuIsUnity(candidates[0]))
|
||||||
|
{
|
||||||
|
return true;
|
||||||
|
}
|
||||||
|
|
||||||
|
for (int j = 1; j < actual; j++)
|
||||||
{
|
{
|
||||||
if (stepper.ComputeNextIsUnity(candidates[j]))
|
if (stepper.ComputeNextIsUnity(candidates[j]))
|
||||||
{
|
{
|
||||||
@@ -535,6 +548,7 @@ internal static partial class PrimeOrderCalculator
|
|||||||
bool allowGpuBatch = true;
|
bool allowGpuBatch = true;
|
||||||
Span<ulong> candidateSpan = CollectionsMarshal.AsSpan(candidates);
|
Span<ulong> candidateSpan = CollectionsMarshal.AsSpan(candidates);
|
||||||
ExponentRemainderStepper powStepper = ThreadStaticPools.RentExponentStepper(divisorData);
|
ExponentRemainderStepper powStepper = ThreadStaticPools.RentExponentStepper(divisorData);
|
||||||
|
bool powStepperInitialized = false;
|
||||||
|
|
||||||
// DebugLog(() => $"Checking candidates ({candidateCount} candidates, {powBudget} pow budget)");
|
// DebugLog(() => $"Checking candidates ({candidateCount} candidates, {powBudget} pow budget)");
|
||||||
int index = 0;
|
int index = 0;
|
||||||
@@ -594,7 +608,33 @@ internal static partial class PrimeOrderCalculator
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
for (int i = 0; i < batchSize && powUsed < powBudget; i++)
|
int startIndex = 0;
|
||||||
|
|
||||||
|
if (!gpuSuccess && !powStepperInitialized)
|
||||||
|
{
|
||||||
|
ulong candidate = batch[0];
|
||||||
|
powUsed++;
|
||||||
|
|
||||||
|
bool equalsOne = powStepper.InitializeCpuIsUnity(candidate);
|
||||||
|
powStepperInitialized = true;
|
||||||
|
if (equalsOne && TryConfirmCandidateCpu(prime, candidate, divisorData, config, ref powUsed, powBudget))
|
||||||
|
{
|
||||||
|
if (gpuPool is not null)
|
||||||
|
{
|
||||||
|
pool.Return(gpuPool, clearArray: false);
|
||||||
|
}
|
||||||
|
|
||||||
|
candidates.Clear();
|
||||||
|
ThreadStaticPools.ReturnUlongList(candidates);
|
||||||
|
ThreadStaticPools.ReturnExponentStepper(powStepper);
|
||||||
|
result = candidate;
|
||||||
|
return true;
|
||||||
|
}
|
||||||
|
|
||||||
|
startIndex = 1;
|
||||||
|
}
|
||||||
|
|
||||||
|
for (int i = startIndex; i < batchSize && powUsed < powBudget; i++)
|
||||||
{
|
{
|
||||||
ulong candidate = batch[i];
|
ulong candidate = batch[i];
|
||||||
powUsed++;
|
powUsed++;
|
||||||
@@ -936,7 +976,18 @@ internal static partial class PrimeOrderCalculator
|
|||||||
|
|
||||||
stepper.Reset();
|
stepper.Reset();
|
||||||
|
|
||||||
for (int j = 0; j < actual; j++)
|
if (powUsed >= powBudget && powBudget > 0)
|
||||||
|
{
|
||||||
|
return true;
|
||||||
|
}
|
||||||
|
|
||||||
|
powUsed++;
|
||||||
|
if (stepper.InitializeCpuIsUnity(candidates[0]))
|
||||||
|
{
|
||||||
|
return true;
|
||||||
|
}
|
||||||
|
|
||||||
|
for (int j = 1; j < actual; j++)
|
||||||
{
|
{
|
||||||
if (powUsed >= powBudget && powBudget > 0)
|
if (powUsed >= powBudget && powBudget > 0)
|
||||||
{
|
{
|
||||||
|
|||||||
Reference in new issue
Block a user