diff --git a/EvenPerfectBitScanner.Benchmarks/Mod10UInt128Benchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod10UInt128Benchmarks.cs index ff414ea..866965b 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod10UInt128Benchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod10UInt128Benchmarks.cs @@ -1,5 +1,4 @@ using BenchmarkDotNet.Attributes; -using PerfectNumbers.Core; namespace EvenPerfectBitScanner.Benchmarks; @@ -15,38 +14,45 @@ public class Mod10UInt128Benchmarks return Value % 10UL; } + /// + /// Fastest + /// [Benchmark] public ulong ModWithHighLowModulo() { - // Split and fold under mod 10: 2^64 ≡ 6 (mod 10) - ulong value = (ulong)Value; - ulong high = value >> 64; - ulong low = value; - ulong highRem = high % 10UL; - ulong lowRem = low % 10UL; - ulong combined = lowRem + highRem * 6UL; - return combined % 10UL; + return ModWithHighLowModuloHelper(Value); } - [Benchmark] - public ulong ModWithLoopModulo() - { - UInt128 value = Value; + private static ulong ModWithHighLowModuloHelper(UInt128 value128) + { + // Split and fold under mod 10: 2^64 ≡ 6 (mod 10) + ulong value = (ulong)value128; + return ((value % 10UL) + ((value >> 64) % 10UL) * 6UL) % 10UL; + } + + [Benchmark] + public ulong ModWithLoopModulo() + { + return ModWithLoopModuloHelper(Value); + } + + private static ulong ModWithLoopModuloHelper(UInt128 value) + { UInt128 zero = UInt128.Zero; if (value == zero) return 0UL; ulong result = (ulong)value; - UInt128 high = value >> 64; + value >>= 64; - while (high != zero) + while (value != zero) { // 2^64 ≡ 6 (mod 10) - result = (result + (ulong)high * 6UL) % 10UL; - high >>= 64; - } + result = (result + (ulong)value * 6UL) % 10UL; + value >>= 64; + } - return result % 10UL; - } + return result % 10UL; + } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod10UIntBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod10UIntBenchmarks.cs index ccfb1e9..14bf8c6 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod10UIntBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod10UIntBenchmarks.cs @@ -1,3 +1,4 @@ +using System.Runtime.CompilerServices; using BenchmarkDotNet.Attributes; using BenchmarkDotNet.Jobs; using PerfectNumbers.Core; @@ -8,10 +9,16 @@ namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod10UIntBenchmarks { + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static uint Mod10(uint value) => value - (uint)((value * UIntExtensions.Mod5Mask) >> 35) * 10U; + [Params(10U, 2047U, 65535U, 2147483647U)] public uint Value { get; set; } - [Benchmark(Baseline = true)] + /// + /// Fastest + /// + [Benchmark(Baseline = true)] public uint ModuloOperator() { return Value % 10U; @@ -20,7 +27,7 @@ public class Mod10UIntBenchmarks [Benchmark] public uint ExtensionMethod() { - return Value.Mod10(); + return Mod10(Value); } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod10ULongBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod10ULongBenchmarks.cs index 04e9cde..0b4a120 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod10ULongBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod10ULongBenchmarks.cs @@ -1,4 +1,3 @@ -using System.Numerics; using BenchmarkDotNet.Attributes; using PerfectNumbers.Core; @@ -16,7 +15,10 @@ public class Mod10ULongBenchmarks return Value % 10UL; } - [Benchmark] + /// + /// Fastest + /// + [Benchmark] public ulong ModMethodModulo() { return Value.Mod10(); diff --git a/EvenPerfectBitScanner.Benchmarks/Mod10_8_5_3Benchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod10_8_5_3Benchmarks.cs index 99aa38d..9f9154c 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod10_8_5_3Benchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod10_8_5_3Benchmarks.cs @@ -1,4 +1,3 @@ -using System.Numerics; using BenchmarkDotNet.Attributes; using PerfectNumbers.Core; @@ -10,7 +9,10 @@ public class Mod10_8_5_3Benchmarks [Params(3UL, 8191UL, 131071UL, 2147483647UL)] public ulong Value { get; set; } - [Benchmark(Baseline = true)] + /// + /// Fastest + /// + [Benchmark(Baseline = true)] public (ulong Mod10, ulong Mod8, ulong Mod3, ulong Mod5) ModuloOperator() { ulong value = Value; diff --git a/EvenPerfectBitScanner.Benchmarks/Mod10_8_5_3PropertiesBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod10_8_5_3PropertiesBenchmarks.cs index b3924e6..a3d0f2e 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod10_8_5_3PropertiesBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod10_8_5_3PropertiesBenchmarks.cs @@ -1,4 +1,3 @@ -using System.Numerics; using System.Runtime.CompilerServices; using BenchmarkDotNet.Attributes; using PerfectNumbers.Core; @@ -36,10 +35,13 @@ public class Mod10_8_5_3PropertiesBenchmarks return rem >= 5UL ? rem - 5UL : rem; } - [Params(3UL, 8191UL, 131071UL, 2147483647UL)] + [Params(3UL, 6UL, 128UL, 256UL, 2_047UL, 8_191UL, 65_535UL, 131_071UL, 2_147_483_647UL, ulong.MaxValue - 1_024UL, ulong.MaxValue - 511UL, ulong.MaxValue - 255UL, ulong.MaxValue - 127UL, ulong.MaxValue - 63UL, ulong.MaxValue - 31UL, ulong.MaxValue - 17UL)] public ulong Value { get; set; } - [Benchmark(Baseline = true)] + /// + /// Fastest + /// + [Benchmark(Baseline = true)] public (ulong Mod10, ulong Mod8, ulong Mod5, ulong Mod3) LegacyModulo() { UInt128 value = Value; @@ -53,12 +55,7 @@ public class Mod10_8_5_3PropertiesBenchmarks [Benchmark] public (ulong Mod10, ulong Mod8, ulong Mod3, ulong Mod5) ModMethodModulo() { - Value.Mod10_8_5_3(out var mod10, out ulong mod8, out ulong mod5, out ulong mod3); - Mod10R = mod10; - Mod8R = mod8; - Mod5R = mod5; - Mod3R = mod3; - + Value.Mod10_8_5_3(out Mod10R, out Mod8R, out Mod5R, out Mod3R); return (Mod10R, Mod8R, Mod3R, Mod5R); } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod10_8_5_3StepsBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod10_8_5_3StepsBenchmarks.cs index da08b13..7b01796 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod10_8_5_3StepsBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod10_8_5_3StepsBenchmarks.cs @@ -1,4 +1,3 @@ -using System.Numerics; using BenchmarkDotNet.Attributes; using PerfectNumbers.Core; @@ -14,19 +13,21 @@ public class Mod10_8_5_3StepsBenchmarks public (ulong Mod10, ulong Mod8, ulong Mod3, ulong Mod5) LegacyModulo() { UInt128 value = Value; - ulong mod10 = (ulong)(value % 10UL) % 10UL; - ulong mod8 = (ulong)(value % 8UL); - ulong mod3 = (ulong)(value % 3UL); - ulong mod5 = (ulong)(value % 5UL); + ulong mod10 = (((ulong)(value % 10UL)) << 1) % 10UL; + ulong mod8 = (((ulong)(value & 7UL)) << 1) & 7UL; + ulong mod5 = (((ulong)(value % 5UL)) << 1) % 5UL; + ulong mod3 = (((ulong)(value % 3UL)) << 1) % 3UL; return (mod10, mod8, mod3, mod5); } - [Benchmark] + /// + /// Fastest + /// + [Benchmark] public (ulong Mod10, ulong Mod8, ulong Mod3, ulong Mod5) ModMethodModulo() { - Value.Mod10_8_5_3(out ulong mod10, out ulong mod8, out ulong mod5, out ulong mod3); + Value.Mod10_8_5_3Steps(out ulong mod10, out ulong mod8, out ulong mod5, out ulong mod3); // Step residues for q += 2*p - mod10 %= 10UL; return (mod10, mod8, mod3, mod5); } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod11UIntBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod11UIntBenchmarks.cs index f16e36b..f34b2c9 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod11UIntBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod11UIntBenchmarks.cs @@ -1,6 +1,6 @@ +using System.Runtime.CompilerServices; using BenchmarkDotNet.Attributes; using BenchmarkDotNet.Jobs; -using PerfectNumbers.Core; namespace EvenPerfectBitScanner.Benchmarks; @@ -8,10 +8,27 @@ namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod11UIntBenchmarks { + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static uint Mod11(uint value) + { + uint remainder = 0U; + while (value != 0U) + { + remainder += value & 1023U; + value >>= 10; + remainder -= 11U * (remainder / 11U); + } + + return remainder; + } + [Params(11U, 2047U, 65535U, 2147483647U)] public uint Value { get; set; } - [Benchmark(Baseline = true)] + /// + /// Fastest + /// + [Benchmark(Baseline = true)] public uint ModuloOperator() { return Value % 11U; @@ -20,7 +37,7 @@ public class Mod11UIntBenchmarks [Benchmark] public uint ExtensionMethod() { - return Value.Mod11(); + return Mod11(Value); } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod11ULongBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod11ULongBenchmarks.cs index 1c450ce..93f43bc 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod11ULongBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod11ULongBenchmarks.cs @@ -1,6 +1,6 @@ +using System.Runtime.CompilerServices; using BenchmarkDotNet.Attributes; using BenchmarkDotNet.Jobs; -using PerfectNumbers.Core; namespace EvenPerfectBitScanner.Benchmarks; @@ -8,10 +8,25 @@ namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod11ULongBenchmarks { + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong Mod11(ulong value) + { + ulong remainder = ((uint)value % 11U) + ((uint)(value >> 32) % 11U) << 2; + while (remainder >= 11UL) + { + remainder -= 11UL; + } + + return remainder; + } + [Params(11UL, 8191UL, 131071UL, 2147483647UL, ulong.MaxValue - 127UL)] public ulong Value { get; set; } - [Benchmark(Baseline = true)] + /// + /// Fastest + /// + [Benchmark(Baseline = true)] public ulong ModuloOperator() { return Value % 11UL; @@ -20,7 +35,7 @@ public class Mod11ULongBenchmarks [Benchmark] public ulong ExtensionMethod() { - return Value.Mod11(); + return Mod11(Value); } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod11_7_5_3Benchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod11_7_5_3Benchmarks.cs index 599d1fa..c5fed35 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod11_7_5_3Benchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod11_7_5_3Benchmarks.cs @@ -1,32 +1,92 @@ +using System.Runtime.CompilerServices; using BenchmarkDotNet.Attributes; -using PerfectNumbers.Core; namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod11_7_5_3Benchmarks { + private static readonly byte[] Mod7ByteCoefficients = [1, 4, 2, 1, 4, 2, 1, 4]; + private static readonly byte[] Mod11ByteCoefficients = [1, 3, 9, 5, 4, 1, 3, 9]; + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private static ulong ReduceMod7(uint value) + { + while (value >= 7U) + { + value = (value >> 3) + (value & 7U); + if (value >= 7U && value < 14U) + { + value -= 7U; + } + } + + return value; + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private static ulong ReduceMod11(uint value) + { + while (value >= 11U) + { + value = (value & 15U) + ((value >> 4) * 5U); + if (value >= 11U && value < 22U) + { + value -= 11U; + } + } + + return value; + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private static void Mod11_7_5_3(ulong value, out ulong mod11, out ulong mod7, out ulong mod5, out ulong mod3) + { + uint mod3And5Accumulator = 0U; + uint mod7Accumulator = 0U; + uint mod11Accumulator = 0U; + + int index = 0; + byte current; + ulong temp = value; + + byte[] mod7ByteCoefficients = Mod7ByteCoefficients; + byte[] mod11ByteCoefficients = Mod11ByteCoefficients; + do + { + current = (byte)temp; + mod3And5Accumulator += current; + mod7Accumulator += (uint)(current * mod7ByteCoefficients[index]); + mod11Accumulator += (uint)(current * mod11ByteCoefficients[index]); + + temp >>= 8; + index++; + } + while (temp != 0UL); + + mod3 = mod3And5Accumulator % 3; + mod5 = mod3And5Accumulator % 5; + mod7 = (uint)ReduceMod7(mod7Accumulator); + mod11 = (uint)ReduceMod11(mod11Accumulator); + } + [Params(3UL, 8191UL, 131071UL, 2147483647UL)] public ulong Value { get; set; } - [Benchmark(Baseline = true)] + /// + /// Fastest + /// + [Benchmark(Baseline = true)] public (ulong Mod11, ulong Mod7, ulong Mod5, ulong Mod3) ModuloOperator() { ulong value = Value; return (value % 11UL, value % 7UL, value % 5UL, value % 3UL); } - [Benchmark] - public (ulong Mod11, ulong Mod7, ulong Mod5, ulong Mod3) ExtensionMethods() - { - ulong value = Value; - return (value.Mod11(), value.Mod7(), value.Mod5(), value.Mod3()); - } - [Benchmark] public (ulong Mod11, ulong Mod7, ulong Mod5, ulong Mod3) CombinedMethod() { - Value.Mod11_7_5_3(out ulong mod11, out ulong mod7, out ulong mod5, out ulong mod3); + Mod11_7_5_3(Value, out ulong mod11, out ulong mod7, out ulong mod5, out ulong mod3); return (mod11, mod7, mod5, mod3); } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod128UInt128Benchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod128UInt128Benchmarks.cs index f50e7d0..1b805bd 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod128UInt128Benchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod128UInt128Benchmarks.cs @@ -1,4 +1,3 @@ -using System.Numerics; using BenchmarkDotNet.Attributes; using BenchmarkDotNet.Jobs; using PerfectNumbers.Core; @@ -12,12 +11,16 @@ public class Mod128UInt128Benchmarks [Params(128UL, 8191UL, 131071UL, 2147483647UL, ulong.MaxValue - 511UL)] public UInt128 Value { get; set; } - [Benchmark(Baseline = true)] + [Benchmark(Baseline = true)] public ulong ModuloOperator() { return (ulong)(Value % 128UL); } + /// + /// Fastest + /// + /// [Benchmark] public ulong ExtensionMethod() { diff --git a/EvenPerfectBitScanner.Benchmarks/Mod128UIntBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod128UIntBenchmarks.cs index a98d7da..5bef7ce 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod128UIntBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod128UIntBenchmarks.cs @@ -1,6 +1,6 @@ +using System.Runtime.CompilerServices; using BenchmarkDotNet.Attributes; using BenchmarkDotNet.Jobs; -using PerfectNumbers.Core; namespace EvenPerfectBitScanner.Benchmarks; @@ -8,10 +8,17 @@ namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod128UIntBenchmarks { - [Params(128U, 2047U, 65535U, 2147483647U)] + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static uint Mod128(uint value) => value & 127U; + + [Params(3UL, 5UL, 6UL, 8UL, 10UL, 11UL, 32UL, 64UL, 128UL, 256UL, 2_047UL, 8_191UL, 65_535UL, 131_071UL, 2_147_483_647UL, uint.MaxValue - 1_024UL, uint.MaxValue - 511UL, uint.MaxValue - 255UL, uint.MaxValue - 127UL, uint.MaxValue - 63UL, uint.MaxValue - 31UL, uint.MaxValue - 17UL, uint.MaxValue)] public uint Value { get; set; } - [Benchmark(Baseline = true)] + /// + /// Fastest most of the times. It's usually much faster than it's slower. + /// + /// + [Benchmark(Baseline = true)] public uint ModuloOperator() { return Value % 128U; @@ -20,7 +27,7 @@ public class Mod128UIntBenchmarks [Benchmark] public uint ExtensionMethod() { - return Value.Mod128(); + return Mod128(Value); } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod128ULongBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod128ULongBenchmarks.cs index 414245b..9eeb806 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod128ULongBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod128ULongBenchmarks.cs @@ -8,15 +8,19 @@ namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod128ULongBenchmarks { - [Params(128UL, 8191UL, 131071UL, 2147483647UL, ulong.MaxValue - 255UL)] + [Params(3UL, 5UL, 6UL, 8UL, 10UL, 11UL, 32UL, 64UL, 128UL, 256UL, 2_047UL, 8_191UL, 65_535UL, 131_071UL, 2_147_483_647UL, uint.MaxValue - 1_024UL, uint.MaxValue - 511UL, uint.MaxValue - 255UL, uint.MaxValue - 127UL, uint.MaxValue - 63UL, uint.MaxValue - 31UL, uint.MaxValue - 17UL, uint.MaxValue)] public ulong Value { get; set; } - [Benchmark(Baseline = true)] + [Benchmark(Baseline = true)] public ulong ModuloOperator() { return Value % 128UL; } + /// + /// Usually fastest for high values. It's slower for some values and faster for other. + /// + /// [Benchmark] public ulong ExtensionMethod() { diff --git a/EvenPerfectBitScanner.Benchmarks/Mod3UIntBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod3UIntBenchmarks.cs index bd81c8d..dca806d 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod3UIntBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod3UIntBenchmarks.cs @@ -20,7 +20,7 @@ public class Mod3UIntBenchmarks [Benchmark] public uint ExtensionMethod() { - return Value.Mod3(); + return UIntExtensions.Mod3(Value); } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod3ULongBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod3ULongBenchmarks.cs index b96b41d..3535545 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod3ULongBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod3ULongBenchmarks.cs @@ -1,3 +1,4 @@ +using System.Runtime.CompilerServices; using BenchmarkDotNet.Attributes; using BenchmarkDotNet.Jobs; using PerfectNumbers.Core; @@ -8,10 +9,22 @@ namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod3ULongBenchmarks { + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private static ulong Mod3(ulong value) => + ( + (uint)(value & ULongExtensions.WordBitMask) + + (uint)((value >> 16) & ULongExtensions.WordBitMask) + + (uint)((value >> 32) & ULongExtensions.WordBitMask) + + (uint)((value >> 48) & ULongExtensions.WordBitMask) + ) % 3; + [Params(3UL, 8191UL, 131071UL, 2147483647UL, ulong.MaxValue - 17UL)] public ulong Value { get; set; } - [Benchmark(Baseline = true)] + /// + /// Fastest + /// + [Benchmark(Baseline = true)] public ulong ModuloOperator() { return Value % 3UL; @@ -20,7 +33,7 @@ public class Mod3ULongBenchmarks [Benchmark] public ulong ExtensionMethod() { - return Value.Mod3(); + return Mod3(Value); } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod5UIntBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod5UIntBenchmarks.cs index 758ba40..fe334fc 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod5UIntBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod5UIntBenchmarks.cs @@ -8,7 +8,7 @@ namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod5UIntBenchmarks { - [Params(5U, 2047U, 65535U, 2147483647U)] + [Params(3UL, 5UL, 6UL, 8UL, 10UL, 11UL, 32UL, 64UL, 128UL, 256UL, 2_047UL, 8_191UL, 65_535UL, 131_071UL, 2_147_483_647UL, uint.MaxValue - 1_024UL, uint.MaxValue - 511UL, uint.MaxValue - 255UL, uint.MaxValue - 127UL, uint.MaxValue - 63UL, uint.MaxValue - 31UL, uint.MaxValue - 17UL, uint.MaxValue)] public uint Value { get; set; } [Benchmark(Baseline = true)] @@ -17,7 +17,10 @@ public class Mod5UIntBenchmarks return Value % 5U; } - [Benchmark] + /// + /// Mixed results. Sometimes it's faster, in higher ranges, sometimes it's slower. + /// + [Benchmark] public uint ExtensionMethod() { return Value.Mod5(); diff --git a/EvenPerfectBitScanner.Benchmarks/Mod5ULongBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod5ULongBenchmarks.cs index d156306..0b49a43 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod5ULongBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod5ULongBenchmarks.cs @@ -1,3 +1,4 @@ +using System.Runtime.CompilerServices; using BenchmarkDotNet.Attributes; using BenchmarkDotNet.Jobs; using PerfectNumbers.Core; @@ -8,10 +9,22 @@ namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod5ULongBenchmarks { + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong Mod5(ulong value) => + ( + (uint)(value & ULongExtensions.WordBitMask) + + (uint)((value >> 16) & ULongExtensions.WordBitMask) + + (uint)((value >> 32) & ULongExtensions.WordBitMask) + + (uint)((value >> 48) & ULongExtensions.WordBitMask) + ) % 5; + [Params(5UL, 8191UL, 131071UL, 2147483647UL, ulong.MaxValue - 31UL)] public ulong Value { get; set; } - [Benchmark(Baseline = true)] + /// + /// Fastest + /// + [Benchmark(Baseline = true)] public ulong ModuloOperator() { return Value % 5UL; @@ -20,7 +33,7 @@ public class Mod5ULongBenchmarks [Benchmark] public ulong ExtensionMethod() { - return Value.Mod5(); + return Mod5(Value); } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod6ComparisonBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod6ComparisonBenchmarks.cs index fc9ea39..7f00b43 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod6ComparisonBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod6ComparisonBenchmarks.cs @@ -1,4 +1,4 @@ -using System.Numerics; +using System.Runtime.CompilerServices; using BenchmarkDotNet.Attributes; using BenchmarkDotNet.Jobs; using PerfectNumbers.Core; @@ -9,6 +9,11 @@ namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod6ComparisonBenchmarks { + private static readonly byte[] Mod6Lookup = [0, 3, 4, 1, 2, 5]; + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong Mod6(ulong value) => Mod6Lookup[(int)(((value % 3UL) << 1) | (value & 1UL))]; + private const int Iterations = 10_000_000; private const ulong Divisor = 6UL; private static readonly ulong FastDivMul = (ulong)(((UInt128)1 << 64) / Divisor); @@ -71,7 +76,7 @@ public class Mod6ComparisonBenchmarks for (int i = 0; i < Iterations; i++) { - checksum += value.Mod6(); + checksum += Mod6(value); value += (ulong)((i & 15) + 1); } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod6UIntBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod6UIntBenchmarks.cs index 765f8e1..b0ac9c6 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod6UIntBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod6UIntBenchmarks.cs @@ -1,6 +1,6 @@ +using System.Runtime.CompilerServices; using BenchmarkDotNet.Attributes; using BenchmarkDotNet.Jobs; -using PerfectNumbers.Core; namespace EvenPerfectBitScanner.Benchmarks; @@ -8,10 +8,18 @@ namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod6UIntBenchmarks { + private static readonly byte[] Mod6Lookup = [0, 3, 4, 1, 2, 5]; + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static uint Mod6(uint value) => Mod6Lookup[(int)(((value % 3U) << 1) | (value & 1U))]; + [Params(6U, 2047U, 65535U, 2147483647U)] public uint Value { get; set; } - [Benchmark(Baseline = true)] + /// + /// Fastest + /// + [Benchmark(Baseline = true)] public uint ModuloOperator() { return Value % 6U; @@ -20,7 +28,7 @@ public class Mod6UIntBenchmarks [Benchmark] public uint ExtensionMethod() { - return Value.Mod6(); + return Mod6(Value); } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod6ULongBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod6ULongBenchmarks.cs index b8af8e7..257b85a 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod6ULongBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod6ULongBenchmarks.cs @@ -1,4 +1,4 @@ -using System.Numerics; +using System.Runtime.CompilerServices; using BenchmarkDotNet.Attributes; using BenchmarkDotNet.Jobs; using PerfectNumbers.Core; @@ -9,13 +9,20 @@ namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod6ULongBenchmarks { + private static readonly byte[] Mod6Lookup = [0, 3, 4, 1, 2, 5]; private static readonly ulong Divisor = 6UL; private static readonly ulong FastDivMul = (ulong)(((UInt128)1 << 64) / Divisor); + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong Mod6(ulong value) => Mod6Lookup[(int)(((value % 3UL) << 1) | (value & 1UL))]; + [Params(3UL, 131071UL, ulong.MaxValue)] public ulong Value { get; set; } - [Benchmark(Baseline = true)] + /// + /// Fastest + /// + [Benchmark(Baseline = true)] public ulong ModuloOperator() { return Value % Divisor; @@ -38,6 +45,6 @@ public class Mod6ULongBenchmarks [Benchmark] public ulong ExtensionMethod() { - return Value.Mod6(); + return Mod6(Value); } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod7UInt128Benchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod7UInt128Benchmarks.cs index 7dacd31..0f50d3d 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod7UInt128Benchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod7UInt128Benchmarks.cs @@ -17,7 +17,11 @@ public class Mod7UInt128Benchmarks return (ulong)(Value % 7UL); } - [Benchmark] + /// + /// Fastest + /// + /// + [Benchmark] public ulong ExtensionMethod() { return Value.Mod7(); diff --git a/EvenPerfectBitScanner.Benchmarks/Mod7UIntBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod7UIntBenchmarks.cs index 4d8faac..d7202cb 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod7UIntBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod7UIntBenchmarks.cs @@ -1,6 +1,6 @@ +using System.Runtime.CompilerServices; using BenchmarkDotNet.Attributes; using BenchmarkDotNet.Jobs; -using PerfectNumbers.Core; namespace EvenPerfectBitScanner.Benchmarks; @@ -8,10 +8,30 @@ namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod7UIntBenchmarks { + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static uint Mod7(uint value) + { + uint remainder = 0U; + while (value != 0U) + { + remainder += value & 7U; + value >>= 3; + if (remainder >= 7U) + { + remainder -= 7U; + } + } + + return remainder >= 7U ? remainder - 7U : remainder; + } + [Params(7U, 2047U, 65535U, 2147483647U)] public uint Value { get; set; } - [Benchmark(Baseline = true)] + /// + /// Fastest + /// + [Benchmark(Baseline = true)] public uint ModuloOperator() { return Value % 7U; @@ -20,7 +40,7 @@ public class Mod7UIntBenchmarks [Benchmark] public uint ExtensionMethod() { - return Value.Mod7(); + return Mod7(Value); } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod7ULongBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod7ULongBenchmarks.cs index b59bdd4..92d1c14 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod7ULongBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod7ULongBenchmarks.cs @@ -1,6 +1,6 @@ +using System.Runtime.CompilerServices; using BenchmarkDotNet.Attributes; using BenchmarkDotNet.Jobs; -using PerfectNumbers.Core; namespace EvenPerfectBitScanner.Benchmarks; @@ -8,10 +8,25 @@ namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod7ULongBenchmarks { + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong Mod7(ulong value) + { + ulong remainder = ((uint)value % 7U) + (((uint)(value >> 32) % 7U) << 2); + while (remainder >= 7UL) + { + remainder -= 7UL; + } + + return remainder; + } + [Params(7UL, 8191UL, 131071UL, 2147483647UL, ulong.MaxValue - 63UL)] public ulong Value { get; set; } - [Benchmark(Baseline = true)] + /// + /// Fastest + /// + [Benchmark(Baseline = true)] public ulong ModuloOperator() { return Value % 7UL; @@ -20,7 +35,7 @@ public class Mod7ULongBenchmarks [Benchmark] public ulong ExtensionMethod() { - return Value.Mod7(); + return Mod7(Value); } } diff --git a/EvenPerfectBitScanner.Benchmarks/Mod8UInt128Benchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod8UInt128Benchmarks.cs index ee2c7cf..3d6dec0 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod8UInt128Benchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod8UInt128Benchmarks.cs @@ -1,4 +1,3 @@ -using System.Numerics; using BenchmarkDotNet.Attributes; using BenchmarkDotNet.Jobs; using PerfectNumbers.Core; @@ -18,7 +17,10 @@ public class Mod8UInt128Benchmarks return (ulong)(Value % 8UL); } - [Benchmark] + /// + /// Fastest + /// + [Benchmark] public ulong ExtensionMethod() { return Value.Mod8(); diff --git a/EvenPerfectBitScanner.Benchmarks/Mod8UIntBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod8UIntBenchmarks.cs index 72e0e4a..189ac59 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod8UIntBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod8UIntBenchmarks.cs @@ -17,6 +17,9 @@ public class Mod8UIntBenchmarks return Value % 8U; } + /// + /// Fastest + /// [Benchmark] public uint ExtensionMethod() { diff --git a/EvenPerfectBitScanner.Benchmarks/Mod8ULongBenchmarks.cs b/EvenPerfectBitScanner.Benchmarks/Mod8ULongBenchmarks.cs index 6af8ff5..17e849f 100644 --- a/EvenPerfectBitScanner.Benchmarks/Mod8ULongBenchmarks.cs +++ b/EvenPerfectBitScanner.Benchmarks/Mod8ULongBenchmarks.cs @@ -1,6 +1,6 @@ +using System.Runtime.CompilerServices; using BenchmarkDotNet.Attributes; using BenchmarkDotNet.Jobs; -using PerfectNumbers.Core; namespace EvenPerfectBitScanner.Benchmarks; @@ -8,19 +8,25 @@ namespace EvenPerfectBitScanner.Benchmarks; [MemoryDiagnoser] public class Mod8ULongBenchmarks { + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong Mod8(ulong value) => value & 7UL; + [Params(8UL, 8191UL, 131071UL, 2147483647UL, ulong.MaxValue - 127UL)] public ulong Value { get; set; } + /// + /// Fastest + /// [Benchmark(Baseline = true)] public ulong ModuloOperator() { return Value % 8UL; } - [Benchmark] + [Benchmark] public ulong ExtensionMethod() { - return Value.Mod8(); + return Mod8(Value); } } diff --git a/EvenPerfectBitScanner/Program.cs b/EvenPerfectBitScanner/Program.cs index 2f3cf7d..e08bb18 100644 --- a/EvenPerfectBitScanner/Program.cs +++ b/EvenPerfectBitScanner/Program.cs @@ -3,7 +3,6 @@ using System.Globalization; using System.Numerics; using System.Runtime.CompilerServices; using System.Text; -using System.Threading; using Open.Collections; using Open.Numeric.Primes; using PerfectNumbers.Core; @@ -1674,7 +1673,7 @@ internal static class Program // Optional: RLE blacklist and binary-threshold filters on p (safe only when configured) if (p <= _rleHardMaxP) { - if (!_rleOnlyLast7 || p % 10UL == 7UL) + if (!_rleOnlyLast7 || p.Mod10() == 7UL) { if (RleBlacklist.IsLoaded() && RleBlacklist.Matches(p)) { diff --git a/PerfectNumbers.Core.Tests/UIntExtensionsTests.cs b/PerfectNumbers.Core.Tests/UIntExtensionsTests.cs index 7a221a1..5f2db93 100644 --- a/PerfectNumbers.Core.Tests/UIntExtensionsTests.cs +++ b/PerfectNumbers.Core.Tests/UIntExtensionsTests.cs @@ -13,7 +13,7 @@ public class UIntExtensionsTests [InlineData(uint.MaxValue)] public void Mod3_matches_operator(uint value) { - value.Mod3().Should().Be(value % 3U); + UIntExtensions.Mod3(value).Should().Be(value % 3U); } [Theory] @@ -27,28 +27,6 @@ public class UIntExtensionsTests value.Mod5().Should().Be(value % 5U); } - [Theory] - [Trait("Category", "Fast")] - [InlineData(0U)] - [InlineData(1U)] - [InlineData(1234567890U)] - [InlineData(uint.MaxValue)] - public void Mod6_matches_operator(uint value) - { - value.Mod6().Should().Be(value % 6U); - } - - [Theory] - [Trait("Category", "Fast")] - [InlineData(0U)] - [InlineData(1U)] - [InlineData(1234567890U)] - [InlineData(uint.MaxValue)] - public void Mod7_matches_operator(uint value) - { - value.Mod7().Should().Be(value % 7U); - } - [Theory] [Trait("Category", "Fast")] [InlineData(0U)] @@ -59,38 +37,5 @@ public class UIntExtensionsTests { value.Mod8().Should().Be(value % 8U); } - - [Theory] - [Trait("Category", "Fast")] - [InlineData(0U)] - [InlineData(1U)] - [InlineData(1234567890U)] - [InlineData(uint.MaxValue)] - public void Mod10_matches_operator(uint value) - { - value.Mod10().Should().Be(value % 10U); - } - - [Theory] - [Trait("Category", "Fast")] - [InlineData(0U)] - [InlineData(1U)] - [InlineData(1234567890U)] - [InlineData(uint.MaxValue)] - public void Mod11_matches_operator(uint value) - { - value.Mod11().Should().Be(value % 11U); - } - - [Theory] - [Trait("Category", "Fast")] - [InlineData(0U)] - [InlineData(1U)] - [InlineData(1234567890U)] - [InlineData(uint.MaxValue)] - public void Mod128_matches_operator(uint value) - { - value.Mod128().Should().Be(value % 128U); - } } diff --git a/PerfectNumbers.Core.Tests/ULongExtensionsTests.cs b/PerfectNumbers.Core.Tests/ULongExtensionsTests.cs index 7bf7e60..b192dab 100644 --- a/PerfectNumbers.Core.Tests/ULongExtensionsTests.cs +++ b/PerfectNumbers.Core.Tests/ULongExtensionsTests.cs @@ -21,66 +21,6 @@ public class ULongExtensionsTests ULongExtensions.Mod10(value).Should().Be(expected); } - [Theory] - [Trait("Category", "Fast")] - [InlineData(0UL)] - [InlineData(ulong.MaxValue)] - [InlineData(1234567890123456789UL)] - public void Mod3_matches_operator(ulong value) - { - value.Mod3().Should().Be(value % 3UL); - } - - [Theory] - [Trait("Category", "Fast")] - [InlineData(0UL)] - [InlineData(ulong.MaxValue)] - [InlineData(1234567890123456789UL)] - public void Mod5_matches_operator(ulong value) - { - value.Mod5().Should().Be(value % 5UL); - } - - [Theory] - [Trait("Category", "Fast")] - [InlineData(0UL)] - [InlineData(ulong.MaxValue)] - [InlineData(1234567890123456789UL)] - public void Mod6_matches_operator(ulong value) - { - value.Mod6().Should().Be(value % 6UL); - } - - [Theory] - [Trait("Category", "Fast")] - [InlineData(0UL)] - [InlineData(ulong.MaxValue)] - [InlineData(1234567890123456789UL)] - public void Mod7_matches_operator(ulong value) - { - value.Mod7().Should().Be(value % 7UL); - } - - [Theory] - [Trait("Category", "Fast")] - [InlineData(0UL)] - [InlineData(ulong.MaxValue)] - [InlineData(1234567890123456789UL)] - public void Mod8_matches_operator(ulong value) - { - value.Mod8().Should().Be(value % 8UL); - } - - [Theory] - [Trait("Category", "Fast")] - [InlineData(0UL)] - [InlineData(ulong.MaxValue)] - [InlineData(1234567890123456789UL)] - public void Mod11_matches_operator(ulong value) - { - value.Mod11().Should().Be(value % 11UL); - } - [Theory] [Trait("Category", "Fast")] [InlineData(0UL)] @@ -111,22 +51,6 @@ public class ULongExtensionsTests step3.Should().Be(expected3); } - [Theory] - [Trait("Category", "Fast")] - [InlineData(0UL)] - [InlineData(1UL)] - [InlineData(1234567890123456789UL)] - [InlineData(ulong.MaxValue)] - public void Mod11_7_5_3_matches_operator(ulong value) - { - value.Mod11_7_5_3(out ulong mod11, out ulong mod7, out ulong mod5, out ulong mod3); - - mod11.Should().Be(value % 11UL); - mod7.Should().Be(value % 7UL); - mod5.Should().Be(value % 5UL); - mod3.Should().Be(value % 3UL); - } - [Theory] [Trait("Category", "Fast")] [InlineData(2UL, true)] diff --git a/PerfectNumbers.Core/Cpu/MersenneNumberLucasLehmerCpuTester.cs b/PerfectNumbers.Core/Cpu/MersenneNumberLucasLehmerCpuTester.cs index 93dddb5..93aafb1 100644 --- a/PerfectNumbers.Core/Cpu/MersenneNumberLucasLehmerCpuTester.cs +++ b/PerfectNumbers.Core/Cpu/MersenneNumberLucasLehmerCpuTester.cs @@ -9,10 +9,10 @@ public class MersenneNumberLucasLehmerCpuTester // - If p ≡ 1 (mod 4) and p shares a factor with (p-1), reject fast. // - If p is divisible by divisors specific to numbers ending with 1 or 7, reject fast. if ( - (exponent.Mod3() == 0UL && exponent != 3UL) || - (exponent.Mod5() == 0UL && exponent != 5UL) || - (exponent.Mod7() == 0UL && exponent != 7UL) || - (exponent.Mod11() == 0UL && exponent != 11UL) + ((exponent % 3UL) == 0UL && exponent != 3UL) || + ((exponent % 5UL) == 0UL && exponent != 5UL) || + ((exponent % 7UL) == 0UL && exponent != 7UL) || + ((exponent % 11UL) == 0UL && exponent != 11UL) ) { return false; diff --git a/PerfectNumbers.Core/Gpu/GpuKernelPool.cs b/PerfectNumbers.Core/Gpu/GpuKernelPool.cs index 81a3cb4..9ee67ef 100644 --- a/PerfectNumbers.Core/Gpu/GpuKernelPool.cs +++ b/PerfectNumbers.Core/Gpu/GpuKernelPool.cs @@ -1,6 +1,5 @@ using System.Collections.Concurrent; using System.Runtime.CompilerServices; -using System.Threading; using ILGPU; using ILGPU.Runtime; using static PerfectNumbers.Core.Gpu.GpuContextPool; @@ -200,8 +199,8 @@ public ref struct GpuKernelLease(IDisposable limiter, GpuContextLease gpu, Kerne ArrayView1D smallCycles) { ulong idx = (ulong)index.X; - ulong idxMod3 = idx.Mod3(); - ulong idxMod5 = idx.Mod5(); + ulong idxMod3 = idx % 3UL; + ulong idxMod5 = idx % 5UL; // residue automaton ulong step10 = (exponent.Mod10() << 1).Mod10(); ulong r10 = q0m10 + (step10 * idx).Mod10(); @@ -217,7 +216,7 @@ public ref struct GpuKernelLease(IDisposable limiter, GpuContextLease gpu, Kerne } else { - ulong step3 = (exponent.Mod3() << 1).Mod3(); + ulong step3 = ((exponent % 3UL) << 1) % 3UL; ulong r3 = q0m3 + (step3 * idxMod3); r3 -= (r3 >= 3UL) ? 3UL : 0UL; if (r3 == 0UL) { @@ -225,7 +224,7 @@ public ref struct GpuKernelLease(IDisposable limiter, GpuContextLease gpu, Kerne } else { - ulong step5Loc = (exponent.Mod5() << 1).Mod5(); + ulong step5Loc = ((exponent % 5UL) << 1) % 5UL; ulong r5 = q0m5 + (step5Loc * idxMod5); if (r5 >= 5UL) r5 -= 5UL; if (r5 == 0UL) { @@ -360,8 +359,8 @@ public ref struct GpuKernelLease(IDisposable limiter, GpuContextLease gpu, Kerne ArrayView1D smallPrimesPow2LastSeven) { ulong idx = (ulong)index.X; - ulong idxMod3 = idx.Mod3(); - ulong idxMod5 = idx.Mod5(); + ulong idxMod3 = idx % 3UL; + ulong idxMod5 = idx % 5UL; ulong r10 = ra.Q0M10 + (ra.Step10 * idx).Mod10(); r10 -= (r10 >= 10UL) ? 10UL : 0UL; bool shouldCheck = r10 != 5UL; if (shouldCheck) @@ -493,8 +492,8 @@ public ref struct GpuKernelLease(IDisposable limiter, GpuContextLease gpu, Kerne ResidueAutomatonArgs ra, ArrayView found, ArrayView1D smallCycles) { ulong idx = (ulong)index.X; - ulong idxMod3 = idx.Mod3(); - ulong idxMod5 = idx.Mod5(); + ulong idxMod3 = idx % 3UL; + ulong idxMod5 = idx % 5UL; ulong r10 = ra.Q0M10 + (ra.Step10 * idx).Mod10(); r10 -= (r10 >= 10UL) ? 10UL : 0UL; bool shouldCheck = r10 != 5UL; if (shouldCheck) @@ -579,8 +578,8 @@ public ref struct GpuKernelLease(IDisposable limiter, GpuContextLease gpu, Kerne ResidueAutomatonArgs ra, ArrayView found, ArrayView1D smallCycles) { ulong idx = (ulong)index.X; - ulong idxMod3 = idx.Mod3(); - ulong idxMod5 = idx.Mod5(); + ulong idxMod3 = idx % 3UL; + ulong idxMod5 = idx % 5UL; ulong r10 = ra.Q0M10 + (ra.Step10 * idx).Mod10(); r10 -= (r10 >= 10UL) ? 10UL : 0UL; bool shouldCheck = r10 != 5UL; if (shouldCheck) diff --git a/PerfectNumbers.Core/Gpu/MersenneNumberIncrementalGpuTester.cs b/PerfectNumbers.Core/Gpu/MersenneNumberIncrementalGpuTester.cs index 076f6f8..6eca1a6 100644 --- a/PerfectNumbers.Core/Gpu/MersenneNumberIncrementalGpuTester.cs +++ b/PerfectNumbers.Core/Gpu/MersenneNumberIncrementalGpuTester.cs @@ -1,6 +1,5 @@ using System.Buffers; using System.Runtime.InteropServices; -using ILGPU; using ILGPU.Runtime; namespace PerfectNumbers.Core.Gpu; @@ -10,98 +9,95 @@ public class MersenneNumberIncrementalGpuTester(GpuKernelType kernelType, bool u private readonly GpuKernelType _kernelType = kernelType; private readonly bool _useGpuOrder = useGpuOrder; - public void Scan(ulong exponent, UInt128 twoP, bool lastIsSeven, UInt128 maxK, ref bool isPrime) - { - var gpuLease = GpuKernelPool.GetKernel(_useGpuOrder); - var execution = gpuLease.EnterExecutionScope(); - var accelerator = gpuLease.Accelerator; - var stream = gpuLease.Stream; - int batchSize = GpuConstants.ScanBatchSize; // large batch improves GPU occupancy - UInt128 kStart = 1UL; - ulong divMul = (ulong)((((UInt128)1 << 64) - UInt128.One) / exponent) + 1UL; - byte last = lastIsSeven ? (byte)1 : (byte)0; // ILGPU kernels do not support bool parameters + public void Scan(ulong exponent, UInt128 twoP, bool lastIsSeven, UInt128 maxK, ref bool isPrime) + { + var gpuLease = GpuKernelPool.GetKernel(_useGpuOrder); + var execution = gpuLease.EnterExecutionScope(); + var accelerator = gpuLease.Accelerator; + var stream = gpuLease.Stream; + int batchSize = GpuConstants.ScanBatchSize; // large batch improves GPU occupancy + UInt128 kStart = 1UL; + ulong divMul = (ulong)((((UInt128)1 << 64) - UInt128.One) / exponent) + 1UL; + byte last = lastIsSeven ? (byte)1 : (byte)0; // ILGPU kernels do not support bool parameters - var pow2Kernel = gpuLease.Pow2ModKernel; - var incKernel = gpuLease.IncrementalKernel; - ulong step10 = (exponent.Mod10() << 1).Mod10(); - ulong step8 = ((exponent & 7UL) << 1) & 7UL; - ulong step3 = (exponent.Mod3() << 1).Mod3(); - ulong step5 = (exponent.Mod5() << 1).Mod5(); - GpuUInt128 twoPGpu = (GpuUInt128)twoP; - var smallCyclesView = GpuKernelPool.EnsureSmallCyclesOnDevice(accelerator); - ResiduePrimeViews primeViews = default; - if (_kernelType == GpuKernelType.Pow2Mod) - { - primeViews = GpuKernelPool.EnsureSmallPrimesOnDevice(accelerator); - } + var pow2Kernel = gpuLease.Pow2ModKernel; + var incKernel = gpuLease.IncrementalKernel; + exponent.Mod10_8_5_3Steps(out ulong step10, out ulong step8, out ulong step5, out ulong step3); + GpuUInt128 twoPGpu = (GpuUInt128)twoP; + var smallCyclesView = GpuKernelPool.EnsureSmallCyclesOnDevice(accelerator); + ResiduePrimeViews primeViews = default; + if (_kernelType == GpuKernelType.Pow2Mod) + { + primeViews = GpuKernelPool.EnsureSmallPrimesOnDevice(accelerator); + } - var orderBuffer = accelerator.Allocate1D(batchSize); - // Avoid giant stack allocations that can trigger StackOverflow when batchSize is large. - // Rent a reusable array from the shared pool instead. - ulong[] orderArray = ArrayPool.Shared.Rent(batchSize); - UInt128 remaining; - int currentSize; - int i; - UInt128 q = UInt128.Zero; - try - { - while (kStart <= maxK && Volatile.Read(ref isPrime)) - { - remaining = maxK - kStart + UInt128.One; - currentSize = remaining > (UInt128)batchSize ? batchSize : (int)remaining; - Span orders = orderArray.AsSpan(0, currentSize); - // Precompute residue automaton bases for this batch - UInt128 q0 = twoP * kStart + UInt128.One; - if (_kernelType == GpuKernelType.Pow2Mod) - { - q0.Mod10_8_5_3(out ulong q0m10, out ulong q0m8, out ulong q0m5, out ulong q0m3); - var kernelArgs = new ResidueAutomatonArgs(q0m10, step10, q0m8, step8, q0m3, step3, q0m5, step5); - pow2Kernel(stream, currentSize, exponent, twoPGpu, (GpuUInt128)kStart, last, divMul, - kernelArgs, orderBuffer.View, smallCyclesView, primeViews.LastOne, primeViews.LastSeven, - primeViews.LastOnePow2, primeViews.LastSevenPow2); - } - else - { - q0.Mod10_8_5_3(out ulong q0m10, out ulong q0m8, out ulong q0m5, out ulong q0m3); - incKernel(stream, currentSize, exponent, twoPGpu, (GpuUInt128)kStart, last, divMul, - q0m10, q0m8, q0m3, q0m5, orderBuffer.View, smallCyclesView); - } + var orderBuffer = accelerator.Allocate1D(batchSize); + // Avoid giant stack allocations that can trigger StackOverflow when batchSize is large. + // Rent a reusable array from the shared pool instead. + ulong[] orderArray = ArrayPool.Shared.Rent(batchSize); + UInt128 remaining; + int currentSize; + int i; + UInt128 q = UInt128.Zero; + try + { + while (kStart <= maxK && Volatile.Read(ref isPrime)) + { + remaining = maxK - kStart + UInt128.One; + currentSize = remaining > (UInt128)batchSize ? batchSize : (int)remaining; + Span orders = orderArray.AsSpan(0, currentSize); + // Precompute residue automaton bases for this batch + UInt128 q0 = twoP * kStart + UInt128.One; + if (_kernelType == GpuKernelType.Pow2Mod) + { + q0.Mod10_8_5_3(out ulong q0m10, out ulong q0m8, out ulong q0m5, out ulong q0m3); + var kernelArgs = new ResidueAutomatonArgs(q0m10, step10, q0m8, step8, q0m3, step3, q0m5, step5); + pow2Kernel(stream, currentSize, exponent, twoPGpu, (GpuUInt128)kStart, last, divMul, + kernelArgs, orderBuffer.View, smallCyclesView, primeViews.LastOne, primeViews.LastSeven, + primeViews.LastOnePow2, primeViews.LastSevenPow2); + } + else + { + q0.Mod10_8_5_3(out ulong q0m10, out ulong q0m8, out ulong q0m5, out ulong q0m3); + incKernel(stream, currentSize, exponent, twoPGpu, (GpuUInt128)kStart, last, divMul, + q0m10, q0m8, q0m3, q0m5, orderBuffer.View, smallCyclesView); + } - stream.Synchronize(); - orderBuffer.View.CopyToCPU(ref MemoryMarshal.GetReference(orders), currentSize); - if (_kernelType == GpuKernelType.Pow2Mod) - { - for (i = 0; i < currentSize; i++) - { - if (orders[i] != 0UL) - { - Volatile.Write(ref isPrime, false); - break; - } - } - } - else - { - q = twoP * kStart + UInt128.One; - for (i = 0; i < currentSize && Volatile.Read(ref isPrime); i++, q += twoP) - { - if (orders[i] != 0UL && q.IsPrimeCandidate()) - { - Volatile.Write(ref isPrime, false); - } - } - } + stream.Synchronize(); + orderBuffer.View.CopyToCPU(ref MemoryMarshal.GetReference(orders), currentSize); + if (_kernelType == GpuKernelType.Pow2Mod) + { + for (i = 0; i < currentSize; i++) + { + if (orders[i] != 0UL) + { + Volatile.Write(ref isPrime, false); + break; + } + } + } + else + { + q = twoP * kStart + UInt128.One; + for (i = 0; i < currentSize && Volatile.Read(ref isPrime); i++, q += twoP) + { + if (orders[i] != 0UL && q.IsPrimeCandidate()) + { + Volatile.Write(ref isPrime, false); + } + } + } - kStart += (UInt128)currentSize; - } - } - finally - { - ArrayPool.Shared.Return(orderArray); - orderBuffer.Dispose(); - execution.Dispose(); - gpuLease.Dispose(); - } + kStart += (UInt128)currentSize; + } + } + finally + { + ArrayPool.Shared.Return(orderArray); + orderBuffer.Dispose(); + execution.Dispose(); + gpuLease.Dispose(); + } - } + } } diff --git a/PerfectNumbers.Core/Gpu/MersenneNumberLucasLehmerGpuTester.cs b/PerfectNumbers.Core/Gpu/MersenneNumberLucasLehmerGpuTester.cs index 0bdca4c..3afcf72 100644 --- a/PerfectNumbers.Core/Gpu/MersenneNumberLucasLehmerGpuTester.cs +++ b/PerfectNumbers.Core/Gpu/MersenneNumberLucasLehmerGpuTester.cs @@ -51,7 +51,7 @@ public class MersenneNumberLucasLehmerGpuTester // Early rejections aligned with incremental/order sieves, but safe for small p: // - If 3 | p and p != 3, then 7 | M_p -> composite. // - If p ≡ 1 (mod 4) and p shares a factor with (p-1), reject fast. - if (exponent.Mod3() == 0UL && exponent != 3UL) + if ((exponent % 3UL) == 0UL && exponent != 3UL) { return false; } diff --git a/PerfectNumbers.Core/Gpu/MersenneNumberOrderGpuTester.cs b/PerfectNumbers.Core/Gpu/MersenneNumberOrderGpuTester.cs index 0d1e738..3708378 100644 --- a/PerfectNumbers.Core/Gpu/MersenneNumberOrderGpuTester.cs +++ b/PerfectNumbers.Core/Gpu/MersenneNumberOrderGpuTester.cs @@ -25,12 +25,7 @@ public class MersenneNumberOrderGpuTester(GpuKernelType kernelType, bool useGpuO }; var foundBuffer = accelerator.Allocate1D(1); - UInt128 exponent128 = exponent; - exponent128.Mod10_8_5_3(out ulong exponentMod10, out ulong exponentMod8, out ulong exponentMod5, out ulong exponentMod3); - ulong step10 = (exponentMod10 << 1).Mod10(); - ulong step8 = (exponentMod8 << 1) & 7UL; - ulong step3 = (exponentMod3 << 1).Mod3(); - ulong step5 = (exponentMod5 << 1).Mod5(); + exponent.Mod10_8_5_3Steps(out ulong step10, out ulong step8, out ulong step5, out ulong step3); try { UInt128 remaining; diff --git a/PerfectNumbers.Core/Gpu/MersenneNumberResidueGpuTester.cs b/PerfectNumbers.Core/Gpu/MersenneNumberResidueGpuTester.cs index a5fb8ae..207c4f1 100644 --- a/PerfectNumbers.Core/Gpu/MersenneNumberResidueGpuTester.cs +++ b/PerfectNumbers.Core/Gpu/MersenneNumberResidueGpuTester.cs @@ -24,7 +24,11 @@ public class MersenneNumberResidueGpuTester(bool useGpuOrder) UInt128 kStart = 1UL; byte last = lastIsSeven ? (byte)1 : (byte)0; var kernel = gpuLease.Pow2ModKernel; - exponent.Mod10_8_5_3Steps(out ulong step10, out ulong step8, out ulong step5, out ulong step3); + ulong step3 = exponent % 3UL, + step5 = exponent % 5UL, + step8 = exponent % 8UL, + step10 = exponent.Mod10(); + GpuUInt128 twoPGpu = (GpuUInt128)twoP; var orderBuffer = accelerator.Allocate1D(batchSize); diff --git a/PerfectNumbers.Core/MersenneDivisorCycles.cs b/PerfectNumbers.Core/MersenneDivisorCycles.cs index edf7930..96c51a6 100644 --- a/PerfectNumbers.Core/MersenneDivisorCycles.cs +++ b/PerfectNumbers.Core/MersenneDivisorCycles.cs @@ -85,7 +85,7 @@ public class MersenneDivisorCycles } } } - + // Binary search for divisor in sorted _table int mid, right, left = 0; ulong d; @@ -189,7 +189,7 @@ public class MersenneDivisorCycles for (int i = 0; i < localCycleIndex; i++) { (d, cycle) = localCycles[i]; - if (d.Mod3() == 0UL || d.Mod5() == 0UL || d.Mod7() == 0UL || d.Mod11() == 0UL) + if ((d % 3UL) == 0UL || (d % 5UL) == 0UL || (d % 7UL) == 0UL || (d % 11UL) == 0UL) continue; writer.Write(d); @@ -222,15 +222,15 @@ public class MersenneDivisorCycles int count = (int)Math.Min(batchSizeUL, maxDivisor), i, idx; ulong[] divisors, outCycles, validDivisors; - divisors = pool.Rent(count); - outCycles = pool.Rent(count); - var lease = GpuKernelPool.GetKernel(useGpuOrder: true); - var accelerator = lease.Accelerator; - var stream = accelerator.CreateStream(); - var kernel = accelerator.LoadAutoGroupedStreamKernel< - Index1D, - ArrayView1D, - ArrayView1D>(GpuDivisorCycleKernel); + divisors = pool.Rent(count); + outCycles = pool.Rent(count); + var lease = GpuKernelPool.GetKernel(useGpuOrder: true); + var accelerator = lease.Accelerator; + var stream = accelerator.CreateStream(); + var kernel = accelerator.LoadAutoGroupedStreamKernel< + Index1D, + ArrayView1D, + ArrayView1D>(GpuDivisorCycleKernel); using Stream outputStream = new FileStream(path, FileMode.OpenOrCreate, FileAccess.ReadWrite, FileShare.Read, BufferSize10M, useAsync: true); @@ -249,7 +249,7 @@ public class MersenneDivisorCycles { while (skipCount > 0L && start <= maxDivisor) { - if ((start & 1UL) == 0UL || start.Mod3() == 0UL || start.Mod5() == 0UL || start.Mod7() == 0UL || start.Mod11() == 0UL) + if ((start & 1UL) == 0UL || (start % 3UL) == 0UL || (start % 5UL) == 0UL || (start % 7UL) == 0UL || (start % 11UL) == 0UL) { start++; continue; @@ -280,10 +280,10 @@ public class MersenneDivisorCycles var bufferDiv = accelerator.Allocate1D(validDivisors); var bufferCycle = accelerator.Allocate1D(idx); - kernel( - idx, - bufferDiv.View, - bufferCycle.View); + kernel( + idx, + bufferDiv.View, + bufferCycle.View); accelerator.Synchronize(); @@ -296,7 +296,7 @@ public class MersenneDivisorCycles { d = divisors[i]; - if (d.Mod3() == 0UL || d.Mod5() == 0UL || d.Mod7() == 0UL || d.Mod11() == 0UL) + if ((d % 3UL) == 0UL || (d % 5UL) == 0UL || (d % 7UL) == 0UL || (d % 11UL) == 0UL) continue; writer.Write(d); @@ -308,16 +308,16 @@ public class MersenneDivisorCycles start = end + 1UL; } } - finally - { - pool.Return(divisors, clearArray: false); - pool.Return(outCycles, clearArray: false); - } + finally + { + pool.Return(divisors, clearArray: false); + pool.Return(outCycles, clearArray: false); + } - stream.Dispose(); - lease.Dispose(); + stream.Dispose(); + lease.Dispose(); - } + } const byte ByteZero = 0; const byte ByteOne = 1; diff --git a/PerfectNumbers.Core/MersenneNumberTester.cs b/PerfectNumbers.Core/MersenneNumberTester.cs index 8f9fe6b..9a24f46 100644 --- a/PerfectNumbers.Core/MersenneNumberTester.cs +++ b/PerfectNumbers.Core/MersenneNumberTester.cs @@ -203,7 +203,7 @@ public sealed class MersenneNumberTester( { // Safe early rejections using known orders for tiny primes: // 7 | M_p iff 3 | p. Avoid rejecting p == 3 where M_3 == 7 is prime. - if (exponent.Mod3() == 0UL && exponent != 3UL) + if ((exponent % 3UL) == 0UL && exponent != 3UL) { return false; } diff --git a/PerfectNumbers.Core/PrimeTester.cs b/PerfectNumbers.Core/PrimeTester.cs index fa962d8..8e9fde4 100644 --- a/PerfectNumbers.Core/PrimeTester.cs +++ b/PerfectNumbers.Core/PrimeTester.cs @@ -53,7 +53,7 @@ public sealed class PrimeTester(bool useInternal = false) { result = false; } - else if (n > 5UL && n.Mod5() == 0UL) + else if (n > 5UL && (n % 5UL) == 0UL) { result = false; } @@ -260,7 +260,7 @@ public sealed class PrimeTester(bool useInternal = false) return; } - if ((n & 1UL) == 0UL || (n > 5UL && n.Mod5() == 0UL)) + if ((n & 1UL) == 0UL || (n > 5UL && (n % 5UL) == 0UL)) { results[index] = 0; return; diff --git a/PerfectNumbers.Core/UInt128Extensions.cs b/PerfectNumbers.Core/UInt128Extensions.cs index ee50019..0474b48 100644 --- a/PerfectNumbers.Core/UInt128Extensions.cs +++ b/PerfectNumbers.Core/UInt128Extensions.cs @@ -96,109 +96,147 @@ public static class UInt128Extensions } [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static bool IsPrimeCandidate(this UInt128 n) - { - UInt128 p; + public static bool IsPrimeCandidate(this UInt128 n) + { + UInt128 p, zero = UInt128.Zero; - uint[] smallPrimes = PrimesGenerator.SmallPrimes; - ulong[] smallPrimesPow2 = PrimesGenerator.SmallPrimesPow2; - ulong i, smallPrimesCount = (ulong)smallPrimes.Length; - for (i = 0UL; i < smallPrimesCount; i++) - { - p = smallPrimes[i]; - if (smallPrimesPow2[i] > n) - { - break; - } + uint[] smallPrimes = PrimesGenerator.SmallPrimes; + ulong[] smallPrimesPow2 = PrimesGenerator.SmallPrimesPow2; + ulong i, smallPrimesCount = (ulong)smallPrimes.Length; + for (i = 0UL; i < smallPrimesCount; i++) + { + if (smallPrimesPow2[i] > n) + { + break; + } - if (n % p == UInt128.Zero) - { - return n == p; - } - } + p = smallPrimes[i]; + if (n % p == zero) + { + return n == p; + } + } return true; } [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static void Mod10_8_5_3(this UInt128 value, out ulong mod10, out ulong mod8, out ulong mod5, out ulong mod3) - { - UInt128 temp = UInt128.Zero; - uint byteSum = 0U; - byte current = 0; - uint accumulator = 0U; - - temp = value; - do - { - current = (byte)temp; - byteSum += current; - temp >>= 8; - } - while (temp != UInt128.Zero); - - mod8 = (ulong)value & 7UL; - - accumulator = byteSum.Mod5(); - mod5 = accumulator; - mod10 = (mod8 & 1UL) == 0UL - ? accumulator switch - { - 0U => 0UL, - 1U => 6UL, - 2U => 2UL, - 3U => 8UL, - _ => 4UL, - } - : accumulator switch - { - 0U => 5UL, - 1U => 1UL, - 2U => 7UL, - 3U => 3UL, - _ => 9UL, - }; - - // Reusing accumulator to hold the mod 3 reduction now that the mod 5 value has been stored. - accumulator = byteSum.Mod3(); - mod3 = accumulator; - } + public static ulong Mod3(this UInt128 value) + { + ulong remainder = ((ulong)value) % 3UL + ((ulong)(value >> 64)) % 3UL; + return remainder >= 3UL ? remainder - 3UL : remainder; + } [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static UInt128 Pow(this UInt128 value, ulong exponent) + public static ulong Mod5(this UInt128 value) { - UInt128 result = UInt128.One; - UInt128 baseValue = value; + ulong remainder = (((ulong)value) % 5UL) + (((ulong)(value >> 64)) % 5UL); + return remainder >= 5UL ? remainder - 5UL : remainder; + } - while (exponent != 0UL) + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong Mod6(this UInt128 value) => ((value.Mod3() << 1) | ((ulong)value & 1UL)) switch + { + 0UL => 0UL, + 1UL => 3UL, + 2UL => 4UL, + 3UL => 1UL, + 4UL => 2UL, + _ => 5UL, + }; + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong Mod7(this UInt128 value) + { + ulong remainder = (((ulong)value) % 7UL) + (((ulong)(value >> 64)) % 7UL) << 1; + while (remainder >= 7UL) { - if ((exponent & 1UL) != 0UL) - { - result *= baseValue; - } - - baseValue *= baseValue; - exponent >>= 1; + remainder -= 7UL; } - return result; + return remainder; } + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong Mod8(this UInt128 value) => (ulong)value & 7UL; + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong Mod10(this UInt128 value128) + { + // Split and fold under mod 10: 2^64 ≡ 6 (mod 10) + ulong value = (ulong)value128; + return ((value % 10UL) + ((value >> 64) % 10UL) * 6UL) % 10UL; + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static void Mod10_8_5_3(this UInt128 value, out ulong mod10, out ulong mod8, out ulong mod5, out ulong mod3) + { + byte current; + uint byteSum = 0U; + mod8 = (ulong)value & 7UL; + UInt128 zero = UInt128.Zero; + do + { + current = (byte)value; + byteSum += current; + value >>= 8; + } + while (value != zero); + + mod5 = byteSum.Mod5(); + mod10 = (mod8 & 1UL) == 0UL + ? mod5 switch + { + 0U => 0UL, + 1U => 6UL, + 2U => 2UL, + 3U => 8UL, + _ => 4UL, + } + : mod5 switch + { + 0U => 5UL, + 1U => 1UL, + 2U => 7UL, + 3U => 3UL, + _ => 9UL, + }; + + mod3 = byteSum % 3U; + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong Mod11(this UInt128 value) + { + ulong remainder = (((ulong)value) % 11UL) + (((ulong)(value >> 64)) % 11UL) * 5UL; + while (remainder >= 11UL) + { + remainder -= 11UL; + } + + return remainder; + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong Mod128(this UInt128 value) => (ulong)value & 127UL; + + [MethodImpl(MethodImplOptions.AggressiveInlining)] public static UInt128 ModPow(this UInt128 value, UInt128 exponent, UInt128 modulus) { - UInt128 zero = UInt128.Zero; - UInt128 result = UInt128.One; - UInt128 baseValue = value % modulus; + UInt128 zero = UInt128.Zero; + UInt128 one = UInt128.One; + UInt128 result = UInt128.One; + value %= modulus; - while (exponent != zero) - { - if ((exponent & UInt128.One) != zero) - { - result = MulMod(result, baseValue, modulus); - } + while (exponent != zero) + { + if ((exponent & one) != zero) + { + result = MulMod(result, value, modulus); + } - baseValue = MulMod(baseValue, baseValue, modulus); + value = MulMod(value, value, modulus); exponent >>= 1; } @@ -206,87 +244,8 @@ public static class UInt128Extensions } [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod10(this UInt128 value) - { - UInt128 zero = UInt128.Zero; - if (value == zero) - return 0UL; - - ulong result = (ulong)value; - UInt128 high = value >> 64; - - while (high != zero) - { - // 2^64 ≡ 6 (mod 10) - result = (result + (ulong)high * 6UL) % 10UL; - high >>= 64; - } - - return result % 10UL; - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod8(this UInt128 value) => (ulong)value & 7UL; - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod3(this UInt128 value) - { - ulong remainder = ((ulong)value % 3UL) + ((ulong)(value >> 64) % 3UL); - return remainder >= 3UL ? remainder - 3UL : remainder; - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod5(this UInt128 value) - { - ulong remainder = ((ulong)value % 5UL) + ((ulong)(value >> 64) % 5UL); - return remainder >= 5UL ? remainder - 5UL : remainder; - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod6(this UInt128 value) - { - return ((value.Mod3() << 1) | ((ulong)value & 1UL)) switch - { - 0UL => 0UL, - 1UL => 3UL, - 2UL => 4UL, - 3UL => 1UL, - 4UL => 2UL, - _ => 5UL, - }; - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod7(this UInt128 value) - { - ulong remainder = ((ulong)value % 7UL) + ((ulong)(value >> 64) % 7UL) * 2UL; - while (remainder >= 7UL) - { - remainder -= 7UL; - } - - return remainder; - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod11(this UInt128 value) - { - ulong remainder = ((ulong)value % 11UL) + ((ulong)(value >> 64) % 11UL) * 5UL; - while (remainder >= 11UL) - { - remainder -= 11UL; - } - - return remainder; - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod128(this UInt128 value) => (ulong)value & 127UL; - - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static UInt128 Mul64(this UInt128 a, UInt128 b) - { + public static UInt128 Mul64(this UInt128 a, UInt128 b) + { ulong aLow = (ulong)a; ulong bLow = (ulong)b; return ((UInt128)(aLow * (ulong)(b >> 64) + aLow.MulHigh(bLow)) << 64) | (aLow * bLow); @@ -295,15 +254,16 @@ public static class UInt128Extensions [MethodImpl(MethodImplOptions.AggressiveInlining)] public static UInt128 MulMod(this UInt128 a, UInt128 b, UInt128 modulus) { - UInt128 zero = UInt128.Zero; - UInt128 result = UInt128.Zero; + UInt128 one = UInt128.One, + result = UInt128.Zero, + zero = UInt128.Zero; - while (b != zero) - { - if ((b & UInt128.One) != zero) - { - result += a; - if (result >= modulus) + while (b != zero) + { + if ((b & one) != zero) + { + result += a; + if (result >= modulus) { result -= modulus; } @@ -320,4 +280,22 @@ public static class UInt128Extensions return result; } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static UInt128 Pow(this UInt128 value, ulong exponent) + { + UInt128 result = UInt128.One; + while (exponent != 0UL) + { + if ((exponent & 1UL) != 0UL) + { + result *= value; + } + + value *= value; + exponent >>= 1; + } + + return result; + } } diff --git a/PerfectNumbers.Core/UIntExtensions.cs b/PerfectNumbers.Core/UIntExtensions.cs index ba23e63..e057cfe 100644 --- a/PerfectNumbers.Core/UIntExtensions.cs +++ b/PerfectNumbers.Core/UIntExtensions.cs @@ -4,28 +4,16 @@ namespace PerfectNumbers.Core; public static class UIntExtensions { - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static uint Mod3(this uint value) => value % 3U; + private const ulong Mod3Mask = 0xAAAAAAABUL; + public const ulong Mod5Mask = 0xCCCCCCCDUL; - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static uint Mod5(this uint value) => value % 5U; + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static uint Mod3(uint value) => value - (uint)((value * Mod3Mask) >> 33) * 3U; - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static uint Mod6(this uint value) => value % 6U; + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static uint Mod5(this uint value) => value - (uint)((value * Mod5Mask) >> 34) * 5U; - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static uint Mod7(this uint value) => value % 7U; - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static uint Mod8(this uint value) => value & 7U; - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static uint Mod10(this uint value) => value % 10U; - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static uint Mod11(this uint value) => value % 11U; - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static uint Mod128(this uint value) => value & 127U; + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static uint Mod8(this uint value) => value & 7U; } diff --git a/PerfectNumbers.Core/ULongExtensions.cs b/PerfectNumbers.Core/ULongExtensions.cs index a311acb..bc43055 100644 --- a/PerfectNumbers.Core/ULongExtensions.cs +++ b/PerfectNumbers.Core/ULongExtensions.cs @@ -48,196 +48,110 @@ public static class ULongExtensions [MethodImpl(MethodImplOptions.AggressiveInlining)] public static ulong FastDiv64(this ulong value, ulong divisor, ulong mul) { - ulong q = value.MulHigh(mul); - if ((UInt128)value - q.Mul64(divisor) >= divisor) + mul = value.MulHigh(mul); + if ((UInt128)value - mul.Mul64(divisor) >= divisor) { - q++; + mul++; } - return q; + return mul; } - private static readonly byte[] Mod7ByteCoefficients = { 1, 4, 2, 1, 4, 2, 1, 4 }; - private static readonly byte[] Mod11ByteCoefficients = { 1, 3, 9, 5, 4, 1, 3, 9 }; + public const ulong WordBitMask = 0xFFFFUL; - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static void Mod10_8_5_3(this ulong value, out ulong mod10, out ulong mod8, out ulong mod5, out ulong mod3) - { - ulong temp = 0UL; - uint byteSum = 0U; - byte current = 0; - uint accumulator = 0U; + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static bool IsPrimeCandidate(this ulong n) + { + int i = 0; + uint[] smallPrimes = PrimesGenerator.SmallPrimes; + ulong[] smallPrimesPow2 = PrimesGenerator.SmallPrimesPow2; + int len = smallPrimes.Length; + ulong p; + for (; i < len; i++) + { + if (smallPrimesPow2[i] > n) + { + break; + } - temp = value; - do - { - current = (byte)temp; - byteSum += current; - temp >>= 8; - } - while (temp != 0UL); + p = smallPrimes[i]; + if ((n % p) == 0UL) + { + return n == p; + } + } - mod8 = value & 7UL; + return true; + } - accumulator = byteSum.Mod5(); - mod5 = accumulator; - mod10 = (mod8 & 1UL) == 0UL - ? accumulator switch - { - 0U => 0UL, - 1U => 6UL, - 2U => 2UL, - 3U => 8UL, - _ => 4UL, - } - : accumulator switch - { - 0U => 5UL, - 1U => 1UL, - 2U => 7UL, - 3U => 3UL, - _ => 9UL, - }; + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong Mod10(this ulong value) => (value & 1UL) == 0UL + ? value % 5UL switch + { + 0UL => 0UL, + 1UL => 6UL, + 2UL => 2UL, + 3UL => 8UL, + _ => 4UL, + } + : value % 5UL switch + { + 0UL => 5UL, + 1UL => 1UL, + 2UL => 7UL, + 3UL => 3UL, + _ => 9UL, + }; - // Reusing accumulator to hold the mod 3 reduction now that the mod 5 value has been stored. - accumulator = byteSum.Mod3(); - mod3 = accumulator; - } + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong Mod128(this ulong value) => value & 127UL; - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static void Mod11_7_5_3(this ulong value, out ulong mod11, out ulong mod7, out ulong mod5, out ulong mod3) - { - uint mod3Accumulator = 0U; - uint mod5Accumulator = 0U; - uint mod7Accumulator = 0U; - uint mod11Accumulator = 0U; - ulong temp = 0UL; - int index = 0; - byte current = 0; + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static void Mod10_8_5_3(this ulong value, out ulong mod10, out ulong mod8, out ulong mod5, out ulong mod3) + { + mod8 = value & 7UL; - temp = value; - do - { - current = (byte)temp; - mod3Accumulator += current; - mod5Accumulator += current; - mod7Accumulator += (uint)(current * Mod7ByteCoefficients[index]); - mod11Accumulator += (uint)(current * Mod11ByteCoefficients[index]); + byte current; + uint byteSum = 0U; + do + { + current = (byte)value; + byteSum += current; + value >>= 8; + } + while (value != 0UL); - temp >>= 8; - index++; - } - while (temp != 0UL); + mod3 = byteSum % 3U; + mod5 = byteSum % 5U; + mod10 = (mod8 & 1UL) == 0UL + ? mod5 switch + { + 0UL => 0UL, + 1UL => 6UL, + 2UL => 2UL, + 3UL => 8UL, + _ => 4UL, + } + : mod5 switch + { + 0UL => 5UL, + 1UL => 1UL, + 2UL => 7UL, + 3UL => 3UL, + _ => 9UL, + }; + } - // Reusing mod3Accumulator to store the reduced mod 3 remainder before exposing it. - mod3Accumulator = mod3Accumulator.Mod3(); - // Reusing mod5Accumulator to store the reduced mod 5 remainder before exposing it. - mod5Accumulator = mod5Accumulator.Mod5(); - // Reusing mod7Accumulator to store the reduced mod 7 remainder before exposing it. - mod7Accumulator = (uint)ReduceMod7(mod7Accumulator); - // Reusing mod11Accumulator to store the reduced mod 11 remainder before exposing it. - mod11Accumulator = (uint)ReduceMod11(mod11Accumulator); + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static void Mod10_8_5_3Steps(this ulong value, out ulong step10, out ulong step8, out ulong step5, out ulong step3) + { + value.Mod10_8_5_3(out step10, out step8, out step5, out step3); - mod3 = mod3Accumulator; - mod5 = mod5Accumulator; - mod7 = mod7Accumulator; - mod11 = mod11Accumulator; - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static void Mod10_8_5_3Steps(this ulong value, out ulong step10, out ulong step8, out ulong step5, out ulong step3) - { - ulong mod10 = 0UL; - ulong mod8 = 0UL; - ulong mod5 = 0UL; - ulong mod3 = 0UL; - - value.Mod10_8_5_3(out mod10, out mod8, out mod5, out mod3); - - step10 = (mod10 << 1).Mod10(); - step8 = ((mod8 << 1) & 7UL); - step5 = (mod5 << 1).Mod5(); - step3 = (mod3 << 1).Mod3(); - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static bool IsPrimeCandidate(this ulong n) - { - int i = 0; - ulong p = 0UL; - uint[] smallPrimes = PrimesGenerator.SmallPrimes; - ulong[] smallPrimesPow2 = PrimesGenerator.SmallPrimesPow2; - int len = smallPrimes.Length; - for (; i < len; i++) - { - if (smallPrimesPow2[i] > n) - { - break; - } - - p = smallPrimes[i]; - if ((n % p) == 0UL) - { - return n == p; - } - } - - return true; - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod10(this ulong value) => value % 10UL; - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod8(this ulong value) => value & 7UL; - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod3(this ulong value) => value % 3UL; - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod5(this ulong value) => value % 5UL; - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod6(this ulong value) => value % 6UL; - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod7(this ulong value) => value % 7UL; - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod11(this ulong value) => value % 11UL; - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong Mod128(this ulong value) => value & 127UL; - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - private static ulong ReduceMod7(uint value) - { - while (value >= 7U) - { - value = (value >> 3) + (value & 7U); - if (value >= 7U && value < 14U) - { - value -= 7U; - } - } - - return value; - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - private static ulong ReduceMod11(uint value) - { - while (value >= 11U) - { - value = (value & 15U) + ((value >> 4) * 5U); - if (value >= 11U && value < 22U) - { - value -= 11U; - } - } - - return value; - } + step10 = (step10 << 1).Mod10(); + step8 = (step8 << 1) & 7UL; + step5 = (step5 << 1) % 5UL; + step3 = (step3 << 1) % 3UL; + } [MethodImpl(MethodImplOptions.AggressiveInlining)] public static UInt128 Mul64(this ulong a, ulong b) => ((UInt128)a.MulHigh(b) << 64) | (UInt128)(a * b); @@ -252,384 +166,326 @@ public static class ULongExtensions ulong w1 = xLow * yHigh; ulong w2 = xHigh * yLow; - ulong carry = - (((xLow * yLow) >> 32) + - (uint)w1 + - (uint)w2) >> 32; - return (xHigh * yHigh) + (w1 >> 32) + (w2 >> 32) + carry; + return (xHigh * yHigh) + (w1 >> 32) + (w2 >> 32) + + ( + ((xLow * yLow) >> 32) + + (uint)w1 + + (uint)w2 + ) >> 32; } [MethodImpl(MethodImplOptions.AggressiveInlining)] public static ulong ModPow64(this ulong value, ulong exponent, ulong modulus) { ulong result = 1UL; - ulong baseValue = value % modulus; + value %= modulus; while (exponent != 0UL) { if ((exponent & 1UL) != 0UL) { - result = MulMod64(result, baseValue, modulus); + result = MulMod64(result, value, modulus); } - baseValue = MulMod64(baseValue, baseValue, modulus); + value = MulMod64(value, value, modulus); exponent >>= 1; } return result; } - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static ulong MulMod64(this ulong a, ulong b, ulong modulus) - { - if (modulus <= 1UL) - { - return 0UL; - } - - if (TryGetMersenneExponent(modulus, out int exponent)) - { - return MulModMersenne(a, b, modulus, exponent); - } - - return (ulong)(((UInt128)a * b) % modulus); - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - private static bool TryGetMersenneExponent(ulong modulus, out int exponent) - { - if (modulus == 0UL) - { - exponent = 0; - return false; - } - - ulong plusOne = modulus + 1UL; - if (plusOne == 0UL) - { - exponent = 64; - return true; - } - - if ((plusOne & (plusOne - 1UL)) != 0UL) - { - exponent = 0; - return false; - } - - exponent = BitOperations.TrailingZeroCount(plusOne); - return true; - } - - [MethodImpl(MethodImplOptions.AggressiveInlining)] - private static ulong MulModMersenne(ulong a, ulong b, ulong modulus, int exponent) - { - UInt128 mask = (UInt128)modulus; - UInt128 product = (UInt128)a * b; - UInt128 folded = (product & mask) + (product >> exponent); - - while (folded > mask) - { - folded = (folded & mask) + (folded >> exponent); - } - - ulong result = (ulong)folded; - if (result >= modulus) - { - result -= modulus; - } - - return result; - } + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static ulong MulMod64(this ulong a, ulong b, ulong modulus) => (ulong)(UInt128)(((a % modulus) * (b % modulus)) % modulus); [MethodImpl(MethodImplOptions.AggressiveInlining)] - public static UInt128 PowMod(this ulong exponent, UInt128 modulus) - { - UInt128 result = UInt128.One; - ulong exponentLoopIndex = 0UL; + public static UInt128 PowMod(this ulong exponent, UInt128 modulus) + { + UInt128 result = UInt128.One; + ulong exponentLoopIndex = 0UL; - // Return 1 because 2^0 = 1 - if (exponent == 0UL) - return result; + // Return 1 because 2^0 = 1 + if (exponent == 0UL) + return result; - // Any number mod 1 is 0 - if (modulus == UInt128.One) - return UInt128.Zero; + // Any number mod 1 is 0 + if (modulus == UInt128.One) + return UInt128.Zero; - // For small exponents, do classic method - if (exponent < 64 || modulus < 4) - { - for (; exponentLoopIndex < exponent; exponentLoopIndex++) - { - result <<= 1; - if (result >= modulus) - result -= modulus; - } + // For small exponents, do classic method + if (exponent < 64 || modulus < 4) + { + for (; exponentLoopIndex < exponent; exponentLoopIndex++) + { + result <<= 1; + if (result >= modulus) + result -= modulus; + } - return result; - } + return result; + } - // Special case: if modulus is a power of two, use bitmasking for efficiency and correctness - if ((modulus & (modulus - 1)) == 0) - { - result = UInt128.One << (int)(exponent & 127); - return result & (modulus - 1); - } + // Special case: if modulus is a power of two, use bitmasking for efficiency and correctness + if ((modulus & (modulus - 1)) == 0) + { + result = UInt128.One << (int)(exponent & 127); + return result & (modulus - 1); + } - // Reusing exponentLoopIndex to iterate again over the exponent range for the general-case accumulation. - exponentLoopIndex = 0UL; - // Reusing result after resetting it for the general modulus accumulation phase. - result = UInt128.One; - for (; exponentLoopIndex < exponent; exponentLoopIndex++) - { - result <<= 1; - if (result >= modulus) - result -= modulus; - } + // Reusing exponentLoopIndex to iterate again over the exponent range for the general-case accumulation. + exponentLoopIndex = 0UL; + // Reusing result after resetting it for the general modulus accumulation phase. + result = UInt128.One; + for (; exponentLoopIndex < exponent; exponentLoopIndex++) + { + result <<= 1; + if (result >= modulus) + result -= modulus; + } - return result; - } + return result; + } /// /// Computes 2^exponent mod modulus using a known cycle length. /// - public static UInt128 PowModWithCycle(this ulong exponent, UInt128 modulus, ulong cycleLength) - { - UInt128 result = UInt128.One; - ulong exponentLoopIndex = 0UL; - ulong rotationCount = 0UL; + public static UInt128 PowModWithCycle(this ulong exponent, UInt128 modulus, ulong cycleLength) + { + UInt128 result = UInt128.One; + ulong exponentLoopIndex = 0UL; - // Return 1 because 2^0 = 1 - if (exponent == 0UL) - return result; + // Return 1 because 2^0 = 1 + if (exponent == 0UL) + return result; - // Any number mod 1 is 0 - if (modulus == UInt128.One) - return UInt128.Zero; + // Any number mod 1 is 0 + if (modulus == UInt128.One) + return UInt128.Zero; - // For small exponents, do classic method - if (exponent < 64 || modulus < 4) - { - for (; exponentLoopIndex < exponent; exponentLoopIndex++) - { - result <<= 1; - if (result >= modulus) - result -= modulus; - } + // For small exponents, do classic method + if (exponent < 64 || modulus < 4) + { + for (; exponentLoopIndex < exponent; exponentLoopIndex++) + { + result <<= 1; + if (result >= modulus) + result -= modulus; + } - return result; - } + return result; + } - // Special case: if modulus is a power of two, use bitmasking for efficiency and correctness - if ((modulus & (modulus - 1)) == 0) - { - result = UInt128.One << (int)(exponent & 127); - return result & (modulus - 1); - } + // Special case: if modulus is a power of two, use bitmasking for efficiency and correctness + if ((modulus & (modulus - 1)) == 0) + { + result = UInt128.One << (int)(exponent & 127); + return result & (modulus - 1); + } - // Reusing exponentLoopIndex to iterate over the rotation count for the cycle-aware path. - exponentLoopIndex = 0UL; - // Reusing result after resetting it for the rotation-accumulation pass. - result = UInt128.One; - rotationCount = exponent % cycleLength; - for (; exponentLoopIndex < rotationCount; exponentLoopIndex++) - { - result <<= 1; - if (result >= modulus) - result -= modulus; - } + // Reusing exponentLoopIndex to iterate over the rotation count for the cycle-aware path. + exponentLoopIndex = 0UL; + // Reusing result after resetting it for the rotation-accumulation pass. + result = UInt128.One; + ulong rotationCount = exponent % cycleLength; + for (; exponentLoopIndex < rotationCount; exponentLoopIndex++) + { + result <<= 1; + if (result >= modulus) + result -= modulus; + } - return result; - } + return result; + } /// /// Computes 2^exponent mod modulus using a known cycle length. /// - public static UInt128 PowModWithCycle(this ulong exponent, UInt128 modulus, UInt128 cycleLength) - { - UInt128 result = UInt128.One; - ulong exponentLoopIndex = 0UL; - UInt128 rotationIndex = UInt128.Zero; - UInt128 rotationCount = UInt128.Zero; + public static UInt128 PowModWithCycle(this ulong exponent, UInt128 modulus, UInt128 cycleLength) + { + UInt128 result = UInt128.One; + ulong exponentLoopIndex = 0UL; - // Return 1 because 2^0 = 1 - if (exponent == UInt128.Zero) - return result; + // Return 1 because 2^0 = 1 + if (exponent == UInt128.Zero) + return result; - // Any number mod 1 is 0 - if (modulus == UInt128.One) - return UInt128.Zero; + // Any number mod 1 is 0 + if (modulus == UInt128.One) + return UInt128.Zero; - // For small exponents, do classic method - if (exponent < UInt128Numbers.SixtyFour || modulus < UInt128Numbers.Four) - { - for (; exponentLoopIndex < exponent; exponentLoopIndex++) - { - result <<= 1; - if (result >= modulus) - result -= modulus; - } + // For small exponents, do classic method + if (exponent < UInt128Numbers.SixtyFour || modulus < UInt128Numbers.Four) + { + for (; exponentLoopIndex < exponent; exponentLoopIndex++) + { + result <<= 1; + if (result >= modulus) + result -= modulus; + } - return result; - } + return result; + } - // Special case: if modulus is a power of two, use bitmasking for efficiency and correctness - if ((modulus & (modulus - UInt128.One)) == UInt128.Zero) - { - result = UInt128.One << (int)(exponent & UInt128Numbers.OneHundredTwentySeven); - return result & (modulus - 1); - } + // Special case: if modulus is a power of two, use bitmasking for efficiency and correctness + if ((modulus & (modulus - UInt128.One)) == UInt128.Zero) + { + result = UInt128.One << (int)(exponent & UInt128Numbers.OneHundredTwentySeven); + return result & (modulus - 1); + } - // Reusing result after resetting it for the rotation-driven accumulation phase. - result = UInt128.One; - rotationCount = exponent % cycleLength; - rotationIndex = UInt128.Zero; - while (rotationIndex < rotationCount) - { - result <<= 1; - if (result >= modulus) - result -= modulus; + // Reusing result after resetting it for the rotation-driven accumulation phase. + result = UInt128.One; + UInt128 rotationCount = exponent % cycleLength; + UInt128 rotationIndex = UInt128.Zero; + while (rotationIndex < rotationCount) + { + result <<= 1; + if (result >= modulus) + result -= modulus; - rotationIndex += UInt128.One; - } + rotationIndex += UInt128.One; + } - return result; - } + return result; + } /// /// Computes 2^exponent mod modulus using a known cycle length. /// - public static UInt128 PowModWithCycle(this UInt128 exponent, UInt128 modulus, UInt128 cycleLength) - { - UInt128 result = UInt128.One; - ulong exponentLoopIndex = 0UL; - UInt128 rotationIndex = UInt128.Zero; - UInt128 rotationCount = UInt128.Zero; + public static UInt128 PowModWithCycle(this UInt128 exponent, UInt128 modulus, UInt128 cycleLength) + { + UInt128 one = UInt128.One, + result = one, + zero = UInt128.Zero; + ulong exponentLoopIndex = 0UL; - // Return 1 because 2^0 = 1 - if (exponent == UInt128.Zero) - return result; + // Return 1 because 2^0 = 1 + if (exponent == zero) + return result; - // Any number mod 1 is 0 - if (modulus == UInt128.One) - return UInt128.Zero; + // Any number mod 1 is 0 + if (modulus == one) + return zero; - // For small exponents, do classic method - if (exponent < UInt128Numbers.SixtyFour || modulus < UInt128Numbers.Four) - { - for (; exponentLoopIndex < exponent; exponentLoopIndex++) - { - result <<= 1; - if (result >= modulus) - result -= modulus; - } + // For small exponents, do classic method + if (exponent < UInt128Numbers.SixtyFour || modulus < UInt128Numbers.Four) + { + for (; exponentLoopIndex < exponent; exponentLoopIndex++) + { + result <<= 1; + if (result >= modulus) + result -= modulus; + } - return result; - } + return result; + } - // Special case: if modulus is a power of two, use bitmasking for efficiency and correctness - if ((modulus & (modulus - UInt128.One)) == UInt128.Zero) - { - result = UInt128.One << (int)(exponent & UInt128Numbers.OneHundredTwentySeven); - return result & (modulus - 1); - } + // Special case: if modulus is a power of two, use bitmasking for efficiency and correctness + if ((modulus & (modulus - one)) == zero) + { + result = one << (int)(exponent & UInt128Numbers.OneHundredTwentySeven); + return result & (modulus - 1); + } - // Reusing result after resetting it for the rotation-driven accumulation phase. - result = UInt128.One; - rotationCount = exponent % cycleLength; - rotationIndex = UInt128.Zero; - while (rotationIndex < rotationCount) - { - result <<= 1; - if (result >= modulus) - result -= modulus; + // Reusing result after resetting it for the rotation-driven accumulation phase. + result = one; + UInt128 rotationCount = exponent % cycleLength; + + // We're reusing "zero" as rotation index for just a little better performance + while (zero < rotationCount) + { + result <<= 1; + if (result >= modulus) + result -= modulus; - rotationIndex += UInt128.One; - } + zero += one; + } - return result; - } + return result; + } /// /// Computes 2^exponent mod modulus using iterative CRT composition from mod 10 up to modulus. /// Only for modulus >= 10 and reasonable size. /// - public static UInt128 PowModCrt(this ulong exponent, UInt128 modulus, MersenneDivisorCycles cycles) - { - if (modulus < 10) - return PowMod(exponent, modulus); // fallback to classic + public static UInt128 PowModCrt(this ulong exponent, UInt128 modulus, MersenneDivisorCycles cycles) + { + if (modulus < 10) + return PowMod(exponent, modulus); // fallback to classic - // Use cycle length 4 for mod 10 - UInt128 result = PowModWithCycle(exponent, 10, 4); - UInt128 currentModulus = 10; - UInt128 modulusCandidate = 11; - UInt128 cycle = UInt128.Zero; - UInt128 remainderForCandidate = UInt128.Zero; + // Use cycle length 4 for mod 10 + UInt128 currentModulus = 10, + cycle, + modulusCandidate = 11, + remainderForCandidate, + result = PowModWithCycle(exponent, 10, 4), + zero = UInt128.Zero; - for (; modulusCandidate <= modulus; modulusCandidate++) - { - cycle = MersenneDivisorCycles.GetCycle(modulusCandidate); - remainderForCandidate = cycle > UInt128.Zero - ? PowModWithCycle(exponent, modulusCandidate, cycle) - : PowMod(exponent, modulusCandidate); + for (; modulusCandidate <= modulus; modulusCandidate++) + { + cycle = MersenneDivisorCycles.GetCycle(modulusCandidate); + remainderForCandidate = cycle > zero + ? PowModWithCycle(exponent, modulusCandidate, cycle) + : PowMod(exponent, modulusCandidate); - // Solve x ≡ result mod currentModulus - // x ≡ remM mod m - // Find x mod (currentModulus * m) - // Since currentModulus and m are coprime, use CRT: - // x = result + currentModulus * t, where t ≡ (remM - result) * inv(currentModulus, m) mod m + // Solve x ≡ result mod currentModulus + // x ≡ remM mod m + // Find x mod (currentModulus * m) + // Since currentModulus and m are coprime, use CRT: + // x = result + currentModulus * t, where t ≡ (remM - result) * inv(currentModulus, m) mod m - result += currentModulus * ((remainderForCandidate + modulusCandidate - (result % modulusCandidate)) * ModInverse(currentModulus, modulusCandidate) % modulusCandidate); - currentModulus *= modulusCandidate; + result += currentModulus * ((remainderForCandidate + modulusCandidate - (result % modulusCandidate)) * ModInverse(currentModulus, modulusCandidate) % modulusCandidate); + currentModulus *= modulusCandidate; - if (currentModulus >= modulus) - break; - } + if (currentModulus >= modulus) + break; + } - return result % modulus; - } + return result % modulus; + } - // Helper: modular inverse (extended Euclidean algorithm) - private static UInt128 ModInverse(UInt128 a, UInt128 m) - { - UInt128 m0 = m, temp; - UInt128 x0 = 0, x1 = 1; - UInt128 originalA = UInt128.Zero; - UInt128 originalM = UInt128.Zero; - if (m == 1) - { - return 0; - } + // Helper: modular inverse (extended Euclidean algorithm) + private static UInt128 ModInverse(UInt128 a, UInt128 m) + { + UInt128 m0 = m, + originalA, + originalM, + temp, + x0 = 0, + x1 = 1; - while (a > 1) - { - originalA = a; - originalM = m; - m = originalA % originalM; - a = originalM; - temp = x0; - x0 = x1 - (originalA / originalM) * x0; - x1 = temp; - } + if (m == 1) + { + return 0; + } - if (x1 < 0) - { - x1 += m0; - } + while (a > 1) + { + originalA = a; + originalM = m; + m = originalA % originalM; + a = originalM; + temp = x0; + x0 = x1 - (originalA / originalM) * x0; + x1 = temp; + } - return x1; - } + if (x1 < 0) + { + x1 += m0; + } + + return x1; + } [MethodImpl(MethodImplOptions.AggressiveInlining)] public static bool SharesFactorWithExponentMinusOne(this ulong exponent) { - ulong value = exponent - 1UL; + ulong prime, value = exponent - 1UL; value >>= BitOperations.TrailingZeroCount(value); uint[] smallPrimes = PrimesGenerator.SmallPrimes; ulong[] smallPrimesPow2 = PrimesGenerator.SmallPrimesPow2; int i = 1, smallPrimesLength = smallPrimes.Length; - ulong prime; for (; i < smallPrimesLength && smallPrimesPow2[i] <= value; i++) {