mlemanczyk/even-perfect-numbers-scanner
1
0 Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Heuristic tables split into 2 groups

Browse Source
This commit is contained in:
mlemanczyk committed 2025-11-17 16:04:27 +01:00
1 parent c413ee1162
commit 183b147135
26 files changed
+1172 -940
Download Patch File Download Diff File Expand all files Collapse all files

No files matched your search

PerfectNumbers.Core/Cpu/MersenneCpuDivisorScanSession.cs
+1 -1
View File
@@ -62,7 +62,7 @@ internal sealed class MersenneCpuDivisorScanSession : IMersenneNumberDivisorByDi
}
}
public void Dispose()
public void Return()
{
ThreadStaticPools.ReturnMersenneCpuDivisorSession(this);
}
PerfectNumbers.Core/Cpu/MersenneNumberDivisorByDivisorCpuTester.cs
+487 -475
View File
@@ -2,486 +2,498 @@ using System.Collections.Generic;
using System.Runtime.CompilerServices;
using PerfectNumbers.Core;
using PerfectNumbers.Core.Gpu;
using PerfectNumbers.Core.Gpu.Accelerators;
namespace PerfectNumbers.Core.Cpu;
public sealed class MersenneNumberDivisorByDivisorCpuTester : IMersenneNumberDivisorByDivisorTester
{
private ulong _divisorLimit;
private int _batchSize = 1_024;
[ThreadStatic]
private static GpuUInt128WorkSet _divisorScanGpuWorkSet;
private struct GpuUInt128WorkSet
{
public GpuUInt128 Step;
public GpuUInt128 Divisor;
public GpuUInt128 Limit;
}
public int BatchSize
{
get => _batchSize;
set => _batchSize = Math.Max(1, value);
}
public void ConfigureFromMaxPrime(ulong maxPrime)
{
_divisorLimit = ComputeDivisorLimitFromMaxPrime(maxPrime);
}
public ulong DivisorLimit
{
[MethodImpl(MethodImplOptions.AggressiveInlining)]
get
{
return _divisorLimit;
}
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public ulong GetAllowedMaxDivisor(ulong prime) => ComputeAllowedMaxDivisor(prime, _divisorLimit);
public bool IsPrime(ulong prime, out bool divisorsExhausted)
{
ulong allowedMax = ComputeAllowedMaxDivisor(prime, _divisorLimit);
// The CPU by-divisor run always hands us primes with enormous divisor limits, so the fallback below never executes.
// if (allowedMax < 3UL)
// {
// // EvenPerfectBitScanner routes primes below the small-divisor cutoff to the GPU path, so the CPU path still sees
// // trivial candidates during targeted tests. Short-circuit here to keep those runs aligned with the production flow.
// divisorsExhausted = true;
// return true;
// }
bool processedAll;
bool composite = CheckDivisors(
prime,
allowedMax,
out processedAll);
if (composite)
{
divisorsExhausted = true;
return false;
}
divisorsExhausted = processedAll || composite;
return true;
}
public void PrepareCandidates(in ReadOnlySpan<ulong> primes, Span<ulong> allowedMaxValues)
{
ulong divisorLimit = _divisorLimit;
int length = primes.Length;
for (int index = 0; index < length; index++)
{
allowedMaxValues[index] = ComputeAllowedMaxDivisor(primes[index], divisorLimit);
}
}
public IMersenneNumberDivisorByDivisorTester.IDivisorScanSession CreateDivisorSession()
{
MersenneCpuDivisorScanSession? session = ThreadStaticPools.RentMersenneCpuDivisorSession();
if (session is not null)
{
return session;
}
return new MersenneCpuDivisorScanSession();
}
private bool CheckDivisors(
ulong prime,
ulong allowedMax,
out bool processedAll)
{
processedAll = false;
// The EvenPerfectBitScanner feeds primes >= 138,000,000 here, so allowedMax >= 3 in production runs.
// Keeping the guard commented out documents the reasoning for benchmarks and tests.
// if (allowedMax < 3UL)
// {
// return false;
// }
ref GpuUInt128WorkSet workSet = ref _divisorScanGpuWorkSet;
ref GpuUInt128 step = ref workSet.Step;
step.High = 0UL;
step.Low = prime;
step.ShiftLeft(1);
ref GpuUInt128 limit = ref workSet.Limit;
limit.High = 0UL;
limit.Low = allowedMax;
ref GpuUInt128 divisor = ref workSet.Divisor;
divisor.High = step.High;
divisor.Low = step.Low;
divisor.Add(1UL);
if (divisor.CompareTo(limit) > 0)
{
processedAll = true;
return false;
}
// Intentionally recomputes factorizations without a per-thread cache.
// The previous factor cache recorded virtually no hits and only slowed down the scan.
DivisorCycleCache cycleCache = DivisorCycleCache.Shared;
ulong stepHigh = step.High;
ulong stepLow = step.Low;
ulong limitHigh = limit.High;
ulong divisorLow = divisor.Low;
// Keep the divisibility filters aligned with the divisor-cycle generator so the
// CPU path never requests cycles that were skipped during cache creation.
LastDigit lastDigit = (prime & 3UL) == 3UL ? LastDigit.Seven : LastDigit.One;
ushort decimalMask = DivisorGenerator.GetDecimalMask(lastDigit);
byte step10;
byte step8;
byte step3;
byte step7;
byte step11;
byte remainder10;
byte remainder8;
byte remainder3;
byte remainder7;
byte remainder11;
if (stepHigh == 0UL && limitHigh == 0UL)
{
step10 = (byte)(stepLow % 10UL);
step8 = (byte)(stepLow & 7UL);
step3 = (byte)(stepLow % 3UL);
step7 = (byte)(stepLow % 7UL);
step11 = (byte)(stepLow % 11UL);
remainder10 = (byte)(divisorLow % 10UL);
remainder8 = (byte)(divisorLow & 7UL);
remainder3 = (byte)(divisorLow % 3UL);
remainder7 = (byte)(divisorLow % 7UL);
remainder11 = (byte)(divisorLow % 11UL);
return CheckDivisors64(
prime,
stepLow,
limit.Low,
divisorLow,
decimalMask,
step10,
step8,
step3,
step7,
step11,
remainder10,
remainder8,
remainder3,
remainder7,
remainder11,
out processedAll);
}
ulong divisorHigh = divisor.High;
step10 = (byte)((((stepHigh % 10UL) * 6UL) + (stepLow % 10UL)) % 10UL);
step8 = (byte)(stepLow % 8UL);
step3 = (byte)(((stepHigh % 3UL) + (stepLow % 3UL)) % 3UL);
step7 = (byte)((((stepHigh % 7UL) * 2UL) + (stepLow % 7UL)) % 7UL);
step11 = (byte)((((stepHigh % 11UL) * 5UL) + (stepLow % 11UL)) % 11UL);
remainder10 = (byte)((((divisorHigh % 10UL) * 6UL) + (divisorLow % 10UL)) % 10UL);
remainder8 = (byte)(divisorLow % 8UL);
remainder3 = (byte)(((divisorHigh % 3UL) + (divisorLow % 3UL)) % 3UL);
remainder7 = (byte)((((divisorHigh % 7UL) * 2UL) + (divisorLow % 7UL)) % 7UL);
remainder11 = (byte)((((divisorHigh % 11UL) * 5UL) + (divisorLow % 11UL)) % 11UL);
while (divisor.CompareTo(limit) <= 0)
{
bool passesSmallModuli = remainder3 != 0 && remainder7 != 0 && remainder11 != 0;
if (passesSmallModuli && (remainder8 == 1 || remainder8 == 7) && ((decimalMask >> remainder10) & 1) != 0)
{
ulong candidate = divisor.Low;
MontgomeryDivisorData divisorData = MontgomeryDivisorData.FromModulus(candidate);
ulong divisorCycle;
// Divisors generated from 2 * k * p + 1 exceed the small-cycle snapshot when p >= 138,000,000, so the short path below never runs.
// if (candidate <= PerfectNumberConstants.MaxQForDivisorCycles)
// {
// divisorCycle = cycleCache.GetCycleLength(candidate);
// }
// else
{
if (!MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(
candidate,
prime,
divisorData,
out ulong computedCycle,
out bool primeOrderFailed) || computedCycle == 0UL)
{
// Divisors produced by 2 * k * p + 1 always exceed PerfectNumberConstants.MaxQForDivisorCycles
// for the exponents scanned here, so skip the unused cache fallback and compute directly.
divisorCycle = MersenneDivisorCycles.CalculateCycleLength(
candidate,
divisorData,
skipPrimeOrderHeuristic: primeOrderFailed);
}
else
{
divisorCycle = computedCycle;
}
}
if (divisorCycle == prime)
{
// A cycle equal to the tested exponent (which is prime in this path) guarantees that the candidate divides
// the corresponding Mersenne number because the order of 2 modulo the divisor is exactly p.
processedAll = true;
return true;
}
if (divisorCycle == 0UL)
{
Console.WriteLine($"Divisor cycle was not calculated for {prime}");
}
}
divisor.Add(step);
remainder10 = AddMod10(remainder10, step10);
remainder8 = AddMod8(remainder8, step8);
remainder3 = AddMod3(remainder3, step3);
remainder7 = AddMod7(remainder7, step7);
remainder11 = AddMod11(remainder11, step11);
}
processedAll = true;
return false;
}
private static bool CheckDivisors64(
ulong prime,
ulong step,
ulong limit,
ulong divisor,
ushort decimalMask,
byte step10,
byte step8,
byte step3,
byte step7,
byte step11,
byte remainder10,
byte remainder8,
byte remainder3,
byte remainder7,
byte remainder11,
out bool processedAll)
{
processedAll = false;
bool canAdvance = step <= limit;
if (!canAdvance)
{
if (divisor <= limit)
{
if (remainder3 != 0 && remainder7 != 0 && remainder11 != 0 && (remainder8 == 1 || remainder8 == 7) && ((decimalMask >> remainder10) & 1) != 0)
{
MontgomeryDivisorData divisorData = MontgomeryDivisorData.FromModulus(divisor);
ulong divisorCycle;
if (!MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(
divisor,
prime,
divisorData,
out ulong computedCycle,
out bool primeOrderFailed) || computedCycle == 0UL)
{
divisorCycle = MersenneDivisorCycles.CalculateCycleLength(
divisor,
divisorData,
skipPrimeOrderHeuristic: primeOrderFailed);
}
else
{
divisorCycle = computedCycle;
}
if (divisorCycle == prime)
{
processedAll = true;
return true;
}
if (divisorCycle == 0UL)
{
Console.WriteLine($"Divisor cycle was not calculated for {prime}");
}
}
}
processedAll = true;
return false;
}
if (divisor > limit)
{
processedAll = true;
return false;
}
ulong remainingIterations = ((limit - divisor) / step) + 1UL;
while (true)
{
bool passesSmallModuli = remainder3 != 0 && remainder7 != 0 && remainder11 != 0;
if (passesSmallModuli && (remainder8 == 1 || remainder8 == 7) && ((decimalMask >> remainder10) & 1) != 0)
{
MontgomeryDivisorData divisorData = MontgomeryDivisorData.FromModulus(divisor);
ulong divisorCycle;
if (!MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(
divisor,
prime,
divisorData,
out ulong computedCycle,
out bool primeOrderFailed) || computedCycle == 0UL)
{
divisorCycle = MersenneDivisorCycles.CalculateCycleLength(
divisor,
divisorData,
skipPrimeOrderHeuristic: primeOrderFailed);
}
else
{
divisorCycle = computedCycle;
}
if (divisorCycle == prime)
{
processedAll = true;
return true;
}
if (divisorCycle == 0UL)
{
Console.WriteLine($"Divisor cycle was not calculated for {prime}");
}
}
remainingIterations--;
if (remainingIterations == 0UL)
{
processedAll = true;
return false;
}
divisor += step;
remainder10 = AddMod10(remainder10, step10);
remainder8 = AddMod8(remainder8, step8);
remainder3 = AddMod3(remainder3, step3);
remainder7 = AddMod7(remainder7, step7);
remainder11 = AddMod11(remainder11, step11);
}
}
private static byte CheckDivisor(ulong prime, ulong divisorCycle, in MontgomeryDivisorData divisorData)
{
ulong residue = prime.Pow2MontgomeryModWithCycleCpu(divisorCycle, divisorData);
return residue == 1UL ? (byte)1 : (byte)0;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static byte AddMod3(byte value, byte delta)
{
const int Modulus = 3;
int sum = value + delta;
if (sum >= Modulus)
{
sum -= Modulus;
}
return (byte)sum;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static byte AddMod7(byte value, byte delta)
{
const int Modulus = 7;
int sum = value + delta;
if (sum >= Modulus)
{
sum -= Modulus;
}
return (byte)sum;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static byte AddMod10(byte value, byte delta)
{
const int Modulus = 10;
int sum = value + delta;
if (sum >= Modulus)
{
sum -= Modulus;
}
return (byte)sum;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static byte AddMod11(byte value, byte delta)
{
const int Modulus = 11;
int sum = value + delta;
if (sum >= Modulus)
{
sum -= Modulus;
}
return (byte)sum;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static byte AddMod8(byte value, byte delta)
{
return (byte)((value + delta) & 7);
}
private static ulong ComputeDivisorLimitFromMaxPrime(ulong maxPrime)
{
// The by-divisor CPU configuration only supplies primes greater than 1, so the guard below never trips.
// if (maxPrime <= 1UL)
// {
// return 0UL;
// }
if (maxPrime - 1UL >= 64UL)
{
return ulong.MaxValue;
}
return (1UL << (int)(maxPrime - 1UL)) - 1UL;
}
private static ulong ComputeAllowedMaxDivisor(ulong prime, ulong divisorLimit)
{
// Production by-divisor scans only handle primes, so inputs never fall below 2.
// if (prime <= 1UL)
// {
// return 0UL;
// }
if (prime - 1UL >= 64UL)
{
return divisorLimit;
}
return Math.Min((1UL << (int)(prime - 1UL)) - 1UL, divisorLimit);
}
private ulong _divisorLimit;
private int _batchSize = 1_024;
[ThreadStatic]
private static GpuUInt128WorkSet _divisorScanGpuWorkSet;
private struct GpuUInt128WorkSet
{
public GpuUInt128 Step;
public GpuUInt128 Divisor;
public GpuUInt128 Limit;
}
public int BatchSize
{
get => _batchSize;
set => _batchSize = Math.Max(1, value);
}
public void ConfigureFromMaxPrime(ulong maxPrime)
{
_divisorLimit = ComputeDivisorLimitFromMaxPrime(maxPrime);
}
public ulong DivisorLimit
{
[MethodImpl(MethodImplOptions.AggressiveInlining)]
get
{
return _divisorLimit;
}
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public ulong GetAllowedMaxDivisor(ulong prime) => ComputeAllowedMaxDivisor(prime, _divisorLimit);
public bool IsPrime(ulong prime, out bool divisorsExhausted)
{
ulong allowedMax = ComputeAllowedMaxDivisor(prime, _divisorLimit);
// The CPU by-divisor run always hands us primes with enormous divisor limits, so the fallback below never executes.
// if (allowedMax < 3UL)
// {
// // EvenPerfectBitScanner routes primes below the small-divisor cutoff to the GPU path, so the CPU path still sees
// // trivial candidates during targeted tests. Short-circuit here to keep those runs aligned with the production flow.
// divisorsExhausted = true;
// return true;
// }
bool processedAll;
bool composite = CheckDivisors(
prime,
allowedMax,
out processedAll);
if (composite)
{
divisorsExhausted = true;
return false;
}
divisorsExhausted = processedAll || composite;
return true;
}
public void PrepareCandidates(in ReadOnlySpan<ulong> primes, Span<ulong> allowedMaxValues)
{
ulong divisorLimit = _divisorLimit;
int length = primes.Length;
for (int index = 0; index < length; index++)
{
allowedMaxValues[index] = ComputeAllowedMaxDivisor(primes[index], divisorLimit);
}
}
public IMersenneNumberDivisorByDivisorTester.IDivisorScanSession CreateDivisorSession()
{
MersenneCpuDivisorScanSession? session = ThreadStaticPools.RentMersenneCpuDivisorSession();
if (session is not null)
{
return session;
}
return new MersenneCpuDivisorScanSession();
}
private bool CheckDivisors(
ulong prime,
ulong allowedMax,
out bool processedAll)
{
processedAll = false;
// The EvenPerfectBitScanner feeds primes >= 138,000,000 here, so allowedMax >= 3 in production runs.
// Keeping the guard commented out documents the reasoning for benchmarks and tests.
// if (allowedMax < 3UL)
// {
// return false;
// }
ref GpuUInt128WorkSet workSet = ref _divisorScanGpuWorkSet;
ref GpuUInt128 step = ref workSet.Step;
step.High = 0UL;
step.Low = prime;
step.ShiftLeft(1);
ref GpuUInt128 limit = ref workSet.Limit;
limit.High = 0UL;
limit.Low = allowedMax;
ref GpuUInt128 divisor = ref workSet.Divisor;
divisor.High = step.High;
divisor.Low = step.Low;
divisor.Add(1UL);
if (divisor.CompareTo(limit) > 0)
{
processedAll = true;
return false;
}
// Intentionally recomputes factorizations without a per-thread cache.
// The previous factor cache recorded virtually no hits and only slowed down the scan.
_ = DivisorCycleCache.Shared;
ulong stepHigh = step.High;
ulong stepLow = step.Low;
ulong limitHigh = limit.High;
ulong divisorLow = divisor.Low;
// Keep the divisibility filters aligned with the divisor-cycle generator so the
// CPU path never requests cycles that were skipped during cache creation.
LastDigit lastDigit = (prime & 3UL) == 3UL ? LastDigit.Seven : LastDigit.One;
ushort decimalMask = DivisorGenerator.GetDecimalMask(lastDigit);
byte step10;
byte step8;
byte step3;
byte step7;
byte step11;
byte remainder10;
byte remainder8;
byte remainder3;
byte remainder7;
byte remainder11;
if (stepHigh == 0UL && limitHigh == 0UL)
{
step10 = (byte)(stepLow % 10UL);
step8 = (byte)(stepLow & 7UL);
step3 = (byte)(stepLow % 3UL);
step7 = (byte)(stepLow % 7UL);
step11 = (byte)(stepLow % 11UL);
remainder10 = (byte)(divisorLow % 10UL);
remainder8 = (byte)(divisorLow & 7UL);
remainder3 = (byte)(divisorLow % 3UL);
remainder7 = (byte)(divisorLow % 7UL);
remainder11 = (byte)(divisorLow % 11UL);
return CheckDivisors64(
prime,
stepLow,
limit.Low,
divisorLow,
decimalMask,
step10,
step8,
step3,
step7,
step11,
remainder10,
remainder8,
remainder3,
remainder7,
remainder11,
out processedAll);
}
ulong divisorHigh = divisor.High;
step10 = (byte)((((stepHigh % 10UL) * 6UL) + (stepLow % 10UL)) % 10UL);
step8 = (byte)(stepLow % 8UL);
step3 = (byte)(((stepHigh % 3UL) + (stepLow % 3UL)) % 3UL);
step7 = (byte)((((stepHigh % 7UL) * 2UL) + (stepLow % 7UL)) % 7UL);
step11 = (byte)((((stepHigh % 11UL) * 5UL) + (stepLow % 11UL)) % 11UL);
remainder10 = (byte)((((divisorHigh % 10UL) * 6UL) + (divisorLow % 10UL)) % 10UL);
remainder8 = (byte)(divisorLow % 8UL);
remainder3 = (byte)(((divisorHigh % 3UL) + (divisorLow % 3UL)) % 3UL);
remainder7 = (byte)((((divisorHigh % 7UL) * 2UL) + (divisorLow % 7UL)) % 7UL);
remainder11 = (byte)((((divisorHigh % 11UL) * 5UL) + (divisorLow % 11UL)) % 11UL);
var gpu = HeuristicCombinedPrimeTesterAccelerator.Rent(1);
while (divisor.CompareTo(limit) <= 0)
{
bool passesSmallModuli = remainder3 != 0 && remainder7 != 0 && remainder11 != 0;
if (passesSmallModuli && (remainder8 == 1 || remainder8 == 7) && ((decimalMask >> remainder10) & 1) != 0)
{
ulong candidate = divisor.Low;
MontgomeryDivisorData divisorData = MontgomeryDivisorData.FromModulus(candidate);
ulong divisorCycle;
// Divisors generated from 2 * k * p + 1 exceed the small-cycle snapshot when p >= 138,000,000, so the short path below never runs.
// if (candidate <= PerfectNumberConstants.MaxQForDivisorCycles)
// {
// divisorCycle = cycleCache.GetCycleLength(candidate);
// }
// else
{
if (!MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(
gpu,
candidate,
prime,
divisorData,
out ulong computedCycle,
out bool primeOrderFailed) || computedCycle == 0UL)
{
// Divisors produced by 2 * k * p + 1 always exceed PerfectNumberConstants.MaxQForDivisorCycles
// for the exponents scanned here, so skip the unused cache fallback and compute directly.
divisorCycle = MersenneDivisorCycles.CalculateCycleLength(
candidate,
divisorData,
skipPrimeOrderHeuristic: primeOrderFailed);
}
else
{
divisorCycle = computedCycle;
}
}
if (divisorCycle == prime)
{
// A cycle equal to the tested exponent (which is prime in this path) guarantees that the candidate divides
// the corresponding Mersenne number because the order of 2 modulo the divisor is exactly p.
processedAll = true;
HeuristicCombinedPrimeTesterAccelerator.Return(gpu);
return true;
}
if (divisorCycle == 0UL)
{
Console.WriteLine($"Divisor cycle was not calculated for {prime}");
}
}
divisor.Add(step);
remainder10 = AddMod10(remainder10, step10);
remainder8 = AddMod8(remainder8, step8);
remainder3 = AddMod3(remainder3, step3);
remainder7 = AddMod7(remainder7, step7);
remainder11 = AddMod11(remainder11, step11);
}
HeuristicCombinedPrimeTesterAccelerator.Return(gpu);
processedAll = true;
return false;
}
private static bool CheckDivisors64(
ulong prime,
ulong step,
ulong limit,
ulong divisor,
ushort decimalMask,
byte step10,
byte step8,
byte step3,
byte step7,
byte step11,
byte remainder10,
byte remainder8,
byte remainder3,
byte remainder7,
byte remainder11,
out bool processedAll)
{
processedAll = false;
bool canAdvance = step <= limit;
var gpu = HeuristicCombinedPrimeTesterAccelerator.Rent(1);
if (!canAdvance)
{
if (divisor <= limit)
{
if (remainder3 != 0 && remainder7 != 0 && remainder11 != 0 && (remainder8 == 1 || remainder8 == 7) && ((decimalMask >> remainder10) & 1) != 0)
{
MontgomeryDivisorData divisorData = MontgomeryDivisorData.FromModulus(divisor);
ulong divisorCycle;
if (!MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(
gpu,
divisor,
prime,
divisorData,
out ulong computedCycle,
out bool primeOrderFailed) || computedCycle == 0UL)
{
divisorCycle = MersenneDivisorCycles.CalculateCycleLength(
divisor,
divisorData,
skipPrimeOrderHeuristic: primeOrderFailed);
}
else
{
divisorCycle = computedCycle;
}
if (divisorCycle == prime)
{
processedAll = true;
HeuristicCombinedPrimeTesterAccelerator.Return(gpu);
return true;
}
if (divisorCycle == 0UL)
{
Console.WriteLine($"Divisor cycle was not calculated for {prime}");
}
}
}
processedAll = true;
return false;
}
if (divisor > limit)
{
processedAll = true;
return false;
}
ulong remainingIterations = ((limit - divisor) / step) + 1UL;
while (true)
{
bool passesSmallModuli = remainder3 != 0 && remainder7 != 0 && remainder11 != 0;
if (passesSmallModuli && (remainder8 == 1 || remainder8 == 7) && ((decimalMask >> remainder10) & 1) != 0)
{
MontgomeryDivisorData divisorData = MontgomeryDivisorData.FromModulus(divisor);
ulong divisorCycle;
if (!MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(
gpu,
divisor,
prime,
divisorData,
out ulong computedCycle,
out bool primeOrderFailed) || computedCycle == 0UL)
{
divisorCycle = MersenneDivisorCycles.CalculateCycleLength(
divisor,
divisorData,
skipPrimeOrderHeuristic: primeOrderFailed);
}
else
{
divisorCycle = computedCycle;
}
if (divisorCycle == prime)
{
processedAll = true;
HeuristicCombinedPrimeTesterAccelerator.Return(gpu);
return true;
}
if (divisorCycle == 0UL)
{
Console.WriteLine($"Divisor cycle was not calculated for {prime}");
}
}
remainingIterations--;
if (remainingIterations == 0UL)
{
processedAll = true;
return false;
}
divisor += step;
remainder10 = AddMod10(remainder10, step10);
remainder8 = AddMod8(remainder8, step8);
remainder3 = AddMod3(remainder3, step3);
remainder7 = AddMod7(remainder7, step7);
remainder11 = AddMod11(remainder11, step11);
}
}
private static byte CheckDivisor(ulong prime, ulong divisorCycle, in MontgomeryDivisorData divisorData)
{
ulong residue = prime.Pow2MontgomeryModWithCycleCpu(divisorCycle, divisorData);
return residue == 1UL ? (byte)1 : (byte)0;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static byte AddMod3(byte value, byte delta)
{
const int Modulus = 3;
int sum = value + delta;
if (sum >= Modulus)
{
sum -= Modulus;
}
return (byte)sum;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static byte AddMod7(byte value, byte delta)
{
const int Modulus = 7;
int sum = value + delta;
if (sum >= Modulus)
{
sum -= Modulus;
}
return (byte)sum;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static byte AddMod10(byte value, byte delta)
{
const int Modulus = 10;
int sum = value + delta;
if (sum >= Modulus)
{
sum -= Modulus;
}
return (byte)sum;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static byte AddMod11(byte value, byte delta)
{
const int Modulus = 11;
int sum = value + delta;
if (sum >= Modulus)
{
sum -= Modulus;
}
return (byte)sum;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static byte AddMod8(byte value, byte delta)
{
return (byte)((value + delta) & 7);
}
private static ulong ComputeDivisorLimitFromMaxPrime(ulong maxPrime)
{
// The by-divisor CPU configuration only supplies primes greater than 1, so the guard below never trips.
// if (maxPrime <= 1UL)
// {
// return 0UL;
// }
if (maxPrime - 1UL >= 64UL)
{
return ulong.MaxValue;
}
return (1UL << (int)(maxPrime - 1UL)) - 1UL;
}
private static ulong ComputeAllowedMaxDivisor(ulong prime, ulong divisorLimit)
{
// Production by-divisor scans only handle primes, so inputs never fall below 2.
// if (prime <= 1UL)
// {
// return 0UL;
// }
if (prime - 1UL >= 64UL)
{
return divisorLimit;
}
return Math.Min((1UL << (int)(prime - 1UL)) - 1UL, divisorLimit);
}
}
PerfectNumbers.Core/Gpu/Accelerators/AcceleratorExtensions.cs
+34
View File
@@ -0,0 +1,34 @@
using ILGPU;
using ILGPU.Runtime;
namespace PerfectNumbers.Core.Gpu.Accelerators;
public static class AcceleratorExtensions
{
public static MemoryBuffer1D<ulong, Stride1D.Dense> CopySpanToDevice(this Accelerator accelerator, AcceleratorStream stream, in ulong[] span)
{
MemoryBuffer1D<ulong, Stride1D.Dense>? buffer;
buffer = accelerator.Allocate1D(stream, span);
return buffer;
}
public static MemoryBuffer1D<ulong, Stride1D.Dense> CopyUintSpanToDevice(this Accelerator accelerator, AcceleratorStream stream, ReadOnlySpan<uint> span)
{
MemoryBuffer1D<ulong, Stride1D.Dense>? buffer;
buffer = accelerator.Allocate1D<ulong>(span.Length);
if (!span.IsEmpty)
{
var converted = new ulong[span.Length];
for (int i = 0; i < span.Length; i++)
{
converted[i] = span[i];
}
buffer.View.CopyFromCPU(stream, converted);
}
return buffer;
}
}
PerfectNumbers.Core/Gpu/Accelerators/HeuristicCombinedPrimeTesterAccelerator.cs
+147
View File
@@ -0,0 +1,147 @@
using System.Runtime.CompilerServices;
using ILGPU;
using ILGPU.Runtime;
namespace PerfectNumbers.Core.Gpu.Accelerators;
public sealed class HeuristicCombinedPrimeTesterAccelerator
{
#region Pool
[ThreadStatic]
private static Queue<HeuristicCombinedPrimeTesterAccelerator>? Pool;
public static void WarmUp(int acceleratorIndex, AcceleratorStream stream)
{
Console.WriteLine($"Preparing heuristic accelerator {acceleratorIndex}...");
HeuristicGpuCombinedTables.WarmUp(acceleratorIndex, stream);
// _ = GpuKernelPool.GetOrAddKernels(accelerator, stream, KernelType.None);
// KernelContainer kernels = GpuKernelPool.GetOrAddKernels(accelerator, stream);
// GpuStaticTableInitializer.EnsureStaticTables(accelerator, kernels, stream);
Console.WriteLine($"Heuristic accelerator {acceleratorIndex} is ready");
}
public static HeuristicCombinedPrimeTesterAccelerator Rent(int minBufferCapacity = 1)
{
var pool = Pool ??= new();
if (!pool.TryDequeue(out var lease))
{
lease = new HeuristicCombinedPrimeTesterAccelerator(minBufferCapacity);
}
else
{
lease.EnsureCapacity(minBufferCapacity);
}
return lease;
}
public static void Return(HeuristicCombinedPrimeTesterAccelerator lease)
{
var pool = Pool ??= new();
pool.Enqueue(lease);
}
public static void DisposeAll()
{
if (Pool != null)
{
while (Pool.TryDequeue(out var lease))
{
lease.DisposeResources();
}
}
}
public static void Clear()
{
if (Pool != null)
{
var retained = new List<HeuristicCombinedPrimeTesterAccelerator>();
while (Pool.TryDequeue(out var lease))
{
lease.DisposeResources();
}
}
}
#endregion
public readonly HeuristicGpuCombinedDivisorTables DivisorTables;
internal readonly Context Context;
public readonly Accelerator Accelerator;
public readonly Kernel HeuristicCombinedTrialDivisionKernel;
public MemoryBuffer1D<ulong, Stride1D.Dense>? Input = null;
public MemoryBuffer1D<byte, Stride1D.Dense>? Output = null;
public MemoryBuffer1D<int, Stride1D.Dense>? HeuristicFlag = null!;
public int BufferCapacity;
private HeuristicGpuCombinedTables _sharedTables;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding1 => _sharedTables.HeuristicCombinedDivisorsEnding1;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding1 => _sharedTables.HeuristicCombinedDivisorSquaresEnding1;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding3 => _sharedTables.HeuristicCombinedDivisorsEnding3;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding3 => _sharedTables.HeuristicCombinedDivisorSquaresEnding3;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding7 => _sharedTables.HeuristicCombinedDivisorsEnding7;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding7 => _sharedTables.HeuristicCombinedDivisorSquaresEnding7;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding9 => _sharedTables.HeuristicCombinedDivisorsEnding9;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding9 => _sharedTables.HeuristicCombinedDivisorSquaresEnding9;
private static readonly Accelerator[] _accelerators = AcceleratorPool.Shared.Accelerators;
internal readonly int AcceleratorIndex;
public HeuristicCombinedPrimeTesterAccelerator(int minBufferCapacity)
{
var acceleratorIndex = AcceleratorPool.Shared.Rent();
AcceleratorIndex = acceleratorIndex;
Accelerator accelerator = _accelerators[acceleratorIndex];
Context = accelerator.Context;
Accelerator = accelerator;
EnsureCapacity(minBufferCapacity);
HeuristicCombinedTrialDivisionKernel = KernelUtil.GetKernel(accelerator.LoadAutoGroupedStreamKernel<Index1D, byte, ArrayView<int>, ulong, ulong, HeuristicGpuCombinedDivisorTables>(PrimeTesterKernels.HeuristicTrialCombinedDivisionKernel));
AcceleratorStream stream = AcceleratorStreamPool.Rent(acceleratorIndex);
// GpuStaticTableInitializer.EnsureStaticTables(accelerator, kernels, stream);
_sharedTables = HeuristicGpuCombinedTables.GetStaticTables(acceleratorIndex);
DivisorTables = _sharedTables.CreateHeuristicDivisorTables();
stream.Synchronize();
AcceleratorStreamPool.Return(acceleratorIndex, stream);
}
public void EnsureCapacity(int minCapacity)
{
if (minCapacity <= BufferCapacity)
{
return;
}
Input?.Dispose();
Output?.Dispose();
HeuristicFlag?.Dispose();
Input = Accelerator.Allocate1D<ulong>(minCapacity);
Output = Accelerator.Allocate1D<byte>(minCapacity);
HeuristicFlag = Accelerator.Allocate1D<int>(minCapacity);
BufferCapacity = minCapacity;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public void Return() => Return(this);
public void DisposeResources()
{
Input?.Dispose();
Output?.Dispose();
HeuristicFlag?.Dispose();
// ReleaseSharedTables(_sharedTables);
// These resources are shared between GPU leases
// Accelerator.Dispose();
// AcceleratorContext.Dispose();
}
}
PerfectNumbers.Core/Gpu/Accelerators/HeuristicGpuCombinedTables.cs
+69
View File
@@ -0,0 +1,69 @@
using ILGPU;
using ILGPU.Runtime;
namespace PerfectNumbers.Core.Gpu.Accelerators;
internal sealed class HeuristicGpuCombinedTables
{
private static readonly Accelerator[] _accelerators = AcceleratorPool.Shared.Accelerators;
private static readonly HeuristicGpuCombinedTables[] _sharedTables = new HeuristicGpuCombinedTables[PerfectNumberConstants.RollingAccelerators];
internal static void WarmUp(int acceleratorIndex, AcceleratorStream stream) => _sharedTables[acceleratorIndex] = new(_accelerators[acceleratorIndex], stream);
internal static HeuristicGpuCombinedTables GetStaticTables(int acceleratorIndex) => _sharedTables[acceleratorIndex];
internal HeuristicGpuCombinedTables(Accelerator accelerator, AcceleratorStream stream)
{
var combinedEnding1 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding1;
HeuristicCombinedDivisorsEnding1 = accelerator.CopySpanToDevice(stream, combinedEnding1);
var combinedSquaresEnding1 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding1Squares;
HeuristicCombinedDivisorSquaresEnding1 = accelerator.CopySpanToDevice(stream, combinedSquaresEnding1);
var combinedEnding3 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding3;
HeuristicCombinedDivisorsEnding3 = accelerator.CopySpanToDevice(stream, combinedEnding3);
var combinedSquaresEnding3 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding3Squares;
HeuristicCombinedDivisorSquaresEnding3 = accelerator.CopySpanToDevice(stream, combinedSquaresEnding3);
var combinedEnding7 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding7;
HeuristicCombinedDivisorsEnding7 = accelerator.CopySpanToDevice(stream, combinedEnding7);
var combinedSquaresEnding7 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding7Squares;
HeuristicCombinedDivisorSquaresEnding7 = accelerator.CopySpanToDevice(stream, combinedSquaresEnding7);
var combinedEnding9 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding9;
HeuristicCombinedDivisorsEnding9 = accelerator.CopySpanToDevice(stream, combinedEnding9);
var combinedSquaresEnding9 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding9Squares;
HeuristicCombinedDivisorSquaresEnding9 = accelerator.CopySpanToDevice(stream, combinedSquaresEnding9);
}
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding1;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding1;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding3;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding3;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding7;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding7;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding9;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding9;
internal HeuristicGpuCombinedDivisorTables CreateHeuristicDivisorTables() => new(
HeuristicCombinedDivisorsEnding1.View,
HeuristicCombinedDivisorsEnding3.View,
HeuristicCombinedDivisorsEnding7.View,
HeuristicCombinedDivisorsEnding9.View,
HeuristicCombinedDivisorSquaresEnding1.View,
HeuristicCombinedDivisorSquaresEnding3.View,
HeuristicCombinedDivisorSquaresEnding7.View,
HeuristicCombinedDivisorSquaresEnding9.View);
internal void Dispose()
{
HeuristicCombinedDivisorsEnding1.Dispose();
HeuristicCombinedDivisorSquaresEnding1.Dispose();
HeuristicCombinedDivisorsEnding3.Dispose();
HeuristicCombinedDivisorSquaresEnding3.Dispose();
HeuristicCombinedDivisorsEnding7.Dispose();
HeuristicCombinedDivisorSquaresEnding7.Dispose();
HeuristicCombinedDivisorsEnding9.Dispose();
HeuristicCombinedDivisorSquaresEnding9.Dispose();
}
}
PerfectNumbers.Core/Gpu/Accelerators/HeuristicGpuGroupABTables.cs
+79
View File
@@ -0,0 +1,79 @@
using ILGPU;
using ILGPU.Runtime;
namespace PerfectNumbers.Core.Gpu.Accelerators;
internal sealed class HeuristicGpuGroupABTables
{
private static readonly Accelerator[] _accelerators = AcceleratorPool.Shared.Accelerators;
private static readonly HeuristicGpuGroupABTables[] _sharedTables = new HeuristicGpuGroupABTables[PerfectNumberConstants.RollingAccelerators];
private static readonly SemaphoreSlim[] _locks = [..Enumerable.Range(1, PerfectNumberConstants.RollingAccelerators).Select(_ => new SemaphoreSlim(1))];
internal static HeuristicGpuGroupABTables EnsureStaticTables(int acceleratorIndex, AcceleratorStream stream)
{
var @lock = _locks[acceleratorIndex];
@lock.Wait();
if (_sharedTables[acceleratorIndex] is {} existing)
{
@lock.Release();
return existing;
}
var accelerator = _accelerators[acceleratorIndex];
existing = new(accelerator, stream);
_sharedTables[acceleratorIndex] = existing;
@lock.Release();
return existing;
}
internal HeuristicGpuGroupABTables(Accelerator accelerator, AcceleratorStream stream)
{
var heuristicGroupA = HeuristicPrimeSieves.GroupADivisorsStorage;
HeuristicGroupADivisors = accelerator.CopySpanToDevice(stream, heuristicGroupA);
var heuristicGroupASquares = HeuristicPrimeSieves.GroupADivisorSquaresStorage;
HeuristicGroupADivisorSquares = accelerator.CopySpanToDevice(stream, heuristicGroupASquares);
HeuristicGroupBDivisorsEnding1 = accelerator.CopyUintSpanToDevice(stream, DivisorGenerator.SmallPrimesLastOneWithoutLastThree);
HeuristicGroupBDivisorSquaresEnding1 = accelerator.CopySpanToDevice(stream, DivisorGenerator.SmallPrimesPow2LastOneWithoutLastThree);
HeuristicGroupBDivisorsEnding7 = accelerator.CopyUintSpanToDevice(stream, DivisorGenerator.SmallPrimesLastSevenWithoutLastThree);
HeuristicGroupBDivisorSquaresEnding7 = accelerator.CopySpanToDevice(stream, DivisorGenerator.SmallPrimesPow2LastSevenWithoutLastThree);
HeuristicGroupBDivisorsEnding9 = accelerator.CopyUintSpanToDevice(stream, DivisorGenerator.SmallPrimesLastNineWithoutLastThree);
HeuristicGroupBDivisorSquaresEnding9 = accelerator.CopySpanToDevice(stream, DivisorGenerator.SmallPrimesPow2LastNineWithoutLastThree);
}
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupADivisors;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupADivisorSquares;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorsEnding1;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorSquaresEnding1;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorsEnding7;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorSquaresEnding7;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorsEnding9;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorSquaresEnding9;
internal HeuristicGpuGroupABDivisorTables CreateHeuristicDivisorTables() => new(
HeuristicGroupADivisors.View,
HeuristicGroupADivisorSquares.View,
HeuristicGroupBDivisorsEnding1.View,
HeuristicGroupBDivisorSquaresEnding1.View,
HeuristicGroupBDivisorsEnding7.View,
HeuristicGroupBDivisorSquaresEnding7.View,
HeuristicGroupBDivisorsEnding9.View,
HeuristicGroupBDivisorSquaresEnding9.View);
internal void Dispose()
{
HeuristicGroupADivisors.Dispose();
HeuristicGroupADivisorSquares.Dispose();
HeuristicGroupBDivisorsEnding1.Dispose();
HeuristicGroupBDivisorSquaresEnding1.Dispose();
HeuristicGroupBDivisorsEnding7.Dispose();
HeuristicGroupBDivisorSquaresEnding7.Dispose();
HeuristicGroupBDivisorsEnding9.Dispose();
HeuristicGroupBDivisorSquaresEnding9.Dispose();
}
}
PerfectNumbers.Core/Gpu/Accelerators/HeuristicGpuTables.cs
-150
View File
@@ -1,150 +0,0 @@
using ILGPU;
using ILGPU.Runtime;
namespace PerfectNumbers.Core.Gpu.Accelerators;
internal sealed class HeuristicGpuTables
{
private static readonly Accelerator[] _accelerators = AcceleratorPool.Shared.Accelerators;
private static readonly HeuristicGpuTables[] _sharedTables = new HeuristicGpuTables[PerfectNumberConstants.RollingAccelerators];
private static readonly SemaphoreSlim[] _locks = [..Enumerable.Range(1, PerfectNumberConstants.RollingAccelerators).Select(_ => new SemaphoreSlim(1))];
internal static HeuristicGpuTables EnsureStaticTables(int acceleratorIndex, AcceleratorStream stream)
{
var @lock = _locks[acceleratorIndex];
@lock.Wait();
if (_sharedTables[acceleratorIndex] is {} existing)
{
@lock.Release();
return existing;
}
var accelerator = _accelerators[acceleratorIndex];
existing = new(accelerator, stream);
_sharedTables[acceleratorIndex] = existing;
@lock.Release();
return existing;
}
internal HeuristicGpuTables(Accelerator accelerator, AcceleratorStream stream)
{
var heuristicGroupA = HeuristicPrimeSieves.GroupADivisorsStorage;
HeuristicGroupADivisors = CopySpanToDevice(accelerator, stream, heuristicGroupA);
var heuristicGroupASquares = HeuristicPrimeSieves.GroupADivisorSquaresStorage;
HeuristicGroupADivisorSquares = CopySpanToDevice(accelerator, stream, heuristicGroupASquares);
HeuristicGroupBDivisorsEnding1 = CopyUintSpanToDevice(accelerator, stream, DivisorGenerator.SmallPrimesLastOneWithoutLastThree);
HeuristicGroupBDivisorSquaresEnding1 = CopySpanToDevice(accelerator, stream, DivisorGenerator.SmallPrimesPow2LastOneWithoutLastThree);
HeuristicGroupBDivisorsEnding7 = CopyUintSpanToDevice(accelerator, stream, DivisorGenerator.SmallPrimesLastSevenWithoutLastThree);
HeuristicGroupBDivisorSquaresEnding7 = CopySpanToDevice(accelerator, stream, DivisorGenerator.SmallPrimesPow2LastSevenWithoutLastThree);
HeuristicGroupBDivisorsEnding9 = CopyUintSpanToDevice(accelerator, stream, DivisorGenerator.SmallPrimesLastNineWithoutLastThree);
HeuristicGroupBDivisorSquaresEnding9 = CopySpanToDevice(accelerator, stream, DivisorGenerator.SmallPrimesPow2LastNineWithoutLastThree);
var combinedEnding1 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding1;
HeuristicCombinedDivisorsEnding1 = CopySpanToDevice(accelerator, stream, combinedEnding1);
var combinedSquaresEnding1 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding1Squares;
HeuristicCombinedDivisorSquaresEnding1 = CopySpanToDevice(accelerator, stream, combinedSquaresEnding1);
var combinedEnding3 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding3;
HeuristicCombinedDivisorsEnding3 = CopySpanToDevice(accelerator, stream, combinedEnding3);
var combinedSquaresEnding3 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding3Squares;
HeuristicCombinedDivisorSquaresEnding3 = CopySpanToDevice(accelerator, stream, combinedSquaresEnding3);
var combinedEnding7 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding7;
HeuristicCombinedDivisorsEnding7 = CopySpanToDevice(accelerator, stream, combinedEnding7);
var combinedSquaresEnding7 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding7Squares;
HeuristicCombinedDivisorSquaresEnding7 = CopySpanToDevice(accelerator, stream, combinedSquaresEnding7);
var combinedEnding9 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding9;
HeuristicCombinedDivisorsEnding9 = CopySpanToDevice(accelerator, stream, combinedEnding9);
var combinedSquaresEnding9 = HeuristicCombinedPrimeTester.CombinedDivisorsEnding9Squares;
HeuristicCombinedDivisorSquaresEnding9 = CopySpanToDevice(accelerator, stream, combinedSquaresEnding9);
}
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupADivisors;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupADivisorSquares;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorsEnding1;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorSquaresEnding1;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorsEnding7;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorSquaresEnding7;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorsEnding9;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorSquaresEnding9;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding1;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding1;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding3;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding3;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding7;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding7;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding9;
internal readonly MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding9;
internal HeuristicGpuDivisorTables CreateHeuristicDivisorTables() => new(
HeuristicCombinedDivisorsEnding1.View,
HeuristicCombinedDivisorsEnding3.View,
HeuristicCombinedDivisorsEnding7.View,
HeuristicCombinedDivisorsEnding9.View,
HeuristicCombinedDivisorSquaresEnding1.View,
HeuristicCombinedDivisorSquaresEnding3.View,
HeuristicCombinedDivisorSquaresEnding7.View,
HeuristicCombinedDivisorSquaresEnding9.View,
HeuristicGroupADivisors.View,
HeuristicGroupADivisorSquares.View,
HeuristicGroupBDivisorsEnding1.View,
HeuristicGroupBDivisorSquaresEnding1.View,
HeuristicGroupBDivisorsEnding7.View,
HeuristicGroupBDivisorSquaresEnding7.View,
HeuristicGroupBDivisorsEnding9.View,
HeuristicGroupBDivisorSquaresEnding9.View);
private static MemoryBuffer1D<ulong, Stride1D.Dense> CopySpanToDevice(Accelerator accelerator, AcceleratorStream stream, in ulong[] span)
{
MemoryBuffer1D<ulong, Stride1D.Dense>? buffer;
buffer = accelerator.Allocate1D(stream, span);
return buffer;
}
private static MemoryBuffer1D<ulong, Stride1D.Dense> CopyUintSpanToDevice(Accelerator accelerator, AcceleratorStream stream, ReadOnlySpan<uint> span)
{
MemoryBuffer1D<ulong, Stride1D.Dense>? buffer;
buffer = accelerator.Allocate1D<ulong>(span.Length);
if (!span.IsEmpty)
{
var converted = new ulong[span.Length];
for (int i = 0; i < span.Length; i++)
{
converted[i] = span[i];
}
buffer.View.CopyFromCPU(stream, converted);
}
return buffer;
}
internal void Dispose()
{
HeuristicGroupADivisors.Dispose();
HeuristicGroupADivisorSquares.Dispose();
HeuristicGroupBDivisorsEnding1.Dispose();
HeuristicGroupBDivisorSquaresEnding1.Dispose();
HeuristicGroupBDivisorsEnding7.Dispose();
HeuristicGroupBDivisorSquaresEnding7.Dispose();
HeuristicGroupBDivisorsEnding9.Dispose();
HeuristicGroupBDivisorSquaresEnding9.Dispose();
HeuristicCombinedDivisorsEnding1.Dispose();
HeuristicCombinedDivisorSquaresEnding1.Dispose();
HeuristicCombinedDivisorsEnding3.Dispose();
HeuristicCombinedDivisorSquaresEnding3.Dispose();
HeuristicCombinedDivisorsEnding7.Dispose();
HeuristicCombinedDivisorSquaresEnding7.Dispose();
HeuristicCombinedDivisorsEnding9.Dispose();
HeuristicCombinedDivisorSquaresEnding9.Dispose();
}
}
PerfectNumbers.Core/Gpu/Accelerators/HeuristicPrimeTesterAccelerator.cs → PerfectNumbers.Core/Gpu/Accelerators/HeuristicGroupABPrimeTesterAccelerator.cs
+20 -35
View File
@@ -1,15 +1,14 @@
using System.Runtime.CompilerServices;
using ILGPU;
using ILGPU.Runtime;
namespace PerfectNumbers.Core.Gpu.Accelerators;
internal sealed class HeuristicPrimeTesterAccelerator
internal sealed class HeuristicGroupABPrimeTesterAccelerator
{
#region Pool
[ThreadStatic]
private static Queue<HeuristicPrimeTesterAccelerator>? Pool;
private static Queue<HeuristicGroupABPrimeTesterAccelerator>? Pool;
internal static void WarmUp()
{
@@ -20,23 +19,22 @@ internal sealed class HeuristicPrimeTesterAccelerator
Console.WriteLine($"Preparing heuristic accelerator {i}...");
var accelerator = accelerators[i];
AcceleratorStream stream = accelerator.CreateStream();
Accelerators.HeuristicGpuTables.EnsureStaticTables(i, stream);
_ = HeuristicGpuGroupABTables.EnsureStaticTables(i, stream);
// _ = GpuKernelPool.GetOrAddKernels(accelerator, stream, KernelType.None);
// KernelContainer kernels = GpuKernelPool.GetOrAddKernels(accelerator, stream);
// GpuStaticTableInitializer.EnsureStaticTables(accelerator, kernels, stream);
stream.Synchronize();
stream.Dispose();
accelerator.Synchronize();
Console.WriteLine($"Heuristic accelerator {i} is ready");
}
}
internal static HeuristicPrimeTesterAccelerator Rent(int minBufferCapacity = 1)
internal static HeuristicGroupABPrimeTesterAccelerator Rent(int minBufferCapacity = 1)
{
var pool = Pool ??= new();
if (!pool.TryDequeue(out var lease))
{
lease = new HeuristicPrimeTesterAccelerator(minBufferCapacity);
lease = new HeuristicGroupABPrimeTesterAccelerator(minBufferCapacity);
}
else
{
@@ -46,7 +44,7 @@ internal sealed class HeuristicPrimeTesterAccelerator
return lease;
}
internal static void Return(HeuristicPrimeTesterAccelerator lease)
internal static void Return(HeuristicGroupABPrimeTesterAccelerator lease)
{
var pool = Pool ??= new();
pool.Enqueue(lease);
@@ -67,7 +65,7 @@ internal sealed class HeuristicPrimeTesterAccelerator
{
if (Pool != null)
{
var retained = new List<HeuristicPrimeTesterAccelerator>();
var retained = new List<HeuristicGroupABPrimeTesterAccelerator>();
while (Pool.TryDequeue(out var lease))
{
lease.DisposeResources();
@@ -77,10 +75,9 @@ internal sealed class HeuristicPrimeTesterAccelerator
#endregion
private readonly HeuristicGpuDivisorTables _heuristicDivisorTables;
public readonly HeuristicGpuGroupABDivisorTables DivisorTables;
internal readonly Context Context;
public readonly Accelerator Accelerator;
public AcceleratorStream? Stream;
[ThreadStatic]
public static Dictionary<Accelerator, Action<AcceleratorStream, Index1D, ArrayView<ulong>, ArrayView<uint>, ArrayView<uint>, ArrayView<uint>, ArrayView<uint>, ArrayView<ulong>, ArrayView<ulong>, ArrayView<ulong>, ArrayView<ulong>, ArrayView<byte>>>? SmallPrimeSieveKernel;
@@ -100,12 +97,12 @@ internal sealed class HeuristicPrimeTesterAccelerator
public readonly Kernel HeuristicTrialDivisionKernel;
public MemoryBuffer1D<ulong, Stride1D.Dense> Input = null!;
public MemoryBuffer1D<byte, Stride1D.Dense> Output = null!;
public MemoryBuffer1D<int, Stride1D.Dense> HeuristicFlag = null!;
public MemoryBuffer1D<ulong, Stride1D.Dense>? Input = null;
public MemoryBuffer1D<byte, Stride1D.Dense>? Output = null;
public MemoryBuffer1D<int, Stride1D.Dense>? HeuristicFlag = null!;
public int BufferCapacity;
private HeuristicGpuTables _sharedTables;
private HeuristicGpuGroupABTables _sharedTables;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupADivisors => _sharedTables.HeuristicGroupADivisors;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupADivisorSquares => _sharedTables.HeuristicGroupADivisorSquares;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorsEnding1 => _sharedTables.HeuristicGroupBDivisorsEnding1;
@@ -114,36 +111,28 @@ internal sealed class HeuristicPrimeTesterAccelerator
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorSquaresEnding7 => _sharedTables.HeuristicGroupBDivisorSquaresEnding7;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorsEnding9 => _sharedTables.HeuristicGroupBDivisorsEnding9;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicGroupBDivisorSquaresEnding9 => _sharedTables.HeuristicGroupBDivisorSquaresEnding9;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding1 => _sharedTables.HeuristicCombinedDivisorsEnding1;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding1 => _sharedTables.HeuristicCombinedDivisorSquaresEnding1;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding3 => _sharedTables.HeuristicCombinedDivisorsEnding3;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding3 => _sharedTables.HeuristicCombinedDivisorSquaresEnding3;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding7 => _sharedTables.HeuristicCombinedDivisorsEnding7;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding7 => _sharedTables.HeuristicCombinedDivisorSquaresEnding7;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorsEnding9 => _sharedTables.HeuristicCombinedDivisorsEnding9;
public MemoryBuffer1D<ulong, Stride1D.Dense> HeuristicCombinedDivisorSquaresEnding9 => _sharedTables.HeuristicCombinedDivisorSquaresEnding9;
public HeuristicGpuDivisorTables HeuristicGpuTables => _heuristicDivisorTables;
private static readonly Accelerator[] _accelerators = AcceleratorPool.Shared.Accelerators;
internal readonly int AcceleratorIndex;
internal HeuristicPrimeTesterAccelerator(int minBufferCapacity)
internal HeuristicGroupABPrimeTesterAccelerator(int minBufferCapacity)
{
var acceleratorIndex = AcceleratorPool.Shared.Rent();
AcceleratorIndex = acceleratorIndex;
Accelerator accelerator = _accelerators[acceleratorIndex];
Context = accelerator.Context;
Accelerator = accelerator;
EnsureCapacity(minBufferCapacity);
HeuristicTrialDivisionKernel = KernelUtil.GetKernel(accelerator.LoadAutoGroupedStreamKernel<Index1D, ArrayView<int>, ulong, ulong, HeuristicGpuDivisorTableKind, HeuristicGpuDivisorTables>(PrimeTesterKernels.HeuristicTrialDivisionKernel));
HeuristicTrialDivisionKernel = KernelUtil.GetKernel(accelerator.LoadAutoGroupedStreamKernel<Index1D, ArrayView<int>, ulong, ulong, HeuristicGpuDivisorTableKind, HeuristicGpuGroupABDivisorTables>(PrimeTesterKernels.HeuristicTrialGroupABDivisionKernel));
AcceleratorStream stream = AcceleratorStreamPool.Rent(acceleratorIndex);
// GpuStaticTableInitializer.EnsureStaticTables(accelerator, kernels, stream);
_sharedTables = Accelerators.HeuristicGpuTables.EnsureStaticTables(acceleratorIndex, stream);
_heuristicDivisorTables = _sharedTables.CreateHeuristicDivisorTables();
_sharedTables = HeuristicGpuGroupABTables.EnsureStaticTables(acceleratorIndex, stream);
DivisorTables = _sharedTables.CreateHeuristicDivisorTables();
stream.Synchronize();
AcceleratorStreamPool.Return(acceleratorIndex, stream);
@@ -156,8 +145,8 @@ internal sealed class HeuristicPrimeTesterAccelerator
return;
}
Input.Dispose();
Output.Dispose();
Input?.Dispose();
Output?.Dispose();
HeuristicFlag?.Dispose();
Input = Accelerator.Allocate1D<ulong>(minCapacity);
@@ -169,8 +158,6 @@ internal sealed class HeuristicPrimeTesterAccelerator
public void Dispose()
{
Stream?.Dispose();
Stream = null;
Return(this);
}
@@ -179,8 +166,6 @@ internal sealed class HeuristicPrimeTesterAccelerator
Input?.Dispose();
Output?.Dispose();
HeuristicFlag?.Dispose();
Stream?.Dispose();
Stream = null;
// ReleaseSharedTables(_sharedTables);
// These resources are shared between GPU leases
// Accelerator.Dispose();
PerfectNumbers.Core/Gpu/Accelerators/PrimeOrderCalculatorAccelerator.cs
+2
View File
@@ -68,6 +68,7 @@ public sealed class PrimeOrderCalculatorAccelerator
// Don't take this from the pool as quick uploads of data to the accelerator consumes much of GPU's memory and throws.
AcceleratorStream stream = accelerator.CreateStream();
LastDigitGpuTables.WarmUp(i, stream);
HeuristicCombinedPrimeTesterAccelerator.WarmUp(i, stream);
// SharedHeuristicGpuTables.EnsureStaticTables(accelerator, stream);
// _ = GpuKernelPool.GetOrAddKernels(accelerator, stream, KernelType.None);
// KernelContainer kernels = GpuKernelPool.GetOrAddKernels(accelerator, stream);
@@ -76,6 +77,7 @@ public sealed class PrimeOrderCalculatorAccelerator
stream.Dispose();
accelerator.Synchronize();
}
}
#endregion
PerfectNumbers.Core/Gpu/HeuristicGpuCombinedDivisorTables.cs
+45
View File
@@ -0,0 +1,45 @@
using System.Runtime.CompilerServices;
using ILGPU;
using ILGPU.Runtime;
namespace PerfectNumbers.Core.Gpu;
public readonly struct HeuristicGpuCombinedDivisorTables(
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorsEnding1,
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorsEnding3,
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorsEnding7,
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorsEnding9,
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorSquaresEnding1,
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorSquaresEnding3,
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorSquaresEnding7,
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorSquaresEnding9)
{
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorsEnding1 = combinedDivisorsEnding1;
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorsEnding3 = combinedDivisorsEnding3;
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorsEnding7 = combinedDivisorsEnding7;
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorsEnding9 = combinedDivisorsEnding9;
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorSquaresEnding1 = combinedDivisorSquaresEnding1;
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorSquaresEnding3 = combinedDivisorSquaresEnding3;
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorSquaresEnding7 = combinedDivisorSquaresEnding7;
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorSquaresEnding9 = combinedDivisorSquaresEnding9;
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public ArrayView1D<ulong, Stride1D.Dense> SelectDivisors(byte ending) => ending switch
{
1 => CombinedDivisorsEnding1,
3 => CombinedDivisorsEnding3,
7 => CombinedDivisorsEnding7,
9 => CombinedDivisorsEnding9,
_ => ArrayView1D<ulong, Stride1D.Dense>.Empty,
};
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public ArrayView1D<ulong, Stride1D.Dense> SelectDivisorSquares(byte ending) => ending switch
{
1 => CombinedDivisorSquaresEnding1,
3 => CombinedDivisorSquaresEnding3,
7 => CombinedDivisorSquaresEnding7,
9 => CombinedDivisorSquaresEnding9,
_ => ArrayView1D<ulong, Stride1D.Dense>.Empty,
};
}
PerfectNumbers.Core/Gpu/HeuristicGpuDivisorTables.cs
-139
View File
@@ -1,139 +0,0 @@
using System.Runtime.CompilerServices;
using ILGPU;
using ILGPU.Runtime;
namespace PerfectNumbers.Core.Gpu;
public readonly struct HeuristicGpuDivisorTables
{
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorsEnding1;
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorsEnding3;
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorsEnding7;
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorsEnding9;
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorSquaresEnding1;
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorSquaresEnding3;
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorSquaresEnding7;
public readonly ArrayView1D<ulong, Stride1D.Dense> CombinedDivisorSquaresEnding9;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupADivisors;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupADivisorSquares;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupBDivisorsEnding1;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupBDivisorSquaresEnding1;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupBDivisorsEnding7;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupBDivisorSquaresEnding7;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupBDivisorsEnding9;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupBDivisorSquaresEnding9;
public HeuristicGpuDivisorTables(
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorsEnding1,
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorsEnding3,
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorsEnding7,
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorsEnding9,
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorSquaresEnding1,
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorSquaresEnding3,
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorSquaresEnding7,
ArrayView1D<ulong, Stride1D.Dense> combinedDivisorSquaresEnding9,
ArrayView1D<ulong, Stride1D.Dense> groupADivisors,
ArrayView1D<ulong, Stride1D.Dense> groupADivisorSquares,
ArrayView1D<ulong, Stride1D.Dense> groupBDivisorsEnding1,
ArrayView1D<ulong, Stride1D.Dense> groupBDivisorSquaresEnding1,
ArrayView1D<ulong, Stride1D.Dense> groupBDivisorsEnding7,
ArrayView1D<ulong, Stride1D.Dense> groupBDivisorSquaresEnding7,
ArrayView1D<ulong, Stride1D.Dense> groupBDivisorsEnding9,
ArrayView1D<ulong, Stride1D.Dense> groupBDivisorSquaresEnding9)
{
CombinedDivisorsEnding1 = combinedDivisorsEnding1;
CombinedDivisorsEnding3 = combinedDivisorsEnding3;
CombinedDivisorsEnding7 = combinedDivisorsEnding7;
CombinedDivisorsEnding9 = combinedDivisorsEnding9;
CombinedDivisorSquaresEnding1 = combinedDivisorSquaresEnding1;
CombinedDivisorSquaresEnding3 = combinedDivisorSquaresEnding3;
CombinedDivisorSquaresEnding7 = combinedDivisorSquaresEnding7;
CombinedDivisorSquaresEnding9 = combinedDivisorSquaresEnding9;
GroupADivisors = groupADivisors;
GroupADivisorSquares = groupADivisorSquares;
GroupBDivisorsEnding1 = groupBDivisorsEnding1;
GroupBDivisorSquaresEnding1 = groupBDivisorSquaresEnding1;
GroupBDivisorsEnding7 = groupBDivisorsEnding7;
GroupBDivisorSquaresEnding7 = groupBDivisorSquaresEnding7;
GroupBDivisorsEnding9 = groupBDivisorsEnding9;
GroupBDivisorSquaresEnding9 = groupBDivisorSquaresEnding9;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public ArrayView1D<ulong, Stride1D.Dense> SelectCombined(byte ending)
{
return ending switch
{
1 => CombinedDivisorsEnding1,
3 => CombinedDivisorsEnding3,
7 => CombinedDivisorsEnding7,
9 => CombinedDivisorsEnding9,
_ => ArrayView1D<ulong, Stride1D.Dense>.Empty,
};
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public ArrayView1D<ulong, Stride1D.Dense> SelectCombinedSquares(byte ending)
{
return ending switch
{
1 => CombinedDivisorSquaresEnding1,
3 => CombinedDivisorSquaresEnding3,
7 => CombinedDivisorSquaresEnding7,
9 => CombinedDivisorSquaresEnding9,
_ => ArrayView1D<ulong, Stride1D.Dense>.Empty,
};
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public ArrayView1D<ulong, Stride1D.Dense> SelectGroupB(byte ending)
{
return ending switch
{
1 => GroupBDivisorsEnding1,
7 => GroupBDivisorsEnding7,
9 => GroupBDivisorsEnding9,
_ => ArrayView1D<ulong, Stride1D.Dense>.Empty,
};
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public ArrayView1D<ulong, Stride1D.Dense> SelectGroupBSquares(byte ending)
{
return ending switch
{
1 => GroupBDivisorSquaresEnding1,
7 => GroupBDivisorSquaresEnding7,
9 => GroupBDivisorSquaresEnding9,
_ => ArrayView1D<ulong, Stride1D.Dense>.Empty,
};
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public ArrayView1D<ulong, Stride1D.Dense> SelectDivisors(HeuristicGpuDivisorTableKind kind, byte ending)
{
return kind switch
{
HeuristicGpuDivisorTableKind.GroupA => GroupADivisors,
HeuristicGpuDivisorTableKind.GroupBEnding1 => GroupBDivisorsEnding1,
HeuristicGpuDivisorTableKind.GroupBEnding7 => GroupBDivisorsEnding7,
HeuristicGpuDivisorTableKind.GroupBEnding9 => GroupBDivisorsEnding9,
HeuristicGpuDivisorTableKind.Combined => SelectCombined(ending),
_ => ArrayView1D<ulong, Stride1D.Dense>.Empty,
};
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public ArrayView1D<ulong, Stride1D.Dense> SelectDivisorSquares(HeuristicGpuDivisorTableKind kind, byte ending)
{
return kind switch
{
HeuristicGpuDivisorTableKind.GroupA => GroupADivisorSquares,
HeuristicGpuDivisorTableKind.GroupBEnding1 => GroupBDivisorSquaresEnding1,
HeuristicGpuDivisorTableKind.GroupBEnding7 => GroupBDivisorSquaresEnding7,
HeuristicGpuDivisorTableKind.GroupBEnding9 => GroupBDivisorSquaresEnding9,
HeuristicGpuDivisorTableKind.Combined => SelectCombinedSquares(ending),
_ => ArrayView1D<ulong, Stride1D.Dense>.Empty,
};
}
}
PerfectNumbers.Core/Gpu/HeuristicGpuGroupABDivisorTables.cs
+63
View File
@@ -0,0 +1,63 @@
using System.Runtime.CompilerServices;
using ILGPU;
using ILGPU.Runtime;
namespace PerfectNumbers.Core.Gpu;
public readonly struct HeuristicGpuGroupABDivisorTables(
ArrayView1D<ulong, Stride1D.Dense> groupADivisors,
ArrayView1D<ulong, Stride1D.Dense> groupADivisorSquares,
ArrayView1D<ulong, Stride1D.Dense> groupBDivisorsEnding1,
ArrayView1D<ulong, Stride1D.Dense> groupBDivisorSquaresEnding1,
ArrayView1D<ulong, Stride1D.Dense> groupBDivisorsEnding7,
ArrayView1D<ulong, Stride1D.Dense> groupBDivisorSquaresEnding7,
ArrayView1D<ulong, Stride1D.Dense> groupBDivisorsEnding9,
ArrayView1D<ulong, Stride1D.Dense> groupBDivisorSquaresEnding9)
{
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupADivisors = groupADivisors;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupADivisorSquares = groupADivisorSquares;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupBDivisorsEnding1 = groupBDivisorsEnding1;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupBDivisorSquaresEnding1 = groupBDivisorSquaresEnding1;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupBDivisorsEnding7 = groupBDivisorsEnding7;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupBDivisorSquaresEnding7 = groupBDivisorSquaresEnding7;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupBDivisorsEnding9 = groupBDivisorsEnding9;
public readonly ArrayView1D<ulong, Stride1D.Dense> GroupBDivisorSquaresEnding9 = groupBDivisorSquaresEnding9;
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public ArrayView1D<ulong, Stride1D.Dense> SelectGroupB(byte ending) => ending switch
{
1 => GroupBDivisorsEnding1,
7 => GroupBDivisorsEnding7,
9 => GroupBDivisorsEnding9,
_ => ArrayView1D<ulong, Stride1D.Dense>.Empty,
};
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public ArrayView1D<ulong, Stride1D.Dense> SelectGroupBSquares(byte ending) => ending switch
{
1 => GroupBDivisorSquaresEnding1,
7 => GroupBDivisorSquaresEnding7,
9 => GroupBDivisorSquaresEnding9,
_ => ArrayView1D<ulong, Stride1D.Dense>.Empty,
};
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public ArrayView1D<ulong, Stride1D.Dense> SelectDivisors(HeuristicGpuDivisorTableKind kind) => kind switch
{
HeuristicGpuDivisorTableKind.GroupA => GroupADivisors,
HeuristicGpuDivisorTableKind.GroupBEnding1 => GroupBDivisorsEnding1,
HeuristicGpuDivisorTableKind.GroupBEnding7 => GroupBDivisorsEnding7,
HeuristicGpuDivisorTableKind.GroupBEnding9 => GroupBDivisorsEnding9,
_ => ArrayView1D<ulong, Stride1D.Dense>.Empty,
};
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public ArrayView1D<ulong, Stride1D.Dense> SelectDivisorSquares(HeuristicGpuDivisorTableKind kind) => kind switch
{
HeuristicGpuDivisorTableKind.GroupA => GroupADivisorSquares,
HeuristicGpuDivisorTableKind.GroupBEnding1 => GroupBDivisorSquaresEnding1,
HeuristicGpuDivisorTableKind.GroupBEnding7 => GroupBDivisorSquaresEnding7,
HeuristicGpuDivisorTableKind.GroupBEnding9 => GroupBDivisorSquaresEnding9,
_ => ArrayView1D<ulong, Stride1D.Dense>.Empty,
};
}
PerfectNumbers.Core/Gpu/Kernels/PrimeTesterKernels.cs
+37 -6
View File
@@ -11,7 +11,6 @@ public enum HeuristicGpuDivisorTableKind : byte
GroupBEnding1 = 1,
GroupBEnding7 = 7,
GroupBEnding9 = 9,
Combined = 255,
}
internal static class PrimeTesterKernels
@@ -78,17 +77,49 @@ internal static class PrimeTesterKernels
results[index] = result;
}
public static void HeuristicTrialDivisionKernel(
public static void HeuristicTrialCombinedDivisionKernel(
Index1D index,
byte nMod10,
ArrayView<int> resultFlag,
ulong n,
ulong maxDivisorSquare,
HeuristicGpuCombinedDivisorTables tables)
{
// byte nMod10 = (byte)(n % 10UL);
ArrayView1D<ulong, Stride1D.Dense> divisors = tables.SelectDivisors(nMod10);
ArrayView1D<ulong, Stride1D.Dense> divisorSquares = tables.SelectDivisorSquares(nMod10);
int divisorLength = (int)divisorSquares.Length;
int threadIndex = index;
if (threadIndex >= divisorLength)
{
return;
}
ulong divisorSquare = divisorSquares[threadIndex];
if (divisorSquare > maxDivisorSquare)
{
return;
}
ulong divisor = divisors[threadIndex];
divisor = n % divisor;
if (divisor == 0UL)
{
resultFlag[0] = 1;
}
}
public static void HeuristicTrialGroupABDivisionKernel(
Index1D index,
ArrayView<int> resultFlag,
ulong n,
ulong maxDivisorSquare,
HeuristicGpuDivisorTableKind tableKind,
HeuristicGpuDivisorTables tables)
HeuristicGpuGroupABDivisorTables tables)
{
byte nMod10 = (byte)(n % 10UL);
ArrayView1D<ulong, Stride1D.Dense> divisors = tables.SelectDivisors(tableKind, nMod10);
ArrayView1D<ulong, Stride1D.Dense> divisorSquares = tables.SelectDivisorSquares(tableKind, nMod10);
ArrayView1D<ulong, Stride1D.Dense> divisors = tables.SelectDivisors(tableKind);
ArrayView1D<ulong, Stride1D.Dense> divisorSquares = tables.SelectDivisorSquares(tableKind);
int divisorLength = (int)divisorSquares.Length;
int threadIndex = index;
PerfectNumbers.Core/Gpu/MersenneNumberDivisorByDivisorGpuTester.cs
+22 -12
View File
@@ -8,6 +8,7 @@ using System.Threading;
using PerfectNumbers.Core;
using ILGPU;
using ILGPU.Runtime;
using PerfectNumbers.Core.Gpu.Accelerators;
namespace PerfectNumbers.Core.Gpu;
@@ -106,9 +107,12 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
ulong lastProcessed;
// GpuPrimeWorkLimiter();
int acceleratorIndex = AcceleratorPool.Shared.Rent();
var accelerator = _accelerators[acceleratorIndex];
var stream = accelerator.CreateStream();
// int acceleratorIndex = AcceleratorPool.Shared.Rent();
// var accelerator = _accelerators[acceleratorIndex];
var gpu = HeuristicCombinedPrimeTesterAccelerator.Rent(batchCapacity);
var acceleratorIndex = gpu.AcceleratorIndex;
var accelerator = gpu.Accelerator;
var stream = AcceleratorStreamPool.Rent(acceleratorIndex);
// Monitor.Enter(gpuLease.ExecutionLock);
@@ -116,6 +120,7 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
var resources = RentBatchResources(accelerator, batchCapacity);
composite = CheckDivisors(
gpu,
prime,
allowedMax,
kernel,
@@ -137,8 +142,9 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
out coveredRange,
out processedCount);
stream.Dispose();
AcceleratorStreamPool.Return(acceleratorIndex, stream);
ReturnBatchResources(resources);
HeuristicCombinedPrimeTesterAccelerator.Return(gpu);
// Monitor.Exit(gpuLease.ExecutionLock);
// GpuPrimeWorkLimiter.Release();
@@ -182,6 +188,7 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
}
private static bool CheckDivisors(
HeuristicCombinedPrimeTesterAccelerator gpu,
ulong prime,
ulong allowedMax,
Action<AcceleratorStream, Index1D, ArrayView<GpuDivisorPartialData>, ArrayView<int>, ArrayView<int>, ArrayView<ulong>, ArrayView<ulong>, ArrayView<byte>, ArrayView<int>> kernel,
@@ -205,7 +212,7 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
{
int batchCapacity = (int)divisorDataBuffer.Length;
bool composite = false;
bool processedAll = false;
bool processedAll;
processedCount = 0UL;
lastProcessed = 0UL;
@@ -303,7 +310,7 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
{
ulong divisorValue = (ulong)nextDivisor128;
MontgomeryDivisorData montgomeryData = MontgomeryDivisorData.FromModulus(divisorValue);
ulong divisorCycle = ResolveDivisorCycle(divisorValue, prime, in montgomeryData);
ulong divisorCycle = ResolveDivisorCycle(gpu, divisorValue, prime, in montgomeryData);
if (divisorCycle == prime)
{
processedCount += (ulong)(i + 1);
@@ -375,9 +382,9 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static ulong ResolveDivisorCycle(ulong divisor, ulong prime, in MontgomeryDivisorData divisorData)
private static ulong ResolveDivisorCycle(HeuristicCombinedPrimeTesterAccelerator gpu, ulong divisor, ulong prime, in MontgomeryDivisorData divisorData)
{
if (!MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(divisor, prime, divisorData, out ulong computedCycle, out bool primeOrderFailed) || computedCycle == 0UL)
if (!MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(gpu, divisor, prime, divisorData, out ulong computedCycle, out bool primeOrderFailed) || computedCycle == 0UL)
{
return MersenneDivisorCycles.CalculateCycleLength(divisor, divisorData, skipPrimeOrderHeuristic: primeOrderFailed);
}
@@ -498,6 +505,7 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
private static readonly Accelerator[] _accelerators = AcceleratorPool.Shared.Accelerators;
private readonly MersenneNumberDivisorByDivisorGpuTester _owner;
private readonly HeuristicCombinedPrimeTesterAccelerator _primeTesterAccelerator;
private readonly Accelerator _accelerator;
private readonly AcceleratorStream _stream;
private Action<AcceleratorStream, Index1D, ArrayView<GpuDivisorPartialData>, ArrayView<int>, ArrayView<int>, ArrayView<ulong>, ArrayView<ulong>, ArrayView<byte>, ArrayView<int>> _kernel = null!;
@@ -516,8 +524,8 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
{
_owner = owner;
// GpuPrimeWorkLimiter.Acquire();
int acceleratorIndex = AcceleratorPool.Shared.Rent();
Accelerator accelerator = _accelerators[acceleratorIndex];
_primeTesterAccelerator = HeuristicCombinedPrimeTesterAccelerator.Rent();
Accelerator accelerator = _primeTesterAccelerator.Accelerator;
_accelerator = accelerator;
_stream = accelerator.CreateStream();
_capacity = Math.Max(1, owner._gpuBatchSize);
@@ -600,7 +608,8 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
if (cycle == 0UL)
{
if (!MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(divisor, firstPrime, divisorData, out ulong computedCycle, out bool primeOrderFailed) || computedCycle == 0UL)
var gpu = _primeTesterAccelerator;
if (!MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(gpu, divisor, firstPrime, divisorData, out ulong computedCycle, out bool primeOrderFailed) || computedCycle == 0UL)
{
cycle = MersenneDivisorCycles.CalculateCycleLength(divisor, divisorData, skipPrimeOrderHeuristic: primeOrderFailed);
}
@@ -657,10 +666,11 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
ref byte hitRef = ref MemoryMarshal.GetReference(hitSlice);
hitView.CopyToCPU(stream, ref hitRef, length);
stream.Synchronize();
// Monitor.Exit(_lease.ExecutionLock);
}
public void Dispose()
public void Return()
{
_owner._sessionPool.Add(this);
}
PerfectNumbers.Core/HeuristicCombinedPrimeTester.cs
+29 -18
View File
@@ -3,6 +3,7 @@ using ILGPU.Runtime;
using System.Runtime.CompilerServices;
using PerfectNumbers.Core.Gpu;
using PerfectNumbers.Core.Gpu.Accelerators;
using ILGPU;
namespace PerfectNumbers.Core;
@@ -172,10 +173,20 @@ public sealed class HeuristicCombinedPrimeTester
{
byte nMod10 = (byte)n.Mod10();
ulong maxDivisorSquare = ComputeHeuristicDivisorSquareLimit(n);
return IsPrimeGpu(n, maxDivisorSquare, nMod10);
var gpu = HeuristicCombinedPrimeTesterAccelerator.Rent(1);
bool isPrime = IsPrimeGpu(gpu, n, maxDivisorSquare, nMod10);
HeuristicCombinedPrimeTesterAccelerator.Return(gpu);
return isPrime;
}
public static bool IsPrimeGpu(ulong n, ulong maxDivisorSquare, byte nMod10)
public static bool IsPrimeGpu(HeuristicCombinedPrimeTesterAccelerator gpu, ulong n)
{
byte nMod10 = (byte)n.Mod10();
ulong maxDivisorSquare = ComputeHeuristicDivisorSquareLimit(n);
return IsPrimeGpu(gpu, n, maxDivisorSquare, nMod10);
}
public static bool IsPrimeGpu(HeuristicCombinedPrimeTesterAccelerator gpu, ulong n, ulong maxDivisorSquare, byte nMod10)
{
// TODO: Is this condition ever met on the execution path in EvenPerfectBitScanner?
if (maxDivisorSquare < 9UL)
@@ -183,7 +194,7 @@ public sealed class HeuristicCombinedPrimeTester
return EvaluateWithOpenNumericFallback(n);
}
bool compositeDetected = HeuristicTrialDivisionGpuDetectsDivisor(n, maxDivisorSquare, nMod10);
bool compositeDetected = HeuristicTrialDivisionGpuDetectsDivisor(gpu, n, maxDivisorSquare, nMod10);
return !compositeDetected;
}
@@ -222,40 +233,38 @@ public sealed class HeuristicCombinedPrimeTester
};
private static bool HeuristicTrialDivisionGpuDetectsDivisor(ulong n, ulong maxDivisorSquare, byte nMod10)
private static bool HeuristicTrialDivisionGpuDetectsDivisor(HeuristicCombinedPrimeTesterAccelerator gpu, ulong n, ulong maxDivisorSquare, byte nMod10)
{
// GpuPrimeWorkLimiter.Acquire();
var gpu = HeuristicPrimeTesterAccelerator.Rent(1);
var accelerator = gpu.Accelerator;
var stream = gpu.Stream;
var kernel = gpu.HeuristicTrialDivisionKernel;
var flagView1D = gpu.HeuristicFlag.View;
var flagView = flagView1D.AsContiguous();
// var gpu = HeuristicCombinedPrimeTesterAccelerator.Rent(1);
var acceleratorIndex = gpu.AcceleratorIndex;
var stream = AcceleratorStreamPool.Rent(acceleratorIndex);
var kernel = gpu.HeuristicCombinedTrialDivisionKernel;
var flagView1D = gpu.HeuristicFlag!.View;
bool compositeDetected = false;
int divisorLength = GetCombinedDivisors(nMod10).Length;
int compositeFlag = 0;
flagView1D.CopyFromCPU(stream, ref compositeFlag, 1);
stream!.Synchronize();
kernel.Launch(
var kernelLauncher = kernel.CreateLauncherDelegate<Action<AcceleratorStream, Index1D, byte, ArrayView<int>, ulong, ulong, HeuristicGpuCombinedDivisorTables>>();
kernelLauncher(
stream,
1,
divisorLength,
flagView,
nMod10,
flagView1D,
n,
maxDivisorSquare,
HeuristicGpuDivisorTableKind.Combined,
gpu.HeuristicGpuTables);
gpu.DivisorTables);
flagView1D.CopyToCPU(stream, ref compositeFlag, 1);
stream.Synchronize();
AcceleratorStreamPool.Return(acceleratorIndex, stream);
// GpuPrimeWorkLimiter.Release();
compositeDetected = compositeFlag != 0;
HeuristicPrimeTesterAccelerator.Return(gpu);
HeuristicCombinedPrimeTesterAccelerator.Return(gpu);
return compositeDetected;
}
@@ -490,6 +499,7 @@ public sealed class HeuristicCombinedPrimeTester
[MethodImpl(MethodImplOptions.AggressiveInlining)]
internal static ulong ResolveHeuristicCycleLength(
HeuristicCombinedPrimeTesterAccelerator gpu,
ulong exponent,
in HeuristicDivisorPreparation preparation,
out bool cycleFromHint,
@@ -508,6 +518,7 @@ public sealed class HeuristicCombinedPrimeTester
MontgomeryDivisorData divisorData = preparation.DivisorData;
bool trySuccess = MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(
gpu,
divisor,
exponent,
divisorData,
PerfectNumbers.Core/HeuristicPrimeTester.cs
+59 -48
View File
@@ -3,6 +3,7 @@ using System.Runtime.CompilerServices;
using ILGPU.Runtime;
using PerfectNumbers.Core.Gpu;
using PerfectNumbers.Core.Gpu.Accelerators;
using ILGPU;
namespace PerfectNumbers.Core;
@@ -242,72 +243,80 @@ public sealed class HeuristicPrimeTester
private bool HeuristicTrialDivisionGpuDetectsDivisor(ulong n, ulong maxDivisorSquare, byte nMod10)
{
// GpuPrimeWorkLimiter.Acquire();
var gpu = HeuristicPrimeTesterAccelerator.Rent(1);
var stream = gpu.Stream;
var kernel = gpu.HeuristicTrialDivisionKernel;
var flagView1D = gpu.HeuristicFlag.View;
var flagView = flagView1D.AsContiguous();
var gpu = HeuristicCombinedPrimeTesterAccelerator.Rent(1);
int acceleratorIndex = gpu.AcceleratorIndex;
var stream = AcceleratorStreamPool.Rent(acceleratorIndex);
var kernel = gpu.HeuristicCombinedTrialDivisionKernel;
var flagView1D = gpu.HeuristicFlag!.View;
bool compositeDetected = false;
bool compositeDetected;
int compositeFlag = 0;
int groupALength = HeuristicPrimeSieves.GroupADivisors.Length;
flagView1D.CopyFromCPU(stream, ref compositeFlag, 1);
kernel.Launch(
stream!,
// Index1D index,
// ArrayView<int> resultFlag,
// ulong n,
// ulong maxDivisorSquare,
// HeuristicGpuDivisorTableKind tableKind,
// HeuristicGpuDivisorTables tables)
var kernelLauncher = kernel.CreateLauncherDelegate<Action<AcceleratorStream, Index1D, byte, ArrayView<int>, ulong, ulong, HeuristicGpuDivisorTableKind, HeuristicGpuCombinedDivisorTables>>();
kernelLauncher(
stream,
1,
groupALength,
flagView,
nMod10,
flagView1D,
n,
maxDivisorSquare,
HeuristicGpuDivisorTableKind.GroupA,
gpu.HeuristicGpuTables);
gpu.DivisorTables);
flagView1D.CopyToCPU(stream, ref compositeFlag, 1);
stream!.Synchronize();
compositeDetected = compositeFlag != 0;
// compositeDetected = compositeFlag != 0;
if (!compositeDetected)
{
ReadOnlySpan<byte> endingOrder = GetGroupBEndingOrder(nMod10);
for (int i = 0; i < endingOrder.Length && !compositeDetected; i++)
{
byte ending = endingOrder[i];
HeuristicGpuDivisorTableKind tableKind = ending switch
{
1 => HeuristicGpuDivisorTableKind.GroupBEnding1,
7 => HeuristicGpuDivisorTableKind.GroupBEnding7,
9 => HeuristicGpuDivisorTableKind.GroupBEnding9,
_ => HeuristicGpuDivisorTableKind.GroupA,
};
// if (!compositeDetected)
// {
// ReadOnlySpan<byte> endingOrder = GetGroupBEndingOrder(nMod10);
// for (int i = 0; i < endingOrder.Length && !compositeDetected; i++)
// {
// byte ending = endingOrder[i];
// HeuristicGpuDivisorTableKind tableKind = ending switch
// {
// 1 => HeuristicGpuDivisorTableKind.GroupBEnding1,
// 7 => HeuristicGpuDivisorTableKind.GroupBEnding7,
// 9 => HeuristicGpuDivisorTableKind.GroupBEnding9,
// _ => HeuristicGpuDivisorTableKind.GroupA,
// };
if (tableKind == HeuristicGpuDivisorTableKind.GroupA)
{
continue;
}
// if (tableKind == HeuristicGpuDivisorTableKind.GroupA)
// {
// continue;
// }
int divisorLength = GetGroupBDivisors(ending).Length;
// int divisorLength = GetGroupBDivisors(ending).Length;
compositeFlag = 0;
flagView1D.CopyFromCPU(stream, ref compositeFlag, 1);
kernel.Launch(
stream,
1,
divisorLength,
flagView,
n,
maxDivisorSquare,
tableKind,
gpu.HeuristicGpuTables);
flagView1D.CopyToCPU(stream, ref compositeFlag, 1);
stream.Synchronize();
compositeDetected = compositeFlag != 0;
}
}
// compositeFlag = 0;
// flagView1D.CopyFromCPU(stream, ref compositeFlag, 1);
// kernelLauncher(
// stream,
// 1,
// flagView1D,
// n,
// maxDivisorSquare,
// tableKind,
// gpu.HeuristicGpuTables);
// flagView1D.CopyToCPU(stream, ref compositeFlag, 1);
// stream.Synchronize();
// compositeDetected = compositeFlag != 0;
// }
// }
gpu.Dispose();
AcceleratorStreamPool.Return(acceleratorIndex, stream);
gpu.Return();
// GpuPrimeWorkLimiter.Release();
return compositeDetected;
}
@@ -457,6 +466,7 @@ public sealed class HeuristicPrimeTester
[MethodImpl(MethodImplOptions.AggressiveInlining)]
internal static ulong ResolveHeuristicCycleLength(
HeuristicCombinedPrimeTesterAccelerator gpu,
ulong exponent,
in HeuristicDivisorPreparation preparation,
out bool cycleFromHint,
@@ -475,6 +485,7 @@ public sealed class HeuristicPrimeTester
MontgomeryDivisorData divisorData = preparation.DivisorData;
bool trySuccess = MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(
gpu,
divisor,
exponent,
divisorData,
PerfectNumbers.Core/IMersenneNumberDivisorByDivisorTester.cs
+1 -1
View File
@@ -20,6 +20,6 @@ public interface IMersenneNumberDivisorByDivisorTester
{
void CheckDivisor(ulong divisor, in MontgomeryDivisorData divisorData, ulong divisorCycle, in ReadOnlySpan<ulong> primes, Span<byte> hits);
void Dispose();
void Return();
}
}
PerfectNumbers.Core/MersenneDivisorCycles.cs
+12 -10
View File
@@ -52,7 +52,7 @@ public class MersenneDivisorCycles
return _smallCycles!;
}
public static bool CycleEqualsExponent(ulong divisor, in MontgomeryDivisorData divisorData, ulong exponent)
public static bool CycleEqualsExponent(HeuristicCombinedPrimeTesterAccelerator gpu, ulong divisor, in MontgomeryDivisorData divisorData, ulong exponent)
{
if (CycleEqualsExponentForMersenneCandidate(divisor, divisorData, exponent))
{
@@ -82,7 +82,7 @@ public class MersenneDivisorCycles
return heuristicCycle == exponent;
}
if (TryCalculateCycleLengthForExponentCpu(divisor, exponent, divisorData, out ulong cycleLength, out bool _) && cycleLength != 0UL)
if (TryCalculateCycleLengthForExponentCpu(gpu, divisor, exponent, divisorData, out ulong cycleLength, out bool _) && cycleLength != 0UL)
{
return cycleLength == exponent;
}
@@ -260,6 +260,7 @@ public class MersenneDivisorCycles
}
public static bool TryCalculateCycleLengthForExponentCpu(
HeuristicCombinedPrimeTesterAccelerator gpu,
ulong divisor,
ulong exponent,
in MontgomeryDivisorData divisorData,
@@ -328,8 +329,8 @@ public class MersenneDivisorCycles
factorCounts[2UL] = twoCount;
}
bool exponentAccumulated = AccumulateFactors(exponent, factorCounts);
bool factorsAccumulated = exponentAccumulated && AccumulateFactors(k, factorCounts);
bool exponentAccumulated = AccumulateFactors(gpu, exponent, factorCounts);
bool factorsAccumulated = exponentAccumulated && AccumulateFactors(gpu, k, factorCounts);
if (factorsAccumulated)
{
@@ -344,6 +345,7 @@ public class MersenneDivisorCycles
}
private static bool AccumulateFactors(
HeuristicCombinedPrimeTesterAccelerator gpu,
ulong value,
Dictionary<ulong, int> counts)
{
@@ -353,7 +355,7 @@ public class MersenneDivisorCycles
}
Dictionary<ulong, int> scratch = ThreadStaticPools.RentFactorScratchDictionary();
bool success = TryFactorIntoCountsInternal(value, scratch);
bool success = TryFactorIntoCountsInternal(gpu, value, scratch);
if (!success)
{
scratch.Clear();
@@ -375,7 +377,7 @@ public class MersenneDivisorCycles
// private static ulong _tryFactorIntoCountsInternalHits;
// private static ulong _tryCalculateCycleLengthHeuristicHits;
private static bool TryFactorIntoCountsInternal(ulong value, Dictionary<ulong, int> counts)
private static bool TryFactorIntoCountsInternal(HeuristicCombinedPrimeTesterAccelerator gpu, ulong value, Dictionary<ulong, int> counts)
{
// The EvenPerfectBitScanner only requests factors for values >= 2 here. Leave the guard documented without executing it.
// if (value <= 1UL)
@@ -412,7 +414,7 @@ public class MersenneDivisorCycles
// Console.WriteLine($"MersenneDivisorCycles.TryFactorIntoCountsInternal hits {Volatile.Read(ref _tryFactorIntoCountsInternalHits)}");
// HeuristicPrimeTester primeTester = _primeTester ??= new();
if (PrimeTester.IsPrimeGpu(remaining))
if (HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, remaining))
// if (PrimeTester.IsPrimeCpu(remaining, CancellationToken.None))
{
AddFactor(counts, remaining);
@@ -437,13 +439,13 @@ public class MersenneDivisorCycles
// return false;
// }
if (!TryFactorIntoCountsInternal(factor, counts))
if (!TryFactorIntoCountsInternal(gpu, factor, counts))
{
return false;
}
// Guard the recursive calls to terminate once Pollard-Rho reduces either branch to 1.
if (quotient > 1UL && !TryFactorIntoCountsInternal(quotient, counts))
if (quotient > 1UL && !TryFactorIntoCountsInternal(gpu, quotient, counts))
{
return false;
}
@@ -958,7 +960,7 @@ public class MersenneDivisorCycles
// HeuristicPrimeTester primeTester = _primeTester ??= new();
// if (!skipPrimeOrderHeuristic && HeuristicPrimeTester.Exclusive.IsPrimeCpu(divisor, CancellationToken.None))
// if (!skipPrimeOrderHeuristic && primeTester.IsPrimeGpu(divisor))
if (!skipPrimeOrderHeuristic && PrimeTester.IsPrimeGpu(divisor))
if (!skipPrimeOrderHeuristic && HeuristicCombinedPrimeTester.IsPrimeGpu(divisor))
{
ulong computedOrder = PrimeOrderCalculator.Calculate(
divisor,
PerfectNumbers.Core/PerfectNumberConstants.cs
+1
View File
@@ -21,3 +21,4 @@ public static class PerfectNumberConstants
// TODO: Promote these magic numbers into a runtime profile derived from EvenPerfectBitScanner.Benchmarks so we can retune
// them automatically when the optimal ranges for 2kp+1 >= 138M change.
}
PerfectNumbers.Core/PrimeOrderCalculator.Cpu.cs
+44 -28
View File
@@ -3,6 +3,7 @@ using System.Diagnostics;
using System.Numerics;
using System.Runtime.CompilerServices;
using System.Runtime.InteropServices;
using ILGPU;
using PerfectNumbers.Core.Cpu;
using PerfectNumbers.Core.Gpu;
using PerfectNumbers.Core.Gpu.Accelerators;
@@ -54,34 +55,39 @@ internal static partial class PrimeOrderCalculator
ulong phi = prime - 1UL;
PrimeOrderCalculatorAccelerator gpu = PrimeOrderCalculatorAccelerator.Rent(1);
HeuristicCombinedPrimeTesterAccelerator gpu2 = HeuristicCombinedPrimeTesterAccelerator.Rent(1);
if (IsGpuHeuristicDevice && PrimeOrderGpuHeuristics.TryCalculateOrder(prime, previousOrder, config, divisorData, out ulong gpuOrder))
{
PrimeOrderCalculatorAccelerator.Return(gpu);
HeuristicCombinedPrimeTesterAccelerator.Return(gpu2);
return gpuOrder;
}
PartialFactorResult phiFactors = PartialFactor(gpu, phi, config);
PartialFactorResult phiFactors = PartialFactor(gpu, gpu2, phi, config);
ulong result;
if (phiFactors.Factors is null)
{
result = prime <= PerfectNumberConstants.MaxQForDivisorCycles
? CalculateByFactorizationGpu(gpu, prime, divisorData)
: CalculateByFactorizationCpu(gpu, prime, divisorData);
? CalculateByFactorizationGpu(gpu, gpu2, prime, divisorData)
: CalculateByFactorizationCpu(gpu, gpu2, prime, divisorData);
phiFactors.Dispose();
PrimeOrderCalculatorAccelerator.Return(gpu);
HeuristicCombinedPrimeTesterAccelerator.Return(gpu2);
return result;
}
result = RunHeuristicPipelineCpu(gpu, prime, previousOrder, config, divisorData, phi, phiFactors);
result = RunHeuristicPipelineCpu(gpu, gpu2, prime, previousOrder, config, divisorData, phi, phiFactors);
phiFactors.Dispose();
PrimeOrderCalculatorAccelerator.Return(gpu);
HeuristicCombinedPrimeTesterAccelerator.Return(gpu2);
return result;
}
private static ulong RunHeuristicPipelineCpu(
PrimeOrderCalculatorAccelerator gpu,
HeuristicCombinedPrimeTesterAccelerator gpu2,
ulong prime,
ulong? previousOrder,
in PrimeOrderSearchConfig config,
@@ -98,7 +104,7 @@ internal static partial class PrimeOrderCalculator
ulong candidateOrder = InitializeStartingOrderCpu(gpu, prime, phi, divisorData);
candidateOrder = ExponentLoweringCpu(candidateOrder, prime, phiFactors, divisorData);
if (TryConfirmOrderCpu(gpu, prime, candidateOrder, divisorData, config, out PartialFactorResult? orderFactors))
if (TryConfirmOrderCpu(gpu, gpu2, prime, candidateOrder, divisorData, config, out PartialFactorResult? orderFactors))
{
// GpuPrimeWorkLimiter.Release();
return candidateOrder;
@@ -108,10 +114,10 @@ internal static partial class PrimeOrderCalculator
{
orderFactors?.Dispose();
// GpuPrimeWorkLimiter.Release();
return CalculateByFactorizationGpu(gpu, prime, divisorData);
return CalculateByFactorizationGpu(gpu, gpu2, prime, divisorData);
}
if (TryHeuristicFinishCpu(gpu, prime, candidateOrder, previousOrder, divisorData, config, phiFactors, orderFactors, out ulong order))
if (TryHeuristicFinishCpu(gpu, gpu2, prime, candidateOrder, previousOrder, divisorData, config, phiFactors, orderFactors, out ulong order))
{
// GpuPrimeWorkLimiter.Release();
return order;
@@ -340,6 +346,7 @@ internal static partial class PrimeOrderCalculator
private static bool TryConfirmOrderCpu(
PrimeOrderCalculatorAccelerator gpu,
HeuristicCombinedPrimeTesterAccelerator gpu2,
ulong prime,
ulong order,
in MontgomeryDivisorData divisorData,
@@ -358,7 +365,7 @@ internal static partial class PrimeOrderCalculator
return false;
}
PartialFactorResult factorization = PartialFactor(gpu, order, config);
PartialFactorResult factorization = PartialFactor(gpu, gpu2, order, config);
if (factorization.Factors is null)
{
factorization.Dispose();
@@ -503,6 +510,7 @@ internal static partial class PrimeOrderCalculator
private static bool TryHeuristicFinishCpu(
PrimeOrderCalculatorAccelerator gpu,
HeuristicCombinedPrimeTesterAccelerator gpu2,
ulong prime,
ulong order,
ulong? previousOrder,
@@ -520,7 +528,7 @@ internal static partial class PrimeOrderCalculator
}
// Reuse the partial factorization from TryConfirmOrderCpu when available.
PartialFactorResult orderFactors = cachedOrderFactors ?? PartialFactor(gpu, order, config);
PartialFactorResult orderFactors = cachedOrderFactors ?? PartialFactor(gpu, gpu2, order, config);
try
{
if (orderFactors.Factors is null)
@@ -636,7 +644,7 @@ internal static partial class PrimeOrderCalculator
bool equalsOne = powStepper.InitializeCpuIsUnity(candidate);
powStepperInitialized = true;
if (equalsOne && TryConfirmCandidateCpu(gpu, prime, candidate, divisorData, config, ref powUsed, powBudget))
if (equalsOne && TryConfirmCandidateCpu(gpu, gpu2, prime, candidate, divisorData, config, ref powUsed, powBudget))
{
if (gpuPool is not null)
{
@@ -673,7 +681,7 @@ internal static partial class PrimeOrderCalculator
continue;
}
if (!TryConfirmCandidateCpu(gpu, prime, candidate, divisorData, config, ref powUsed, powBudget))
if (!TryConfirmCandidateCpu(gpu, gpu2, prime, candidate, divisorData, config, ref powUsed, powBudget))
{
continue;
}
@@ -884,9 +892,9 @@ internal static partial class PrimeOrderCalculator
}
}
private static bool TryConfirmCandidateCpu(PrimeOrderCalculatorAccelerator gpu, ulong prime, ulong candidate, in MontgomeryDivisorData divisorData, in PrimeOrderSearchConfig config, ref int powUsed, int powBudget)
private static bool TryConfirmCandidateCpu(PrimeOrderCalculatorAccelerator gpu, HeuristicCombinedPrimeTesterAccelerator gpu2, ulong prime, ulong candidate, in MontgomeryDivisorData divisorData, in PrimeOrderSearchConfig config, ref int powUsed, int powBudget)
{
PartialFactorResult factorization = PartialFactor(gpu, candidate, config);
PartialFactorResult factorization = PartialFactor(gpu, gpu2, candidate, config);
try
{
if (factorization.Factors is null)
@@ -1082,7 +1090,9 @@ internal static partial class PrimeOrderCalculator
// private static ulong _partialFactorPendingHits;
// private static ulong _partialFactorCofactorHits;
private static PartialFactorResult PartialFactor(PrimeOrderCalculatorAccelerator gpu, ulong value, in PrimeOrderSearchConfig config)
private static bool _goToGpu = true;
private static PartialFactorResult PartialFactor(PrimeOrderCalculatorAccelerator gpu, HeuristicCombinedPrimeTesterAccelerator gpu2, ulong value, in PrimeOrderSearchConfig config)
{
if (value <= 1UL)
{
@@ -1169,7 +1179,11 @@ internal static partial class PrimeOrderCalculator
// Console.WriteLine($"Partial factor hits {Volatile.Read(ref _partialFactorHits)}");
// bool isPrime = PrimeTester.IsPrimeCpu(composite, CancellationToken.None);
bool isPrime = PrimeTester.IsPrimeGpu(gpu, composite);
bool goToGpu = Volatile.Read(ref _goToGpu);
Volatile.Write(ref _goToGpu, !goToGpu);
bool isPrime = goToGpu
? HeuristicCombinedPrimeTester.IsPrimeGpu(gpu2, composite)
: PrimeTester.IsPrime(composite);
// bool isPrime = Open.Numeric.Primes.Prime.Numbers.IsPrime(composite);
if (isPrime)
@@ -1241,7 +1255,7 @@ internal static partial class PrimeOrderCalculator
// Console.WriteLine($"Partial factor pending hits {Volatile.Read(ref _partialFactorPendingHits)}");
// bool isPrime = PrimeTester.IsPrimeCpu(composite, CancellationToken.None);
bool isPrime = PrimeTester.IsPrimeGpu(gpu, composite);
bool isPrime = HeuristicCombinedPrimeTester.IsPrimeGpu(gpu2, composite);
// bool isPrime = Open.Numeric.Primes.Prime.Numbers.IsPrime(composite);
entry = entry.WithPrimality(isPrime);
@@ -1278,7 +1292,8 @@ internal static partial class PrimeOrderCalculator
// Atomic.Add(ref _partialFactorCofactorHits, 1UL);
// Console.WriteLine($"Partial factor cofactor hits {Volatile.Read(ref _partialFactorCofactorHits)}");
cofactorIsPrime = PrimeTester.IsPrimeGpu(gpu, cofactor);
cofactorIsPrime = HeuristicCombinedPrimeTester.IsPrimeGpu(gpu2, cofactor);
// cofactorIsPrime = PrimeTester.IsPrimeGpu(gpu, cofactor);
// cofactorIsPrime = PrimeTester.IsPrimeGpu(cofactor);
// cofactorIsPrime = Open.Numeric.Primes.Prime.Numbers.IsPrime(cofactor);
}
@@ -1698,11 +1713,11 @@ internal static partial class PrimeOrderCalculator
exponents = newExponentSpan;
}
private static ulong CalculateByFactorizationCpu(PrimeOrderCalculatorAccelerator gpu, ulong prime, in MontgomeryDivisorData divisorData)
private static ulong CalculateByFactorizationCpu(PrimeOrderCalculatorAccelerator gpu, HeuristicCombinedPrimeTesterAccelerator gpu2, ulong prime, in MontgomeryDivisorData divisorData)
{
ulong phi = prime - 1UL;
Dictionary<ulong, int> counts = new(capacity: 8);
FactorCompletelyCpu(phi, counts);
FactorCompletelyCpu(gpu2, phi, counts);
if (counts.Count == 0)
{
return phi;
@@ -1746,9 +1761,9 @@ internal static partial class PrimeOrderCalculator
return order;
}
private static void FactorCompletelyCpu(ulong value, Dictionary<ulong, int> counts)
private static void FactorCompletelyCpu(HeuristicCombinedPrimeTesterAccelerator gpu, ulong value, Dictionary<ulong, int> counts)
{
FactorCompletely(value, counts, knownComposite: false);
FactorCompletely(gpu, value, counts, knownComposite: false);
}
// private static ulong _factorCompletelyHits;
@@ -1756,7 +1771,7 @@ internal static partial class PrimeOrderCalculator
// [ThreadStatic]
// private static HeuristicPrimeTester? _tester;
private static void FactorCompletely(ulong value, Dictionary<ulong, int> counts, bool knownComposite)
private static void FactorCompletely(HeuristicCombinedPrimeTesterAccelerator gpu, ulong value, Dictionary<ulong, int> counts, bool knownComposite)
{
if (value <= 1UL)
{
@@ -1768,7 +1783,8 @@ internal static partial class PrimeOrderCalculator
// HeuristicPrimeTester tester = _tester ??= new();
// bool isPrime = primeTester.HeuristicIsPrimeGpu(value);
bool isPrime = PrimeTester.IsPrimeGpu(value);
bool isPrime = HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, value);
// bool isPrime = PrimeTester.IsPrimeGpu(value);
// bool isPrime = PrimeTester.IsPrimeGpu(value);
// bool isPrime = Open.Numeric.Primes.Prime.Numbers.IsPrime(value);
if (!knownComposite && isPrime)
@@ -1784,7 +1800,7 @@ internal static partial class PrimeOrderCalculator
// Console.WriteLine($"Factor completely after PollardRho hits {Volatile.Read(ref _factorCompletelyPollardRhoHits)}");
// isPrime = primeTester.HeuristicIsPrimeGpu(factor);
isPrime = PrimeTester.IsPrimeGpu(factor);
isPrime = HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, factor);
// isPrime = PrimeTester.IsPrimeGpu(factor);
// isPrime = Open.Numeric.Primes.Prime.Numbers.IsPrime(factor);
if (isPrime)
@@ -1801,13 +1817,13 @@ internal static partial class PrimeOrderCalculator
if (remaining > 1UL)
{
FactorCompletely(remaining, counts, knownComposite: false);
FactorCompletely(gpu, remaining, counts, knownComposite: false);
}
return;
}
FactorCompletely(factor, counts, knownComposite: true);
FactorCompletely(gpu, factor, counts, knownComposite: true);
if (quotient == factor)
{
@@ -1821,7 +1837,7 @@ internal static partial class PrimeOrderCalculator
// Atomic.Add(ref _partialFactorCofactorHits, 1UL);
// Console.WriteLine($"Partial factor cofactor hits {Volatile.Read(ref _partialFactorCofactorHits)}");
isPrime = PrimeTester.IsPrimeGpu(quotient);
isPrime = HeuristicCombinedPrimeTester.IsPrimeGpu(gpu, quotient);
// isPrime = PrimeTester.IsPrimeGpu(quotient);
// isPrime = Open.Numeric.Primes.Prime.Numbers.IsPrime(quotient);
}
@@ -1832,7 +1848,7 @@ internal static partial class PrimeOrderCalculator
}
else
{
FactorCompletely(quotient, counts, knownComposite: true);
FactorCompletely(gpu, quotient, counts, knownComposite: true);
}
}
PerfectNumbers.Core/PrimeOrderCalculator.Gpu.cs
+2 -2
View File
@@ -12,13 +12,13 @@ internal static partial class PrimeOrderCalculator
private const int GpuSmallPrimeFactorSlots = 64;
private static ulong CalculateByFactorizationGpu(PrimeOrderCalculatorAccelerator gpu, ulong prime, in MontgomeryDivisorData divisorData)
private static ulong CalculateByFactorizationGpu(PrimeOrderCalculatorAccelerator gpu, HeuristicCombinedPrimeTesterAccelerator gpu2, ulong prime, in MontgomeryDivisorData divisorData)
{
ulong phi = prime - 1UL;
Dictionary<ulong, int> counts = gpu.Pow2ModEntriesToTestOnHost;
counts.Clear();
FactorCompletelyCpu(phi, counts);
FactorCompletelyCpu(gpu2, phi, counts);
if (counts.Count == 0)
{
return phi;