diff --git a/PerfectNumbers.Core/MersenneDivisorCycles.cs b/PerfectNumbers.Core/MersenneDivisorCycles.cs index 10e86cc..e660c2a 100644 --- a/PerfectNumbers.Core/MersenneDivisorCycles.cs +++ b/PerfectNumbers.Core/MersenneDivisorCycles.cs @@ -304,8 +304,7 @@ public class MersenneDivisorCycles // return false; // } - // Each divisor candidate is evaluated once along the scan, so a fresh factor map per call avoids keeping an unused cache. - Dictionary factorCounts = new Dictionary(); + Dictionary factorCounts = ThreadStaticPools.RentFactorCountDictionary(); bool success = false; if (twoCount > 0) @@ -323,6 +322,7 @@ public class MersenneDivisorCycles primeOrderFailed = false; } + ThreadStaticPools.ReturnFactorCountDictionary(factorCounts); return success; } @@ -335,11 +335,11 @@ public class MersenneDivisorCycles return true; } - // Factorizations are local to a single reduction attempt; there is no reuse across calls. - Dictionary scratch = new Dictionary(); + Dictionary scratch = ThreadStaticPools.RentFactorScratchDictionary(); bool success = TryFactorIntoCountsInternal(value, scratch); if (!success) { + ThreadStaticPools.ReturnFactorScratchDictionary(scratch); return false; } @@ -349,6 +349,7 @@ public class MersenneDivisorCycles AddFactor(counts, pair.Key, pair.Value); } + ThreadStaticPools.ReturnFactorScratchDictionary(scratch); return true; } @@ -468,21 +469,19 @@ public class MersenneDivisorCycles Array.Sort(primes, 0, count); + const int StackMultiplicityCapacity = 16; + const int HeapMultiplicityCapacity = 256; + + // Scanner workloads keep almost all multiplicities within the stack-backed span, so the stack path handles the hot loop. + Span stackCandidates = stackalloc ulong[StackMultiplicityCapacity]; + Span stackEvaluations = stackalloc bool[StackMultiplicityCapacity]; + + ulong[]? heapCandidateArray = null; + bool[]? heapEvaluationArray = null; + + ExponentRemainderStepper stepper = ThreadStaticPools.RentExponentStepper(divisorData); + ulong order = initialOrder; - - int stackMultiplicity = 0; - for (int i = 0; i < count; i++) - { - int multiplicity = factorCounts[primes[i]]; - if (multiplicity <= 16 && multiplicity > stackMultiplicity) - { - stackMultiplicity = multiplicity; - } - } - - Span stackCandidates = stackMultiplicity > 0 ? stackalloc ulong[stackMultiplicity] : Span.Empty; - Span stackEvaluations = stackMultiplicity > 0 ? stackalloc bool[stackMultiplicity] : Span.Empty; - for (int i = 0; i < count; i++) { ulong prime = primes[i]; @@ -494,25 +493,52 @@ public class MersenneDivisorCycles // continue; // } - if (multiplicity <= 16) + if (multiplicity <= StackMultiplicityCapacity) { - ProcessReduceOrderPrime(stackCandidates.Slice(0, multiplicity), stackEvaluations.Slice(0, multiplicity), ref order, prime, multiplicity, divisorData); + ProcessReduceOrderPrime(stackCandidates.Slice(0, multiplicity), stackEvaluations.Slice(0, multiplicity), ref order, prime, multiplicity, ref stepper); + continue; } - else + + if (multiplicity > HeapMultiplicityCapacity) { - ulong[] rentedCandidates = ThreadStaticPools.UlongPool.Rent(multiplicity); - bool[] rentedEvaluations = ThreadStaticPools.BoolPool.Rent(multiplicity); - ProcessReduceOrderPrime(rentedCandidates.AsSpan(0, multiplicity), rentedEvaluations.AsSpan(0, multiplicity), ref order, prime, multiplicity, divisorData); - ThreadStaticPools.UlongPool.Return(rentedCandidates, clearArray: false); - ThreadStaticPools.BoolPool.Return(rentedEvaluations, clearArray: false); + if (heapCandidateArray is not null) + { + ThreadStaticPools.UlongPool.Return(heapCandidateArray); + } + + if (heapEvaluationArray is not null) + { + ThreadStaticPools.BoolPool.Return(heapEvaluationArray); + } + + ThreadStaticPools.ReturnExponentStepper(stepper); + ThreadStaticPools.UlongPool.Return(primes); + throw new InvalidOperationException($"Factor multiplicity {multiplicity} exceeds the supported limit of {HeapMultiplicityCapacity}."); } + + // Larger multiplicities appear during stress tests; reuse a shared heap buffer in those cases. + heapCandidateArray ??= ThreadStaticPools.UlongPool.Rent(HeapMultiplicityCapacity); + heapEvaluationArray ??= ThreadStaticPools.BoolPool.Rent(HeapMultiplicityCapacity); + ProcessReduceOrderPrime(heapCandidateArray.AsSpan(0, multiplicity), heapEvaluationArray.AsSpan(0, multiplicity), ref order, prime, multiplicity, ref stepper); + } + + ThreadStaticPools.ReturnExponentStepper(stepper); + ThreadStaticPools.UlongPool.Return(primes); + + if (heapCandidateArray is not null) + { + ThreadStaticPools.UlongPool.Return(heapCandidateArray); + } + + if (heapEvaluationArray is not null) + { + ThreadStaticPools.BoolPool.Return(heapEvaluationArray); } - ThreadStaticPools.UlongPool.Return(primes, clearArray: false); return order; } - private static void ProcessReduceOrderPrime(Span candidateBuffer, Span evaluationBuffer, ref ulong order, ulong prime, int multiplicity, in MontgomeryDivisorData divisorData) + private static void ProcessReduceOrderPrime(Span candidateBuffer, Span evaluationBuffer, ref ulong order, ulong prime, int multiplicity, ref ExponentRemainderStepper stepper) { ulong working = order; int actual = 0; @@ -549,7 +575,8 @@ public class MersenneDivisorCycles candidates[right] = tmp; } - ExponentRemainderStepper stepper = ThreadStaticPools.RentExponentStepper(divisorData); + stepper.Reset(); + for (int j = 0; j < actual; j++) { evaluations[j] = stepper.ComputeNextIsUnity(candidates[j]); @@ -564,8 +591,6 @@ public class MersenneDivisorCycles order = candidates[j]; } - - ThreadStaticPools.ReturnExponentStepper(stepper); } private static ulong PollardRho64(ulong n) diff --git a/PerfectNumbers.Core/PrimeOrderCalculator.Cpu.cs b/PerfectNumbers.Core/PrimeOrderCalculator.Cpu.cs index 5cc4ed4..296279f 100644 --- a/PerfectNumbers.Core/PrimeOrderCalculator.Cpu.cs +++ b/PerfectNumbers.Core/PrimeOrderCalculator.Cpu.cs @@ -220,6 +220,9 @@ internal static partial class PrimeOrderCalculator Span stackCandidates = stackalloc ulong[StackExponentCapacity]; Span stackEvaluations = stackalloc bool[StackExponentCapacity]; + ulong[]? heapCandidateArray = null; + bool[]? heapEvaluationArray = null; + for (int i = 0; i < length; i++) { ulong primeFactor = buffer[i].Value; @@ -227,12 +230,23 @@ internal static partial class PrimeOrderCalculator // Factor exponents produced by Pollard-Rho and trial division are always positive in the scanner flow. // if (exponent <= 0) // { - // continue; + // continue; // } if (exponent > ExponentHardLimit) { ThreadStaticPools.ReturnExponentStepper(stepper); + + if (heapCandidateArray is not null) + { + ThreadStaticPools.UlongPool.Return(heapCandidateArray); + } + + if (heapEvaluationArray is not null) + { + ThreadStaticPools.BoolPool.Return(heapEvaluationArray); + } + pool.Return(tempArray, clearArray: false); throw new InvalidOperationException($"Prime factor exponent {exponent} exceeds the supported limit of {ExponentHardLimit}."); } @@ -240,18 +254,26 @@ internal static partial class PrimeOrderCalculator if (exponent <= StackExponentCapacity) { ProcessExponentLoweringPrime(stackCandidates.Slice(0, exponent), stackEvaluations.Slice(0, exponent), ref order, primeFactor, exponent, ref stepper); + continue; } - else - { - ulong[] rentedCandidates = ThreadStaticPools.UlongPool.Rent(exponent); - bool[] rentedEvaluations = ThreadStaticPools.BoolPool.Rent(exponent); - ProcessExponentLoweringPrime(rentedCandidates.AsSpan(0, exponent), rentedEvaluations.AsSpan(0, exponent), ref order, primeFactor, exponent, ref stepper); - ThreadStaticPools.UlongPool.Return(rentedCandidates, clearArray: false); - ThreadStaticPools.BoolPool.Return(rentedEvaluations, clearArray: false); - } + + heapCandidateArray ??= ThreadStaticPools.UlongPool.Rent(ExponentHardLimit); + heapEvaluationArray ??= ThreadStaticPools.BoolPool.Rent(ExponentHardLimit); + ProcessExponentLoweringPrime(heapCandidateArray.AsSpan(0, exponent), heapEvaluationArray.AsSpan(0, exponent), ref order, primeFactor, exponent, ref stepper); } ThreadStaticPools.ReturnExponentStepper(stepper); + + if (heapCandidateArray is not null) + { + ThreadStaticPools.UlongPool.Return(heapCandidateArray); + } + + if (heapEvaluationArray is not null) + { + ThreadStaticPools.BoolPool.Return(heapEvaluationArray); + } + pool.Return(tempArray, clearArray: false); return order; } @@ -374,6 +396,10 @@ internal static partial class PrimeOrderCalculator Span stackBuffer = stackalloc ulong[StackExponentCapacity]; + ulong[]? heapCandidateArray = null; + ExponentRemainderStepper stepper = ThreadStaticPools.RentExponentStepper(divisorData); + bool violates = false; + for (int i = 0; i < length; i++) { ulong primeFactor = span[i].Value; @@ -381,37 +407,51 @@ internal static partial class PrimeOrderCalculator // Factor exponents emitted by TryConfirmCandidateCpu stay positive for production workloads. // if (exponent <= 0) // { - // continue; + // continue; // } if (exponent > ExponentHardLimit) { + ThreadStaticPools.ReturnExponentStepper(stepper); + + if (heapCandidateArray is not null) + { + ThreadStaticPools.UlongPool.Return(heapCandidateArray); + } + throw new InvalidOperationException($"Factor exponent {exponent} exceeds the supported limit of {ExponentHardLimit}."); } if (exponent <= StackExponentCapacity) { - if (ValidateOrderForFactor(stackBuffer.Slice(0, exponent), primeFactor, exponent, order, divisorData)) + if (ValidateOrderForFactor(stackBuffer.Slice(0, exponent), primeFactor, exponent, order, ref stepper)) { - return false; + violates = true; + break; } + + continue; } - else + + heapCandidateArray ??= ThreadStaticPools.UlongPool.Rent(ExponentHardLimit); + if (ValidateOrderForFactor(heapCandidateArray.AsSpan(0, exponent), primeFactor, exponent, order, ref stepper)) { - ulong[] rentedCandidates = ThreadStaticPools.UlongPool.Rent(exponent); - if (ValidateOrderForFactor(rentedCandidates.AsSpan(0, exponent), primeFactor, exponent, order, divisorData)) - { - ThreadStaticPools.UlongPool.Return(rentedCandidates, clearArray: false); - return false; - } - ThreadStaticPools.UlongPool.Return(rentedCandidates, clearArray: false); + violates = true; + break; } } - return true; + ThreadStaticPools.ReturnExponentStepper(stepper); + + if (heapCandidateArray is not null) + { + ThreadStaticPools.UlongPool.Return(heapCandidateArray); + } + + return !violates; } - private static bool ValidateOrderForFactor(Span buffer, ulong primeFactor, int exponent, ulong order, in MontgomeryDivisorData divisorData) + private static bool ValidateOrderForFactor(Span buffer, ulong primeFactor, int exponent, ulong order, ref ExponentRemainderStepper stepper) { ulong working = order; int actual = 0; @@ -441,17 +481,16 @@ internal static partial class PrimeOrderCalculator Span candidates = buffer.Slice(0, actual); candidates.Reverse(); - ExponentRemainderStepper stepper = ThreadStaticPools.RentExponentStepper(divisorData); + stepper.Reset(); + for (int j = 0; j < actual; j++) { if (stepper.ComputeNextIsUnity(candidates[j])) { - ThreadStaticPools.ReturnExponentStepper(stepper); return true; } } - ThreadStaticPools.ReturnExponentStepper(stepper); return false; } @@ -833,6 +872,9 @@ internal static partial class PrimeOrderCalculator const int ExponentHardLimit = 256; Span stackBuffer = stackalloc ulong[StackExponentCapacity]; + ulong[]? heapCandidateArray = null; + ExponentRemainderStepper stepper = ThreadStaticPools.RentExponentStepper(divisorData); + bool violates = false; for (int i = 0; i < length; i++) { @@ -841,34 +883,48 @@ internal static partial class PrimeOrderCalculator // Factorization entries for candidate confirmation are always positive here. // if (exponent <= 0) // { - // continue; + // continue; // } if (exponent > ExponentHardLimit) { + ThreadStaticPools.ReturnExponentStepper(stepper); + + if (heapCandidateArray is not null) + { + ThreadStaticPools.UlongPool.Return(heapCandidateArray); + } + throw new InvalidOperationException($"Candidate factor exponent {exponent} exceeds the supported limit of {ExponentHardLimit}."); } if (exponent <= StackExponentCapacity) { - if (CheckCandidateViolation(stackBuffer.Slice(0, exponent), primeFactor, exponent, candidate, prime, divisorData, ref powUsed, powBudget)) + if (CheckCandidateViolation(stackBuffer.Slice(0, exponent), primeFactor, exponent, candidate, prime, ref powUsed, powBudget, ref stepper)) { - return false; + violates = true; + break; } + + continue; } - else + + heapCandidateArray ??= ThreadStaticPools.UlongPool.Rent(ExponentHardLimit); + if (CheckCandidateViolation(heapCandidateArray.AsSpan(0, exponent), primeFactor, exponent, candidate, prime, ref powUsed, powBudget, ref stepper)) { - ulong[] rentedCandidates = ThreadStaticPools.UlongPool.Rent(exponent); - if (CheckCandidateViolation(rentedCandidates.AsSpan(0, exponent), primeFactor, exponent, candidate, prime, divisorData, ref powUsed, powBudget)) - { - ThreadStaticPools.UlongPool.Return(rentedCandidates, clearArray: false); - return false; - } - ThreadStaticPools.UlongPool.Return(rentedCandidates, clearArray: false); + violates = true; + break; } } - return true; + ThreadStaticPools.ReturnExponentStepper(stepper); + + if (heapCandidateArray is not null) + { + ThreadStaticPools.UlongPool.Return(heapCandidateArray); + } + + return !violates; } finally { @@ -876,7 +932,7 @@ internal static partial class PrimeOrderCalculator } } - private static bool CheckCandidateViolation(Span buffer, ulong primeFactor, int exponent, ulong candidate, ulong prime, in MontgomeryDivisorData divisorData, ref int powUsed, int powBudget) + private static bool CheckCandidateViolation(Span buffer, ulong primeFactor, int exponent, ulong candidate, ulong prime, ref int powUsed, int powBudget, ref ExponentRemainderStepper stepper) { ulong working = candidate; int actual = 0; @@ -906,28 +962,25 @@ internal static partial class PrimeOrderCalculator Span candidates = buffer.Slice(0, actual); candidates.Reverse(); - ExponentRemainderStepper stepper = ThreadStaticPools.RentExponentStepper(divisorData); + stepper.Reset(); + for (int j = 0; j < actual; j++) { if (powUsed >= powBudget && powBudget > 0) { - ThreadStaticPools.ReturnExponentStepper(stepper); return true; } powUsed++; if (stepper.ComputeNextIsUnity(candidates[j])) { - ThreadStaticPools.ReturnExponentStepper(stepper); return true; } } - ThreadStaticPools.ReturnExponentStepper(stepper); return false; } - private static bool Pow2EqualsOneCpu(ulong exponent, ulong prime, in MontgomeryDivisorData divisorData) { if (IsGpuPow2Allowed) diff --git a/PerfectNumbers.Core/ThreadStaticPools.cs b/PerfectNumbers.Core/ThreadStaticPools.cs index 8ba6845..8191676 100644 --- a/PerfectNumbers.Core/ThreadStaticPools.cs +++ b/PerfectNumbers.Core/ThreadStaticPools.cs @@ -251,7 +251,6 @@ namespace PerfectNumbers.Core pool.Add(stack); } - // The factor-count scratch caches were removed because divisor factorizations never repeat across scan iterations. [ThreadStatic] private static List>? _ulongIntDictionaryPool; @@ -294,6 +293,54 @@ namespace PerfectNumbers.Core pool.Add(dictionary); } + [ThreadStatic] + private static Dictionary? _factorCountDictionary; + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static Dictionary RentFactorCountDictionary() + { + Dictionary? dictionary = _factorCountDictionary; + if (dictionary is null) + { + return new Dictionary(8); + } + + _factorCountDictionary = null; + dictionary.Clear(); + return dictionary; + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static void ReturnFactorCountDictionary(Dictionary dictionary) + { + dictionary.Clear(); + _factorCountDictionary = dictionary; + } + + [ThreadStatic] + private static Dictionary? _factorScratchDictionary; + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static Dictionary RentFactorScratchDictionary() + { + Dictionary? dictionary = _factorScratchDictionary; + if (dictionary is null) + { + return new Dictionary(8); + } + + _factorScratchDictionary = null; + dictionary.Clear(); + return dictionary; + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static void ReturnFactorScratchDictionary(Dictionary dictionary) + { + dictionary.Clear(); + _factorScratchDictionary = dictionary; + } + [ThreadStatic] private static ExponentRemainderStepper _exponentRemainderStepper;