diff --git a/EvenPerfectBitScanner/CliArguments.cs b/EvenPerfectBitScanner/CliArguments.cs index c5701fd..56dda07 100644 --- a/EvenPerfectBitScanner/CliArguments.cs +++ b/EvenPerfectBitScanner/CliArguments.cs @@ -25,17 +25,19 @@ namespace EvenPerfectBitScanner; internal readonly bool UseResidue; internal readonly bool UseDivisor; internal readonly bool UseByDivisor; - internal readonly bool UseDivisorPow2Increment; - internal readonly int ByDivisorSpecialRange; + internal readonly bool UseDivisorPow2Increment; + internal readonly bool UseDivisorPredictiveIncrement; + internal readonly int ByDivisorSpecialRange; internal readonly bool UseGcdFilter; - internal readonly bool TestMode; - internal readonly bool UseGpuCycles; - internal readonly ComputationDevice MersenneDevice; - internal readonly ComputationDevice OrderDevice; - internal readonly int ScanBatchSize; - internal readonly int SliceSize; - internal readonly ulong ResidueKMax; - internal readonly string FilterFile; + internal readonly bool TestMode; + internal readonly bool UseGpuCycles; + internal readonly ComputationDevice MersenneDevice; + internal readonly ComputationDevice OrderDevice; + internal readonly int ScanBatchSize; + internal readonly int SliceSize; + internal readonly ulong MaxKLimit; + internal readonly bool MaxKConfigured; + internal readonly string FilterFile; internal readonly bool ForcePrimeKernelsOnCpu; internal readonly ModReductionMode ModReductionMode; internal readonly NttBackend NttBackend; @@ -75,6 +77,7 @@ namespace EvenPerfectBitScanner; bool useByDivisor, bool useGcdFilter, bool useDivisorPow2Increment, + bool useDivisorPredictiveIncrement, int byDivisorSpecialRange, bool testMode, bool useGpuCycles, @@ -82,7 +85,8 @@ namespace EvenPerfectBitScanner; ComputationDevice orderDevice, int scanBatchSize, int sliceSize, - ulong residueKMax, + ulong maxKLimit, + bool maxKConfigured, string filterFile, bool forcePrimeKernelsOnCpu, ModReductionMode modReductionMode, @@ -121,6 +125,7 @@ namespace EvenPerfectBitScanner; UseDivisor = useDivisor; UseByDivisor = useByDivisor; UseDivisorPow2Increment = useDivisorPow2Increment; + UseDivisorPredictiveIncrement = useDivisorPredictiveIncrement; ByDivisorSpecialRange = byDivisorSpecialRange; UseGcdFilter = useGcdFilter; TestMode = testMode; @@ -129,7 +134,8 @@ namespace EvenPerfectBitScanner; OrderDevice = orderDevice; ScanBatchSize = scanBatchSize; SliceSize = sliceSize; - ResidueKMax = residueKMax; + MaxKLimit = maxKLimit; + MaxKConfigured = maxKConfigured; FilterFile = filterFile; ForcePrimeKernelsOnCpu = forcePrimeKernelsOnCpu; ModReductionMode = modReductionMode; @@ -173,6 +179,7 @@ namespace EvenPerfectBitScanner; bool useDivisor = false; bool useByDivisor = false; bool useDivisorPow2Increment = false; + bool useDivisorPredictiveIncrement = false; int byDivisorSpecialRange = 0; bool useGcdFilter = false; bool testMode = false; @@ -181,7 +188,8 @@ namespace EvenPerfectBitScanner; ComputationDevice orderDevice = ComputationDevice.Gpu; int scanBatchSize = 2_097_152; int sliceSize = 32; - ulong residueKMax = 5_000_000UL; + ulong maxKLimit = 5_000_000UL; + bool maxKConfigured = false; string filterFile = string.Empty; bool forcePrimeKernelsOnCpu = false; ModReductionMode modReductionMode = ModReductionMode.Auto; @@ -189,7 +197,7 @@ namespace EvenPerfectBitScanner; string? resultsDirectory = null; string? resultsPrefix = null; int writeBatchSize = CalculationResultHandler.DefaultWriteBatchSize; - string cyclesPath = DefaultCyclesPath; + string cyclesPath = DefaultCyclesPath; bool continueCyclesGeneration = false; int cyclesBatchSize = 512; ulong divisorCyclesSearchLimit = PerfectNumberConstants.ExtraDivisorCycleSearchLimit; @@ -307,14 +315,21 @@ namespace EvenPerfectBitScanner; if (value.Equals("pow2groups", StringComparison.OrdinalIgnoreCase)) { useDivisorPow2Increment = true; + useDivisorPredictiveIncrement = false; } else if (value.Equals("sequential", StringComparison.OrdinalIgnoreCase)) { useDivisorPow2Increment = false; + useDivisorPredictiveIncrement = false; + } + else if (value.Equals("predictive", StringComparison.OrdinalIgnoreCase)) + { + useDivisorPow2Increment = false; + useDivisorPredictiveIncrement = true; } else { - Console.WriteLine("Invalid value for --bydivisor-k-increment. Use onebyone or pow2groups."); + Console.WriteLine("Invalid value for --bydivisor-k-increment. Use sequential, pow2groups, or predictive."); return default; } @@ -346,11 +361,12 @@ namespace EvenPerfectBitScanner; continue; } - if (argument.StartsWith("--residue-max-k=", StringComparison.OrdinalIgnoreCase)) + if (argument.StartsWith("--max-k=", StringComparison.OrdinalIgnoreCase)) { - if (Utf8CliParser.TryParseUInt64(argument.AsSpan("--residue-max-k=".Length), out ulong parsedResidueMax)) + if (Utf8CliParser.TryParseUInt64(argument.AsSpan("--max-k=".Length), out ulong parsedMaxK)) { - residueKMax = parsedResidueMax; + maxKLimit = parsedMaxK; + maxKConfigured = true; } continue; @@ -637,6 +653,7 @@ namespace EvenPerfectBitScanner; useByDivisor, useGcdFilter, useDivisorPow2Increment, + useDivisorPredictiveIncrement, byDivisorSpecialRange, testMode, useGpuCycles, @@ -644,7 +661,8 @@ namespace EvenPerfectBitScanner; orderDevice, scanBatchSize, sliceSize, - residueKMax, + maxKLimit, + maxKConfigured, filterFile, forcePrimeKernelsOnCpu, modReductionMode, @@ -678,9 +696,9 @@ namespace EvenPerfectBitScanner; Console.WriteLine(" --increment=bit|add exponent increment method (default add)"); Console.WriteLine(" --threads= number of worker threads"); Console.WriteLine(" --mersenne=pow2mod|incremental|lucas|residue|divisor|bydivisor Mersenne test method"); - Console.WriteLine(" --bydivisor-k-increment=sequential|pow2groups k iteration strategy for --mersenne=bydivisor (default sequential)"); + Console.WriteLine(" --bydivisor-k-increment=sequential|pow2groups|predictive k iteration strategy for --mersenne=bydivisor (default sequential)"); Console.WriteLine(" --bydivisor-special-range= widen pow2groups special values to +/- n percent (default 0 = points)"); - Console.WriteLine(" --residue-max-k= max k for residue Mersenne test (q = 2*p*k + 1)"); + Console.WriteLine(" --max-k= max k for Mersenne tests (q = 2*p*k + 1), including --mersenne=bydivisor"); Console.WriteLine(" --mersenne-device=cpu|gpu|hybrid Device for Mersenne method (default gpu)"); Console.WriteLine(" --primes-device=cpu|gpu device for prime-scan kernels (default gpu)"); Console.WriteLine(" --gpu-prime-batch= batch size for GPU primality sieve (default 262144)"); diff --git a/EvenPerfectBitScanner/Program.cs b/EvenPerfectBitScanner/Program.cs index 568bae3..fc9559f 100644 --- a/EvenPerfectBitScanner/Program.cs +++ b/EvenPerfectBitScanner/Program.cs @@ -3,9 +3,11 @@ using System.Diagnostics; using System.Runtime.CompilerServices; using EvenPerfectBitScanner.Candidates; using EvenPerfectBitScanner.IO; +using System.Numerics; using Open.Collections; using PerfectNumbers.Core; using PerfectNumbers.Core.Cpu; +using PerfectNumbers.Core.ByDivisor; using PerfectNumbers.Core.Gpu; using PerfectNumbers.Core.Gpu.Accelerators; using MersenneNumberDivisorByDivisorCpuTesterWithCpuOrder = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForOneByOneDivisorSetForCpuOrder; @@ -14,6 +16,9 @@ using MersenneNumberDivisorByDivisorCpuTesterWithGpuOrder = PerfectNumbers.Core. using ByDivisorCpuOrderPow2 = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPow2GroupsDivisorSetForCpuOrder; using ByDivisorHybridOrderPow2 = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPow2GroupsDivisorSetForHybridOrder; using ByDivisorGpuOrderPow2 = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPow2GroupsDivisorSetForGpuOrder; +using ByDivisorCpuOrderPredictive = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPredictiveDivisorSetForCpuOrder; +using ByDivisorHybridOrderPredictive = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPredictiveDivisorSetForHybridOrder; +using ByDivisorGpuOrderPredictive = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPredictiveDivisorSetForGpuOrder; using System.Globalization; namespace EvenPerfectBitScanner; @@ -31,12 +36,20 @@ private static ThreadLocal? MersenneTesters; private static ByDivisorCpuOrderPow2 _byDivisorCpuOrderTesterPow2; private static ByDivisorHybridOrderPow2 _byDivisorHybridOrderTesterPow2; private static ByDivisorGpuOrderPow2 _byDivisorGpuOrderTesterPow2; + private static ByDivisorCpuOrderPredictive _byDivisorCpuOrderTesterPredictive; + private static ByDivisorHybridOrderPredictive _byDivisorHybridOrderTesterPredictive; + private static ByDivisorGpuOrderPredictive _byDivisorGpuOrderTesterPredictive; private static Dictionary? _byDivisorPreviousResults; + private static string? _byDivisorModelPath; private static bool _limitReached; private static ulong _byDivisorStartPrime; private static CliArguments _cliArguments; private static int _testTargetPrimeCount; private static int _testProcessedPrimeCount; + private static readonly object _byDivisorModelRefreshSync = new(); + private static int _byDivisorPendingRebuilds; + private static DateTime _byDivisorLastRebuildUtc = DateTime.MinValue; + private static readonly TimeSpan _byDivisorRebuildMinInterval = TimeSpan.FromSeconds(30); private static ComputationDevice _mersenneDevice; private static OrderDevice _orderCalculationDevice; private static CacheStatus _cacheStatus; @@ -125,7 +138,8 @@ private static ThreadLocal? MersenneTesters; _orderCalculationDevice = orderCalculationDevice; _cacheStatus = CacheStatus.Disabled; int sliceSize = Math.Max(1, _cliArguments.SliceSize); - ulong residueKMax = _cliArguments.ResidueKMax; + ulong maxKLimit = _cliArguments.MaxKLimit; + bool maxKConfigured = _cliArguments.MaxKConfigured; string filterFile = _cliArguments.FilterFile; UInt128 maxPrimeLimit = _cliArguments.MaxPrimeLimit; bool maxPrimeConfigured = _cliArguments.MaxPrimeConfigured; @@ -161,6 +175,12 @@ private static ThreadLocal? MersenneTesters; Directory.CreateDirectory(PerfectNumberConstants.ByDivisorStateDirectory); + if (useByDivisor && !maxKConfigured) + { + // For by-divisor mode without explicit --max-k, rely on per-prime limits. + maxKLimit = ulong.MaxValue; + } + if (!File.Exists(cyclesPath) || continueCyclesGeneration) { @@ -246,7 +266,7 @@ private static ThreadLocal? MersenneTesters; _orderCalculationDevice, mersenneDevice, kernelType, - residueKMax); + maxKLimit); if (!useLucas) { Console.WriteLine("Warming up orders"); @@ -270,12 +290,17 @@ private static ThreadLocal? MersenneTesters; } bool usePow2Divisors = _cliArguments.UseDivisorPow2Increment; + bool usePredictiveDivisors = _cliArguments.UseDivisorPredictiveIncrement; if (EnvironmentConfiguration.UseCpuOrder) { if (usePow2Divisors) { _byDivisorCpuOrderTesterPow2 = new ByDivisorCpuOrderPow2(); } + else if (usePredictiveDivisors) + { + _byDivisorCpuOrderTesterPredictive = new ByDivisorCpuOrderPredictive(); + } else { _byDivisorCpuOrderTester = new MersenneNumberDivisorByDivisorCpuTesterWithCpuOrder(); @@ -287,6 +312,10 @@ private static ThreadLocal? MersenneTesters; { _byDivisorHybridOrderTesterPow2 = new ByDivisorHybridOrderPow2(); } + else if (usePredictiveDivisors) + { + _byDivisorHybridOrderTesterPredictive = new ByDivisorHybridOrderPredictive(); + } else { _byDivisorHybridOrderTester = new MersenneNumberDivisorByDivisorCpuTesterWithHybridOrder(); @@ -298,6 +327,10 @@ private static ThreadLocal? MersenneTesters; { _byDivisorGpuOrderTesterPow2 = new ByDivisorGpuOrderPow2(); } + else if (usePredictiveDivisors) + { + _byDivisorGpuOrderTesterPredictive = new ByDivisorGpuOrderPredictive(); + } else { _byDivisorGpuOrderTester = new MersenneNumberDivisorByDivisorCpuTesterWithGpuOrder(); @@ -402,6 +435,29 @@ private static ThreadLocal? MersenneTesters; CalculationResultHandler.CreateResultsFileWithHeader(); } + if (useByDivisor) + { + _byDivisorModelPath = Path.Combine(PerfectNumberConstants.ByDivisorStateDirectory, "bydivisor_k10_model.json"); + var byDivisorTuning = new ByDivisorScanTuning + { + CheapLimit = EnvironmentConfiguration.ByDivisorCheapKLimit, + }; + EnvironmentConfiguration.ByDivisorScanTuning = byDivisorTuning; + bool modelExists = File.Exists(_byDivisorModelPath); + bool resultsExist = File.Exists(resultsFileName); + bool needsBuild = !modelExists || (resultsExist && File.GetLastWriteTimeUtc(resultsFileName) > File.GetLastWriteTimeUtc(_byDivisorModelPath)); + Console.WriteLine(needsBuild + ? $"Building by-divisor class model from {resultsFileName}..." + : $"Loading by-divisor class model from {_byDivisorModelPath}..."); + EnvironmentConfiguration.ByDivisorClassModel = ByDivisorClassModelLoader.LoadOrBuild( + _byDivisorModelPath, + resultsFileName, + byDivisorTuning, + byDivisorTuning.CheapLimit, + rebuildWhenResultsNewer: true); + Console.WriteLine($"By-divisor class model ready (cheap<={byDivisorTuning.CheapLimit:F0})."); + } + bool useFilter = !testMode && !string.IsNullOrEmpty(filterFile); HashSet filter = []; List byDivisorCandidates = []; @@ -518,6 +574,7 @@ private static ThreadLocal? MersenneTesters; { int byDivisorPrimeRange = Math.Max(1, _cliArguments.BlockSize); bool usePow2Divisors = _cliArguments.UseDivisorPow2Increment; + bool usePredictiveDivisors = _cliArguments.UseDivisorPredictiveIncrement; if (EnvironmentConfiguration.UseCpuOrder) { if (usePow2Divisors) @@ -530,7 +587,21 @@ private static ThreadLocal? MersenneTesters; _byDivisorStartPrime, static () => _lastCompositeP = true, static () => _lastCompositeP = false, - static (candidate, searchedMersenne, detailedCheck, passedAllTests, divisor) => CalculationResultHandler.HandleResult(candidate, divisor, searchedMersenne, detailedCheck, passedAllTests, _lastCompositeP), + HandleByDivisorResult, + threadCount, + byDivisorPrimeRange); + } + else if (usePredictiveDivisors) + { + Console.WriteLine("Using by-divisor predictive increment (CPU order)."); + MersenneNumberDivisorByDivisorTester.Run( + byDivisorCandidates, + ref _byDivisorCpuOrderTesterPredictive, + _byDivisorPreviousResults, + _byDivisorStartPrime, + static () => _lastCompositeP = true, + static () => _lastCompositeP = false, + HandleByDivisorResult, threadCount, byDivisorPrimeRange); } @@ -543,7 +614,7 @@ private static ThreadLocal? MersenneTesters; _byDivisorStartPrime, static () => _lastCompositeP = true, static () => _lastCompositeP = false, - static (candidate, searchedMersenne, detailedCheck, passedAllTests, divisor) => CalculationResultHandler.HandleResult(candidate, divisor, searchedMersenne, detailedCheck, passedAllTests, _lastCompositeP), + HandleByDivisorResult, threadCount, byDivisorPrimeRange); } @@ -560,7 +631,21 @@ private static ThreadLocal? MersenneTesters; _byDivisorStartPrime, static () => _lastCompositeP = true, static () => _lastCompositeP = false, - static (candidate, searchedMersenne, detailedCheck, passedAllTests, divisor) => CalculationResultHandler.HandleResult(candidate, divisor, searchedMersenne, detailedCheck, passedAllTests, _lastCompositeP), + HandleByDivisorResult, + threadCount, + byDivisorPrimeRange); + } + else if (usePredictiveDivisors) + { + Console.WriteLine("Using by-divisor predictive increment (Hybrid order)."); + MersenneNumberDivisorByDivisorTester.Run( + byDivisorCandidates, + ref _byDivisorHybridOrderTesterPredictive, + _byDivisorPreviousResults, + _byDivisorStartPrime, + static () => _lastCompositeP = true, + static () => _lastCompositeP = false, + HandleByDivisorResult, threadCount, byDivisorPrimeRange); } @@ -573,7 +658,7 @@ private static ThreadLocal? MersenneTesters; _byDivisorStartPrime, static () => _lastCompositeP = true, static () => _lastCompositeP = false, - static (candidate, searchedMersenne, detailedCheck, passedAllTests, divisor) => CalculationResultHandler.HandleResult(candidate, divisor, searchedMersenne, detailedCheck, passedAllTests, _lastCompositeP), + HandleByDivisorResult, threadCount, byDivisorPrimeRange); } @@ -590,7 +675,21 @@ private static ThreadLocal? MersenneTesters; _byDivisorStartPrime, static () => _lastCompositeP = true, static () => _lastCompositeP = false, - static (candidate, searchedMersenne, detailedCheck, passedAllTests, divisor) => CalculationResultHandler.HandleResult(candidate, divisor, searchedMersenne, detailedCheck, passedAllTests, _lastCompositeP), + HandleByDivisorResult, + threadCount, + byDivisorPrimeRange); + } + else if (usePredictiveDivisors) + { + Console.WriteLine("Using by-divisor predictive increment (GPU order)."); + MersenneNumberDivisorByDivisorTester.Run( + byDivisorCandidates, + ref _byDivisorGpuOrderTesterPredictive, + _byDivisorPreviousResults, + _byDivisorStartPrime, + static () => _lastCompositeP = true, + static () => _lastCompositeP = false, + HandleByDivisorResult, threadCount, byDivisorPrimeRange); } @@ -603,7 +702,7 @@ private static ThreadLocal? MersenneTesters; _byDivisorStartPrime, static () => _lastCompositeP = true, static () => _lastCompositeP = false, - static (candidate, searchedMersenne, detailedCheck, passedAllTests, divisor) => CalculationResultHandler.HandleResult(candidate, divisor, searchedMersenne, detailedCheck, passedAllTests, _lastCompositeP), + HandleByDivisorResult, threadCount, byDivisorPrimeRange); } @@ -1016,5 +1115,70 @@ private static ThreadLocal? MersenneTesters; detailedCheck = IsMersennePrime(gpu, p); return detailedCheck; } + + private static void HandleByDivisorResult(ulong candidate, bool searchedMersenne, bool detailedCheck, bool passedAllTests, BigInteger divisor) + { + CalculationResultHandler.HandleResult(candidate, divisor, searchedMersenne, detailedCheck, passedAllTests, _lastCompositeP); + if (divisor.IsZero || _byDivisorModelPath is null || EnvironmentConfiguration.ByDivisorClassModel is null) + { + return; + } + + TryRefreshByDivisorModel(); + } + + private static void TryRefreshByDivisorModel() + { + string? modelPath = _byDivisorModelPath; + if (modelPath is null) + { + return; + } + + string resultsPath = CalculationResultHandler.ResultsFileName; + if (string.IsNullOrEmpty(resultsPath)) + { + return; + } + + bool shouldRefresh = false; + lock (_byDivisorModelRefreshSync) + { + _byDivisorPendingRebuilds++; + DateTime now = DateTime.UtcNow; + if (_byDivisorPendingRebuilds >= 3 || now - _byDivisorLastRebuildUtc >= _byDivisorRebuildMinInterval) + { + _byDivisorPendingRebuilds = 0; + _byDivisorLastRebuildUtc = now; + shouldRefresh = true; + } + } + + if (!shouldRefresh) + { + return; + } + + lock (_byDivisorModelRefreshSync) + { + try + { + double cheapLimit = EnvironmentConfiguration.ByDivisorCheapKLimit; + Console.WriteLine("Refreshing by-divisor class model from results..."); + ByDivisorClassModelLoader.RebuildAndSave(modelPath, resultsPath, cheapLimit); + EnvironmentConfiguration.ByDivisorClassModel = ByDivisorClassModelLoader.LoadOrBuild( + modelPath, + resultsPath, + EnvironmentConfiguration.ByDivisorScanTuning, + cheapLimit, + rebuildWhenResultsNewer: false); + Console.WriteLine($"By-divisor class model refreshed (cheap<={cheapLimit:F0})."); + } + catch (Exception ex) + { + Console.WriteLine($"By-divisor model refresh failed: {ex.Message}"); + } + } + } } diff --git a/PerfectNumbers.Core.Tests/ByDivisorClassModelTests.cs b/PerfectNumbers.Core.Tests/ByDivisorClassModelTests.cs new file mode 100644 index 0000000..6d1855c --- /dev/null +++ b/PerfectNumbers.Core.Tests/ByDivisorClassModelTests.cs @@ -0,0 +1,148 @@ +using System; +using System.IO; +using System.Numerics; +using FluentAssertions; +using PerfectNumbers.Core.ByDivisor; +using Xunit; +using System.Text.Json; + +namespace PerfectNumbers.Core.Tests; + +public sealed class ByDivisorClassModelTests +{ + [Fact] + public void LoadOrBuild_builds_model_from_results_file() + { + string tempPath = Path.GetTempFileName(); + try + { + File.WriteAllLines( + tempPath, + new[] + { + "p,divisor,searchedMersenne,detailedCheck,passedAllTests", + "5,31,true,true,false", + "5,0,true,true,true", + "13,8191,true,true,false", + }); + + string modelPath = tempPath + ".model"; + ByDivisorScanTuning tuning = new() + { + CheapLimit = 32, + }; + + ByDivisorClassModel model = ByDivisorClassModelLoader.LoadOrBuild(modelPath, tempPath, tuning, tuning.CheapLimit); + + ByDivisorClassEntry classFive = model.Entries[5]; + classFive.K50.Should().Be(3d); + classFive.K75.Should().Be(3d); + classFive.GapProbability.Should().BeApproximately(0.5d, 1e-12); + classFive.SampleCount.Should().Be(2); + classFive.GapCount.Should().Be(1); + + ByDivisorClassEntry classThirteen = model.Entries[13]; + classThirteen.K50.Should().Be(315d); + classThirteen.K75.Should().Be(315d); + classThirteen.GapProbability.Should().Be(0d); + } + finally + { + TryDelete(tempPath); + TryDelete(tempPath + ".model"); + } + } + + [Fact] + public void BuildPlan_uses_aggressive_tuning_for_hard_class() + { + var model = new ByDivisorClassModel + { + CheapLimit = 50, + }; + model.Entries[0] = new ByDivisorClassEntry + { + K50 = 1000d, + K75 = 2000d, + GapProbability = 0.8d, + Stability = "escaping", + }; + + ByDivisorScanTuning tuning = new() + { + CheapLimit = 50, + Gamma50 = 0.8d, + Gamma75 = 1.2d, + AggressiveGamma50 = 1.0d, + AggressiveGamma75 = 1.4d, + Rho1 = 1.35d, + Rho2 = 1.6d, + AggressiveRho2 = 1.8d, + HardGapThreshold = 0.75d, + }; + + ByDivisorScanPlan plan = model.BuildPlan(2UL, BigInteger.Zero, tuning); + + plan.CheapLimit.Should().Be(50d); + plan.Start.Should().Be(1000d); + plan.Target.Should().Be(2800d); + plan.Rho1.Should().Be(1.35d); + plan.Rho2.Should().Be(1.8d); + } + + [Fact] + public void LoadOrBuild_rebuilds_when_results_newer_than_model() + { + string tempResults = Path.GetTempFileName(); + string tempModel = tempResults + ".model"; + try + { + File.WriteAllLines( + tempResults, + new[] + { + "p,divisor,searchedMersenne,detailedCheck,passedAllTests", + "5,31,true,true,false", + }); + + // Write a stale model + var staleModel = new ByDivisorClassModel + { + CheapLimit = 10, + }; + File.WriteAllText(tempModel, JsonSerializer.Serialize(staleModel)); + File.SetLastWriteTimeUtc(tempModel, DateTime.UtcNow.AddMinutes(-10)); + File.SetLastWriteTimeUtc(tempResults, DateTime.UtcNow); + + ByDivisorScanTuning tuning = new() + { + CheapLimit = 32, + }; + + ByDivisorClassModel model = ByDivisorClassModelLoader.LoadOrBuild(tempModel, tempResults, tuning, tuning.CheapLimit); + + model.Entries[5].K50.Should().Be(3d); + File.GetLastWriteTimeUtc(tempModel).Should().BeAfter(File.GetLastWriteTimeUtc(tempResults).AddSeconds(-2)); + } + finally + { + TryDelete(tempResults); + TryDelete(tempModel); + } + } + + private static void TryDelete(string path) + { + try + { + if (File.Exists(path)) + { + File.Delete(path); + } + } + catch + { + // ignore cleanup failure in tests + } + } +} diff --git a/PerfectNumbers.Core.Tests/MersenneNumberDivisorCpuTesterTests.cs b/PerfectNumbers.Core.Tests/MersenneNumberDivisorCpuTesterTests.cs index 4836180..c07a042 100644 --- a/PerfectNumbers.Core.Tests/MersenneNumberDivisorCpuTesterTests.cs +++ b/PerfectNumbers.Core.Tests/MersenneNumberDivisorCpuTesterTests.cs @@ -1,4 +1,5 @@ using System.Numerics; +using System.Reflection; using FluentAssertions; using PerfectNumbers.Core.Cpu; using MersenneNumberDivisorByDivisorCpuTesterWithCpuOrder = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForOneByOneDivisorSetForCpuOrder; @@ -99,9 +100,38 @@ public class MersenneNumberDivisorCpuTesterTests tester.ResumeFromState(Path.Combine(Path.GetTempPath(), "bydivisor-selfdivisor.bin"), BigInteger.Zero, kSelf); tester.ResetStateTracking(); - tester.IsPrime(31UL, out bool divisorsExhausted, out BigInteger divisor).Should().BeTrue(); - divisorsExhausted.Should().BeTrue(); - divisor.Should().Be(BigInteger.Zero); + tester.IsPrime(31UL, out bool divisorsExhausted, out BigInteger divisor).Should().BeTrue(); + divisorsExhausted.Should().BeTrue(); + divisor.Should().Be(BigInteger.Zero); + } + + [Theory] + [Trait("Category", "Fast")] + [InlineData(127UL, 127, true)] + [InlineData(127UL, 89, false)] + [InlineData(89UL, 61, false)] + public void IsMersenneValue_flags_only_matching_exponent_as_trivial(ulong prime, int exponent, bool expected) + { + BigInteger divisor = (BigInteger.One << exponent) - BigInteger.One; + InvokeIsMersenneValue(prime, divisor).Should().Be(expected); + } + + [Fact] + [Trait("Category", "Fast")] + public void IsMersenneValue_ignores_mersenne_shaped_divisors_above_prime() + { + ulong prime = 89UL; + BigInteger divisor = (BigInteger.One << 101) - BigInteger.One; + InvokeIsMersenneValue(prime, divisor).Should().BeFalse(); + } + + private static bool InvokeIsMersenneValue(ulong prime, BigInteger divisor) + { + MethodInfo? method = typeof(MersenneNumberDivisorByDivisorCpuTesterWithCpuOrder) + .GetMethod("IsMersenneValue", BindingFlags.NonPublic | BindingFlags.Static); + + method.Should().NotBeNull("IsMersenneValue should exist as a private static helper on the tester type."); + return (bool)method!.Invoke(null, new object[] { prime, divisor })!; } } diff --git a/PerfectNumbers.Core/ByDivisor/ByDivisorClassModel.cs b/PerfectNumbers.Core/ByDivisor/ByDivisorClassModel.cs new file mode 100644 index 0000000..49099bc --- /dev/null +++ b/PerfectNumbers.Core/ByDivisor/ByDivisorClassModel.cs @@ -0,0 +1,344 @@ +using System; +using System.Collections.Generic; +using System.Globalization; +using System.IO; +using System.Linq; +using System.Numerics; +using System.Text.Json; +using System.Text.Json.Serialization; + +namespace PerfectNumbers.Core.ByDivisor; + +public sealed class ByDivisorClassEntry +{ + public double K50 { get; set; } + public double K75 { get; set; } + public double GapProbability { get; set; } + public int SampleCount { get; set; } + public int GapCount { get; set; } + public string Stability { get; set; } = "unknown"; +} + +public sealed class ByDivisorScanPlan +{ + public double CheapLimit { get; init; } + public double Start { get; init; } + public double Target { get; init; } + public double Rho1 { get; init; } + public double Rho2 { get; init; } +} + +public sealed class ByDivisorScanTuning +{ + public double CheapLimit { get; init; } = 10_000d; + public double Gamma50 { get; init; } = 0.8d; + public double Gamma75 { get; init; } = 1.2d; + public double AggressiveGamma50 { get; init; } = 1.0d; + public double AggressiveGamma75 { get; init; } = 1.4d; + public double Rho1 { get; init; } = 1.35d; + public double Rho2 { get; init; } = 1.6d; + public double AggressiveRho2 { get; init; } = 1.8d; + public double HardGapThreshold { get; init; } = 0.75d; + + [JsonIgnore] + public static ByDivisorScanTuning Default => new(); +} + +public sealed class ByDivisorClassModel +{ + private const int ClassCount = 1024; + public int Version { get; set; } = 1; + public double CheapLimit { get; set; } = 10_000d; + public ByDivisorClassEntry[] Entries { get; set; } = Enumerable.Range(0, ClassCount) + .Select(_ => new ByDivisorClassEntry { K50 = 1d, K75 = 1d }) + .ToArray(); + + public ByDivisorScanPlan BuildPlan(ulong prime, BigInteger allowedMax, in ByDivisorScanTuning tuning) + { + int suffix = (int)(prime & 1023UL); + ByDivisorClassEntry entry = Entries[suffix]; + + double k50 = entry.K50 > 1d ? entry.K50 : 1d; + double k75 = entry.K75 > 1d ? entry.K75 : Math.Max(k50, 1d); + bool isHard = entry.GapProbability >= tuning.HardGapThreshold || string.Equals(entry.Stability, "escaping", StringComparison.OrdinalIgnoreCase); + + double gamma50 = isHard ? tuning.AggressiveGamma50 : tuning.Gamma50; + double gamma75 = isHard ? tuning.AggressiveGamma75 : tuning.Gamma75; + double start = Math.Ceiling(Math.Max(1d, gamma50 * k50)); + double target = Math.Ceiling(Math.Max(start + 1d, gamma75 * k75)); + + double cheap = tuning.CheapLimit > 1d ? tuning.CheapLimit : CheapLimit; + double rho2 = isHard ? tuning.AggressiveRho2 : tuning.Rho2; + + double maxKAllowed = ComputeMaxKAllowed(prime, allowedMax); + if (maxKAllowed > 0d) + { + cheap = Math.Min(cheap, maxKAllowed); + start = Math.Min(start, maxKAllowed); + target = Math.Min(target, maxKAllowed); + } + + return new ByDivisorScanPlan + { + CheapLimit = cheap, + Start = start, + Target = target, + Rho1 = tuning.Rho1, + Rho2 = rho2, + }; + } + + private static double ComputeMaxKAllowed(ulong prime, in BigInteger allowedMax) + { + if (prime == 0UL || allowedMax.IsZero) + { + return 0d; + } + + BigInteger step = ((BigInteger)prime) << 1; + if (step.IsZero) + { + return 0d; + } + + BigInteger maxK = (allowedMax - BigInteger.One) / step; + if (maxK <= BigInteger.Zero) + { + return 0d; + } + + double value = (double)maxK; + return double.IsInfinity(value) || value <= 0d ? 0d : value; + } +} + +public static class ByDivisorClassModelLoader +{ + public static ByDivisorClassModel LoadOrBuild(string modelPath, string resultsPath, in ByDivisorScanTuning tuning, double cheapLimit, bool rebuildWhenResultsNewer = true) + { + DateTime modelTimestamp = File.Exists(modelPath) ? File.GetLastWriteTimeUtc(modelPath) : DateTime.MinValue; + DateTime resultsTimestamp = File.Exists(resultsPath) ? File.GetLastWriteTimeUtc(resultsPath) : DateTime.MinValue; + + if (File.Exists(modelPath) && (!rebuildWhenResultsNewer || resultsTimestamp <= modelTimestamp)) + { + try + { + ByDivisorClassModel? loaded = JsonSerializer.Deserialize(File.ReadAllText(modelPath)); + if (loaded is not null) + { + loaded.CheapLimit = cheapLimit; + return loaded; + } + } + catch + { + // fall back to rebuild + } + } + + ByDivisorClassModel model = BuildFromResults(resultsPath, cheapLimit); + Save(modelPath, model); + return model; + } + + public static void RebuildAndSave(string modelPath, string resultsPath, double cheapLimit) + { + ByDivisorClassModel model = BuildFromResults(resultsPath, cheapLimit); + Save(modelPath, model); + } + + private static void Save(string modelPath, ByDivisorClassModel model) + { + string? directory = Path.GetDirectoryName(modelPath); + if (!string.IsNullOrEmpty(directory)) + { + Directory.CreateDirectory(directory); + } + + var options = new JsonSerializerOptions + { + WriteIndented = true, + }; + string json = JsonSerializer.Serialize(model, options); + string tempPath = modelPath + ".tmp"; + File.WriteAllText(tempPath, json); + File.Copy(tempPath, modelPath, overwrite: true); + File.Delete(tempPath); + } + + private static ByDivisorClassModel BuildFromResults(string resultsPath, double cheapLimit) + { + var model = new ByDivisorClassModel + { + CheapLimit = cheapLimit, + }; + + if (!File.Exists(resultsPath)) + { + return model; + } + + var samples = new List[1024]; + Span sampleCounts = stackalloc int[1024]; + Span gapCounts = stackalloc int[1024]; + + using FileStream stream = new(resultsPath, FileMode.Open, FileAccess.Read, FileShare.ReadWrite); + using StreamReader reader = new(stream); + string? line; + bool headerSkipped = false; + + while ((line = reader.ReadLine()) is not null) + { + if (!headerSkipped) + { + headerSkipped = true; + continue; + } + + if (string.IsNullOrWhiteSpace(line)) + { + continue; + } + + string[] parts = line.Split(','); + if (parts.Length < 5) + { + continue; + } + + if (!ulong.TryParse(parts[0], NumberStyles.None, CultureInfo.InvariantCulture, out ulong p)) + { + continue; + } + + if (!BigInteger.TryParse(parts[1], NumberStyles.None, CultureInfo.InvariantCulture, out BigInteger divisor)) + { + continue; + } + + if (!bool.TryParse(parts[3], out bool detailedCheck)) + { + detailedCheck = false; + } + + if (!bool.TryParse(parts[4], out bool passedAllTests)) + { + passedAllTests = false; + } + + int suffix = (int)(p & 1023UL); + bool validDivisor = divisor > BigInteger.One && IsValidDivisor(p, divisor); + if (validDivisor) + { + double kValue = ComputeK(divisor, p); + if (!double.IsNaN(kValue) && !double.IsInfinity(kValue)) + { + samples[suffix] ??= new List(); + samples[suffix]!.Add(kValue); + } + } + + if (detailedCheck) + { + sampleCounts[suffix]++; + if (!validDivisor || passedAllTests) + { + gapCounts[suffix]++; + } + } + } + + for (int i = 0; i < 1024; i++) + { + List? classSamples = samples[i]; + double k50 = 1d; + double k75 = 1d; + if (classSamples is not null && classSamples.Count > 0) + { + classSamples.Sort(); + k50 = ComputePercentile(classSamples, 0.5d); + k75 = ComputePercentile(classSamples, 0.75d); + } + + int total = sampleCounts[i]; + int gaps = gapCounts[i]; + double gapProb = total > 0 ? (double)gaps / total : 0d; + + model.Entries[i] = new ByDivisorClassEntry + { + K50 = k50, + K75 = k75, + GapProbability = gapProb, + SampleCount = total, + GapCount = gaps, + Stability = model.Entries[i].Stability, + }; + } + + return model; + } + + private static bool IsValidDivisor(ulong prime, in BigInteger divisor) + { + if (divisor <= BigInteger.One) + { + return false; + } + + BigInteger modulus = ((BigInteger)prime) << 1; + if (modulus.IsZero) + { + return false; + } + + if ((divisor - BigInteger.One) % modulus != BigInteger.Zero) + { + return false; + } + + return BigInteger.ModPow(2, prime, divisor) == BigInteger.One; + } + + private static double ComputeK(in BigInteger divisor, ulong prime) + { + BigInteger numerator = divisor - BigInteger.One; + BigInteger denominator = ((BigInteger)prime) << 1; + if (denominator.IsZero) + { + return 0d; + } + + BigInteger k = numerator / denominator; + if (k <= BigInteger.Zero) + { + return 0d; + } + + double value = (double)k; + if (double.IsInfinity(value) || double.IsNaN(value) || value <= 0d) + { + return double.MaxValue; + } + + return value; + } + + private static double ComputePercentile(List values, double percentile) + { + if (values.Count == 0) + { + return 1d; + } + + double index = (values.Count - 1) * percentile; + int lower = (int)Math.Floor(index); + int upper = (int)Math.Ceiling(index); + if (lower == upper) + { + return values[lower]; + } + + double fraction = index - lower; + return values[lower] + ((values[upper] - values[lower]) * fraction); + } +} diff --git a/PerfectNumbers.Core/Cpu/MersenneNumberDivisorByDivisorCpuTesterWithTemplate.cs b/PerfectNumbers.Core/Cpu/MersenneNumberDivisorByDivisorCpuTesterWithTemplate.cs index 89133a1..54be908 100644 --- a/PerfectNumbers.Core/Cpu/MersenneNumberDivisorByDivisorCpuTesterWithTemplate.cs +++ b/PerfectNumbers.Core/Cpu/MersenneNumberDivisorByDivisorCpuTesterWithTemplate.cs @@ -6,6 +6,7 @@ using System.Collections.Generic; using System.Collections.Concurrent; using PerfectNumbers.Core.Gpu; using PerfectNumbers.Core.Gpu.Accelerators; +using PerfectNumbers.Core.ByDivisor; namespace PerfectNumbers.Core.Cpu; @@ -13,6 +14,7 @@ public enum DivisorSet { OneByOne, Pow2Groups, + Predictive, } [EnumDependentTemplate(typeof(DivisorSet))] @@ -71,19 +73,24 @@ public struct MersenneNumberDivisorByDivisorCpuTesterWithTemplate() : IMersenneN [MethodImpl(MethodImplOptions.AggressiveInlining)] private static bool IsMersenneValue(ulong prime, in BigInteger divisor) { - if (divisor <= BigInteger.One || prime == 0UL || prime > int.MaxValue) - { - return false; - } - - int divisorBitLength = GetBitLengthPortable(divisor); - if (divisorBitLength != (int)prime) + if (divisor <= BigInteger.One || prime == 0UL) { return false; } BigInteger plusOne = divisor + BigInteger.One; - return plusOne > BigInteger.One && IsPowerOfTwoBig(plusOne); + if (plusOne <= BigInteger.One || !IsPowerOfTwoBig(plusOne)) + { + return false; + } + + int exponentBits = GetBitLengthPortable(plusOne) - 1; + if (exponentBits <= 0) + { + return false; + } + + return prime == (ulong)exponentBits; } [MethodImpl(MethodImplOptions.AggressiveInlining)] @@ -331,8 +338,36 @@ public struct MersenneNumberDivisorByDivisorCpuTesterWithTemplate() : IMersenneN step, out processedAll, out foundDivisor); + // Console.WriteLine($"[by-divisor-plan] p={prime} suffix={(prime & 1023UL)} using pow2-groups plan; model skipped."); return compositePow2; #else +#if DivisorSet_Predictive + ByDivisorClassModel? classifier = EnvironmentConfiguration.ByDivisorClassModel; + if (classifier is not null) + { + bool compositePlanned = TryExecutePlannedScan( + classifier, + prime, + allowedMax, + normalizedMinK, + out processedAll, + out foundDivisor); + if (compositePlanned) + { + return true; + } + + if (processedAll) + { + return false; + } + } + else + { + Console.WriteLine($"[by-divisor-plan] p={prime} suffix={(prime & 1023UL)} model unavailable, falling back to sequential scan."); + } +#endif + if (fits64 && composite64) { allowedMax64 = unlimited ? ulong.MaxValue : (ulong)allowedMax; @@ -556,6 +591,400 @@ public struct MersenneNumberDivisorByDivisorCpuTesterWithTemplate() : IMersenneN } #endif + #if !DivisorSet_Pow2Groups && DivisorSet_Predictive + [MethodImpl(MethodImplOptions.AggressiveOptimization)] + private bool TryExecutePlannedScan( + ByDivisorClassModel model, + ulong prime, + BigInteger allowedMax, + BigInteger minK, + out bool processedAll, + out BigInteger foundDivisor) + { + processedAll = false; + foundDivisor = BigInteger.Zero; + + if (minK > ulong.MaxValue) + { + // Console.WriteLine($"[by-divisor-plan] p={prime} using BigInteger plan (minK={minK})."); + return TryExecutePlannedScanBig(model, prime, allowedMax, minK, out processedAll, out foundDivisor); + } + + ulong allowedMax64; + if (allowedMax.IsZero) + { + allowedMax64 = ulong.MaxValue; + } + else if (allowedMax > ulong.MaxValue) + { + // Console.WriteLine($"[by-divisor-plan] p={prime} using BigInteger plan (allowedMax={allowedMax})."); + return TryExecutePlannedScanBig(model, prime, allowedMax, minK, out processedAll, out foundDivisor); + } + else + { + allowedMax64 = (ulong)allowedMax; + } + bool composite = TryExecutePlannedScan64( + model, + prime, + allowedMax64, + (ulong)minK, + out processedAll, + out ulong foundDivisor64); + if (composite) + { + foundDivisor = foundDivisor64; + return true; + } + + if (processedAll) + { + foundDivisor = BigInteger.Zero; + } + + return false; + } + + [MethodImpl(MethodImplOptions.AggressiveOptimization)] + private bool TryExecutePlannedScanBig( + ByDivisorClassModel model, + ulong prime, + in BigInteger allowedMax, + in BigInteger minK, + out bool processedAll, + out BigInteger foundDivisor) + { + processedAll = true; + foundDivisor = BigInteger.Zero; + + if (allowedMax.IsZero) + { + processedAll = false; + return false; + } + + BigInteger step = ((BigInteger)prime) << 1; + BigInteger maxKAllowed = (allowedMax - BigInteger.One) / step; + TryClampMaxKForMersenneTail(prime, step, ref maxKAllowed); + if (maxKAllowed <= BigInteger.Zero || minK > maxKAllowed) + { + processedAll = false; + return false; + } + + ByDivisorScanTuning tuning = EnvironmentConfiguration.ByDivisorScanTuning; + ByDivisorScanPlan plan = model.BuildPlan(prime, allowedMax, tuning); + int suffix = (int)(prime & 1023UL); + // ByDivisorClassEntry entry = model.Entries[suffix]; + // Console.WriteLine($"[by-divisor-plan] p={prime} suffix={suffix} k50={entry.K50:F2} k75={entry.K75:F2} gap={entry.GapProbability:F3} cheap<= {plan.CheapLimit:F0} start={plan.Start:F0} target={plan.Target:F0} rho1={plan.Rho1:F2} rho2={plan.Rho2:F2} maxK={maxKAllowed}"); + + BigInteger cheapLimit = BigInteger.Min((BigInteger)plan.CheapLimit, maxKAllowed); + if (cheapLimit >= minK) + { + BigInteger currentK = minK; + bool compositeCheap = CheckDivisorsBigPlanRange( + prime, + step, + allowedMax, + minK, + cheapLimit, + out bool processedRange, + out BigInteger foundRange, + ref currentK); + if (compositeCheap) + { + foundDivisor = foundRange; + return true; + } + + if (!processedRange) + { + processedAll = false; + return false; + } + } + + BigInteger startPlan = plan.Start <= 0d ? minK : new BigInteger(Math.Ceiling(plan.Start)); + if (startPlan < minK) + { + startPlan = minK; + } + + BigInteger target = plan.Target <= 0d ? startPlan + BigInteger.One : new BigInteger(Math.Ceiling(plan.Target)); + if (target < startPlan) + { + target = startPlan + BigInteger.One; + } + + BigInteger current = cheapLimit + BigInteger.One; + if (current < startPlan) + { + current = startPlan; + } + + bool loggedProgressive = false; + while (current <= maxKAllowed) + { + if (!loggedProgressive) + { + // Console.WriteLine($"[by-divisor-progressive] p={prime} startK={current} target={target} rho={(current < target ? plan.Rho1 : plan.Rho2):F2} maxK={maxKAllowed}"); + loggedProgressive = true; + } + + BigInteger startRange = current; + BigInteger endRange = current; + bool composite = CheckDivisorsBigPlanRange( + prime, + step, + allowedMax, + startRange, + endRange, + out bool processedRange, + out BigInteger foundRange, + ref startRange); + if (composite) + { + foundDivisor = foundRange; + return true; + } + + if (!processedRange) + { + processedAll = false; + return false; + } + + double multiplier = current < target ? plan.Rho1 : plan.Rho2; + BigInteger advanced = new BigInteger(Math.Ceiling((double)current * multiplier)); + if (advanced <= current) + { + advanced = current + BigInteger.One; + } + + if (advanced > maxKAllowed) + { + break; + } + + current = advanced; + } + + return false; + } + + [MethodImpl(MethodImplOptions.AggressiveOptimization)] + private bool CheckDivisorsBigPlanRange( + ulong prime, + in BigInteger step, + in BigInteger allowedMax, + in BigInteger startK, + in BigInteger endK, + out bool processedAll, + out BigInteger foundDivisor, + ref BigInteger currentK) + { + processedAll = true; + foundDivisor = BigInteger.Zero; + + if (startK < BigInteger.One || endK < startK) + { + return false; + } + + BigInteger maxKAllowed = (allowedMax - BigInteger.One) / step; + TryClampMaxKForMersenneTail(prime, step, ref maxKAllowed); + if (startK > maxKAllowed) + { + return false; + } + + BigInteger effectiveEnd = endK > maxKAllowed ? maxKAllowed : endK; + ushort primeDecimalMask = DivisorGenerator.GetDecimalMask((prime & 3UL) == 3UL); + BigInteger divisor = (step * startK) + BigInteger.One; + if (divisor > allowedMax) + { + return false; + } + + BigInteger iterations = (effectiveEnd - startK) + BigInteger.One; + currentK = startK; + + var rem10 = new BigIntegerCycleRemainderStepper(10); + var rem8 = new BigIntegerCycleRemainderStepper(8); + var rem3 = new BigIntegerCycleRemainderStepper(3); + var rem7 = new BigIntegerCycleRemainderStepper(7); + var rem11 = new BigIntegerCycleRemainderStepper(11); + var rem13 = new BigIntegerCycleRemainderStepper(13); + var rem17 = new BigIntegerCycleRemainderStepper(17); + var rem19 = new BigIntegerCycleRemainderStepper(19); + + byte remainder10 = rem10.Initialize(divisor); + BigInteger primeBig = (BigInteger)prime; + + while (iterations.Sign > 0) + { + if (remainder10 == 1 || remainder10 == 7) + { + if (rem3.NextIsNotDivisible(divisor) && + rem7.NextIsNotDivisible(divisor) && + rem11.NextIsNotDivisible(divisor) && + rem13.NextIsNotDivisible(divisor) && + rem17.NextIsNotDivisible(divisor) && + rem19.NextIsNotDivisible(divisor)) + { + byte remainder8 = rem8.ComputeNext(divisor); + if (remainder8 == 1 || remainder8 == 7) + { + if (((primeDecimalMask >> remainder10) & 1) != 0) + { + RecordState(currentK); + BigInteger powResult = BigInteger.ModPow(BigIntegerNumbers.Two, primeBig, divisor); + if (powResult.IsOne && !IsMersenneValue(prime, divisor)) + { + foundDivisor = divisor; + return true; + } + } + } + } + } + + iterations--; + if (iterations.IsZero) + { + break; + } + + currentK++; + divisor += step; + remainder10 = rem10.ComputeNext(divisor); + } + + return false; + } + + [MethodImpl(MethodImplOptions.AggressiveOptimization)] + private bool TryExecutePlannedScan64( + ByDivisorClassModel model, + ulong prime, + ulong allowedMax, + ulong minK, + out bool processedAll, + out ulong foundDivisor) + { + processedAll = true; + foundDivisor = 0UL; + + ulong step64 = prime << 1; + if (step64 == 0UL) + { + processedAll = false; + return false; + } + + ByDivisorScanTuning tuning = EnvironmentConfiguration.ByDivisorScanTuning; + ByDivisorScanPlan plan = model.BuildPlan(prime, allowedMax == ulong.MaxValue ? BigInteger.Zero : new BigInteger(allowedMax), tuning); + int suffix = (int)(prime & 1023UL); + // ByDivisorClassEntry entry = model.Entries[suffix]; + // Console.WriteLine($"[by-divisor-plan] p={prime} suffix={suffix} k50={entry.K50:F2} k75={entry.K75:F2} gap={entry.GapProbability:F3} cheap<= {plan.CheapLimit:F0} start={plan.Start:F0} target={plan.Target:F0} rho1={plan.Rho1:F2} rho2={plan.Rho2:F2}"); + + ulong effectiveAllowed = allowedMax == 0UL ? ulong.MaxValue : allowedMax; + ulong maxKAllowed = effectiveAllowed == ulong.MaxValue ? ulong.MaxValue : ((effectiveAllowed - 1UL) / step64); + if (maxKAllowed == 0UL || maxKAllowed < minK) + { + processedAll = false; + return false; + } + + ulong cheapLimit = (ulong)Math.Min(plan.CheapLimit, maxKAllowed); + if (cheapLimit >= minK) + { + ulong currentK = minK; + bool composite = CheckDivisors64Range( + prime, + step64, + effectiveAllowed, + minK, + cheapLimit, + ref currentK, + out bool processedRange, + out ulong foundRange); + if (composite) + { + foundDivisor = foundRange; + return true; + } + + if (!processedRange) + { + processedAll = false; + return false; + } + } + + double nextK = Math.Max(plan.Start, cheapLimit + 1d); + if (nextK < minK) + { + nextK = minK; + } + + double target = plan.Target <= 0d ? nextK : plan.Target; + if (target < nextK) + { + target = nextK + 1d; + } + + ulong current = (ulong)Math.Ceiling(nextK); + bool loggedProgressive = false; + while (current <= maxKAllowed) + { + if (!loggedProgressive) + { + Console.WriteLine($"[by-divisor-progressive] p={prime} startK={current} target={target:F0} rho={(current < target ? plan.Rho1 : plan.Rho2):F2} maxK={maxKAllowed}"); + loggedProgressive = true; + } + + ulong currentK = current; + bool compositeSingle = CheckDivisors64Range( + prime, + step64, + effectiveAllowed, + currentK, + currentK, + ref currentK, + out bool processedSingle, + out ulong foundSingle); + if (compositeSingle) + { + foundDivisor = foundSingle; + return true; + } + + if (!processedSingle) + { + processedAll = false; + return false; + } + + double multiplier = current < target ? plan.Rho1 : plan.Rho2; + double advanced = Math.Ceiling(current * multiplier); + if (advanced <= current) + { + advanced = current + 1d; + } + + if (advanced > maxKAllowed) + { + break; + } + + current = (ulong)advanced; + } + + return false; + } + #endif + [MethodImpl(MethodImplOptions.AggressiveOptimization)] private bool CheckDivisors64Bit( ulong prime, diff --git a/PerfectNumbers.Core/Cpu/MersenneNumberDivisorByDivisorPrimeScanSessionWithTemplate.cs b/PerfectNumbers.Core/Cpu/MersenneNumberDivisorByDivisorPrimeScanSessionWithTemplate.cs index cd06b14..775051f 100644 --- a/PerfectNumbers.Core/Cpu/MersenneNumberDivisorByDivisorPrimeScanSessionWithTemplate.cs +++ b/PerfectNumbers.Core/Cpu/MersenneNumberDivisorByDivisorPrimeScanSessionWithTemplate.cs @@ -11,6 +11,10 @@ using System.Runtime.CompilerServices; using TesterCpu = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPow2GroupsDivisorSetForCpuOrder; using TesterHybrid = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPow2GroupsDivisorSetForHybridOrder; using TesterGpu = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPow2GroupsDivisorSetForGpuOrder; +#elif DivisorSet_Predictive +using TesterCpu = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPredictiveDivisorSetForCpuOrder; +using TesterHybrid = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPredictiveDivisorSetForHybridOrder; +using TesterGpu = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPredictiveDivisorSetForGpuOrder; #else using TesterCpu = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForOneByOneDivisorSetForCpuOrder; using TesterHybrid = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForOneByOneDivisorSetForHybridOrder; diff --git a/PerfectNumbers.Core/EnvironmentConfiguration.cs b/PerfectNumbers.Core/EnvironmentConfiguration.cs index 624b69e..1be8f67 100644 --- a/PerfectNumbers.Core/EnvironmentConfiguration.cs +++ b/PerfectNumbers.Core/EnvironmentConfiguration.cs @@ -20,6 +20,9 @@ public static class EnvironmentConfiguration public static bool UseCpuOrder; public static bool UsePow2GroupDivisors; public static int ByDivisorSpecialRange; + public static ByDivisor.ByDivisorClassModel? ByDivisorClassModel; + public static ByDivisor.ByDivisorScanTuning ByDivisorScanTuning = ByDivisor.ByDivisorScanTuning.Default; + public static double ByDivisorCheapKLimit = 100_000d; [MethodImpl(MethodImplOptions.AggressiveInlining | MethodImplOptions.AggressiveOptimization)] public static void Initialize() diff --git a/PerfectNumbers.Core/MersenneNumberDivisorByDivisorTester.cs b/PerfectNumbers.Core/MersenneNumberDivisorByDivisorTester.cs index 489a40d..4ea7e0f 100644 --- a/PerfectNumbers.Core/MersenneNumberDivisorByDivisorTester.cs +++ b/PerfectNumbers.Core/MersenneNumberDivisorByDivisorTester.cs @@ -9,6 +9,9 @@ using System.Runtime.CompilerServices; using MersenneNumberDivisorByDivisorCpuTesterWithCpuOrder = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForOneByOneDivisorSetForCpuOrder; using MersenneNumberDivisorByDivisorCpuTesterWithHybridOrder = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForOneByOneDivisorSetForHybridOrder; using MersenneNumberDivisorByDivisorCpuTesterWithGpuOrder = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForOneByOneDivisorSetForGpuOrder; +using MersenneNumberDivisorByDivisorCpuTesterWithPredictiveCpuOrder = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPredictiveDivisorSetForCpuOrder; +using MersenneNumberDivisorByDivisorCpuTesterWithPredictiveHybridOrder = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPredictiveDivisorSetForHybridOrder; +using MersenneNumberDivisorByDivisorCpuTesterWithPredictiveGpuOrder = PerfectNumbers.Core.Cpu.MersenneNumberDivisorByDivisorCpuTesterWithForPredictiveDivisorSetForGpuOrder; namespace PerfectNumbers.Core; @@ -45,6 +48,37 @@ public static class MersenneNumberDivisorByDivisorTester Action pr) => RunCpuOrderFiltered(filteredPrimes, prototype, workerCount, partitionSize, mc, cc, pr)); + [MethodImpl(MethodImplOptions.AggressiveInlining | MethodImplOptions.AggressiveOptimization)] + public static void Run( + List candidates, + ref MersenneNumberDivisorByDivisorCpuTesterWithPredictiveCpuOrder tester, + Dictionary? previousResults, + ulong startPrime, + Action markComposite, + Action clearComposite, + Action printResult, + int threadCount, + int primesPerTask) + => RunStructTester( + candidates, + ref tester, + previousResults, + startPrime, + markComposite, + clearComposite, + printResult, + threadCount, + primesPerTask, + static ( + List filteredPrimes, + in MersenneNumberDivisorByDivisorCpuTesterWithPredictiveCpuOrder prototype, + int workerCount, + int partitionSize, + Action mc, + Action cc, + Action pr) => + RunPredictiveCpuOrderFiltered(filteredPrimes, prototype, workerCount, partitionSize, mc, cc, pr)); + [MethodImpl(MethodImplOptions.AggressiveInlining | MethodImplOptions.AggressiveOptimization)] public static void Run( List candidates, @@ -76,6 +110,37 @@ public static class MersenneNumberDivisorByDivisorTester Action pr) => RunHybridOrderFiltered(filteredPrimes, prototype, workerCount, partitionSize, mc, cc, pr)); + [MethodImpl(MethodImplOptions.AggressiveInlining | MethodImplOptions.AggressiveOptimization)] + public static void Run( + List candidates, + ref MersenneNumberDivisorByDivisorCpuTesterWithPredictiveHybridOrder tester, + Dictionary? previousResults, + ulong startPrime, + Action markComposite, + Action clearComposite, + Action printResult, + int threadCount, + int primesPerTask) + => RunStructTester( + candidates, + ref tester, + previousResults, + startPrime, + markComposite, + clearComposite, + printResult, + threadCount, + primesPerTask, + static ( + List filteredPrimes, + in MersenneNumberDivisorByDivisorCpuTesterWithPredictiveHybridOrder prototype, + int workerCount, + int partitionSize, + Action mc, + Action cc, + Action pr) => + RunPredictiveHybridOrderFiltered(filteredPrimes, prototype, workerCount, partitionSize, mc, cc, pr)); + [MethodImpl(MethodImplOptions.AggressiveInlining | MethodImplOptions.AggressiveOptimization)] public static void Run( List candidates, @@ -107,6 +172,37 @@ public static class MersenneNumberDivisorByDivisorTester Action pr) => RunGpuOrderFiltered(filteredPrimes, prototype, workerCount, partitionSize, mc, cc, pr)); + [MethodImpl(MethodImplOptions.AggressiveInlining | MethodImplOptions.AggressiveOptimization)] + public static void Run( + List candidates, + ref MersenneNumberDivisorByDivisorCpuTesterWithPredictiveGpuOrder tester, + Dictionary? previousResults, + ulong startPrime, + Action markComposite, + Action clearComposite, + Action printResult, + int threadCount, + int primesPerTask) + => RunStructTester( + candidates, + ref tester, + previousResults, + startPrime, + markComposite, + clearComposite, + printResult, + threadCount, + primesPerTask, + static ( + List filteredPrimes, + in MersenneNumberDivisorByDivisorCpuTesterWithPredictiveGpuOrder prototype, + int workerCount, + int partitionSize, + Action mc, + Action cc, + Action pr) => + RunPredictiveGpuOrderFiltered(filteredPrimes, prototype, workerCount, partitionSize, mc, cc, pr)); + [MethodImpl(MethodImplOptions.AggressiveInlining | MethodImplOptions.AggressiveOptimization)] public static void Run( List candidates, @@ -1056,6 +1152,84 @@ public static class MersenneNumberDivisorByDivisorTester startGate.Dispose(); } + private static void RunPredictiveCpuOrderFiltered( + List filteredPrimes, + in MersenneNumberDivisorByDivisorCpuTesterWithPredictiveCpuOrder prototype, + int workerCount, + int partitionSize, + Action markComposite, + Action clearComposite, + Action printResult) + { + MersenneNumberDivisorByDivisorCpuTesterWithPredictiveCpuOrder prototypeCopy = prototype; + int totalCount = filteredPrimes.Count; + if (partitionSize > totalCount) + { + partitionSize = totalCount; + } + + if (workerCount == 1) + { + DeterministicRandomCpu.Initialize(); + var session = new MersenneNumberDivisorByDivisorPrimeScanSessionWithForPredictiveDivisorSetForCpuOrder(prototypeCopy, markComposite, clearComposite, printResult); + foreach (ulong prime in filteredPrimes) + { + session.ProcessPrime(prime); + } + + session.Dispose(); + return; + } + + TaskScheduler scheduler = UnboundedTaskScheduler.Instance; + int nextIndex = 0; + Task[] tasks = new Task[workerCount]; + var startGate = new ManualResetEventSlim(initialState: false); + + for (int taskIndex = 0; taskIndex < workerCount; taskIndex++) + { + int workerIndex = taskIndex; + tasks[taskIndex] = Task.Factory.StartNew( + () => + { + DeterministicRandomCpu.Initialize(); + var session = new MersenneNumberDivisorByDivisorPrimeScanSessionWithForPredictiveDivisorSetForCpuOrder(prototypeCopy, markComposite, clearComposite, printResult); + startGate.Wait(); + Console.WriteLine($"Task started for worker {workerIndex.ToString(CultureInfo.InvariantCulture)}"); + + while (true) + { + int start = Interlocked.Add(ref nextIndex, partitionSize) - partitionSize; + if (start >= totalCount) + { + break; + } + + int end = start + partitionSize; + if (end > totalCount) + { + end = totalCount; + } + + for (int index = start; index < end; index++) + { + session.ProcessPrime(filteredPrimes[index]); + } + } + + session.Dispose(); + Console.WriteLine($"Task finished for worker {workerIndex.ToString(CultureInfo.InvariantCulture)}"); + }, + CancellationToken.None, + TaskCreationOptions.DenyChildAttach, + scheduler); + } + + startGate.Set(); + Task.WaitAll(tasks); + startGate.Dispose(); + } + private static void RunHybridOrderFiltered( List filteredPrimes, in MersenneNumberDivisorByDivisorCpuTesterWithHybridOrder prototype, @@ -1147,6 +1321,84 @@ public static class MersenneNumberDivisorByDivisorTester startGate.Dispose(); } + private static void RunPredictiveHybridOrderFiltered( + List filteredPrimes, + in MersenneNumberDivisorByDivisorCpuTesterWithPredictiveHybridOrder prototype, + int workerCount, + int partitionSize, + Action markComposite, + Action clearComposite, + Action printResult) + { + MersenneNumberDivisorByDivisorCpuTesterWithPredictiveHybridOrder prototypeCopy = prototype; + int totalCount = filteredPrimes.Count; + if (partitionSize > totalCount) + { + partitionSize = totalCount; + } + + if (workerCount == 1) + { + DeterministicRandomCpu.Initialize(); + var session = new MersenneNumberDivisorByDivisorPrimeScanSessionWithForPredictiveDivisorSetForHybridOrder(prototypeCopy, markComposite, clearComposite, printResult); + foreach (ulong prime in filteredPrimes) + { + session.ProcessPrime(prime); + } + + session.Dispose(); + return; + } + + TaskScheduler scheduler = UnboundedTaskScheduler.Instance; + int nextIndex = 0; + Task[] tasks = new Task[workerCount]; + var startGate = new ManualResetEventSlim(initialState: false); + + for (int taskIndex = 0; taskIndex < workerCount; taskIndex++) + { + int workerIndex = taskIndex; + tasks[taskIndex] = Task.Factory.StartNew( + () => + { + DeterministicRandomCpu.Initialize(); + var session = new MersenneNumberDivisorByDivisorPrimeScanSessionWithForPredictiveDivisorSetForHybridOrder(prototypeCopy, markComposite, clearComposite, printResult); + startGate.Wait(); + Console.WriteLine($"Task started for worker {workerIndex.ToString(CultureInfo.InvariantCulture)}"); + + while (true) + { + int start = Interlocked.Add(ref nextIndex, partitionSize) - partitionSize; + if (start >= totalCount) + { + break; + } + + int end = start + partitionSize; + if (end > totalCount) + { + end = totalCount; + } + + for (int index = start; index < end; index++) + { + session.ProcessPrime(filteredPrimes[index]); + } + } + + session.Dispose(); + Console.WriteLine($"Task finished for worker {workerIndex.ToString(CultureInfo.InvariantCulture)}"); + }, + CancellationToken.None, + TaskCreationOptions.DenyChildAttach, + scheduler); + } + + startGate.Set(); + Task.WaitAll(tasks); + startGate.Dispose(); + } + private static void RunGpuOrderFiltered( List filteredPrimes, in MersenneNumberDivisorByDivisorCpuTesterWithGpuOrder prototype, @@ -1238,6 +1490,84 @@ public static class MersenneNumberDivisorByDivisorTester startGate.Dispose(); } + private static void RunPredictiveGpuOrderFiltered( + List filteredPrimes, + in MersenneNumberDivisorByDivisorCpuTesterWithPredictiveGpuOrder prototype, + int workerCount, + int partitionSize, + Action markComposite, + Action clearComposite, + Action printResult) + { + MersenneNumberDivisorByDivisorCpuTesterWithPredictiveGpuOrder prototypeCopy = prototype; + int totalCount = filteredPrimes.Count; + if (partitionSize > totalCount) + { + partitionSize = totalCount; + } + + if (workerCount == 1) + { + DeterministicRandomCpu.Initialize(); + var session = new MersenneNumberDivisorByDivisorPrimeScanSessionWithForPredictiveDivisorSetForGpuOrder(prototypeCopy, markComposite, clearComposite, printResult); + foreach (ulong prime in filteredPrimes) + { + session.ProcessPrime(prime); + } + + session.Dispose(); + return; + } + + TaskScheduler scheduler = UnboundedTaskScheduler.Instance; + int nextIndex = 0; + Task[] tasks = new Task[workerCount]; + var startGate = new ManualResetEventSlim(initialState: false); + + for (int taskIndex = 0; taskIndex < workerCount; taskIndex++) + { + int workerIndex = taskIndex; + tasks[taskIndex] = Task.Factory.StartNew( + () => + { + DeterministicRandomCpu.Initialize(); + var session = new MersenneNumberDivisorByDivisorPrimeScanSessionWithForPredictiveDivisorSetForGpuOrder(prototypeCopy, markComposite, clearComposite, printResult); + startGate.Wait(); + Console.WriteLine($"Task started for worker {workerIndex.ToString(CultureInfo.InvariantCulture)}"); + + while (true) + { + int start = Interlocked.Add(ref nextIndex, partitionSize) - partitionSize; + if (start >= totalCount) + { + break; + } + + int end = start + partitionSize; + if (end > totalCount) + { + end = totalCount; + } + + for (int index = start; index < end; index++) + { + session.ProcessPrime(filteredPrimes[index]); + } + } + + session.Dispose(); + Console.WriteLine($"Task finished for worker {workerIndex.ToString(CultureInfo.InvariantCulture)}"); + }, + CancellationToken.None, + TaskCreationOptions.DenyChildAttach, + scheduler); + } + + startGate.Set(); + Task.WaitAll(tasks); + startGate.Dispose(); + } + private static void RunGpuTesterFiltered( List filteredPrimes, MersenneNumberDivisorByDivisorGpuTester prototype, diff --git a/docs/ByDivisorPredictiveScanPlan.md b/docs/ByDivisorPredictiveScanPlan.md new file mode 100644 index 0000000..41cdb82 --- /dev/null +++ b/docs/ByDivisorPredictiveScanPlan.md @@ -0,0 +1,52 @@ +# By-Divisor Adaptive Scan Plan + +## Goals +- [x] Add deterministic class-based strategy selection for --mersenne=bydivisor (CPU/Hybrid first) using k~q/(2p) stats per 10-bit suffix. +- [x] Persist and reuse class statistics (k10 percentiles, gap probability, stability) to drive k-start/k-step decisions; rebuild from recorded divisors when stats files are missing or stale. +- [x] Skip trivial divisors and previously found divisors for a given p when rebuilding stats or scanning. +- [x] Introduce a classifier component that categorizes p -> class params and exposes recommended scan ranges/steps for k. +- [x] Integrate classifier into by-divisor CPU/hybrid scan path to choose aggressive k progression (cheap low-k sweep + geometric escalation to class targets). +- [x] Update persistence to append new findings and refresh the model for future runs. +- [x] Keep GPU path unchanged for now; guard hooks for future extension. + +## Background / Inputs +- Existing results file (p,divisor,searchedMersenne,detailedCheck,passedAllTests) already loaded in Program.cs (by-divisor flow requires --filter-p and reads prior results). +- Divisor form: q = 2*k*p + 1. We need per-class stats of q/p percentiles -> k percentiles (k10 ~= P10/2, k50 ~= median/2, k75 ~= P75/2) and gap probability (missing rows = gaps). +- Classes keyed by suffix s = p & 1023 (10 LSBs). Stability classification available when rolling envelope data exists; otherwise default to "unknown/assume stable". +- Aggressive strategy (Spec 2): + - Phase A: always scan cheap low k up to K_cheap. + - Phase B: if not found, jump to class start k_start = max(K_cheap+1, ceil(gamma50*k50)), escalate by rho1 until k_target = ceil(gamma75*k75), then rho2. + - Hard/gappy classes (gap_prob >= 0.75 or escaping) use more aggressive multipliers (gamma50=1.0, gamma75=1.4, rho2=1.8). Default gamma50=0.8, gamma75=1.2, rho1=1.35, rho2=1.6. + +## Planned Artifacts +- New docs/stats files to persist model (e.g., docs/stats/bydivisor_k10_model.json or CSV) with suffix->{k50,k75,gap_prob,stability,timestamps}. +- Runtime cache file near results/state directory (e.g., Checks/bydivisor-kmodel.bin or similar) for quick load/update. +- Helper to rebuild model from prior results file when model file missing; updates appended when new divisors found. + +## Implementation Tasks +1) Data model and persistence +- [x] Define DTO/structs for per-suffix stats (k50,k75,gap_prob,stability tier, sample counts, last-updated). +- [x] Choose on-disk format (JSON/CSV) stored under docs/stats and runtime writable copy under Checks/ (or alongside results). Implement load/save with versioning. +- [x] Add builder that reads existing results (p,q,passed) and computes percentiles + gap probability (skip trivial q, enforce q == 1 mod 2p and pow check when needed or reuse flags). +- [x] Handle missing model files by building from results; merge newly found divisors via rebuild after a run. + +2) Classifier +- [x] Implement classifier service (ByDivisorClassModel) returning class params for p: suffix lookup, fallbacks when stats absent. +- [x] Expose recommendation API for scanning: GetScanPlan(p) -> {kCheap, kStart, kTarget, rho1, rho2} with deterministic defaults. +- [x] Flag classes as hard/gappy using thresholds (gap_prob, stability) to pick aggressive multipliers. + +3) Integration with by-divisor CPU/Hybrid flow +- [x] Feed classifier with loaded model at Program start (by-divisor path) using the already loaded results file; share between CPU/Hybrid testers. +- [x] Extend divisor scan session to use Phase A/B plan: cheap linear sweep up to K_cheap, then geometric progression per class (k_start -> k_target with rho1, then rho2) instead of flat increments. +- [x] Ensure pow2-groups path remains unchanged; apply to sequential path first (CPU/Hybrid). Guard GPU path (no-op for now). +- [x] Skip trivial divisors and previously recorded divisors for p when scanning/rebuilding stats. +- [x] Persist updated model after processing primes (additions from new divisors found during run). + +4) Validation +- [ ] Add targeted unit/functional tests (xUnit + FluentAssertions) for classifier defaults, model rebuild from synthetic results, and scan plan generation for easy vs hard classes. +- [ ] Document how to run fast tests (Category=Fast) and keep runtime under 2 minutes. + +## Open Questions / Assumptions +- Results file path reused for stats rebuild; assumes divisor column contains the found q (non-zero only when divisor found). +- Stability classification file presence: if not shipped, treat as stable; allow future plug-in when rolling envelopes added. +- K_cheap initial guess TBD (needs perf check); provisional: small constant (e.g., 10_000) until benchmarked.