Use bit masks for divisor heuristics
This commit is contained in:
1 parent
46a2178d65
commit
bbe0c02069
6 files changed
+83
-43
No files matched your search
@@ -175,17 +175,31 @@ public sealed class MersenneNumberDivisorByDivisorCpuTester : IMersenneNumberDiv
|
||||
// Keep the divisibility filters aligned with the divisor-cycle generator so the
|
||||
// CPU path never requests cycles that were skipped during cache creation.
|
||||
bool lastIsSeven = (prime & 3UL) == 3UL;
|
||||
const ushort DecimalMaskWhenLastIsSeven = (1 << 3) | (1 << 7) | (1 << 9);
|
||||
const ushort DecimalMaskOtherwise = (1 << 1) | (1 << 3) | (1 << 9);
|
||||
const byte AcceptableRemainder8Mask = 0b10000010;
|
||||
const int NonZeroMod3Mask = 0b00000110;
|
||||
const int NonZeroMod5Mask = 0b00011110;
|
||||
const int NonZeroMod7Mask = 0b01111110;
|
||||
const int NonZeroMod11Mask = 0b11111111110;
|
||||
const int RequiredSmallPrimeMask = 0b1111;
|
||||
|
||||
while (divisor.CompareTo(limit) <= 0)
|
||||
{
|
||||
ulong candidate = divisor.Low;
|
||||
processedCount++;
|
||||
|
||||
bool admissible = lastIsSeven
|
||||
? (remainder10 == 3 || remainder10 == 7 || remainder10 == 9)
|
||||
: (remainder10 == 1 || remainder10 == 3 || remainder10 == 9);
|
||||
ushort decimalMask = lastIsSeven ? DecimalMaskWhenLastIsSeven : DecimalMaskOtherwise;
|
||||
bool admissible = ((decimalMask >> remainder10) & 1) != 0;
|
||||
|
||||
if (admissible && (remainder8 == 1 || remainder8 == 7) && remainder3 != 0 && remainder5 != 0 && remainder7 != 0 && remainder11 != 0)
|
||||
bool passesOctalMask = ((AcceptableRemainder8Mask >> remainder8) & 1) != 0;
|
||||
int smallPrimeBits = ((NonZeroMod3Mask >> remainder3) & 1)
|
||||
| (((NonZeroMod5Mask >> remainder5) & 1) << 1)
|
||||
| (((NonZeroMod7Mask >> remainder7) & 1) << 2)
|
||||
| (((NonZeroMod11Mask >> remainder11) & 1) << 3);
|
||||
bool passesSmallPrimeMasks = (smallPrimeBits & RequiredSmallPrimeMask) == RequiredSmallPrimeMask;
|
||||
|
||||
if (admissible && passesOctalMask && passesSmallPrimeMasks)
|
||||
{
|
||||
MontgomeryDivisorData divisorData = MontgomeryDivisorData.FromModulus(candidate);
|
||||
ulong divisorCycle;
|
||||
|
||||
@@ -8,16 +8,16 @@ internal readonly struct GpuDivisorPartialData
|
||||
public readonly ReadOnlyGpuUInt128 ModulusWide;
|
||||
|
||||
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||
private GpuDivisorPartialData(ulong modulus, in ReadOnlyGpuUInt128 modulusWide)
|
||||
public GpuDivisorPartialData(ulong modulus)
|
||||
{
|
||||
Modulus = modulus;
|
||||
ModulusWide = new ReadOnlyGpuUInt128(modulus);
|
||||
}
|
||||
|
||||
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||
public GpuDivisorPartialData(ulong modulus, in ReadOnlyGpuUInt128 modulusWide)
|
||||
{
|
||||
Modulus = modulus;
|
||||
ModulusWide = modulusWide;
|
||||
}
|
||||
|
||||
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||
public static GpuDivisorPartialData Create(ulong modulus)
|
||||
{
|
||||
ReadOnlyGpuUInt128 modulusWide = new(modulus);
|
||||
return new GpuDivisorPartialData(modulus, in modulusWide);
|
||||
}
|
||||
}
|
||||
@@ -258,7 +258,6 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
|
||||
for (int i = 0; i < filteredCount; i++)
|
||||
{
|
||||
ulong divisorValue = filteredDivisorsSpan[i];
|
||||
GpuDivisorPartialData divisorDataValue = GpuDivisorPartialData.Create(divisorValue);
|
||||
MontgomeryDivisorData montgomeryData = MontgomeryDivisorData.FromModulus(divisorValue);
|
||||
ulong divisorCycle = ResolveDivisorCycle(divisorValue, prime, in montgomeryData);
|
||||
if (divisorCycle != prime)
|
||||
@@ -321,27 +320,38 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
|
||||
{
|
||||
// Mirror the CPU-side divisor screen so the GPU path rejects multiples of 3, 5, 7, and 11
|
||||
// before resolving cycles, keeping cache lookups aligned with the persisted snapshot.
|
||||
ulong remainder10 = candidate % 10UL;
|
||||
bool accept10 = lastIsSeven
|
||||
? (remainder10 == 3UL || remainder10 == 7UL || remainder10 == 9UL)
|
||||
: (remainder10 == 1UL || remainder10 == 3UL || remainder10 == 9UL);
|
||||
if (!accept10)
|
||||
const ushort DecimalMaskWhenLastIsSeven = (1 << 3) | (1 << 7) | (1 << 9);
|
||||
const ushort DecimalMaskOtherwise = (1 << 1) | (1 << 3) | (1 << 9);
|
||||
const byte AcceptableRemainder8Mask = 0b10000010;
|
||||
const int NonZeroMod3Mask = 0b00000110;
|
||||
const int NonZeroMod5Mask = 0b00011110;
|
||||
const int NonZeroMod7Mask = 0b01111110;
|
||||
const int NonZeroMod11Mask = 0b11111111110;
|
||||
const int RequiredSmallPrimeMask = 0b1111;
|
||||
|
||||
int remainder10 = (int)(candidate % 10UL);
|
||||
ushort decimalMask = lastIsSeven ? DecimalMaskWhenLastIsSeven : DecimalMaskOtherwise;
|
||||
if (((decimalMask >> remainder10) & 1) == 0)
|
||||
{
|
||||
return false;
|
||||
}
|
||||
|
||||
ulong remainder8 = candidate & 7UL;
|
||||
if (remainder8 != 1UL && remainder8 != 7UL)
|
||||
int remainder8 = (int)(candidate & 7UL);
|
||||
if (((AcceptableRemainder8Mask >> remainder8) & 1) == 0)
|
||||
{
|
||||
return false;
|
||||
}
|
||||
|
||||
if (candidate % 3UL == 0UL || candidate % 5UL == 0UL || candidate % 7UL == 0UL || candidate % 11UL == 0UL)
|
||||
{
|
||||
return false;
|
||||
}
|
||||
int remainder3 = (int)(candidate % 3UL);
|
||||
int remainder5 = (int)(candidate % 5UL);
|
||||
int remainder7 = (int)(candidate % 7UL);
|
||||
int remainder11 = (int)(candidate % 11UL);
|
||||
int smallPrimeBits = (NonZeroMod3Mask >> remainder3) & 1;
|
||||
smallPrimeBits |= ((NonZeroMod5Mask >> remainder5) & 1) << 1;
|
||||
smallPrimeBits |= ((NonZeroMod7Mask >> remainder7) & 1) << 2;
|
||||
smallPrimeBits |= ((NonZeroMod11Mask >> remainder11) & 1) << 3;
|
||||
|
||||
return true;
|
||||
return (smallPrimeBits & RequiredSmallPrimeMask) == RequiredSmallPrimeMask;
|
||||
}
|
||||
|
||||
|
||||
@@ -460,7 +470,7 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
|
||||
// return;
|
||||
// }
|
||||
|
||||
GpuDivisorPartialData cachedData = GpuDivisorPartialData.Create(divisor);
|
||||
GpuDivisorPartialData cachedData = new GpuDivisorPartialData(divisor);
|
||||
|
||||
ulong cycle = divisorCycle;
|
||||
ulong firstPrime = primes[0];
|
||||
|
||||
@@ -58,7 +58,11 @@ public class MersenneDivisorCycles
|
||||
return true;
|
||||
}
|
||||
|
||||
if (exponent <= 1UL || divisor <= 1UL || (divisor & 1UL) == 0UL)
|
||||
const int InvalidParameterMask = 0b111;
|
||||
int invalidBits = (exponent <= 1UL ? 1 : 0)
|
||||
| ((divisor <= 1UL ? 1 : 0) << 1)
|
||||
| ((int)(1UL ^ (divisor & 1UL)) << 2);
|
||||
if ((invalidBits & InvalidParameterMask) != 0)
|
||||
{
|
||||
return false;
|
||||
}
|
||||
@@ -91,10 +95,12 @@ public class MersenneDivisorCycles
|
||||
|
||||
public static bool CycleEqualsExponentForMersenneCandidate(ulong divisor, in MontgomeryDivisorData divisorData, ulong exponent)
|
||||
{
|
||||
if (!IsValidMersenneDivisorCandidate(divisor, exponent))
|
||||
{
|
||||
return false;
|
||||
}
|
||||
// EvenPerfectBitScanner only calls this helper with divisors generated from 2 * k * p + 1 and primes greater than one,
|
||||
// so the validation guard would never trigger in production. Keep it commented out to leave the hot path branch-free.
|
||||
// if (!IsValidMersenneDivisorCandidate(divisor, exponent))
|
||||
// {
|
||||
// return false;
|
||||
// }
|
||||
|
||||
if (divisor <= PerfectNumberConstants.MaxQForDivisorCycles)
|
||||
{
|
||||
@@ -115,7 +121,11 @@ public class MersenneDivisorCycles
|
||||
|
||||
private static bool IsValidMersenneDivisorCandidate(ulong divisor, ulong exponent)
|
||||
{
|
||||
if (exponent <= 1UL || divisor <= 1UL || (divisor & 1UL) == 0UL)
|
||||
const int InvalidParameterMask = 0b111;
|
||||
int invalidBits = (exponent <= 1UL ? 1 : 0)
|
||||
| ((divisor <= 1UL ? 1 : 0) << 1)
|
||||
| ((int)(1UL ^ (divisor & 1UL)) << 2);
|
||||
if ((invalidBits & InvalidParameterMask) != 0)
|
||||
{
|
||||
return false;
|
||||
}
|
||||
@@ -255,17 +265,19 @@ public class MersenneDivisorCycles
|
||||
cycleLength = 0UL;
|
||||
primeOrderFailed = false;
|
||||
|
||||
if ((divisor & (divisor - 1UL)) == 0UL)
|
||||
{
|
||||
cycleLength = 1UL;
|
||||
return true;
|
||||
}
|
||||
// EvenPerfectBitScanner never feeds powers of two into this helper, so leave the branch commented out to avoid redundant work.
|
||||
// if ((divisor & (divisor - 1UL)) == 0UL)
|
||||
// {
|
||||
// cycleLength = 1UL;
|
||||
// return true;
|
||||
// }
|
||||
|
||||
if (divisor <= 3UL || exponent <= 1UL)
|
||||
{
|
||||
cycleLength = CalculateCycleLength(divisor, divisorData);
|
||||
return true;
|
||||
}
|
||||
// The scanner only calls this path with prime exponents greater than one and odd divisors above three, keeping the small-parameter guard unreachable.
|
||||
// if (divisor <= 3UL || exponent <= 1UL)
|
||||
// {
|
||||
// cycleLength = CalculateCycleLength(divisor, divisorData);
|
||||
// return true;
|
||||
// }
|
||||
|
||||
ulong computedOrder = PrimeOrderCalculator.Calculate(
|
||||
divisor,
|
||||
|
||||
@@ -188,7 +188,9 @@ internal static partial class PrimeOrderCalculator
|
||||
private static ulong InitializeStartingOrderCpu(ulong prime, ulong phi, in MontgomeryDivisorData divisorData)
|
||||
{
|
||||
ulong order = phi;
|
||||
if ((prime & 7UL) == 1UL || (prime & 7UL) == 7UL)
|
||||
const int AcceptablePrimeResidueMask = (1 << 1) | (1 << 7);
|
||||
int residue = (int)(prime & 7UL);
|
||||
if (((AcceptablePrimeResidueMask >> residue) & 1) != 0)
|
||||
{
|
||||
ulong half = phi >> 1;
|
||||
if (Pow2EqualsOneCpu(half, prime, divisorData))
|
||||
|
||||
@@ -155,6 +155,8 @@ public static partial class ULongExtensions
|
||||
GpuUInt128 result = new(2UL);
|
||||
if (windowValue == 1UL)
|
||||
{
|
||||
// Single-bit windows show up whenever the bit scan finds an isolated exponent bit,
|
||||
// which happens routinely on the GPU path, so keep the fast return enabled.
|
||||
return 2UL;
|
||||
}
|
||||
|
||||
|
||||
Reference in new issue
Block a user