From 7cbc27d6c3dc8da48f9daedbe55a3d0e1f274bc3 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Marek=20Lema=C5=84czyk?= Date: Mon, 20 Oct 2025 18:36:58 +0200 Subject: [PATCH] Use adaptive pow2 stepping in residue stepper --- PerfectNumbers.Core/MersenneResidueStepper.cs | 23 ++++--------------- 1 file changed, 4 insertions(+), 19 deletions(-) diff --git a/PerfectNumbers.Core/MersenneResidueStepper.cs b/PerfectNumbers.Core/MersenneResidueStepper.cs index a2dbd2a..629b6cf 100644 --- a/PerfectNumbers.Core/MersenneResidueStepper.cs +++ b/PerfectNumbers.Core/MersenneResidueStepper.cs @@ -34,28 +34,12 @@ public class MersenneResidueStepper return State(); } - UInt128 g = UInt128.One; ulong absDelta = (ulong)Math.Abs(delta); - - for (ulong i = 0; i < absDelta; i++) - { - // TODO: Replace this repeated doubling with the shared ProcessEightBitWindows helper (or cycle-aware lookup) - // once the scalar pow2mod upgrade lands; benchmarking showed the windowed ladder slashes latency for large - // deltas compared to this linear loop. - // TODO: Once divisor cycles power the residue stepper, prefer stepping by cycle lengths instead of incrementing - // one bit at a time so large jumps stay on the mandatory cycle-accelerated path. - g <<= 1; - - if (g >= m) - { - g -= m; - } - } + UInt128 g = absDelta.PowMod(m); if (delta > 0) { - // TODO: Once the scalar windowed pow2 helper lands, reroute this `%` path through it to avoid the slower - // multiply-and-mod sequence measured in the legacy benchmarks. + // Reuse the cached 2^p residue and only multiply by the adaptive pow2 step computed above. pow2_mod = (pow2_mod * g) % m; } else @@ -63,7 +47,8 @@ public class MersenneResidueStepper // Modular inverse for negative delta // TODO: Switch to a cached divisor-cycle remainder instead of computing a fresh modular inverse for every // backward jump; the divisor cycle benchmarks demonstrated significant savings once shared cycle data was used. - pow2_mod = (pow2_mod * ModInverse(g, m)) % m; + UInt128 inverse = ModInverse(g, m); + pow2_mod = (pow2_mod * inverse) % m; } p = (ulong)((long)p + delta);