From 8a6fc0c6a5cffc4748acb9f2366ab3de15d20785 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Marek=20Lema=C5=84czyk?= Date: Sun, 21 Sep 2025 02:25:18 +0200 Subject: [PATCH] Revert "Merge pull request #24 from mlemanczyk/codex/modify-evenperfectbitscanner-for-residue-mode" This reverts commit 2d024cc70836b8e3a468f61187f2f0498086fe20, reversing changes made to b1c29ad16e6ab6c308a71db7b094725a9b39ef38. --- EvenPerfectBitScanner.Tests/ProgramTests.cs | 93 +----- EvenPerfectBitScanner/Program.cs | 49 +-- .../MersenneNumberResidueCpuTesterTests.cs | 12 +- .../MersenneNumberResidueGpuTesterTests.cs | 12 +- .../Cpu/MersenneNumberResidueCpuTester.cs | 302 ++++++++---------- .../Gpu/MersenneNumberResidueGpuTester.cs | 136 +++----- PerfectNumbers.Core/MersenneNumberTester.cs | 212 ++++-------- 7 files changed, 265 insertions(+), 551 deletions(-) diff --git a/EvenPerfectBitScanner.Tests/ProgramTests.cs b/EvenPerfectBitScanner.Tests/ProgramTests.cs index 6a0e429..122f26c 100644 --- a/EvenPerfectBitScanner.Tests/ProgramTests.cs +++ b/EvenPerfectBitScanner.Tests/ProgramTests.cs @@ -165,95 +165,6 @@ public class ProgramTests } } - [Fact] - public void Divisor_mode_sets_detailedCheck_when_divisor_is_found() - { - var useDivisorField = typeof(Program).GetField("_useDivisor", BindingFlags.NonPublic | BindingFlags.Static)!; - var divisorField = typeof(Program).GetField("_divisor", BindingFlags.NonPublic | BindingFlags.Static)!; - var testerField = typeof(Program).GetField("_divisorTester", BindingFlags.NonPublic | BindingFlags.Static)!; - var candidatesField = typeof(MersenneNumberDivisorGpuTester).GetField("_divisorCandidates", BindingFlags.NonPublic | BindingFlags.Static)!; - var primeField = typeof(Program).GetField("PrimeTesters", BindingFlags.NonPublic | BindingFlags.Static)!; - var residueField = typeof(Program).GetField("PResidue", BindingFlags.NonPublic | BindingFlags.Static)!; - var forceCpuProp = typeof(GpuContextPool).GetProperty("ForceCpu"); - - useDivisorField.SetValue(null, true); - divisorField.SetValue(null, UInt128.Zero); - testerField.SetValue(null, new MersenneNumberDivisorGpuTester()); - primeField.SetValue(null, new ThreadLocal(() => new PrimeTester(), trackAllValues: true)); - residueField.SetValue(null, new ThreadLocal(() => new ModResidueTracker(ResidueModel.Identity, 2UL, true), trackAllValues: true)); - candidatesField.SetValue(null, new (ulong, uint)[] { (23UL, 11U) }); - forceCpuProp!.SetValue(null, true); - - try - { - Program.IsEvenPerfectCandidate(11UL, 32UL, out bool searched, out bool detailedCheck).Should().BeFalse(); - searched.Should().BeTrue(); - detailedCheck.Should().BeTrue(); - } - finally - { - useDivisorField.SetValue(null, false); - testerField.SetValue(null, null); - candidatesField.SetValue(null, null); - primeField.SetValue(null, null); - residueField.SetValue(null, null); - forceCpuProp!.SetValue(null, false); - } - } - - [Fact] - public void Residue_mode_sets_detailedCheck_based_on_exhaustion() - { - var useDivisorField = typeof(Program).GetField("_useDivisor", BindingFlags.NonPublic | BindingFlags.Static)!; - var residueModeField = typeof(Program).GetField("_useResidueMode", BindingFlags.NonPublic | BindingFlags.Static)!; - var mersenneField = typeof(Program).GetField("MersenneTesters", BindingFlags.NonPublic | BindingFlags.Static)!; - var primeField = typeof(Program).GetField("PrimeTesters", BindingFlags.NonPublic | BindingFlags.Static)!; - var residueField = typeof(Program).GetField("PResidue", BindingFlags.NonPublic | BindingFlags.Static)!; - var forceCpuProp = typeof(GpuContextPool).GetProperty("ForceCpu"); - - const ulong exponent = 7UL; - - useDivisorField.SetValue(null, false); - residueModeField.SetValue(null, true); - primeField.SetValue(null, new ThreadLocal(() => new PrimeTester(), trackAllValues: true)); - residueField.SetValue(null, new ThreadLocal(() => new ModResidueTracker(ResidueModel.Identity, exponent, true), trackAllValues: true)); - forceCpuProp!.SetValue(null, true); - - try - { - mersenneField.SetValue(null, new ThreadLocal(() => new MersenneNumberTester( - useGpuScan: false, - useGpuOrder: false, - useResidue: true, - maxK: 0UL, - residueDivisorSets: 0UL), trackAllValues: true)); - - Program.IsEvenPerfectCandidate(exponent, 0UL, out bool searched, out bool detailedCheck).Should().BeTrue(); - searched.Should().BeTrue(); - detailedCheck.Should().BeFalse(); - - mersenneField.SetValue(null, new ThreadLocal(() => new MersenneNumberTester( - useGpuScan: false, - useGpuOrder: false, - useResidue: true, - maxK: 4UL, - residueDivisorSets: 2UL), trackAllValues: true)); - - Program.IsEvenPerfectCandidate(exponent, 0UL, out searched, out detailedCheck).Should().BeTrue(); - searched.Should().BeTrue(); - detailedCheck.Should().BeTrue(); - } - finally - { - mersenneField.SetValue(null, null); - primeField.SetValue(null, null); - residueField.SetValue(null, null); - residueModeField.SetValue(null, false); - forceCpuProp!.SetValue(null, false); - useDivisorField.SetValue(null, false); - } - } - [Fact] public void Divisor_mode_accepts_large_search_limit() { @@ -335,7 +246,7 @@ public class ProgramTests var forceCpuProp = typeof(GpuContextPool).GetProperty("ForceCpu"); useDivisorField.SetValue(null, false); - mersenneField.SetValue(null, new ThreadLocal(() => new MersenneNumberTester(useIncremental: true, useResidue: true, maxK: 1_000UL, residueDivisorSets: 1UL), trackAllValues: true)); + mersenneField.SetValue(null, new ThreadLocal(() => new MersenneNumberTester(useIncremental: true, useResidue: true, maxK: 1_000UL), trackAllValues: true)); primeField.SetValue(null, new ThreadLocal(() => new PrimeTester(), trackAllValues: true)); residueField.SetValue(null, new ThreadLocal(() => new ModResidueTracker(ResidueModel.Identity, 2UL, true), trackAllValues: true)); forceCpuProp!.SetValue(null, false); @@ -344,7 +255,7 @@ public class ProgramTests { Program.IsEvenPerfectCandidate(29UL, 0UL, out bool searched, out bool detailed).Should().BeFalse(); searched.Should().BeTrue(); - detailed.Should().BeTrue(); + detailed.Should().BeFalse(); Program.IsEvenPerfectCandidate(31UL, 0UL, out searched, out detailed).Should().BeTrue(); searched.Should().BeTrue(); diff --git a/EvenPerfectBitScanner/Program.cs b/EvenPerfectBitScanner/Program.cs index 16e984e..dbd683e 100644 --- a/EvenPerfectBitScanner/Program.cs +++ b/EvenPerfectBitScanner/Program.cs @@ -451,17 +451,16 @@ internal static class Program MersenneTesters = new ThreadLocal(() => { var tester = new MersenneNumberTester( - useIncremental: !useLucas, - useOrderCache: false, - kernelType: kernelType, - useModuloWorkaround: useModuloWorkaround, - useOrder: useOrder, - useGpuLucas: mersenneOnGpu, - useGpuScan: mersenneOnGpu, - useGpuOrder: orderOnGpu, - useResidue: useResidue, - maxK: residueKMax, - residueDivisorSets: divisorCyclesSearchLimit128); + useIncremental: !useLucas, + useOrderCache: false, + kernelType: kernelType, + useModuloWorkaround: useModuloWorkaround, + useOrder: useOrder, + useGpuLucas: mersenneOnGpu, + useGpuScan: mersenneOnGpu, + useGpuOrder: orderOnGpu, + useResidue: useResidue, + maxK: residueKMax); if (!useLucas) { tester.WarmUpOrders(currentP, _orderWarmupLimitOverride ?? 5_000_000UL); @@ -767,7 +766,7 @@ internal static class Program Console.WriteLine(" --divisor-cycles-device=cpu|gpu device for cycles generation (default gpu)"); Console.WriteLine(" --divisor-cycles-batch-size= batch size for cycles generation (default 512)"); Console.WriteLine(" --divisor-cycles-continue continue divisor cycles generation"); - Console.WriteLine(" --divisor-cycles-limit= cycle search iterations when --mersenne=divisor or residue (UInt128 supported; residue path honors full range)"); + Console.WriteLine(" --divisor-cycles-limit= cycle search iterations when --mersenne=divisor"); Console.WriteLine(" --use-order test primality via q order"); Console.WriteLine(" --workaround-mod avoid '%' operator on the GPU"); @@ -1127,30 +1126,8 @@ internal static class Program return _divisorTester!.IsPrime(p, _divisor, divisorCyclesSearchLimit, out detailedCheck); } - bool mersennePrime = MersenneTesters.Value!.IsMersennePrime(p, out bool residueExhausted); - detailedCheck = residueExhausted; - return mersennePrime; - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - internal static UInt128 ResolveResidueMaxK( - bool useResidue, - bool residueMaxExplicit, - bool divisorLimitExplicit, - UInt128 residueMaxK, - UInt128 divisorLimit) - { - if (!useResidue) - { - return residueMaxK; - } - - if (residueMaxExplicit || !divisorLimitExplicit) - { - return residueMaxK; - } - - return divisorLimit; + detailedCheck = MersenneTesters.Value!.IsMersennePrime(p); + return detailedCheck; } // Use ModResidueTracker with a small set of primes to pre-filter composite p. diff --git a/PerfectNumbers.Core.Tests/MersenneNumberResidueCpuTesterTests.cs b/PerfectNumbers.Core.Tests/MersenneNumberResidueCpuTesterTests.cs index c6fb890..883c5b2 100644 --- a/PerfectNumbers.Core.Tests/MersenneNumberResidueCpuTesterTests.cs +++ b/PerfectNumbers.Core.Tests/MersenneNumberResidueCpuTesterTests.cs @@ -69,18 +69,8 @@ public class MersenneNumberResidueCpuTesterTests bool expectedPrime) { bool isPrime = true; - bool exhausted = false; - tester.Scan( - exponent, - (UInt128)exponent << 1, - LastDigitIsSeven(exponent), - (UInt128)perSetLimit, - (UInt128)setCount, - (UInt128)overallLimit, - ref isPrime, - ref exhausted); + tester.Scan(exponent, (UInt128)exponent << 1, LastDigitIsSeven(exponent), (UInt128)maxK, ref isPrime); isPrime.Should().Be(expectedPrime); - exhausted.Should().BeTrue(); } } diff --git a/PerfectNumbers.Core.Tests/MersenneNumberResidueGpuTesterTests.cs b/PerfectNumbers.Core.Tests/MersenneNumberResidueGpuTesterTests.cs index 5b15e5d..1085799 100644 --- a/PerfectNumbers.Core.Tests/MersenneNumberResidueGpuTesterTests.cs +++ b/PerfectNumbers.Core.Tests/MersenneNumberResidueGpuTesterTests.cs @@ -288,18 +288,8 @@ public class MersenneNumberResidueGpuTesterTests bool expectedPrime) { bool isPrime = true; - bool exhausted = false; - tester.Scan( - exponent, - (UInt128)exponent << 1, - LastDigitIsSeven(exponent), - (UInt128)perSetLimit, - (UInt128)setCount, - (UInt128)overallLimit, - ref isPrime, - ref exhausted); + tester.Scan(exponent, (UInt128)exponent << 1, LastDigitIsSeven(exponent), (UInt128)maxK, ref isPrime); isPrime.Should().Be(expectedPrime); - exhausted.Should().BeTrue(); } } diff --git a/PerfectNumbers.Core/Cpu/MersenneNumberResidueCpuTester.cs b/PerfectNumbers.Core/Cpu/MersenneNumberResidueCpuTester.cs index 6bc8966..9d912c5 100644 --- a/PerfectNumbers.Core/Cpu/MersenneNumberResidueCpuTester.cs +++ b/PerfectNumbers.Core/Cpu/MersenneNumberResidueCpuTester.cs @@ -5,32 +5,25 @@ public class MersenneNumberResidueCpuTester private ModResidueTracker? _mersenneResidueTracker; // CPU residue variant using tracker + unrolled residue updates (no method calls in the loop). - public void Scan( - ulong exponent, - UInt128 twoP, - bool lastIsSeven, - UInt128 perSetLimit, - UInt128 setCount, - UInt128 overallLimit, - ref bool isPrime, - ref bool divisorsExhausted) - { - if (setCount == UInt128.Zero || perSetLimit == UInt128.Zero || overallLimit == UInt128.Zero) - { - return; - } + public void Scan(ulong exponent, UInt128 twoP, bool lastIsSeven, UInt128 maxK, ref bool isPrime) + { + // Initialize/update Mersenne residue tracker and start a merge walk for ascending q divisors - // Initialize/update Mersenne residue tracker and start a merge walk for ascending q divisors + _mersenneResidueTracker ??= new(ResidueModel.Mersenne); + // For small divisors q <= 4,000,000 use cached cycle length to fast-path checks + // and avoid per-q powmods when possible when scanning lanes. + _mersenneResidueTracker.BeginMerge(exponent); - _mersenneResidueTracker ??= new(ResidueModel.Mersenne); - // For small divisors q <= 4,000,000 use cached cycle length to fast-path checks - // and avoid per-q powmods when possible when scanning lanes. - _mersenneResidueTracker.BeginMerge(exponent); - - ulong step10 = (ulong)(twoP % 10UL) % 10UL; - ulong step8 = (ulong)(twoP % 8UL); - ulong step3 = (ulong)(twoP % 3UL); - ulong step5 = (ulong)(twoP % 5UL); + // Start at k = 1: q0 = 2*p*1 + 1 + UInt128 k = UInt128.One; + UInt128 q = twoP + UInt128.One; + // Current residues + q.Mod10_8_5_3(out ulong r10, out ulong r8, out ulong r5, out ulong r3); + // Step residues for q += 2*p + ulong step10 = (ulong)(twoP % 10UL) % 10UL; + ulong step8 = (ulong)(twoP % 8UL); + ulong step3 = (ulong)(twoP % 3UL); + ulong step5 = (ulong)(twoP % 5UL); // Allowed last-digit sets for q depending on last digit of M_p // lastIsSeven == true => allow {7,9} @@ -38,162 +31,139 @@ public class MersenneNumberResidueCpuTester int allowMask = lastIsSeven ? CpuConstants.LastSevenMask10 : CpuConstants.LastOneMask10; // Predeclare temps to avoid redeclaration overhead in the loop - UInt128 remaining; - int i, lanes; - // per-lane residues - ulong r10_0, r10_1, r10_2, r10_3, r8_0, r8_1, r8_2, r8_3, r3_0, r3_1, r3_2, r3_3, r5_0, r5_1, r5_2, r5_3; - // q lanes - UInt128 twoP2 = twoP + twoP, - twoP3 = twoP2 + twoP, - twoP4 = twoP3 + twoP, - qIncremental; + UInt128 remaining; + int i, lanes; + // per-lane residues + ulong r10_0, r10_1, r10_2, r10_3, r8_0, r8_1, r8_2, r8_3, r3_0, r3_1, r3_2, r3_3, r5_0, r5_1, r5_2, r5_3; + // q lanes + UInt128 twoP2 = twoP + twoP, + twoP3 = twoP2 + twoP, + twoP4 = twoP3 + twoP, + qIncremental; - bool localIsPrime, isDivider; + bool localIsPrime, isDivider; - ModResidueTracker tracker = _mersenneResidueTracker!; - UInt128 limitInclusive = overallLimit + UInt128.One; - for (UInt128 setIndex = UInt128.Zero; setIndex < setCount && Volatile.Read(ref isPrime); setIndex++) - { - UInt128 setOffset = perSetLimit * setIndex; - UInt128 k = setOffset + UInt128.One; - if (k >= limitInclusive) - { - break; - } + maxK++; - UInt128 setLimitExclusive = setOffset + perSetLimit + UInt128.One; - if (setLimitExclusive > limitInclusive) - { - setLimitExclusive = limitInclusive; - } + ModResidueTracker tracker = _mersenneResidueTracker!; + while (k < maxK && (localIsPrime = Volatile.Read(ref isPrime))) + { + remaining = maxK - k; + if (remaining == UInt128.Zero) + { + break; + } - UInt128 q = twoP * k + UInt128.One; - q.Mod10_8_5_3(out ulong r10, out ulong r8, out ulong r5, out ulong r3); + lanes = remaining > 4UL ? 4 : (int)remaining; - while (k < setLimitExclusive && (localIsPrime = Volatile.Read(ref isPrime))) - { - remaining = setLimitExclusive - k; - if (remaining == UInt128.Zero) - { - break; - } + // Compute lane residues incrementally (avoid multiplications, minimize reductions) + // r10 + r10_0 = r10; + r10_1 = r10_0 + step10; if (r10_1 >= 20UL) r10_1 -= 20UL; if (r10_1 >= 10UL) r10_1 -= 10UL; + r10_2 = r10_1 + step10; if (r10_2 >= 20UL) r10_2 -= 20UL; if (r10_2 >= 10UL) r10_2 -= 10UL; + r10_3 = r10_2 + step10; if (r10_3 >= 20UL) r10_3 -= 20UL; if (r10_3 >= 10UL) r10_3 -= 10UL; - lanes = remaining > 4UL ? 4 : (int)remaining; - if (lanes <= 0) - { - break; - } + // r8 + r8_0 = r8 & 7UL; + r8_1 = (r8_0 + step8) & 7UL; + r8_2 = (r8_1 + step8) & 7UL; + r8_3 = (r8_2 + step8) & 7UL; - // Compute lane residues incrementally (avoid multiplications, minimize reductions) - // r10 - r10_0 = r10; - r10_1 = r10_0 + step10; if (r10_1 >= 20UL) r10_1 -= 20UL; if (r10_1 >= 10UL) r10_1 -= 10UL; - r10_2 = r10_1 + step10; if (r10_2 >= 20UL) r10_2 -= 20UL; if (r10_2 >= 10UL) r10_2 -= 10UL; - r10_3 = r10_2 + step10; if (r10_3 >= 20UL) r10_3 -= 20UL; if (r10_3 >= 10UL) r10_3 -= 10UL; + // r3 (pattern 0,1,2,0) + r3_0 = r3; if (r3_0 >= 3UL) r3_0 -= 3UL; + r3_1 = r3_0 + step3; if (r3_1 >= 3UL) r3_1 -= 3UL; + r3_2 = r3_1 + step3; if (r3_2 >= 3UL) r3_2 -= 3UL; + r3_3 = r3; if (r3_3 >= 3UL) r3_3 -= 3UL; - // r8 - r8_0 = r8 & 7UL; - r8_1 = (r8_0 + step8) & 7UL; - r8_2 = (r8_1 + step8) & 7UL; - r8_3 = (r8_2 + step8) & 7UL; + // r5 (pattern 0,1,2,3) + r5_0 = r5; if (r5_0 >= 5UL) r5_0 -= 5UL; + r5_1 = r5_0 + step5; if (r5_1 >= 10UL) r5_1 -= 10UL; if (r5_1 >= 5UL) r5_1 -= 5UL; + r5_2 = r5_1 + step5; if (r5_2 >= 10UL) r5_2 -= 10UL; if (r5_2 >= 5UL) r5_2 -= 5UL; + r5_3 = r5_2 + step5; if (r5_3 >= 10UL) r5_3 -= 10UL; if (r5_3 >= 5UL) r5_3 -= 5UL; - // r3 (pattern 0,1,2,0) - r3_0 = r3; if (r3_0 >= 3UL) r3_0 -= 3UL; - r3_1 = r3_0 + step3; if (r3_1 >= 3UL) r3_1 -= 3UL; - r3_2 = r3_1 + step3; if (r3_2 >= 3UL) r3_2 -= 3UL; - r3_3 = r3; if (r3_3 >= 3UL) r3_3 -= 3UL; + // lane 0 + if (lanes >= 1 && + (((allowMask >> (int)r10_0) & 1) != 0 && (r8_0 == 1UL || r8_0 == 7UL) && r3_0 != 0UL && r5_0 != 0UL)) + { + // cycle-based quick check for small q + if (localIsPrime && tracker.MergeOrAppend(exponent, q, out isDivider) && isDivider && + (q <= ulong.MaxValue ? ((ulong)q).IsPrimeCandidate() : q.IsPrimeCandidate())) + { + localIsPrime = false; + } + } - // r5 (pattern 0,1,2,3) - r5_0 = r5; if (r5_0 >= 5UL) r5_0 -= 5UL; - r5_1 = r5_0 + step5; if (r5_1 >= 10UL) r5_1 -= 10UL; if (r5_1 >= 5UL) r5_1 -= 5UL; - r5_2 = r5_1 + step5; if (r5_2 >= 10UL) r5_2 -= 10UL; if (r5_2 >= 5UL) r5_2 -= 5UL; - r5_3 = r5_2 + step5; if (r5_3 >= 10UL) r5_3 -= 10UL; if (r5_3 >= 5UL) r5_3 -= 5UL; + // lane 1 + qIncremental = q + twoP; + if (lanes >= 2 && + ((allowMask >> (int)r10_1) & 1) != 0 && (r8_1 == 1UL || r8_1 == 7UL) && r3_1 != 0UL && r5_1 != 0UL && + localIsPrime) + { + if (localIsPrime && tracker.MergeOrAppend(exponent, qIncremental, out isDivider) && isDivider && + (qIncremental <= ulong.MaxValue ? ((ulong)qIncremental).IsPrimeCandidate() : qIncremental.IsPrimeCandidate())) + { + localIsPrime = false; + } + } - // lane 0 - if (lanes >= 1 && - (((allowMask >> (int)r10_0) & 1) != 0 && (r8_0 == 1UL || r8_0 == 7UL) && r3_0 != 0UL && r5_0 != 0UL)) - { - if (localIsPrime && tracker.MergeOrAppend(exponent, q, out isDivider) && isDivider && - (q <= ulong.MaxValue ? ((ulong)q).IsPrimeCandidate() : q.IsPrimeCandidate())) - { - localIsPrime = false; - } - } + // lane 2 + qIncremental += twoP; + if (lanes >= 3 && + ((allowMask >> (int)r10_2) & 1) != 0 && (r8_2 == 1UL || r8_2 == 7UL) && r3_2 != 0UL && r5_2 != 0UL && + localIsPrime) + { + if (localIsPrime && tracker.MergeOrAppend(exponent, qIncremental, out isDivider) && isDivider && + (qIncremental <= ulong.MaxValue ? ((ulong)qIncremental).IsPrimeCandidate() : qIncremental.IsPrimeCandidate())) + { + localIsPrime = false; + } + } - // lane 1 - qIncremental = q + twoP; - if (lanes >= 2 && - ((allowMask >> (int)r10_1) & 1) != 0 && (r8_1 == 1UL || r8_1 == 7UL) && r3_1 != 0UL && r5_1 != 0UL && - localIsPrime) - { - if (localIsPrime && tracker.MergeOrAppend(exponent, qIncremental, out isDivider) && isDivider && - (qIncremental <= ulong.MaxValue ? ((ulong)qIncremental).IsPrimeCandidate() : qIncremental.IsPrimeCandidate())) - { - localIsPrime = false; - } - } + // lane 3 + qIncremental += twoP; + if (lanes >= 4 && + ((allowMask >> (int)r10_3) & 1) != 0 && (r8_3 == 1UL || r8_3 == 7UL) && r3_3 != 0UL && r5_3 != 0UL && + localIsPrime) + { + if (localIsPrime && tracker.MergeOrAppend(exponent, qIncremental, out isDivider) && isDivider && + (qIncremental <= ulong.MaxValue ? ((ulong)qIncremental).IsPrimeCandidate() : qIncremental.IsPrimeCandidate())) + { + localIsPrime = false; + } + } - // lane 2 - qIncremental += twoP; - if (lanes >= 3 && - ((allowMask >> (int)r10_2) & 1) != 0 && (r8_2 == 1UL || r8_2 == 7UL) && r3_2 != 0UL && r5_2 != 0UL && - localIsPrime) - { - if (localIsPrime && tracker.MergeOrAppend(exponent, qIncremental, out isDivider) && isDivider && - (qIncremental <= ulong.MaxValue ? ((ulong)qIncremental).IsPrimeCandidate() : qIncremental.IsPrimeCandidate())) - { - localIsPrime = false; - } - } + if (!localIsPrime) + { + Volatile.Write(ref isPrime, false); + return; + } + else if (!Volatile.Read(ref isPrime)) + { + return; + } - // lane 3 - qIncremental += twoP; - if (lanes >= 4 && - ((allowMask >> (int)r10_3) & 1) != 0 && (r8_3 == 1UL || r8_3 == 7UL) && r3_3 != 0UL && r5_3 != 0UL && - localIsPrime) - { - if (localIsPrime && tracker.MergeOrAppend(exponent, qIncremental, out isDivider) && isDivider && - (qIncremental <= ulong.MaxValue ? ((ulong)qIncremental).IsPrimeCandidate() : qIncremental.IsPrimeCandidate())) - { - localIsPrime = false; - } - } + // advance base by processed lanes + k += (UInt128)lanes; + q += lanes switch + { + 0 => 0UL, + 1 => twoP, + 2 => twoP2, + 3 => twoP3, + 4 => twoP4, + _ => throw new ArgumentException("Unsupproted value", nameof(twoP)) + }; - if (!localIsPrime) - { - Volatile.Write(ref isPrime, false); - divisorsExhausted = true; - return; - } - else if (!Volatile.Read(ref isPrime)) - { - divisorsExhausted = true; - return; - } + // advance residues by 'lanes' steps (consistent with base advance) + for (i = 0; i < lanes; i++) + { + r10 += step10; if (r10 >= 20UL) r10 -= 20UL; if (r10 >= 10UL) r10 -= 10UL; + r8 = (r8 + step8) & 7UL; + r3 += step3; if (r3 >= 3UL) r3 -= 3UL; + r5 += step5; if (r5 >= 10UL) r5 -= 10UL; if (r5 >= 5UL) r5 -= 5UL; + } + } - // advance base by processed lanes - k += (UInt128)lanes; - q += lanes switch - { - 0 => 0UL, - 1 => twoP, - 2 => twoP2, - 3 => twoP3, - 4 => twoP4, - _ => throw new ArgumentException("Unsupproted value", nameof(twoP)) - }; - - // advance residues by 'lanes' steps (consistent with base advance) - for (i = 0; i < lanes; i++) - { - r10 += step10; if (r10 >= 20UL) r10 -= 20UL; if (r10 >= 10UL) r10 -= 10UL; - r8 = (r8 + step8) & 7UL; - r3 += step3; if (r3 >= 3UL) r3 -= 3UL; - r5 += step5; if (r5 >= 10UL) r5 -= 10UL; if (r5 >= 5UL) r5 -= 5UL; - } - - divisorsExhausted = true; - } -} - } + return; + } } diff --git a/PerfectNumbers.Core/Gpu/MersenneNumberResidueGpuTester.cs b/PerfectNumbers.Core/Gpu/MersenneNumberResidueGpuTester.cs index 8c1783a..a9ee58d 100644 --- a/PerfectNumbers.Core/Gpu/MersenneNumberResidueGpuTester.cs +++ b/PerfectNumbers.Core/Gpu/MersenneNumberResidueGpuTester.cs @@ -12,36 +12,22 @@ public class MersenneNumberResidueGpuTester(bool useGpuOrder) // via GpuKernelPool to avoid host round-trips. // GPU residue variant: check 2^p % q == 1 for q = 2*p*k + 1. - public void Scan( - ulong exponent, - UInt128 twoP, - bool lastIsSeven, - UInt128 perSetLimit, - UInt128 setCount, - UInt128 overallLimit, - ref bool isPrime, - ref bool divisorsExhausted) - { - if (setCount == UInt128.Zero || perSetLimit == UInt128.Zero || overallLimit == UInt128.Zero) - { - return; - } - - var gpuLease = GpuKernelPool.GetKernel(_useGpuOrder); - var execution = gpuLease.EnterExecutionScope(); + public void Scan(ulong exponent, UInt128 twoP, bool lastIsSeven, UInt128 maxK, ref bool isPrime) + { + var gpuLease = GpuKernelPool.GetKernel(_useGpuOrder); var accelerator = gpuLease.Accelerator; - var stream = gpuLease.Stream; - // Ensure device has small cycles and primes tables for in-kernel lookup - var smallCyclesView = GpuKernelPool.EnsureSmallCyclesOnDevice(accelerator); - ResiduePrimeViews primeViews = GpuKernelPool.EnsureSmallPrimesOnDevice(accelerator); - int batchSize = GpuConstants.ScanBatchSize; - byte last = lastIsSeven ? (byte)1 : (byte)0; - var kernel = gpuLease.Pow2ModKernel; - ulong step10 = ((exponent % 10UL) << 1) % 10UL; - ulong step8 = ((exponent & 7UL) << 1) & 7UL; - ulong step3 = ((exponent % 3UL) << 1) % 3UL; - ulong step5 = ((exponent % 5UL) << 1) % 5UL; - GpuUInt128 twoPGpu = (GpuUInt128)twoP; + // Ensure device has small cycles and primes tables for in-kernel lookup + var smallCyclesView = GpuKernelPool.EnsureSmallCyclesOnDevice(accelerator); + ResiduePrimeViews primeViews = GpuKernelPool.EnsureSmallPrimesOnDevice(accelerator); + int batchSize = GpuConstants.ScanBatchSize; + UInt128 kStart = 1UL; + byte last = lastIsSeven ? (byte)1 : (byte)0; + var kernel = gpuLease.Pow2ModKernel; + ulong step10 = ((exponent % 10UL) << 1) % 10UL; + ulong step8 = ((exponent & 7UL) << 1) & 7UL; + ulong step3 = ((exponent % 3UL) << 1) % 3UL; + ulong step5 = ((exponent % 5UL) << 1) % 5UL; + GpuUInt128 twoPGpu = (GpuUInt128)twoP; var orderBuffer = accelerator.Allocate1D(batchSize); // Device memory returned by Allocate1D is not guaranteed to be zeroed. @@ -52,79 +38,57 @@ public class MersenneNumberResidueGpuTester(bool useGpuOrder) // testers such as the order kernels. orderBuffer.MemSetToZero(); ulong[] orderArray = ArrayPool.Shared.Rent(batchSize); - UInt128 batchSize128 = (UInt128)batchSize; - UInt128 limitInclusive = overallLimit == UInt128.MaxValue - ? overallLimit - : checked(overallLimit + UInt128.One); + UInt128 remaining; + int currentSize; + int i; + UInt128 batchSize128 = (UInt128)batchSize, q; Span orders = orderArray.AsSpan(0, batchSize); try { - for (UInt128 setIndex = UInt128.Zero; setIndex < setCount && Volatile.Read(ref isPrime); setIndex++) + while (kStart < maxK && Volatile.Read(ref isPrime)) { - UInt128 setOffset = checked(perSetLimit * setIndex); - UInt128 setStart = checked(setOffset + UInt128.One); - if (setStart >= limitInclusive) + remaining = maxK - kStart; + if (remaining > batchSize128) + { + currentSize = batchSize; + } + else + { + currentSize = (int)remaining; + orders = orderArray.AsSpan(0, currentSize); + } + + q = twoP * kStart + UInt128.One; + q.Mod10_8_5_3(out ulong q0m10, out ulong q0m8, out ulong q0m5, out ulong q0m3); + var ra = new ResidueAutomatonArgs(q0m10, step10, q0m8, step8, q0m3, step3, q0m5, step5); + + kernel(currentSize, exponent, twoPGpu, (GpuUInt128)kStart, last, 0UL, + ra, orderBuffer.View, smallCyclesView, primeViews.LastOne, primeViews.LastSeven, primeViews.LastOnePow2, primeViews.LastSevenPow2); + + accelerator.Synchronize(); + orderBuffer.View.CopyToCPU(ref MemoryMarshal.GetReference(orders), currentSize); + if (!Volatile.Read(ref isPrime)) { break; } - UInt128 setLimitExclusive = checked(setOffset + perSetLimit); - setLimitExclusive = setLimitExclusive >= limitInclusive - ? limitInclusive - : checked(setLimitExclusive + UInt128.One); - if (setLimitExclusive > limitInclusive) + for (i = 0; i < currentSize; i++) { - setLimitExclusive = limitInclusive; - } - - UInt128 kStart = setStart; - while (kStart < setLimitExclusive && Volatile.Read(ref isPrime)) - { - UInt128 remaining = checked(setLimitExclusive - kStart); - int currentSize = remaining > batchSize128 ? batchSize : (int)remaining; - if (currentSize <= 0) + if (orders[i] != 0UL) { + Volatile.Write(ref isPrime, false); break; } - orders = orderArray.AsSpan(0, currentSize); - - UInt128 q = checked(twoP * kStart); - q = checked(q + UInt128.One); - q.Mod10_8_5_3(out ulong q0m10, out ulong q0m8, out ulong q0m5, out ulong q0m3); - var ra = new ResidueAutomatonArgs(q0m10, step10, q0m8, step8, q0m3, step3, q0m5, step5); - - kernel(stream, currentSize, exponent, twoPGpu, (GpuUInt128)kStart, last, 0UL, - ra, orderBuffer.View, smallCyclesView, primeViews.LastOne, primeViews.LastSeven, primeViews.LastOnePow2, primeViews.LastSevenPow2); - - stream.Synchronize(); - orderBuffer.View.CopyToCPU(ref MemoryMarshal.GetReference(orders), currentSize); - stream.Synchronize(); - if (!Volatile.Read(ref isPrime)) - { - break; - } - - for (int i = 0; i < currentSize; i++) - { - if (orders[i] != 0UL) - { - Volatile.Write(ref isPrime, false); - break; - } - } - - kStart = checked(kStart + (UInt128)currentSize); } + + kStart += batchSize128; } } finally { ArrayPool.Shared.Return(orderArray); - orderBuffer.Dispose(); - execution.Dispose(); - gpuLease.Dispose(); - } - - divisorsExhausted = true; - } + orderBuffer.Dispose(); + gpuLease.Dispose(); + } + } } diff --git a/PerfectNumbers.Core/MersenneNumberTester.cs b/PerfectNumbers.Core/MersenneNumberTester.cs index 8a30580..15f7920 100644 --- a/PerfectNumbers.Core/MersenneNumberTester.cs +++ b/PerfectNumbers.Core/MersenneNumberTester.cs @@ -16,28 +16,26 @@ public enum GpuKernelType } public sealed class MersenneNumberTester( - bool useIncremental = true, - bool useOrderCache = false, - GpuKernelType kernelType = GpuKernelType.Incremental, - bool useModuloWorkaround = false, - bool useOrder = false, - bool useGpuLucas = true, - bool useGpuScan = true, - bool useGpuOrder = true, - bool useResidue = true, - UInt128? maxK = null, - UInt128? residueDivisorSets = null) + bool useIncremental = true, + bool useOrderCache = false, + GpuKernelType kernelType = GpuKernelType.Incremental, + bool useModuloWorkaround = false, + bool useOrder = false, + bool useGpuLucas = true, + bool useGpuScan = true, + bool useGpuOrder = true, + bool useResidue = true, + ulong maxK = 5_000_000UL) { - private readonly bool _useResidue = useResidue; - private readonly bool _useIncremental = useIncremental && !useResidue; - private readonly UInt128 _maxK = maxK ?? (UInt128)5_000_000UL; - private readonly UInt128 _residueDivisorSetCount = residueDivisorSets ?? (UInt128)PerfectNumberConstants.ExtraDivisorCycleSearchLimit; + private readonly bool _useResidue = useResidue; + private readonly bool _useIncremental = useIncremental && !useResidue; + private readonly ulong _maxK = maxK; private readonly GpuKernelType _kernelType = kernelType; - private readonly bool _useModuloWorkaround = useModuloWorkaround; - private readonly bool _useGpuLucas = useGpuLucas; - private readonly bool _useGpuScan = useGpuScan; // device for pow2mod/incremental scanning - private readonly bool _useGpuOrder = useGpuOrder; // device for order computations - private static readonly ConcurrentDictionary OrderCache = new(); + private readonly bool _useModuloWorkaround = useModuloWorkaround; + private readonly bool _useGpuLucas = useGpuLucas; + private readonly bool _useGpuScan = useGpuScan; // device for pow2mod/incremental scanning + private readonly bool _useGpuOrder = useGpuOrder; // device for order computations + private static readonly ConcurrentDictionary OrderCache = new(); private readonly MersenneNumberIncrementalGpuTester _incrementalGpuTester = new(kernelType, useGpuOrder); private readonly MersenneNumberIncrementalCpuTester _incrementalCpuTester = new(kernelType); private readonly MersenneNumberOrderGpuTester _orderGpuTester = new(kernelType, useGpuOrder); @@ -50,24 +48,7 @@ public sealed class MersenneNumberTester( // TODO: Ensure Program passes useResidue correctly and that residue-vs-incremental-vs-LL selection respects CLI flags. // Also consider injecting a shared cycles cache (MersenneDivisorCycles.Shared) to both CPU and GPU testers. - [MethodImpl(MethodImplOptions.AggressiveInlining)] - private static UInt128 DivideRoundUp(UInt128 value, UInt128 divisor) - { - if (divisor == UInt128.Zero) - { - return UInt128.Zero; - } - - UInt128 result = value / divisor; - if (result * divisor == value) - { - return result; - } - - return result + UInt128.One; - } - - public void WarmUpOrders(ulong exponent, ulong limit = 5_000_000UL) + public void WarmUpOrders(ulong exponent, ulong limit = 5_000_000UL) { bool cacheEnabled = useOrderCache; if (!cacheEnabled) @@ -218,116 +199,51 @@ public sealed class MersenneNumberTester( } } - public bool IsMersennePrime(ulong exponent) => IsMersennePrime(exponent, out _); + public bool IsMersennePrime(ulong exponent) + { + // Safe early rejections using known orders for tiny primes: + // 7 | M_p iff 3 | p. Avoid rejecting p == 3 where M_3 == 7 is prime. + if ((exponent % 3UL) == 0UL && exponent != 3UL) + { + return false; + } - public bool IsMersennePrime(ulong exponent, out bool divisorsExhausted) - { - divisorsExhausted = false; + if ((exponent & 3UL) == 1UL && exponent.SharesFactorWithExponentMinusOne()) + { + return false; + } - // Safe early rejections using known orders for tiny primes: - // 7 | M_p iff 3 | p. Avoid rejecting p == 3 where M_3 == 7 is prime. - if ((exponent % 3UL) == 0UL && exponent != 3UL) - { - divisorsExhausted = true; - return false; - } - - if ((exponent & 3UL) == 1UL && exponent.SharesFactorWithExponentMinusOne()) - { - divisorsExhausted = true; - return false; - } - - // Skip even exponents: for even p > 2, M_p is composite; here p is large. - if ((exponent & 1UL) == 0UL) - { - divisorsExhausted = true; - return false; - } + // Skip even exponents: for even p > 2, M_p is composite; here p is large. + if ((exponent & 1UL) == 0UL) + { + return false; + } bool prePrime = true; - UInt128 twoP = (UInt128)exponent << 1; // 2 * p - // UInt128 maxDivisor = new(ulong.MaxValue, ulong.MaxValue); - UInt128 maxK = UInt128Numbers.Two.Pow(exponent) / 2 + 1; - bool lastIsSeven = (exponent & 3UL) == 3UL; + UInt128 twoP = (UInt128)exponent << 1; // 2 * p + // UInt128 maxDivisor = new(ulong.MaxValue, ulong.MaxValue); + UInt128 maxK = UInt128Numbers.Two.Pow(exponent) / 2 + 1; + bool lastIsSeven = (exponent & 3UL) == 3UL; - // In residue mode the residue scan is the final primality check. - if (_useResidue) - { - if (_maxK < maxK) - { - maxK = _maxK; - } + // In residue mode the residue scan is the final primality check. + if (_useResidue) + { + if ((UInt128)_maxK < maxK) + { + maxK = _maxK; + } - if (_maxK == UInt128.Zero || _residueDivisorSetCount == UInt128.Zero) - { - return prePrime; - } + if (_useGpuScan) + { + _residueGpuTester!.Scan(exponent, twoP, lastIsSeven, maxK, ref prePrime); + } + else + { + _residueCpuTester!.Scan(exponent, twoP, lastIsSeven, maxK, ref prePrime); + } - UInt128 totalLimit = _maxK; - if (totalLimit > maxK) - { - totalLimit = maxK; - } - - if (totalLimit == UInt128.Zero) - { - return true; - } - - UInt128 requestedSets = _residueDivisorSetCount; - if (_useGpuScan) - { - UInt128 maxUsefulSets = DivideRoundUp(totalLimit, (UInt128)GpuConstants.ScanBatchSize); - if (maxUsefulSets == UInt128.Zero) - { - maxUsefulSets = UInt128.One; - } - - if (requestedSets == UInt128.Zero || requestedSets > maxUsefulSets) - { - requestedSets = maxUsefulSets; - } - } - - if (requestedSets == UInt128.Zero) - { - requestedSets = UInt128.One; - } - - UInt128 perSetLimit = DivideRoundUp(totalLimit, requestedSets); - if (perSetLimit == UInt128.Zero) - { - perSetLimit = UInt128.One; - } - - UInt128 setsNeeded = DivideRoundUp(totalLimit, perSetLimit); - UInt128 effectiveSets = requestedSets; - if (setsNeeded < effectiveSets) - { - effectiveSets = setsNeeded; - } - - bool scanCompleted = false; - - if (_useGpuScan) - { - _residueGpuTester!.Scan(exponent, twoP, lastIsSeven, perSetLimit, effectiveSets, totalLimit, ref prePrime, ref scanCompleted); - } - else - { - _residueCpuTester!.Scan(exponent, twoP, lastIsSeven, perSetLimit, effectiveSets, totalLimit, ref prePrime, ref scanCompleted); - } - - if (!prePrime) - { - divisorsExhausted = true; - return false; - } - - divisorsExhausted = scanCompleted; - return true; - } + return prePrime; + } if (!_useIncremental) { @@ -347,15 +263,12 @@ public sealed class MersenneNumberTester( if (!prePrime) { - divisorsExhausted = true; return false; } - bool lucasPrime = _useGpuLucas - ? _lucasLehmerGpuTester.IsPrime(exponent, runOnGpu: true) - : _lucasLehmerCpuTester.IsPrime(exponent); - divisorsExhausted = true; - return lucasPrime; + return _useGpuLucas + ? _lucasLehmerGpuTester.IsPrime(exponent, runOnGpu: true) + : _lucasLehmerCpuTester.IsPrime(exponent); } bool isPrime = true; @@ -383,9 +296,8 @@ public sealed class MersenneNumberTester( } } - divisorsExhausted = true; - return isPrime; - } + return isPrime; + } // Limit of k to scan in the pre-Lucas prescan public const int PreLucasPreScanK = 5_000_000;