From 6a1198c417426ab03e60ce5b35979d7ab6f9d82a Mon Sep 17 00:00:00 2001 From: mlemanczyk Date: Sat, 29 Nov 2025 12:09:33 +0100 Subject: [PATCH] CPU vs GPU work delegation adjusted to better utilize CPU. --- .../PrimeTesterBenchmarks.cs | 63 +++- EvenPerfectBitScanner.Tests/ProgramTests.cs | 2 +- EvenPerfectBitScanner/Program.cs | 2 +- .../PrimeTester/PrimeTesterIsPrimeTests.cs | 6 +- PerfectNumbers.Core/MersenneDivisorCycles.cs | 4 +- .../MersennePrimeFactorTester.cs | 4 +- PerfectNumbers.Core/PerfectNumberConstants.cs | 2 +- .../PrimeOrderCalculator.Cpu.cs | 32 +- .../PrimeOrderCalculator.Shared.cs | 2 +- PerfectNumbers.Core/PrimeTester.cs | 323 +----------------- PerfectNumbers.Core/PrimeTesterByLastDigit.cs | 239 +++++++++++++ .../PrimeTesterSeverOrOthers.cs | 241 +++++++++++++ 12 files changed, 576 insertions(+), 344 deletions(-) create mode 100644 PerfectNumbers.Core/PrimeTesterByLastDigit.cs create mode 100644 PerfectNumbers.Core/PrimeTesterSeverOrOthers.cs diff --git a/EvenPerfectBitScanner.Benchmarks/PrimeTesterBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/PrimeTesterBenchmarks.cs index b27b469..924ea47 100644 --- a/EvenPerfectBitScanner.Benchmarks/PrimeTesterBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/PrimeTesterBenchmarks.cs @@ -136,7 +136,7 @@ public class PrimeTesterBenchmarks for (int i = 0; i < values.Length; i++) { - if (PrimeTester.IsPrime(values[i])) + if (PrimeTester.IsPrimeCpu(values[i])) { primeCount++; } @@ -153,7 +153,24 @@ public class PrimeTesterBenchmarks for (int i = 0; i < values.Length; i++) { - if (PrimeTester.IsPrime(values[i])) + if (PrimeTesterByLastDigit.IsPrimeCpu(values[i])) + { + primeCount++; + } + } + + return primeCount; + } + + [Benchmark] + public int SevenOrOthersCpu() + { + ulong[] values = _candidates; + int primeCount = 0; + + for (int i = 0; i < values.Length; i++) + { + if (PrimeTesterSeverOrOthers.IsPrimeCpu(values[i])) { primeCount++; } @@ -195,12 +212,54 @@ public class PrimeTesterBenchmarks 4_000_000UL, 4_096); + yield return PrimeBenchmarkCase.FromRange( + "≥16_000_000", + "Sampled admissible odd values ≥ 16,000,000 (k ≤ 5, not divisible by 5)", + 16_000_001UL, + 16_500_000UL, + 4_096); + + yield return PrimeBenchmarkCase.FromRange( + "≥32_000_000", + "Sampled admissible odd values ≥ 32,000,000 (k ≤ 5, not divisible by 5)", + 32_000_001UL, + 32_500_000UL, + 4_096); + + yield return PrimeBenchmarkCase.FromRange( + "≥64_000_000", + "Sampled admissible odd values ≥ 64,000,000 (k ≤ 5, not divisible by 5)", + 64_000_001UL, + 64_500_000UL, + 4_096); + + yield return PrimeBenchmarkCase.FromRange( + "≥100_000_000", + "Sampled admissible odd values ≥ 100,000,000 (k ≤ 5, not divisible by 5)", + 100_000_001UL, + 100_500_000UL, + 4_096); + yield return PrimeBenchmarkCase.FromRange( "≥138_000_000", "Sampled admissible odd values ≥ 138,000,000 (k ≤ 5, not divisible by 5)", 138_000_001UL, 138_500_000UL, 4_096); + + yield return PrimeBenchmarkCase.FromRange( + "≥276_000_000", + "Sampled admissible odd values ≥ 276,000,000 (k ≤ 5, not divisible by 5)", + 276_000_001UL, + 276_500_000UL, + 4_096); + + yield return PrimeBenchmarkCase.FromRange( + "≥552_000_000", + "Sampled admissible odd values ≥ 552,000,000 (k ≤ 5, not divisible by 5)", + 552_000_001UL, + 552_500_000UL, + 4_096); } public readonly record struct PrimeBenchmarkCase(string Name, string Description, ulong[] Candidates) diff --git a/EvenPerfectBitScanner.Tests/ProgramTests.cs b/EvenPerfectBitScanner.Tests/ProgramTests.cs index 2007013..9e3e315 100644 --- a/EvenPerfectBitScanner.Tests/ProgramTests.cs +++ b/EvenPerfectBitScanner.Tests/ProgramTests.cs @@ -52,7 +52,7 @@ public class ProgramTests [InlineData(81UL, false)] public void IsPrime_identifies_primes_correctly(ulong n, bool expected) { - PrimeTester.IsPrime(n).Should().Be(expected); + PrimeTester.IsPrimeCpu(n).Should().Be(expected); } [Fact] diff --git a/EvenPerfectBitScanner/Program.cs b/EvenPerfectBitScanner/Program.cs index 559c6cd..3c5f0fc 100644 --- a/EvenPerfectBitScanner/Program.cs +++ b/EvenPerfectBitScanner/Program.cs @@ -631,7 +631,7 @@ internal static class Program if (!_cliArguments.UseByDivisor) { if (!(_runPrimesOnCpu - ? PrimeTester.IsPrime(p) + ? PrimeTester.IsPrimeCpu(p) : HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, p))) { _lastCompositeP = true; diff --git a/PerfectNumbers.Core.Tests/PrimeTester/PrimeTesterIsPrimeTests.cs b/PerfectNumbers.Core.Tests/PrimeTester/PrimeTesterIsPrimeTests.cs index 383fc28..9f9a43b 100644 --- a/PerfectNumbers.Core.Tests/PrimeTester/PrimeTesterIsPrimeTests.cs +++ b/PerfectNumbers.Core.Tests/PrimeTester/PrimeTesterIsPrimeTests.cs @@ -27,7 +27,7 @@ public class PrimeTesterIsPrimeTests { var tester = new PrimeTester(); - PrimeTester.IsPrime(n).Should().Be(expected); + PrimeTester.IsPrimeCpu(n).Should().Be(expected); } [Theory] @@ -56,8 +56,8 @@ public class PrimeTesterIsPrimeTests public void IsPrimeGpu_falls_back_to_cpu_when_forced() { var tester = new PrimeTester(); - PrimeTester.IsPrime(137UL).Should().BeTrue(); - PrimeTester.IsPrime(341UL).Should().BeFalse(); + PrimeTester.IsPrimeCpu(137UL).Should().BeTrue(); + PrimeTester.IsPrimeCpu(341UL).Should().BeFalse(); } [Fact] diff --git a/PerfectNumbers.Core/MersenneDivisorCycles.cs b/PerfectNumbers.Core/MersenneDivisorCycles.cs index c66aece..02bc79e 100644 --- a/PerfectNumbers.Core/MersenneDivisorCycles.cs +++ b/PerfectNumbers.Core/MersenneDivisorCycles.cs @@ -415,7 +415,7 @@ public class MersenneDivisorCycles // Console.WriteLine($"MersenneDivisorCycles.TryFactorIntoCountsInternal hits {Volatile.Read(ref _tryFactorIntoCountsInternalHits)}"); // HeuristicPrimeTester primeTester = _primeTester ??= new(); - if (HeuristicCombinedPrimeTester.IsPrimeCpu(remaining)) + if (PrimeTesterByLastDigit.IsPrimeCpu(remaining)) // if (HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, remaining)) // if (PrimeTester.IsPrime(remaining)) { @@ -966,7 +966,7 @@ public class MersenneDivisorCycles // HeuristicPrimeTester primeTester = _primeTester ??= new(); // if (!skipPrimeOrderHeuristic && HeuristicPrimeTester.Exclusive.IsPrimeCpu(divisor, CancellationToken.None)) - if (!skipPrimeOrderHeuristic && HeuristicCombinedPrimeTester.IsPrimeCpu(divisor)) + if (!skipPrimeOrderHeuristic && PrimeTesterByLastDigit.IsPrimeCpu(divisor)) // if (!skipPrimeOrderHeuristic && HeuristicCombinedPrimeTester.IsPrimeGpu(divisor)) { ulong computedOrder = PrimeOrderCalculator.Calculate( diff --git a/PerfectNumbers.Core/MersennePrimeFactorTester.cs b/PerfectNumbers.Core/MersennePrimeFactorTester.cs index 2393e2b..1ecdce6 100644 --- a/PerfectNumbers.Core/MersennePrimeFactorTester.cs +++ b/PerfectNumbers.Core/MersennePrimeFactorTester.cs @@ -35,7 +35,7 @@ public static class MersennePrimeFactorTester // Atomic.Add(ref _isPrimeFactorHits, 1UL); // Console.WriteLine($"MersennePrimeFactorTester.IsPrimeFactor hits {Volatile.Read(ref _isPrimeFactorHits)}"); - if (!PrimeTester.IsPrime(q)) + if (!PrimeTester.IsPrimeCpu(q)) { return false; } @@ -256,7 +256,7 @@ public static class MersennePrimeFactorTester // Console.WriteLine($"MersennePrimeFactorTester.Factor64 hits {Volatile.Read(ref _factor64Hits)}"); // if (HeuristicPrimeTester.Exclusive.IsPrime(m, ct)) - if (PrimeTester.IsPrime(m)) + if (PrimeTester.IsPrimeCpu(m)) { dict[m] = 1; continue; diff --git a/PerfectNumbers.Core/PerfectNumberConstants.cs b/PerfectNumbers.Core/PerfectNumberConstants.cs index 8c9f271..6b74b6c 100644 --- a/PerfectNumbers.Core/PerfectNumberConstants.cs +++ b/PerfectNumbers.Core/PerfectNumberConstants.cs @@ -8,7 +8,7 @@ public static class PerfectNumberConstants public const int DefaultSmallPrimeFactorSlotCount = 64; public const int DefaultSpecialMaxFactorCapacity = 1024; public const ulong ExtraDivisorCycleSearchLimit = 64UL; - public const int GpuFrequency = 192 - 1; + public const int GpuFrequency = 2048 - 1; public const int MaxQForDivisorCycles = 4_000_000; public static readonly uint PrimesLimit = 1_000_000; //(ulong)Array.MaxLength;// 1_000_000; public const int MaxOddPowersCount = 128; diff --git a/PerfectNumbers.Core/PrimeOrderCalculator.Cpu.cs b/PerfectNumbers.Core/PrimeOrderCalculator.Cpu.cs index fb3c11f..404fc2a 100644 --- a/PerfectNumbers.Core/PrimeOrderCalculator.Cpu.cs +++ b/PerfectNumbers.Core/PrimeOrderCalculator.Cpu.cs @@ -68,8 +68,11 @@ internal static partial class PrimeOrderCalculator // result = prime <= PerfectNumberConstants.MaxQForDivisorCycles // ? CalculateByFactorizationGpu(gpu, prime, divisorData) // : CalculateByFactorizationCpu(gpu, prime, divisorData); - result = CalculateByFactorizationGpu(gpu, prime, divisorData); - // result = CalculateByFactorizationCpu(gpu, prime, divisorData); + // result = RunOnCpu() + // ? CalculateByFactorizationCpu(gpu, prime, divisorData) + // : CalculateByFactorizationGpu(gpu, prime, divisorData); + // result = CalculateByFactorizationGpu(gpu, prime, divisorData); + result = CalculateByFactorizationCpu(gpu, prime, divisorData); phiFactors.Dispose(); return result; @@ -105,8 +108,10 @@ internal static partial class PrimeOrderCalculator if (config.Mode == PrimeOrderMode.Strict) { orderFactors?.Dispose(); - // return CalculateByFactorizationCpu(gpu, prime, divisorData); - return CalculateByFactorizationGpu(gpu, prime, divisorData); + return CalculateByFactorizationCpu(gpu, prime, divisorData); + // return RunOnCpu() + // ? CalculateByFactorizationCpu(gpu, prime, divisorData) + // : CalculateByFactorizationGpu(gpu, prime, divisorData); } if (TryHeuristicFinishCpu(gpu, prime, candidateOrder, previousOrder, divisorData, config, phiFactors, orderFactors, out ulong order)) @@ -132,7 +137,7 @@ internal static partial class PrimeOrderCalculator Span stackBuffer = stackalloc ulong[length]; - if (RunOnCpu()) + if (length <= 8) { return EvaluateSpecialMaxCandidatesCpu(stackBuffer, factorSpan, phi, prime, divisorData); } @@ -1175,7 +1180,7 @@ internal static partial class PrimeOrderCalculator { // bool isPrime = HeuristicCombinedPrimeTester.IsPrimeCpu(composite); // bool isPrime = PrimeTester.IsPrime(composite); - bool isPrime = HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, composite); + bool isPrime = RunOnCpu() ? PrimeTesterByLastDigit.IsPrimeCpu(composite) : HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, composite); // bool isPrime = PrimeTester.IsPrimeInternal(composite, CancellationToken.None); // Atomic.Add(ref _partialFactorPendingHits, 1UL); @@ -1218,7 +1223,8 @@ internal static partial class PrimeOrderCalculator { // cofactorIsPrime = PrimeTester.IsPrime(cofactor); // cofactorIsPrime = HeuristicCombinedPrimeTester.IsPrimeCpu(cofactor); - cofactorIsPrime = HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, cofactor); + cofactorIsPrime = RunOnCpu() ? PrimeTesterByLastDigit.IsPrimeCpu(cofactor) : HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, cofactor); + // cofactorIsPrime = HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, cofactor); // cofactorIsPrime = PrimeTester.IsPrimeInternal(cofactor, CancellationToken.None); // Atomic.Add(ref _partialFactorCofactorHits, 1UL); @@ -1319,7 +1325,7 @@ internal static partial class PrimeOrderCalculator while (compositeStack.Count > 0) { ulong composite = compositeStack.Pop(); - bool isPrime = RunOnCpu() ? PrimeTester.IsPrime(composite) : HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, composite); + bool isPrime = RunOnCpu() ? PrimeTesterByLastDigit.IsPrimeCpu(composite) : HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, composite); if (isPrime) { @@ -1354,7 +1360,7 @@ internal static partial class PrimeOrderCalculator while (compositeStack.Count > 0) { ulong composite = compositeStack.Pop(); - bool isPrime = PrimeTester.IsPrime(composite); + bool isPrime = PrimeTesterByLastDigit.IsPrimeCpu(composite); if (isPrime) { @@ -1644,7 +1650,7 @@ internal static partial class PrimeOrderCalculator x = AdvancePolynomialCpu(x, c, n); y = AdvancePolynomialTwiceCpu(y, c, n); ulong diff = x > y ? x - y : y - x; - d = BinaryGcd(diff, n); + d = diff.BinaryGcd(n); } if (d != n) @@ -1862,7 +1868,7 @@ internal static partial class PrimeOrderCalculator // bool isPrime = primeTester.HeuristicIsPrimeGpu(value); // bool isPrime = HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, value); // bool isPrime = PrimeTester.IsPrimeGpu(value); - bool isPrime = PrimeTester.IsPrime(value); + bool isPrime = PrimeTesterByLastDigit.IsPrimeCpu(value); // bool isPrime = RunOnGpu() // ? HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, value) // : PrimeTester.IsPrime(value); @@ -1885,7 +1891,7 @@ internal static partial class PrimeOrderCalculator // isPrime = RunOnGpu() // ? HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, factor) // : PrimeTester.IsPrime(factor); - isPrime = PrimeTester.IsPrime(factor); + isPrime = PrimeTesterByLastDigit.IsPrimeCpu(factor); // isPrime = Open.Numeric.Primes.Prime.Numbers.IsPrime(factor); if (isPrime) { @@ -1922,7 +1928,7 @@ internal static partial class PrimeOrderCalculator // Console.WriteLine($"Partial factor cofactor hits {Volatile.Read(ref _partialFactorCofactorHits)}"); // isPrime = HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, quotient); - isPrime = PrimeTester.IsPrime(quotient); + isPrime = PrimeTesterByLastDigit.IsPrimeCpu(quotient); // isPrime = PrimeTester.IsPrimeGpu(quotient); // isPrime = Open.Numeric.Primes.Prime.Numbers.IsPrime(quotient); } diff --git a/PerfectNumbers.Core/PrimeOrderCalculator.Shared.cs b/PerfectNumbers.Core/PrimeOrderCalculator.Shared.cs index 465b633..b6d90c9 100644 --- a/PerfectNumbers.Core/PrimeOrderCalculator.Shared.cs +++ b/PerfectNumbers.Core/PrimeOrderCalculator.Shared.cs @@ -953,7 +953,7 @@ internal static partial class PrimeOrderCalculator if (value <= ulong.MaxValue) { // return HeuristicPrimeTester.Exclusive.IsPrimeCpu((ulong)value, CancellationToken.None); - return PrimeTester.IsPrime((ulong)value); + return PrimeTester.IsPrimeCpu((ulong)value); } if ((value & UInt128.One) == UInt128.Zero) diff --git a/PerfectNumbers.Core/PrimeTester.cs b/PerfectNumbers.Core/PrimeTester.cs index 28c30d7..2cf98d9 100644 --- a/PerfectNumbers.Core/PrimeTester.cs +++ b/PerfectNumbers.Core/PrimeTester.cs @@ -6,16 +6,16 @@ using PerfectNumbers.Core.Gpu.Accelerators; namespace PerfectNumbers.Core; -public sealed class PrimeTesterByLastDigit +public sealed class PrimeTester { - public PrimeTesterByLastDigit() + public PrimeTester() { } [ThreadStatic] - private static PrimeTesterByLastDigit? _tester; + private static PrimeTester? _tester; - public static PrimeTesterByLastDigit Exclusive => _tester ??= new(); + public static PrimeTester Exclusive => _tester ??= new(); public static bool IsPrimeCpu(ulong n) { @@ -34,249 +34,6 @@ public sealed class PrimeTesterByLastDigit bool requiresTrialDivision = result && n >= 7UL; - if (requiresTrialDivision) - { - // EvenPerfectBitScanner streams exponents starting at 136,279,841, so the Mod10/GCD guard never fires on the - // production path. Leave the logic commented out as instrumentation for diagnostic builds. - // bool sharesMaxExponentFactor = n.Mod10() == 1UL && SharesFactorWithMaxExponent(n); - // result &= !sharesMaxExponentFactor; - - if (result) - { - uint[] smallPrimeDivisors; - ulong[] smallPrimeDivisorsMul; - - ulong nMod10 = n.Mod10(); - switch (nMod10) - { - case 1UL: - smallPrimeDivisors = PrimesGenerator.SmallPrimesLastOne; - smallPrimeDivisorsMul = PrimesGenerator.SmallPrimesPow2LastOne; - break; - case 3UL: - smallPrimeDivisors = DivisorGenerator.SmallPrimesLastThree; - smallPrimeDivisorsMul = DivisorGenerator.SmallPrimesPow2LastThree; - break; - case 7UL: - smallPrimeDivisors = PrimesGenerator.SmallPrimesLastSeven; - smallPrimeDivisorsMul = PrimesGenerator.SmallPrimesPow2LastSeven; - break; - default: - smallPrimeDivisors = DivisorGenerator.SmallPrimesLastNine; - smallPrimeDivisorsMul = DivisorGenerator.SmallPrimesPow2LastNine; - break; - } - - int smallPrimeDivisorsLength = smallPrimeDivisors.Length; - for (int i = 0; i < smallPrimeDivisorsLength; i++) - { - ulong value = smallPrimeDivisorsMul[i]; - if (value > n) - { - break; - } - - value = n % smallPrimeDivisors[i]; - if (value == 0) - { - result = false; - break; - } - } - } - } - - return result; - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static bool IsPrimeGpu(PrimeOrderCalculatorAccelerator gpu, ulong n) - { - byte flag = 0; - var inputView = gpu.InputView; - var outputView = gpu.OutputByteView; - - int acceleratorIndex = gpu.AcceleratorIndex; - AcceleratorStream stream = AcceleratorStreamPool.Rent(acceleratorIndex); - inputView.CopyFromCPU(stream, ref n, 1); - - var kernelLauncher = gpu.SmallPrimeSieveKernelLauncher; - - kernelLauncher( - stream, - 1, - inputView, - gpu.DevicePrimesLastOne, - gpu.DevicePrimesLastSeven, - gpu.DevicePrimesLastThree, - gpu.DevicePrimesLastNine, - gpu.DevicePrimesPow2LastOne, - gpu.DevicePrimesPow2LastSeven, - gpu.DevicePrimesPow2LastThree, - gpu.DevicePrimesPow2LastNine, - outputView); - - outputView.CopyToCPU(stream, ref flag, 1); - stream.Synchronize(); - AcceleratorStreamPool.Return(acceleratorIndex, stream); - - return flag != 0; - } - - public static int GpuBatchSize { get; set; } = 262_144; - - private static readonly object GpuWarmUpLock = new(); - private static int WarmedGpuLeaseCount; - - public static void WarmUpGpuKernels(int threadCount) - { - int target = threadCount >> 2; - if (target == 0) - { - target = threadCount; - } - - lock (GpuWarmUpLock) - { - if (target <= WarmedGpuLeaseCount) - { - return; - } - - PrimeOrderCalculatorAccelerator.WarmUp(); - WarmedGpuLeaseCount = target; - } - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static void IsPrimeBatchGpu(PrimeOrderCalculatorAccelerator gpu, ReadOnlySpan values, Span results) - { - // GpuPrimeWorkLimiter.Acquire(); - int acceleratorIndex = gpu.AcceleratorIndex; - int totalLength = values.Length; - int batchSize = GpuBatchSize; - - var inputView = gpu.InputView; - var outputView = gpu.OutputByteView; - - int pos = 0; - AcceleratorStream stream = AcceleratorStreamPool.Rent(acceleratorIndex); - var kernelLauncher = gpu.SmallPrimeSieveKernelLauncher; - - while (pos < totalLength) - { - int remaining = totalLength - pos; - int count = remaining > batchSize ? batchSize : remaining; - - var valueSlice = values.Slice(pos, count); - inputView.CopyFromCPU(stream, valueSlice); - - kernelLauncher( - stream, - count, - inputView, - gpu.DevicePrimesLastOne, - gpu.DevicePrimesLastSeven, - gpu.DevicePrimesLastThree, - gpu.DevicePrimesLastNine, - gpu.DevicePrimesPow2LastOne, - gpu.DevicePrimesPow2LastSeven, - gpu.DevicePrimesPow2LastThree, - gpu.DevicePrimesPow2LastNine, - outputView); - - var resultSlice = results.Slice(pos, count); - outputView.CopyToCPU(stream, resultSlice); - - pos += count; - } - - stream.Synchronize(); - AcceleratorStreamPool.Return(acceleratorIndex, stream); - // GpuPrimeWorkLimiter.Release(); - } - - // Expose cache clearing for accelerator disposal coordination - public static void ClearGpuCaches() - { - PrimeOrderCalculatorAccelerator.Clear(); - } - - internal static void DisposeGpuContexts() - { - PrimeOrderCalculatorAccelerator.DisposeAll(); - lock (GpuWarmUpLock) - { - WarmedGpuLeaseCount = 0; - } - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - internal static bool SharesFactorWithMaxExponent(ulong n) - { - // TODO: Replace this on-the-fly GCD probe with the cached factor table derived from - // ResidueComputationBenchmarks so divisor-cycle metadata can short-circuit the test - // instead of recomputing binary GCD for every candidate. - ulong m = (ulong)BitOperations.Log2(n); - return n.BinaryGcd(m) != 1UL; - } - - internal static void SharesFactorWithMaxExponentBatch(PrimeOrderCalculatorAccelerator gpu, ReadOnlySpan values, Span results) - { - // TODO: Route this batch helper through the shared GPU kernel pool from - // GpuUInt128BinaryGcdBenchmarks so we reuse cached kernels, pinned host buffers, - // and divisor-cycle staging instead of allocating new device buffers per call. - // Check in benchmarks, which implementation was the fastest, is compatible with GPU, - // and implement it. - - int length = values.Length; - - int acceleratorIndex = gpu.AcceleratorIndex; - gpu.EnsureCapacity(0, length); - var inputBufferView = gpu.InputView; - var resultBufferView = gpu.OutputByteView; - var kernelLauncher = gpu.SharesFactorKernelLauncher; - - AcceleratorStream stream = AcceleratorStreamPool.Rent(acceleratorIndex); - inputBufferView.CopyFromCPU(stream, values); - - kernelLauncher(stream, length, inputBufferView, resultBufferView); - - resultBufferView.CopyToCPU(stream, in results); - stream.Synchronize(); - - AcceleratorStreamPool.Return(acceleratorIndex, stream); - } -} - -public sealed class PrimeTester -{ - public PrimeTester() - { - } - - [ThreadStatic] - private static PrimeTester? _tester; - - public static PrimeTester Exclusive => _tester ??= new(); - - public static bool IsPrime(ulong n) - { - if (n <= 1UL) - { - return false; - } - - if (n == 2UL) - { - throw new InvalidOperationException("PrimeTester.IsPrime encountered the sentinel input 2."); - } - - bool isOdd = (n & 1UL) != 0UL; - bool result = isOdd; - - bool requiresTrialDivision = result && n >= 7UL; - if (requiresTrialDivision) { // EvenPerfectBitScanner streams exponents starting at 136,279,841, so the Mod10/GCD guard never fires on the @@ -308,76 +65,6 @@ public sealed class PrimeTester return result; } - public static bool IsPrimeLastDigit(ulong n) - { - if (n <= 1UL) - { - return false; - } - - if (n == 2UL) - { - throw new InvalidOperationException("PrimeTester.IsPrime encountered the sentinel input 2."); - } - - bool isOdd = (n & 1UL) != 0UL; - bool result = isOdd; - - bool requiresTrialDivision = result && n >= 7UL; - - if (requiresTrialDivision) - { - // EvenPerfectBitScanner streams exponents starting at 136,279,841, so the Mod10/GCD guard never fires on the - // production path. Leave the logic commented out as instrumentation for diagnostic builds. - // bool sharesMaxExponentFactor = n.Mod10() == 1UL && SharesFactorWithMaxExponent(n); - // result &= !sharesMaxExponentFactor; - - if (result) - { - uint[] smallPrimeDivisors = PrimesGenerator.SmallPrimes; - ulong[] smallPrimeDivisorsMul = PrimesGenerator.SmallPrimesPow2; - - ulong nMod10 = n.Mod10(); - switch (nMod10) - { - case 1UL: - smallPrimeDivisors = PrimesGenerator.SmallPrimesLastOne; - smallPrimeDivisorsMul = PrimesGenerator.SmallPrimesPow2LastOne; - break; - case 3UL: - smallPrimeDivisors = DivisorGenerator.SmallPrimesLastThree; - smallPrimeDivisorsMul = DivisorGenerator.SmallPrimesPow2LastThree; - break; - case 7UL: - smallPrimeDivisors = PrimesGenerator.SmallPrimesLastSeven; - smallPrimeDivisorsMul = PrimesGenerator.SmallPrimesPow2LastSeven; - break; - case 9UL: - smallPrimeDivisors = DivisorGenerator.SmallPrimesLastNine; - smallPrimeDivisorsMul = DivisorGenerator.SmallPrimesPow2LastNine; - break; - } - - int smallPrimeDivisorsLength = smallPrimeDivisors.Length; - for (int i = 0; i < smallPrimeDivisorsLength; i++) - { - if (smallPrimeDivisorsMul[i] > n) - { - break; - } - - if (n % smallPrimeDivisors[i] == 0) - { - result = false; - break; - } - } - } - } - - return result; - } - [MethodImpl(MethodImplOptions.AggressiveInlining)] public static bool IsPrimeGpu(PrimeOrderCalculatorAccelerator gpu, ulong n) { @@ -536,4 +223,4 @@ public sealed class PrimeTester AcceleratorStreamPool.Return(acceleratorIndex, stream); } -} +} \ No newline at end of file diff --git a/PerfectNumbers.Core/PrimeTesterByLastDigit.cs b/PerfectNumbers.Core/PrimeTesterByLastDigit.cs new file mode 100644 index 0000000..df0ac27 --- /dev/null +++ b/PerfectNumbers.Core/PrimeTesterByLastDigit.cs @@ -0,0 +1,239 @@ +using System.Numerics; +using System.Runtime.CompilerServices; +using ILGPU.Runtime; +using PerfectNumbers.Core.Gpu.Accelerators; + +namespace PerfectNumbers.Core; + +public sealed class PrimeTesterByLastDigit +{ + public PrimeTesterByLastDigit() + { + } + + [ThreadStatic] + private static PrimeTesterByLastDigit? _tester; + + public static PrimeTesterByLastDigit Exclusive => _tester ??= new(); + + public static bool IsPrimeCpu(ulong n) + { + // The below IFs never trigger in production code of EvenPerfectBitScanner on --mersenne=bydivisor path. + // if (n <= 1UL) + // { + // throw new InvalidOperationException("PrimeTester.IsPrime encountered the sentinel input 2."); + // // return false; + // } + + // if (n == 2UL) + // { + // throw new InvalidOperationException("PrimeTester.IsPrime encountered the sentinel input 2."); + // } + + // EvenPerfectBitScanner streams exponents starting at 136,279,841, so the Mod10/GCD guard never fires on the + // production path. Leave the logic commented out as instrumentation for diagnostic builds. + // bool sharesMaxExponentFactor = n.Mod10() == 1UL && SharesFactorWithMaxExponent(n); + // result &= !sharesMaxExponentFactor; + + uint[] smallPrimeDivisors; + ulong[] smallPrimeDivisorsMul; + + ulong nMod10 = n.Mod10(); + switch (nMod10) + { + case 1UL: + smallPrimeDivisors = PrimesGenerator.SmallPrimesLastOne; + smallPrimeDivisorsMul = PrimesGenerator.SmallPrimesPow2LastOne; + break; + case 3UL: + smallPrimeDivisors = DivisorGenerator.SmallPrimesLastThree; + smallPrimeDivisorsMul = DivisorGenerator.SmallPrimesPow2LastThree; + break; + case 7UL: + smallPrimeDivisors = PrimesGenerator.SmallPrimesLastSeven; + smallPrimeDivisorsMul = PrimesGenerator.SmallPrimesPow2LastSeven; + break; + default: + smallPrimeDivisors = DivisorGenerator.SmallPrimesLastNine; + smallPrimeDivisorsMul = DivisorGenerator.SmallPrimesPow2LastNine; + break; + } + + int smallPrimeDivisorsLength = smallPrimeDivisors.Length; + for (int i = 0; i < smallPrimeDivisorsLength; i++) + { + ulong value = smallPrimeDivisorsMul[i]; + if (value > n) + { + break; + } + + value = n % smallPrimeDivisors[i]; + if (value == 0) + { + return false; + } + } + + return true; + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static bool IsPrimeGpu(PrimeOrderCalculatorAccelerator gpu, ulong n) + { + byte flag = 0; + var inputView = gpu.InputView; + var outputView = gpu.OutputByteView; + + int acceleratorIndex = gpu.AcceleratorIndex; + AcceleratorStream stream = AcceleratorStreamPool.Rent(acceleratorIndex); + inputView.CopyFromCPU(stream, ref n, 1); + + var kernelLauncher = gpu.SmallPrimeSieveKernelLauncher; + + kernelLauncher( + stream, + 1, + inputView, + gpu.DevicePrimesLastOne, + gpu.DevicePrimesLastSeven, + gpu.DevicePrimesLastThree, + gpu.DevicePrimesLastNine, + gpu.DevicePrimesPow2LastOne, + gpu.DevicePrimesPow2LastSeven, + gpu.DevicePrimesPow2LastThree, + gpu.DevicePrimesPow2LastNine, + outputView); + + outputView.CopyToCPU(stream, ref flag, 1); + stream.Synchronize(); + AcceleratorStreamPool.Return(acceleratorIndex, stream); + + return flag != 0; + } + + public static int GpuBatchSize { get; set; } = 262_144; + + private static readonly object GpuWarmUpLock = new(); + private static int WarmedGpuLeaseCount; + + public static void WarmUpGpuKernels(int threadCount) + { + int target = threadCount >> 2; + if (target == 0) + { + target = threadCount; + } + + lock (GpuWarmUpLock) + { + if (target <= WarmedGpuLeaseCount) + { + return; + } + + PrimeOrderCalculatorAccelerator.WarmUp(); + WarmedGpuLeaseCount = target; + } + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static void IsPrimeBatchGpu(PrimeOrderCalculatorAccelerator gpu, ReadOnlySpan values, Span results) + { + // GpuPrimeWorkLimiter.Acquire(); + int acceleratorIndex = gpu.AcceleratorIndex; + int totalLength = values.Length; + int batchSize = GpuBatchSize; + + var inputView = gpu.InputView; + var outputView = gpu.OutputByteView; + + int pos = 0; + AcceleratorStream stream = AcceleratorStreamPool.Rent(acceleratorIndex); + var kernelLauncher = gpu.SmallPrimeSieveKernelLauncher; + + while (pos < totalLength) + { + int remaining = totalLength - pos; + int count = remaining > batchSize ? batchSize : remaining; + + var valueSlice = values.Slice(pos, count); + inputView.CopyFromCPU(stream, valueSlice); + + kernelLauncher( + stream, + count, + inputView, + gpu.DevicePrimesLastOne, + gpu.DevicePrimesLastSeven, + gpu.DevicePrimesLastThree, + gpu.DevicePrimesLastNine, + gpu.DevicePrimesPow2LastOne, + gpu.DevicePrimesPow2LastSeven, + gpu.DevicePrimesPow2LastThree, + gpu.DevicePrimesPow2LastNine, + outputView); + + var resultSlice = results.Slice(pos, count); + outputView.CopyToCPU(stream, resultSlice); + + pos += count; + } + + stream.Synchronize(); + AcceleratorStreamPool.Return(acceleratorIndex, stream); + // GpuPrimeWorkLimiter.Release(); + } + + // Expose cache clearing for accelerator disposal coordination + public static void ClearGpuCaches() + { + PrimeOrderCalculatorAccelerator.Clear(); + } + + internal static void DisposeGpuContexts() + { + PrimeOrderCalculatorAccelerator.DisposeAll(); + lock (GpuWarmUpLock) + { + WarmedGpuLeaseCount = 0; + } + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + internal static bool SharesFactorWithMaxExponent(ulong n) + { + // TODO: Replace this on-the-fly GCD probe with the cached factor table derived from + // ResidueComputationBenchmarks so divisor-cycle metadata can short-circuit the test + // instead of recomputing binary GCD for every candidate. + ulong m = (ulong)BitOperations.Log2(n); + return n.BinaryGcd(m) != 1UL; + } + + internal static void SharesFactorWithMaxExponentBatch(PrimeOrderCalculatorAccelerator gpu, ReadOnlySpan values, Span results) + { + // TODO: Route this batch helper through the shared GPU kernel pool from + // GpuUInt128BinaryGcdBenchmarks so we reuse cached kernels, pinned host buffers, + // and divisor-cycle staging instead of allocating new device buffers per call. + // Check in benchmarks, which implementation was the fastest, is compatible with GPU, + // and implement it. + + int length = values.Length; + + int acceleratorIndex = gpu.AcceleratorIndex; + gpu.EnsureCapacity(0, length); + var inputBufferView = gpu.InputView; + var resultBufferView = gpu.OutputByteView; + var kernelLauncher = gpu.SharesFactorKernelLauncher; + + AcceleratorStream stream = AcceleratorStreamPool.Rent(acceleratorIndex); + inputBufferView.CopyFromCPU(stream, values); + + kernelLauncher(stream, length, inputBufferView, resultBufferView); + + resultBufferView.CopyToCPU(stream, in results); + stream.Synchronize(); + + AcceleratorStreamPool.Return(acceleratorIndex, stream); + } +} diff --git a/PerfectNumbers.Core/PrimeTesterSeverOrOthers.cs b/PerfectNumbers.Core/PrimeTesterSeverOrOthers.cs new file mode 100644 index 0000000..6cf103a --- /dev/null +++ b/PerfectNumbers.Core/PrimeTesterSeverOrOthers.cs @@ -0,0 +1,241 @@ +using System.Numerics; +using System.Runtime.CompilerServices; +using ILGPU.Runtime; +using PerfectNumbers.Core.Gpu.Accelerators; + +namespace PerfectNumbers.Core; + +public sealed class PrimeTesterSeverOrOthers +{ + public PrimeTesterSeverOrOthers() + { + } + + [ThreadStatic] + private static PrimeTesterSeverOrOthers? _tester; + + public static PrimeTesterSeverOrOthers Exclusive => _tester ??= new(); + + public static bool IsPrimeCpu(ulong n) + { + if (n <= 1UL) + { + throw new InvalidOperationException("PrimeTester.IsPrime encountered the sentinel input 2."); + // return false; + } + + if (n == 2UL) + { + throw new InvalidOperationException("PrimeTester.IsPrime encountered the sentinel input 2."); + } + + bool isOdd = (n & 1UL) != 0UL; + bool result = isOdd; + + bool requiresTrialDivision = result && n >= 7UL; + + if (requiresTrialDivision) + { + // EvenPerfectBitScanner streams exponents starting at 136,279,841, so the Mod10/GCD guard never fires on the + // production path. Leave the logic commented out as instrumentation for diagnostic builds. + // bool sharesMaxExponentFactor = n.Mod10() == 1UL && SharesFactorWithMaxExponent(n); + // result &= !sharesMaxExponentFactor; + + if (result) + { + uint[] smallPrimeDivisors; + ulong[] smallPrimeDivisorsMul; + + ulong nMod10 = n.Mod10(); + if (nMod10 == 7) + { + smallPrimeDivisors = PrimesGenerator.SmallPrimesLastSeven; + smallPrimeDivisorsMul = PrimesGenerator.SmallPrimesPow2LastSeven; + } + else + { + smallPrimeDivisors = PrimesGenerator.SmallPrimesLastOne; + smallPrimeDivisorsMul = PrimesGenerator.SmallPrimesPow2LastOne; + } + + int smallPrimeDivisorsLength = smallPrimeDivisors.Length; + for (int i = 0; i < smallPrimeDivisorsLength; i++) + { + ulong value = smallPrimeDivisorsMul[i]; + if (value > n) + { + break; + } + + value = n % smallPrimeDivisors[i]; + if (value == 0) + { + result = false; + break; + } + } + } + } + + return result; + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static bool IsPrimeGpu(PrimeOrderCalculatorAccelerator gpu, ulong n) + { + byte flag = 0; + var inputView = gpu.InputView; + var outputView = gpu.OutputByteView; + + int acceleratorIndex = gpu.AcceleratorIndex; + AcceleratorStream stream = AcceleratorStreamPool.Rent(acceleratorIndex); + inputView.CopyFromCPU(stream, ref n, 1); + + var kernelLauncher = gpu.SmallPrimeSieveKernelLauncher; + + kernelLauncher( + stream, + 1, + inputView, + gpu.DevicePrimesLastOne, + gpu.DevicePrimesLastSeven, + gpu.DevicePrimesLastThree, + gpu.DevicePrimesLastNine, + gpu.DevicePrimesPow2LastOne, + gpu.DevicePrimesPow2LastSeven, + gpu.DevicePrimesPow2LastThree, + gpu.DevicePrimesPow2LastNine, + outputView); + + outputView.CopyToCPU(stream, ref flag, 1); + stream.Synchronize(); + AcceleratorStreamPool.Return(acceleratorIndex, stream); + + return flag != 0; + } + + public static int GpuBatchSize { get; set; } = 262_144; + + private static readonly object GpuWarmUpLock = new(); + private static int WarmedGpuLeaseCount; + + public static void WarmUpGpuKernels(int threadCount) + { + int target = threadCount >> 2; + if (target == 0) + { + target = threadCount; + } + + lock (GpuWarmUpLock) + { + if (target <= WarmedGpuLeaseCount) + { + return; + } + + PrimeOrderCalculatorAccelerator.WarmUp(); + WarmedGpuLeaseCount = target; + } + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static void IsPrimeBatchGpu(PrimeOrderCalculatorAccelerator gpu, ReadOnlySpan values, Span results) + { + // GpuPrimeWorkLimiter.Acquire(); + int acceleratorIndex = gpu.AcceleratorIndex; + int totalLength = values.Length; + int batchSize = GpuBatchSize; + + var inputView = gpu.InputView; + var outputView = gpu.OutputByteView; + + int pos = 0; + AcceleratorStream stream = AcceleratorStreamPool.Rent(acceleratorIndex); + var kernelLauncher = gpu.SmallPrimeSieveKernelLauncher; + + while (pos < totalLength) + { + int remaining = totalLength - pos; + int count = remaining > batchSize ? batchSize : remaining; + + var valueSlice = values.Slice(pos, count); + inputView.CopyFromCPU(stream, valueSlice); + + kernelLauncher( + stream, + count, + inputView, + gpu.DevicePrimesLastOne, + gpu.DevicePrimesLastSeven, + gpu.DevicePrimesLastThree, + gpu.DevicePrimesLastNine, + gpu.DevicePrimesPow2LastOne, + gpu.DevicePrimesPow2LastSeven, + gpu.DevicePrimesPow2LastThree, + gpu.DevicePrimesPow2LastNine, + outputView); + + var resultSlice = results.Slice(pos, count); + outputView.CopyToCPU(stream, resultSlice); + + pos += count; + } + + stream.Synchronize(); + AcceleratorStreamPool.Return(acceleratorIndex, stream); + // GpuPrimeWorkLimiter.Release(); + } + + // Expose cache clearing for accelerator disposal coordination + public static void ClearGpuCaches() + { + PrimeOrderCalculatorAccelerator.Clear(); + } + + internal static void DisposeGpuContexts() + { + PrimeOrderCalculatorAccelerator.DisposeAll(); + lock (GpuWarmUpLock) + { + WarmedGpuLeaseCount = 0; + } + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + internal static bool SharesFactorWithMaxExponent(ulong n) + { + // TODO: Replace this on-the-fly GCD probe with the cached factor table derived from + // ResidueComputationBenchmarks so divisor-cycle metadata can short-circuit the test + // instead of recomputing binary GCD for every candidate. + ulong m = (ulong)BitOperations.Log2(n); + return n.BinaryGcd(m) != 1UL; + } + + internal static void SharesFactorWithMaxExponentBatch(PrimeOrderCalculatorAccelerator gpu, ReadOnlySpan values, Span results) + { + // TODO: Route this batch helper through the shared GPU kernel pool from + // GpuUInt128BinaryGcdBenchmarks so we reuse cached kernels, pinned host buffers, + // and divisor-cycle staging instead of allocating new device buffers per call. + // Check in benchmarks, which implementation was the fastest, is compatible with GPU, + // and implement it. + + int length = values.Length; + + int acceleratorIndex = gpu.AcceleratorIndex; + gpu.EnsureCapacity(0, length); + var inputBufferView = gpu.InputView; + var resultBufferView = gpu.OutputByteView; + var kernelLauncher = gpu.SharesFactorKernelLauncher; + + AcceleratorStream stream = AcceleratorStreamPool.Rent(acceleratorIndex); + inputBufferView.CopyFromCPU(stream, values); + + kernelLauncher(stream, length, inputBufferView, resultBufferView); + + resultBufferView.CopyToCPU(stream, in results); + stream.Synchronize(); + + AcceleratorStreamPool.Return(acceleratorIndex, stream); + } +}