Remove Montgomery data from GPU divisor flow
This commit is contained in:
1 parent
69acd3dad9
commit
46a2178d65
4 files changed
+11
-152
No files matched your search
@@ -1,5 +1,4 @@
|
||||
using System.Runtime.CompilerServices;
|
||||
using PerfectNumbers.Core;
|
||||
|
||||
namespace PerfectNumbers.Core.Gpu;
|
||||
|
||||
@@ -7,32 +6,18 @@ internal readonly struct GpuDivisorPartialData
|
||||
{
|
||||
public readonly ulong Modulus;
|
||||
public readonly ReadOnlyGpuUInt128 ModulusWide;
|
||||
private readonly MontgomeryDivisorData _montgomeryData;
|
||||
|
||||
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||
private GpuDivisorPartialData(ulong modulus, ReadOnlyGpuUInt128 modulusWide, in MontgomeryDivisorData montgomeryData)
|
||||
private GpuDivisorPartialData(ulong modulus, in ReadOnlyGpuUInt128 modulusWide)
|
||||
{
|
||||
Modulus = modulus;
|
||||
ModulusWide = modulusWide;
|
||||
_montgomeryData = montgomeryData;
|
||||
}
|
||||
|
||||
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||
public static GpuDivisorPartialData Create(ulong modulus)
|
||||
{
|
||||
MontgomeryDivisorData montgomery = MontgomeryDivisorData.FromModulus(modulus);
|
||||
return new GpuDivisorPartialData(modulus, new ReadOnlyGpuUInt128(modulus), in montgomery);
|
||||
}
|
||||
|
||||
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||
public static GpuDivisorPartialData FromMontgomery(in MontgomeryDivisorData divisor)
|
||||
{
|
||||
return new GpuDivisorPartialData(divisor.Modulus, new ReadOnlyGpuUInt128(divisor.Modulus), in divisor);
|
||||
}
|
||||
|
||||
public MontgomeryDivisorData MontgomeryData
|
||||
{
|
||||
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||
get => _montgomeryData;
|
||||
ReadOnlyGpuUInt128 modulusWide = new(modulus);
|
||||
return new GpuDivisorPartialData(modulus, in modulusWide);
|
||||
}
|
||||
}
|
||||
@@ -83,30 +83,4 @@ internal static class DivisorByDivisorKernels
|
||||
}
|
||||
}
|
||||
|
||||
public static void ConvertMontgomeryResultsToHitsKernel(
|
||||
Index1D index,
|
||||
ArrayView<ulong> results,
|
||||
ArrayView<byte> hits)
|
||||
{
|
||||
int globalIndex = index;
|
||||
// EvenPerfectBitScanner launches this kernel with an extent matching the result and hit spans,
|
||||
// so the defensive bounds guard stays commented out to keep the conversion branch-free.
|
||||
// int length = (int)results.Length;
|
||||
// if (globalIndex >= length)
|
||||
// {
|
||||
// return;
|
||||
// }
|
||||
|
||||
ulong montgomeryResult = results[globalIndex];
|
||||
hits[globalIndex] = montgomeryResult == 1UL ? (byte)1 : (byte)0;
|
||||
}
|
||||
|
||||
public static void ComputeMontgomeryExponentKernel(Index1D index, GpuDivisorPartialData divisor, ArrayView<ulong> exponents, ArrayView<ulong> results)
|
||||
{
|
||||
ulong modulus = divisor.Modulus;
|
||||
// Each divisor flowing through this kernel follows q = 2kp + 1 with k >= 1, so modulus is always odd and exceeds one.
|
||||
ulong exponent = exponents[index];
|
||||
ulong montgomeryResult = ULongExtensions.Pow2MontgomeryModWindowedGpuKeepMontgomery(divisor.MontgomeryData, exponent);
|
||||
results[index] = montgomeryResult;
|
||||
}
|
||||
}
|
||||
@@ -259,7 +259,8 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
|
||||
{
|
||||
ulong divisorValue = filteredDivisorsSpan[i];
|
||||
GpuDivisorPartialData divisorDataValue = GpuDivisorPartialData.Create(divisorValue);
|
||||
ulong divisorCycle = ResolveDivisorCycle(divisorValue, prime, divisorDataValue);
|
||||
MontgomeryDivisorData montgomeryData = MontgomeryDivisorData.FromModulus(divisorValue);
|
||||
ulong divisorCycle = ResolveDivisorCycle(divisorValue, prime, in montgomeryData);
|
||||
if (divisorCycle != prime)
|
||||
{
|
||||
continue;
|
||||
@@ -305,11 +306,11 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
|
||||
}
|
||||
|
||||
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||
private static ulong ResolveDivisorCycle(ulong divisor, ulong prime, in GpuDivisorPartialData divisorData)
|
||||
private static ulong ResolveDivisorCycle(ulong divisor, ulong prime, in MontgomeryDivisorData divisorData)
|
||||
{
|
||||
if (!MersenneDivisorCycles.TryCalculateCycleLengthForExponentGpu(divisor, prime, divisorData, out ulong computedCycle, out bool primeOrderFailed) || computedCycle == 0UL)
|
||||
if (!MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(divisor, prime, divisorData, out ulong computedCycle, out bool primeOrderFailed) || computedCycle == 0UL)
|
||||
{
|
||||
return MersenneDivisorCycles.CalculateCycleLengthGpu(divisor, divisorData, skipPrimeOrderHeuristic: primeOrderFailed);
|
||||
return MersenneDivisorCycles.CalculateCycleLength(divisor, divisorData, skipPrimeOrderHeuristic: primeOrderFailed);
|
||||
}
|
||||
|
||||
return computedCycle;
|
||||
@@ -459,20 +460,16 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
|
||||
// return;
|
||||
// }
|
||||
|
||||
GpuDivisorPartialData cachedData = GpuDivisorPartialData.FromMontgomery(divisorData);
|
||||
if (cachedData.Modulus != divisor)
|
||||
{
|
||||
cachedData = GpuDivisorPartialData.Create(divisor);
|
||||
}
|
||||
GpuDivisorPartialData cachedData = GpuDivisorPartialData.Create(divisor);
|
||||
|
||||
ulong cycle = divisorCycle;
|
||||
ulong firstPrime = primes[0];
|
||||
|
||||
if (cycle == 0UL)
|
||||
{
|
||||
if (!MersenneDivisorCycles.TryCalculateCycleLengthForExponentGpu(divisor, firstPrime, cachedData, out ulong computedCycle, out bool primeOrderFailed) || computedCycle == 0UL)
|
||||
if (!MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(divisor, firstPrime, divisorData, out ulong computedCycle, out bool primeOrderFailed) || computedCycle == 0UL)
|
||||
{
|
||||
cycle = MersenneDivisorCycles.CalculateCycleLengthGpu(divisor, cachedData, skipPrimeOrderHeuristic: primeOrderFailed);
|
||||
cycle = MersenneDivisorCycles.CalculateCycleLength(divisor, divisorData, skipPrimeOrderHeuristic: primeOrderFailed);
|
||||
}
|
||||
else
|
||||
{
|
||||
|
||||
@@ -326,103 +326,6 @@ public class MersenneDivisorCycles
|
||||
ThreadStaticPools.ReturnFactorCountDictionary(factorCounts);
|
||||
return success;
|
||||
}
|
||||
public static bool TryCalculateCycleLengthForExponentGpu(
|
||||
ulong divisor,
|
||||
ulong exponent,
|
||||
in MontgomeryDivisorData divisorData,
|
||||
out ulong cycleLength,
|
||||
out bool primeOrderFailed)
|
||||
{
|
||||
cycleLength = 0UL;
|
||||
primeOrderFailed = false;
|
||||
|
||||
if ((divisor & (divisor - 1UL)) == 0UL)
|
||||
{
|
||||
cycleLength = 1UL;
|
||||
return true;
|
||||
}
|
||||
|
||||
if (divisor <= 3UL || exponent <= 1UL)
|
||||
{
|
||||
cycleLength = CalculateCycleLength(divisor, divisorData);
|
||||
return true;
|
||||
}
|
||||
|
||||
ulong computedOrder = PrimeOrderCalculator.Calculate(
|
||||
divisor,
|
||||
previousOrder: null,
|
||||
divisorData,
|
||||
PrimeOrderCalculator.PrimeOrderSearchConfig.HeuristicDefault,
|
||||
PrimeOrderCalculator.PrimeOrderHeuristicDevice.Cpu);
|
||||
if (computedOrder != 0UL)
|
||||
{
|
||||
cycleLength = computedOrder;
|
||||
return true;
|
||||
}
|
||||
|
||||
primeOrderFailed = true;
|
||||
|
||||
ulong phi = divisor - 1UL;
|
||||
// A zero totient would imply divisor == 0, which never occurs on the tested path. Keep the guard documented without
|
||||
// branching in the hot loop.
|
||||
// if (phi == 0UL)
|
||||
// {
|
||||
// return false;
|
||||
// }
|
||||
|
||||
int twoCount = BitOperations.TrailingZeroCount(phi);
|
||||
ulong reducedPhi = phi >> twoCount;
|
||||
if (reducedPhi == 0UL || (reducedPhi % exponent) != 0UL)
|
||||
{
|
||||
return false;
|
||||
}
|
||||
|
||||
ulong k = reducedPhi / exponent;
|
||||
// k == 0 would indicate reducedPhi < exponent, which the EvenPerfectBitScanner call flow prevents. Document the
|
||||
// assumption instead of branching.
|
||||
// if (k == 0UL)
|
||||
// {
|
||||
// return false;
|
||||
// }
|
||||
|
||||
Dictionary<ulong, int> factorCounts = ThreadStaticPools.RentFactorCountDictionary();
|
||||
bool success = false;
|
||||
|
||||
if (twoCount > 0)
|
||||
{
|
||||
factorCounts[2UL] = twoCount;
|
||||
}
|
||||
|
||||
bool exponentAccumulated = AccumulateFactors(exponent, factorCounts);
|
||||
bool factorsAccumulated = exponentAccumulated && AccumulateFactors(k, factorCounts);
|
||||
|
||||
if (factorsAccumulated)
|
||||
{
|
||||
cycleLength = ReduceOrder(divisorData, phi, factorCounts);
|
||||
success = true;
|
||||
primeOrderFailed = false;
|
||||
}
|
||||
|
||||
factorCounts.Clear();
|
||||
ThreadStaticPools.ReturnFactorCountDictionary(factorCounts);
|
||||
return success;
|
||||
}
|
||||
|
||||
internal static bool TryCalculateCycleLengthForExponentGpu(
|
||||
ulong divisor,
|
||||
ulong exponent,
|
||||
in GpuDivisorPartialData divisorData,
|
||||
out ulong cycleLength,
|
||||
out bool primeOrderFailed) =>
|
||||
TryCalculateCycleLengthForExponentGpu(divisor, exponent, divisorData.MontgomeryData, out cycleLength, out primeOrderFailed);
|
||||
|
||||
internal static ulong CalculateCycleLengthGpu(
|
||||
ulong divisor,
|
||||
in GpuDivisorPartialData divisorData,
|
||||
bool skipPrimeOrderHeuristic = false)
|
||||
{
|
||||
return CalculateCycleLength(divisor, divisorData.MontgomeryData, skipPrimeOrderHeuristic);
|
||||
}
|
||||
|
||||
private static bool AccumulateFactors(
|
||||
ulong value,
|
||||
|
||||
Reference in new issue
Block a user