Restore struct residue args and add residue benchmark project

This commit is contained in:
Marek Lemańczyk committed 2025-09-21 05:29:13 +02:00
1 parent 4fe701425e
commit 5ab69f8c01
10 files changed
+326 -165

No files matched your search

+4 -4
View File
@@ -158,7 +158,7 @@ public ref struct GpuKernelLease(IDisposable limiter, GpuContextLease gpu, Kerne
}
}
public readonly Action<AcceleratorStream, Index1D, ulong, GpuUInt128, GpuUInt128, byte, ulong, ResidueAutomatonArgs, ArrayView<int>, ArrayView1D<uint, Stride1D.Dense>> IncrementalOrderKernel
public readonly Action<AcceleratorStream, Index1D, ulong, GpuUInt128, GpuUInt128, byte, ulong, ResidueAutomatonArgs, ArrayView<int>, ArrayView1D<uint, Stride1D.Dense>> IncrementalOrderKernel
{
get
{
@@ -175,7 +175,7 @@ public ref struct GpuKernelLease(IDisposable limiter, GpuContextLease gpu, Kerne
}
}
public readonly Action<AcceleratorStream, Index1D, ulong, GpuUInt128, GpuUInt128, byte, ulong, ResidueAutomatonArgs, ArrayView<int>, ArrayView1D<uint, Stride1D.Dense>> Pow2ModOrderKernel
public readonly Action<AcceleratorStream, Index1D, ulong, GpuUInt128, GpuUInt128, byte, ulong, ResidueAutomatonArgs, ArrayView<int>, ArrayView1D<uint, Stride1D.Dense>> Pow2ModOrderKernel
{
get
{
@@ -486,7 +486,7 @@ public ref struct GpuKernelLease(IDisposable limiter, GpuContextLease gpu, Kerne
}
private static void IncrementalOrderKernelScan(Index1D index, ulong exponent, GpuUInt128 twoP, GpuUInt128 kStart, byte lastIsSeven, ulong divMul,
ResidueAutomatonArgs ra, ArrayView<int> found, ArrayView1D<uint, Stride1D.Dense> smallCycles)
ResidueAutomatonArgs ra, ArrayView<int> found, ArrayView1D<uint, Stride1D.Dense> smallCycles)
{
ulong idx = (ulong)index.X;
ulong r10 = ra.Q0M10 + (ra.Step10 * idx) % 10UL; r10 -= (r10 >= 10UL) ? 10UL : 0UL;
@@ -570,7 +570,7 @@ public ref struct GpuKernelLease(IDisposable limiter, GpuContextLease gpu, Kerne
}
private static void Pow2ModOrderKernelScan(Index1D index, ulong exponent, GpuUInt128 twoP, GpuUInt128 kStart, byte lastIsSeven, ulong _,
ResidueAutomatonArgs ra, ArrayView<int> found, ArrayView1D<uint, Stride1D.Dense> smallCycles)
ResidueAutomatonArgs ra, ArrayView<int> found, ArrayView1D<uint, Stride1D.Dense> smallCycles)
{
ulong idx = (ulong)index.X;
ulong r10 = ra.Q0M10 + (ra.Step10 * idx) % 10UL; r10 -= (r10 >= 10UL) ? 10UL : 0UL;
@@ -23,10 +23,12 @@ public class MersenneNumberIncrementalGpuTester(GpuKernelType kernelType, bool u
var pow2Kernel = gpuLease.Pow2ModKernel;
var incKernel = gpuLease.IncrementalKernel;
ulong step10 = ((exponent % 10UL) << 1) % 10UL;
ulong step8 = ((exponent & 7UL) << 1) & 7UL;
ulong step3 = ((exponent % 3UL) << 1) % 3UL;
ulong step5 = ((exponent % 5UL) << 1) % 5UL;
UInt128 exponent128 = exponent;
exponent128.Mod10_8_5_3(out ulong exponentMod10, out ulong exponentMod8, out ulong exponentMod5, out ulong exponentMod3);
ulong step10 = (exponentMod10 << 1) % 10UL;
ulong step8 = (exponentMod8 << 1) & 7UL;
ulong step3 = (exponentMod3 << 1) % 3UL;
ulong step5 = (exponentMod5 << 1) % 5UL;
GpuUInt128 twoPGpu = (GpuUInt128)twoP;
var smallCyclesView = GpuKernelPool.EnsureSmallCyclesOnDevice(accelerator);
ResiduePrimeViews primeViews = default;
@@ -55,9 +57,9 @@ public class MersenneNumberIncrementalGpuTester(GpuKernelType kernelType, bool u
if (_kernelType == GpuKernelType.Pow2Mod)
{
q0.Mod10_8_5_3(out ulong q0m10, out ulong q0m8, out ulong q0m5, out ulong q0m3);
var ra = new ResidueAutomatonArgs(q0m10, step10, q0m8, step8, q0m3, step3, q0m5, step5);
var kernelArgs = new ResidueAutomatonArgs(q0m10, step10, q0m8, step8, q0m3, step3, q0m5, step5);
pow2Kernel(stream, currentSize, exponent, twoPGpu, (GpuUInt128)kStart, last, divMul,
ra, orderBuffer.View, smallCyclesView, primeViews.LastOne, primeViews.LastSeven,
kernelArgs, orderBuffer.View, smallCyclesView, primeViews.LastOne, primeViews.LastSeven,
primeViews.LastOnePow2, primeViews.LastSevenPow2);
}
else
@@ -25,6 +25,12 @@ public class MersenneNumberOrderGpuTester(GpuKernelType kernelType, bool useGpuO
};
var foundBuffer = accelerator.Allocate1D<int>(1);
UInt128 exponent128 = exponent;
exponent128.Mod10_8_5_3(out ulong exponentMod10, out ulong exponentMod8, out ulong exponentMod5, out ulong exponentMod3);
ulong step10 = (exponentMod10 << 1) % 10UL;
ulong step8 = (exponentMod8 << 1) & 7UL;
ulong step3 = (exponentMod3 << 1) % 3UL;
ulong step5 = (exponentMod5 << 1) % 5UL;
try
{
UInt128 remaining;
@@ -37,21 +43,14 @@ public class MersenneNumberOrderGpuTester(GpuKernelType kernelType, bool useGpuO
foundBuffer.MemSetToZero();
// Precompute residue automaton bases for this batch
UInt128 q0 = twoP * kStart + 1UL;
ulong q0m10 = q0.Mod10();
ulong q0m8 = q0.Mod8();
ulong q0m3 = q0.Mod3();
ulong q0m5 = q0.Mod5();
ulong step10 = ((exponent % 10UL) << 1) % 10UL;
ulong step8 = ((exponent & 7UL) << 1) & 7UL;
ulong step3 = ((exponent % 3UL) << 1) % 3UL;
ulong step5 = ((exponent % 5UL) << 1) % 5UL;
var ra = new ResidueAutomatonArgs(q0m10, step10, q0m8, step8, q0m3, step3, q0m5, step5);
q0.Mod10_8_5_3(out ulong q0m10, out ulong q0m8, out ulong q0m5, out ulong q0m3);
var kernelArgs = new ResidueAutomatonArgs(q0m10, step10, q0m8, step8, q0m3, step3, q0m5, step5);
// Ensure device has small cycles table for early rejections in order kernels
var smallCyclesView = GpuKernelPool.EnsureSmallCyclesOnDevice(accelerator);
kernel(stream, currentSize, exponent, (GpuUInt128)twoP, (GpuUInt128)kStart, last, divMul,
ra, foundBuffer.View, smallCyclesView);
kernelArgs, foundBuffer.View, smallCyclesView);
stream.Synchronize();
foundBuffer.View.CopyToCPU(ref found, 1);
@@ -20,74 +20,89 @@ public class MersenneNumberResidueGpuTester(bool useGpuOrder)
// Ensure device has small cycles and primes tables for in-kernel lookup
var smallCyclesView = GpuKernelPool.EnsureSmallCyclesOnDevice(accelerator);
ResiduePrimeViews primeViews = GpuKernelPool.EnsureSmallPrimesOnDevice(accelerator);
int batchSize = GpuConstants.ScanBatchSize;
UInt128 kStart = 1UL;
byte last = lastIsSeven ? (byte)1 : (byte)0;
var kernel = gpuLease.Pow2ModKernel;
ulong step10 = ((exponent % 10UL) << 1) % 10UL;
ulong step8 = ((exponent & 7UL) << 1) & 7UL;
ulong step3 = ((exponent % 3UL) << 1) % 3UL;
ulong step5 = ((exponent % 5UL) << 1) % 5UL;
GpuUInt128 twoPGpu = (GpuUInt128)twoP;
int batchSize = GpuConstants.ScanBatchSize;
UInt128 kStart = 1UL;
byte last = lastIsSeven ? (byte)1 : (byte)0;
var kernel = gpuLease.Pow2ModKernel;
UInt128 exponent128 = exponent;
exponent128.Mod10_8_5_3(out ulong exponentMod10, out ulong exponentMod8, out ulong exponentMod5, out ulong exponentMod3);
ulong step10 = (exponentMod10 << 1) % 10UL;
ulong step8 = (exponentMod8 << 1) & 7UL;
ulong step3 = (exponentMod3 << 1) % 3UL;
ulong step5 = (exponentMod5 << 1) % 5UL;
GpuUInt128 twoPGpu = (GpuUInt128)twoP;
var orderBuffer = accelerator.Allocate1D<ulong>(batchSize);
// Device memory returned by Allocate1D is not guaranteed to be zeroed.
// Ensure the buffer starts with a known state so that lanes skipped by the
// kernel (for example due to early rejections) do not leave stale garbage
// values that would be interpreted as composite witnesses when copied back
// to the host. This mirrors the explicit zeroing performed by other GPU
// testers such as the order kernels.
orderBuffer.MemSetToZero();
ulong[] orderArray = ArrayPool<ulong>.Shared.Rent(batchSize);
UInt128 remaining;
int currentSize;
int i;
UInt128 batchSize128 = (UInt128)batchSize, q;
Span<ulong> orders = orderArray.AsSpan(0, batchSize);
ref var ordersRef = ref MemoryMarshal.GetReference(orders);
while (kStart < maxK && Volatile.Read(ref isPrime))
{
remaining = maxK - kStart;
if (remaining > batchSize128)
{
currentSize = batchSize;
}
else
{
currentSize = (int)remaining;
orders = orderArray.AsSpan(0, currentSize);
ordersRef = ref MemoryMarshal.GetReference(orders);
}
var orderBuffer = accelerator.Allocate1D<ulong>(batchSize);
// Device memory returned by Allocate1D is not guaranteed to be zeroed.
// Ensure the buffer starts with a known state so that lanes skipped by the
// kernel (for example due to early rejections) do not leave stale garbage
// values that would be interpreted as composite witnesses when copied back
// to the host. This mirrors the explicit zeroing performed by other GPU
// testers such as the order kernels.
orderBuffer.MemSetToZero();
ulong[] orderArray = ArrayPool<ulong>.Shared.Rent(batchSize);
UInt128 batchSize128 = (UInt128)batchSize;
UInt128 fullBatchStep = twoP * batchSize128;
UInt128 q = twoP * kStart + UInt128.One;
q = twoP * kStart + UInt128.One;
q.Mod10_8_5_3(out ulong q0m10, out ulong q0m8, out ulong q0m5, out ulong q0m3);
var ra = new ResidueAutomatonArgs(q0m10, step10, q0m8, step8, q0m3, step3, q0m5, step5);
try
{
while (kStart < maxK && Volatile.Read(ref isPrime))
{
UInt128 remaining = maxK - kStart;
int currentSize = remaining > batchSize128 ? batchSize : (int)remaining;
Span<ulong> orders = orderArray.AsSpan(0, currentSize);
ref ulong ordersRef = ref MemoryMarshal.GetReference(orders);
kernel(stream, currentSize, exponent, twoPGpu, (GpuUInt128)kStart, last, 0UL,
ra, orderBuffer.View, smallCyclesView, primeViews.LastOne, primeViews.LastSeven, primeViews.LastOnePow2, primeViews.LastSevenPow2);
q.Mod10_8_5_3(out ulong q0m10, out ulong q0m8, out ulong q0m5, out ulong q0m3);
var kernelArgs = new ResidueAutomatonArgs(q0m10, step10, q0m8, step8, q0m3, step3, q0m5, step5);
accelerator.Synchronize();
orderBuffer.View.CopyToCPU(ref ordersRef, currentSize);
if (!Volatile.Read(ref isPrime))
{
break;
}
kernel(stream, currentSize, exponent, twoPGpu, (GpuUInt128)kStart, last, 0UL,
kernelArgs, orderBuffer.View, smallCyclesView, primeViews.LastOne, primeViews.LastSeven, primeViews.LastOnePow2, primeViews.LastSevenPow2);
for (i = 0; i < currentSize; i++)
{
if (orders[i] != 0UL)
{
Volatile.Write(ref isPrime, false);
goto cleanup;
}
}
accelerator.Synchronize();
orderBuffer.View.CopyToCPU(ref ordersRef, currentSize);
if (!Volatile.Read(ref isPrime))
{
break;
}
kStart += batchSize128;
}
bool compositeFound = false;
for (int i = 0; i < currentSize; i++)
{
if (orders[i] != 0UL)
{
Volatile.Write(ref isPrime, false);
compositeFound = true;
break;
}
}
cleanup:
ArrayPool<ulong>.Shared.Return(orderArray);
orderBuffer.Dispose();
gpuLease.Dispose();
}
if (compositeFound)
{
break;
}
UInt128 processed = (UInt128)currentSize;
kStart += processed;
if (kStart < maxK)
{
if (processed == batchSize128)
{
q += fullBatchStep;
}
else
{
q += twoP * processed;
}
}
}
}
finally
{
ArrayPool<ulong>.Shared.Return(orderArray);
orderBuffer.Dispose();
gpuLease.Dispose();
}
}
}