Merge pull request #248 from mlemanczyk/codex/optimize-heuristicprimetester-for-divisor-selection
Add filtered group B divisor tables for heuristic GPU checks
This commit is contained in:
13 files changed
+1124
-100
No files matched your search
@@ -12,6 +12,7 @@ namespace EvenPerfectBitScanner.Benchmarks;
|
|||||||
[MemoryDiagnoser]
|
[MemoryDiagnoser]
|
||||||
public class PrimeTesterBenchmarks
|
public class PrimeTesterBenchmarks
|
||||||
{
|
{
|
||||||
|
private static readonly HeuristicCombinedPrimeTester _heuristicCombinedTester = new();
|
||||||
private static readonly HeuristicPrimeTester _heuristicTester = new();
|
private static readonly HeuristicPrimeTester _heuristicTester = new();
|
||||||
private ulong[] _candidates = [];
|
private ulong[] _candidates = [];
|
||||||
|
|
||||||
@@ -71,6 +72,42 @@ public class PrimeTesterBenchmarks
|
|||||||
return primeCount;
|
return primeCount;
|
||||||
}
|
}
|
||||||
|
|
||||||
|
[Benchmark]
|
||||||
|
public int HeuristicCombinedGpu()
|
||||||
|
{
|
||||||
|
HeuristicCombinedPrimeTester tester = _heuristicCombinedTester;
|
||||||
|
ulong[] values = _candidates;
|
||||||
|
int primeCount = 0;
|
||||||
|
|
||||||
|
for (int i = 0; i < values.Length; i++)
|
||||||
|
{
|
||||||
|
if (HeuristicCombinedPrimeTester.IsPrimeGpu(values[i]))
|
||||||
|
{
|
||||||
|
primeCount++;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
return primeCount;
|
||||||
|
}
|
||||||
|
|
||||||
|
[Benchmark]
|
||||||
|
public int HeuristicCombinedCpu()
|
||||||
|
{
|
||||||
|
HeuristicCombinedPrimeTester tester = _heuristicCombinedTester;
|
||||||
|
ulong[] values = _candidates;
|
||||||
|
int primeCount = 0;
|
||||||
|
|
||||||
|
for (int i = 0; i < values.Length; i++)
|
||||||
|
{
|
||||||
|
if (HeuristicCombinedPrimeTester.IsPrimeCpu(values[i]))
|
||||||
|
{
|
||||||
|
primeCount++;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
return primeCount;
|
||||||
|
}
|
||||||
|
|
||||||
private static readonly PrimeTester _nonHeuristicTester = new();
|
private static readonly PrimeTester _nonHeuristicTester = new();
|
||||||
|
|
||||||
[Benchmark]
|
[Benchmark]
|
||||||
|
|||||||
@@ -55,7 +55,7 @@ public class ProgramTests
|
|||||||
[InlineData(81UL, false)]
|
[InlineData(81UL, false)]
|
||||||
public void IsPrime_identifies_primes_correctly(ulong n, bool expected)
|
public void IsPrime_identifies_primes_correctly(ulong n, bool expected)
|
||||||
{
|
{
|
||||||
new PrimeTester().IsPrime(n, CancellationToken.None).Should().Be(expected);
|
PrimeTester.IsPrime(n, CancellationToken.None).Should().Be(expected);
|
||||||
}
|
}
|
||||||
|
|
||||||
[Fact]
|
[Fact]
|
||||||
|
|||||||
@@ -622,7 +622,7 @@ internal static class Program
|
|||||||
if (!_cliArguments.UseByDivisor)
|
if (!_cliArguments.UseByDivisor)
|
||||||
{
|
{
|
||||||
if (!(GpuContextPool.ForceCpu
|
if (!(GpuContextPool.ForceCpu
|
||||||
? PrimeTesters.Value!.IsPrime(p, CancellationToken.None)
|
? PrimeTester.IsPrime(p, CancellationToken.None)
|
||||||
: PrimeTesters.Value!.IsPrimeGpu(p, CancellationToken.None)))
|
: PrimeTesters.Value!.IsPrimeGpu(p, CancellationToken.None)))
|
||||||
{
|
{
|
||||||
_lastCompositeP = true;
|
_lastCompositeP = true;
|
||||||
|
|||||||
@@ -0,0 +1,243 @@
|
|||||||
|
using System;
|
||||||
|
using System.Numerics;
|
||||||
|
using System.Collections.Generic;
|
||||||
|
using PerfectNumbers.Core;
|
||||||
|
using FluentAssertions;
|
||||||
|
using Xunit;
|
||||||
|
|
||||||
|
namespace PerfectNumbers.Core.Tests;
|
||||||
|
|
||||||
|
public class HeuristicCombinedPrimeTesterTests
|
||||||
|
{
|
||||||
|
[Theory]
|
||||||
|
[Trait("Category", "Fast")]
|
||||||
|
[InlineData((byte)1)]
|
||||||
|
[InlineData((byte)7)]
|
||||||
|
[InlineData((byte)9)]
|
||||||
|
public void Enumerator_uses_one_a_one_b_schedule_by_default(byte lastDigit)
|
||||||
|
{
|
||||||
|
Span<HeuristicCombinedPrimeTester.HeuristicGroupBSequenceState> groupBState =
|
||||||
|
stackalloc HeuristicCombinedPrimeTester.HeuristicGroupBSequenceState[4];
|
||||||
|
|
||||||
|
var enumerator = HeuristicCombinedPrimeTester.CreateHeuristicDivisorEnumerator(500, lastDigit, groupBState);
|
||||||
|
var results = new List<HeuristicCombinedPrimeTester.HeuristicDivisorCandidate>(32);
|
||||||
|
|
||||||
|
while (results.Count < 20 && enumerator.TryGetNext(out var candidate))
|
||||||
|
{
|
||||||
|
results.Add(candidate);
|
||||||
|
}
|
||||||
|
|
||||||
|
const int segmentsToCheck = 6;
|
||||||
|
results.Count.Should().BeGreaterThanOrEqualTo(4 + segmentsToCheck * 2);
|
||||||
|
results[0].Value.Should().Be(3UL);
|
||||||
|
results[0].Group.Should().Be(HeuristicCombinedPrimeTester.HeuristicDivisorGroup.GroupAConstant);
|
||||||
|
results[1].Value.Should().Be(7UL);
|
||||||
|
results[1].Group.Should().Be(HeuristicCombinedPrimeTester.HeuristicDivisorGroup.GroupAConstant);
|
||||||
|
results[2].Value.Should().Be(11UL);
|
||||||
|
results[2].Group.Should().Be(HeuristicCombinedPrimeTester.HeuristicDivisorGroup.GroupAConstant);
|
||||||
|
results[3].Value.Should().Be(13UL);
|
||||||
|
results[3].Group.Should().Be(HeuristicCombinedPrimeTester.HeuristicDivisorGroup.GroupAConstant);
|
||||||
|
|
||||||
|
int index = 4;
|
||||||
|
|
||||||
|
for (int segment = 0; segment < segmentsToCheck; segment++)
|
||||||
|
{
|
||||||
|
HeuristicCombinedPrimeTester.HeuristicDivisorCandidate groupACandidate = results[index++];
|
||||||
|
groupACandidate.Value.Should().BeGreaterThan(13UL);
|
||||||
|
(groupACandidate.Value % 10UL).Should().Be(3UL);
|
||||||
|
groupACandidate.Group.Should().NotBe(HeuristicCombinedPrimeTester.HeuristicDivisorGroup.GroupB);
|
||||||
|
|
||||||
|
HeuristicCombinedPrimeTester.HeuristicDivisorCandidate groupBCandidate = results[index++];
|
||||||
|
groupBCandidate.Value.Should().BeGreaterThan(13UL);
|
||||||
|
(groupBCandidate.Value % 10UL).Should().NotBe(3UL);
|
||||||
|
groupBCandidate.Group.Should().Be(HeuristicCombinedPrimeTester.HeuristicDivisorGroup.GroupB);
|
||||||
|
}
|
||||||
|
|
||||||
|
}
|
||||||
|
|
||||||
|
[Theory]
|
||||||
|
[Trait("Category", "Fast")]
|
||||||
|
[InlineData((byte)1)]
|
||||||
|
[InlineData((byte)7)]
|
||||||
|
[InlineData((byte)9)]
|
||||||
|
public void One_a_one_b_pattern_interleaves_group_a_and_group_b_values(byte lastDigit)
|
||||||
|
{
|
||||||
|
ReadOnlySpan<uint> combined = HeuristicCombinedPrimeTester.GetCombinedDivisors(
|
||||||
|
lastDigit,
|
||||||
|
HeuristicCombinedPrimeTester.CombinedDivisorPattern.OneAOneB);
|
||||||
|
|
||||||
|
combined.Length.Should().BeGreaterThan(18);
|
||||||
|
combined[0].Should().Be(3U);
|
||||||
|
combined[1].Should().Be(7U);
|
||||||
|
combined[2].Should().Be(11U);
|
||||||
|
combined[3].Should().Be(13U);
|
||||||
|
|
||||||
|
int index = 4;
|
||||||
|
const int segmentsToCheck = 6;
|
||||||
|
|
||||||
|
for (int segment = 0; segment < segmentsToCheck; segment++)
|
||||||
|
{
|
||||||
|
uint groupAValue = combined[index++];
|
||||||
|
groupAValue.Should().BeGreaterThan(13U);
|
||||||
|
(groupAValue % 10U).Should().Be(3U);
|
||||||
|
|
||||||
|
uint groupBValue = combined[index++];
|
||||||
|
groupBValue.Should().BeGreaterThan(13U);
|
||||||
|
(groupBValue % 10U).Should().NotBe(3U);
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
[Theory]
|
||||||
|
[Trait("Category", "Fast")]
|
||||||
|
[InlineData((byte)1)]
|
||||||
|
[InlineData((byte)7)]
|
||||||
|
[InlineData((byte)9)]
|
||||||
|
public void Two_a_one_b_pattern_interleaves_two_group_a_values_before_group_b(byte lastDigit)
|
||||||
|
{
|
||||||
|
ReadOnlySpan<uint> combined = HeuristicCombinedPrimeTester.GetCombinedDivisors(
|
||||||
|
lastDigit,
|
||||||
|
HeuristicCombinedPrimeTester.CombinedDivisorPattern.TwoAOneB);
|
||||||
|
|
||||||
|
combined.Length.Should().BeGreaterThan(20);
|
||||||
|
combined[0].Should().Be(3U);
|
||||||
|
combined[1].Should().Be(7U);
|
||||||
|
combined[2].Should().Be(11U);
|
||||||
|
combined[3].Should().Be(13U);
|
||||||
|
|
||||||
|
int index = 4;
|
||||||
|
const int segmentsToCheck = 5;
|
||||||
|
|
||||||
|
for (int segment = 0; segment < segmentsToCheck; segment++)
|
||||||
|
{
|
||||||
|
for (int i = 0; i < 2; i++)
|
||||||
|
{
|
||||||
|
uint groupAValue = combined[index++];
|
||||||
|
groupAValue.Should().BeGreaterThan(13U);
|
||||||
|
(groupAValue % 10U).Should().Be(3U);
|
||||||
|
}
|
||||||
|
|
||||||
|
uint groupBValue = combined[index++];
|
||||||
|
groupBValue.Should().BeGreaterThan(13U);
|
||||||
|
(groupBValue % 10U).Should().NotBe(3U);
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
[Theory]
|
||||||
|
[Trait("Category", "Fast")]
|
||||||
|
[InlineData((byte)1)]
|
||||||
|
[InlineData((byte)7)]
|
||||||
|
[InlineData((byte)9)]
|
||||||
|
public void Three_a_two_b_pattern_interleaves_three_group_a_values_before_two_group_b_values(byte lastDigit)
|
||||||
|
{
|
||||||
|
ReadOnlySpan<uint> combined = HeuristicCombinedPrimeTester.GetCombinedDivisors(
|
||||||
|
lastDigit,
|
||||||
|
HeuristicCombinedPrimeTester.CombinedDivisorPattern.ThreeATwoB);
|
||||||
|
|
||||||
|
combined.Length.Should().BeGreaterThan(28);
|
||||||
|
combined[0].Should().Be(3U);
|
||||||
|
combined[1].Should().Be(7U);
|
||||||
|
combined[2].Should().Be(11U);
|
||||||
|
combined[3].Should().Be(13U);
|
||||||
|
|
||||||
|
int index = 4;
|
||||||
|
const int segmentsToCheck = 4;
|
||||||
|
|
||||||
|
for (int segment = 0; segment < segmentsToCheck; segment++)
|
||||||
|
{
|
||||||
|
for (int i = 0; i < 3; i++)
|
||||||
|
{
|
||||||
|
uint groupAValue = combined[index++];
|
||||||
|
groupAValue.Should().BeGreaterThan(13U);
|
||||||
|
(groupAValue % 10U).Should().Be(3U);
|
||||||
|
}
|
||||||
|
|
||||||
|
for (int i = 0; i < 2; i++)
|
||||||
|
{
|
||||||
|
uint groupBValue = combined[index++];
|
||||||
|
groupBValue.Should().BeGreaterThan(13U);
|
||||||
|
(groupBValue % 10U).Should().NotBe(3U);
|
||||||
|
}
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
[Theory]
|
||||||
|
[Trait("Category", "Fast")]
|
||||||
|
[InlineData((byte)1)]
|
||||||
|
[InlineData((byte)7)]
|
||||||
|
[InlineData((byte)9)]
|
||||||
|
public void Three_a_one_b_pattern_interleaves_three_group_a_values_before_group_b(byte lastDigit)
|
||||||
|
{
|
||||||
|
ReadOnlySpan<uint> combined = HeuristicCombinedPrimeTester.GetCombinedDivisors(
|
||||||
|
lastDigit,
|
||||||
|
HeuristicCombinedPrimeTester.CombinedDivisorPattern.ThreeAOneB);
|
||||||
|
|
||||||
|
combined.Length.Should().BeGreaterThan(24);
|
||||||
|
combined[0].Should().Be(3U);
|
||||||
|
combined[1].Should().Be(7U);
|
||||||
|
combined[2].Should().Be(11U);
|
||||||
|
combined[3].Should().Be(13U);
|
||||||
|
|
||||||
|
int index = 4;
|
||||||
|
const int segmentsToCheck = 4;
|
||||||
|
|
||||||
|
for (int segment = 0; segment < segmentsToCheck; segment++)
|
||||||
|
{
|
||||||
|
for (int i = 0; i < 3; i++)
|
||||||
|
{
|
||||||
|
uint groupAValue = combined[index++];
|
||||||
|
groupAValue.Should().BeGreaterThan(13U);
|
||||||
|
(groupAValue % 10U).Should().Be(3U);
|
||||||
|
}
|
||||||
|
|
||||||
|
uint groupBValue = combined[index++];
|
||||||
|
groupBValue.Should().BeGreaterThan(13U);
|
||||||
|
(groupBValue % 10U).Should().NotBe(3U);
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
[Theory]
|
||||||
|
[Trait("Category", "Fast")]
|
||||||
|
[InlineData(0UL)]
|
||||||
|
[InlineData(1UL)]
|
||||||
|
[InlineData(2UL)]
|
||||||
|
[InlineData(3UL)]
|
||||||
|
[InlineData(4UL)]
|
||||||
|
[InlineData(5UL)]
|
||||||
|
[InlineData(8UL)]
|
||||||
|
[InlineData(9UL)]
|
||||||
|
[InlineData(15UL)]
|
||||||
|
[InlineData(16UL)]
|
||||||
|
[InlineData(17UL)]
|
||||||
|
[InlineData(24UL)]
|
||||||
|
[InlineData(25UL)]
|
||||||
|
[InlineData(26UL)]
|
||||||
|
[InlineData(255UL)]
|
||||||
|
[InlineData(256UL)]
|
||||||
|
[InlineData(257UL)]
|
||||||
|
[InlineData(1_000UL)]
|
||||||
|
[InlineData(1_048_575UL)]
|
||||||
|
[InlineData(1_048_576UL)]
|
||||||
|
[InlineData(1_048_577UL)]
|
||||||
|
[InlineData(4_000_000UL)]
|
||||||
|
[InlineData(4_000_001UL)]
|
||||||
|
[InlineData(4_000_004UL)]
|
||||||
|
[InlineData(65_535_999UL)]
|
||||||
|
[InlineData(65_536_000UL)]
|
||||||
|
[InlineData(65_536_001UL)]
|
||||||
|
[InlineData(4_294_967_295UL)]
|
||||||
|
[InlineData(4_294_967_296UL)]
|
||||||
|
[InlineData(4_294_967_297UL)]
|
||||||
|
[InlineData(18_446_744_073_709_551_615UL)]
|
||||||
|
[InlineData(18_446_744_073_709_551_614UL)]
|
||||||
|
[InlineData(18_446_744_073_709_551_613UL)]
|
||||||
|
public void ComputeHeuristicSqrt_returns_floor_square_root(ulong value)
|
||||||
|
{
|
||||||
|
ulong sqrt = HeuristicCombinedPrimeTester.ComputeHeuristicSqrt(value);
|
||||||
|
|
||||||
|
bool squareBelowOrEqual = (UInt128)sqrt * sqrt <= (UInt128)value;
|
||||||
|
bool nextSquareAbove = (UInt128)(sqrt + 1UL) * (UInt128)(sqrt + 1UL) > (UInt128)value;
|
||||||
|
|
||||||
|
squareBelowOrEqual.Should().BeTrue();
|
||||||
|
nextSquareAbove.Should().BeTrue();
|
||||||
|
}
|
||||||
|
}
|
||||||
@@ -26,7 +26,7 @@ public class PrimeTesterIsPrimeTests
|
|||||||
{
|
{
|
||||||
var tester = new PrimeTester();
|
var tester = new PrimeTester();
|
||||||
|
|
||||||
tester.IsPrime(n, CancellationToken.None).Should().Be(expected);
|
PrimeTester.IsPrime(n, CancellationToken.None).Should().Be(expected);
|
||||||
}
|
}
|
||||||
|
|
||||||
[Theory]
|
[Theory]
|
||||||
|
|||||||
@@ -36,15 +36,30 @@ public class PrimesGeneratorTests
|
|||||||
var lastThreePow2 = DivisorGenerator.SmallPrimesPow2LastThree;
|
var lastThreePow2 = DivisorGenerator.SmallPrimesPow2LastThree;
|
||||||
var lastNine = DivisorGenerator.SmallPrimesLastNine;
|
var lastNine = DivisorGenerator.SmallPrimesLastNine;
|
||||||
var lastNinePow2 = DivisorGenerator.SmallPrimesPow2LastNine;
|
var lastNinePow2 = DivisorGenerator.SmallPrimesPow2LastNine;
|
||||||
|
var lastOneWithoutThree = DivisorGenerator.SmallPrimesLastOneWithoutLastThree;
|
||||||
|
var lastOneWithoutThreePow2 = DivisorGenerator.SmallPrimesPow2LastOneWithoutLastThree;
|
||||||
|
var lastSevenWithoutThree = DivisorGenerator.SmallPrimesLastSevenWithoutLastThree;
|
||||||
|
var lastSevenWithoutThreePow2 = DivisorGenerator.SmallPrimesPow2LastSevenWithoutLastThree;
|
||||||
|
var lastNineWithoutThree = DivisorGenerator.SmallPrimesLastNineWithoutLastThree;
|
||||||
|
var lastNineWithoutThreePow2 = DivisorGenerator.SmallPrimesPow2LastNineWithoutLastThree;
|
||||||
|
|
||||||
lastOne.Length.Should().Be((int)PerfectNumberConstants.PrimesLimit);
|
lastOne.Length.Should().Be((int)PerfectNumberConstants.PrimesLimit);
|
||||||
lastSeven.Length.Should().Be((int)PerfectNumberConstants.PrimesLimit);
|
lastSeven.Length.Should().Be((int)PerfectNumberConstants.PrimesLimit);
|
||||||
lastThree.Length.Should().Be((int)PerfectNumberConstants.PrimesLimit);
|
lastThree.Length.Should().Be((int)PerfectNumberConstants.PrimesLimit);
|
||||||
lastNine.Length.Should().Be((int)PerfectNumberConstants.PrimesLimit);
|
lastNine.Length.Should().Be((int)PerfectNumberConstants.PrimesLimit);
|
||||||
|
lastOneWithoutThree.Length.Should().BeGreaterThan(0);
|
||||||
|
lastSevenWithoutThree.Length.Should().BeGreaterThan(0);
|
||||||
|
lastNineWithoutThree.Length.Should().BeGreaterThan(0);
|
||||||
|
lastOneWithoutThree.Length.Should().BeLessThan(lastOne.Length);
|
||||||
|
lastSevenWithoutThree.Length.Should().BeLessThan(lastSeven.Length);
|
||||||
|
lastNineWithoutThree.Length.Should().BeLessThan(lastNine.Length);
|
||||||
lastOnePow2.Length.Should().Be(lastOne.Length);
|
lastOnePow2.Length.Should().Be(lastOne.Length);
|
||||||
lastSevenPow2.Length.Should().Be(lastSeven.Length);
|
lastSevenPow2.Length.Should().Be(lastSeven.Length);
|
||||||
lastThreePow2.Length.Should().Be(lastThree.Length);
|
lastThreePow2.Length.Should().Be(lastThree.Length);
|
||||||
lastNinePow2.Length.Should().Be(lastNine.Length);
|
lastNinePow2.Length.Should().Be(lastNine.Length);
|
||||||
|
lastOneWithoutThreePow2.Length.Should().Be(lastOneWithoutThree.Length);
|
||||||
|
lastSevenWithoutThreePow2.Length.Should().Be(lastSevenWithoutThree.Length);
|
||||||
|
lastNineWithoutThreePow2.Length.Should().Be(lastNineWithoutThree.Length);
|
||||||
|
|
||||||
lastOne[..5].Should().Equal([3U, 7U, 11U, 13U, 19U]);
|
lastOne[..5].Should().Equal([3U, 7U, 11U, 13U, 19U]);
|
||||||
lastSeven[..5].Should().Equal([3U, 7U, 11U, 13U, 17U]);
|
lastSeven[..5].Should().Equal([3U, 7U, 11U, 13U, 17U]);
|
||||||
@@ -78,6 +93,47 @@ public class PrimesGeneratorTests
|
|||||||
lastThreePow2[i].Should().Be(prime * (ulong)prime);
|
lastThreePow2[i].Should().Be(prime * (ulong)prime);
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
|
lastOneWithoutThree.Should().NotBeEmpty();
|
||||||
|
lastSevenWithoutThree.Should().NotBeEmpty();
|
||||||
|
lastNineWithoutThree.Should().NotBeEmpty();
|
||||||
|
|
||||||
|
for (int i = 0; i < lastOneWithoutThree.Length; i++)
|
||||||
|
{
|
||||||
|
uint prime = lastOneWithoutThree[i];
|
||||||
|
(prime % 10U != 3U).Should().BeTrue();
|
||||||
|
if (i > 0)
|
||||||
|
{
|
||||||
|
prime.Should().BeGreaterThan(lastOneWithoutThree[i - 1]);
|
||||||
|
}
|
||||||
|
|
||||||
|
lastOneWithoutThreePow2[i].Should().Be(prime * (ulong)prime);
|
||||||
|
}
|
||||||
|
|
||||||
|
for (int i = 0; i < lastSevenWithoutThree.Length; i++)
|
||||||
|
{
|
||||||
|
uint prime = lastSevenWithoutThree[i];
|
||||||
|
(prime % 10U != 3U).Should().BeTrue();
|
||||||
|
if (i > 0)
|
||||||
|
{
|
||||||
|
prime.Should().BeGreaterThan(lastSevenWithoutThree[i - 1]);
|
||||||
|
}
|
||||||
|
|
||||||
|
lastSevenWithoutThreePow2[i].Should().Be(prime * (ulong)prime);
|
||||||
|
}
|
||||||
|
|
||||||
|
for (int i = 0; i < lastNineWithoutThree.Length; i++)
|
||||||
|
{
|
||||||
|
uint prime = lastNineWithoutThree[i];
|
||||||
|
(prime % 10U != 3U).Should().BeTrue();
|
||||||
|
if (i > 0)
|
||||||
|
{
|
||||||
|
prime.Should().BeGreaterThan(lastNineWithoutThree[i - 1]);
|
||||||
|
}
|
||||||
|
|
||||||
|
lastNineWithoutThreePow2[i].Should().Be(prime * (ulong)prime);
|
||||||
|
}
|
||||||
|
|
||||||
for (int i = 0; i < lastNine.Length; i++)
|
for (int i = 0; i < lastNine.Length; i++)
|
||||||
{
|
{
|
||||||
uint prime = lastNine[i];
|
uint prime = lastNine[i];
|
||||||
|
|||||||
@@ -16,6 +16,12 @@ internal static class DivisorGenerator
|
|||||||
public static readonly ulong[] SmallPrimesPow2LastThree;
|
public static readonly ulong[] SmallPrimesPow2LastThree;
|
||||||
public static readonly uint[] SmallPrimesLastNine;
|
public static readonly uint[] SmallPrimesLastNine;
|
||||||
public static readonly ulong[] SmallPrimesPow2LastNine;
|
public static readonly ulong[] SmallPrimesPow2LastNine;
|
||||||
|
public static readonly uint[] SmallPrimesLastOneWithoutLastThree;
|
||||||
|
public static readonly ulong[] SmallPrimesPow2LastOneWithoutLastThree;
|
||||||
|
public static readonly uint[] SmallPrimesLastSevenWithoutLastThree;
|
||||||
|
public static readonly ulong[] SmallPrimesPow2LastSevenWithoutLastThree;
|
||||||
|
public static readonly uint[] SmallPrimesLastNineWithoutLastThree;
|
||||||
|
public static readonly ulong[] SmallPrimesPow2LastNineWithoutLastThree;
|
||||||
|
|
||||||
private const ushort DecimalMaskWhenLastIsSeven = (1 << 3) | (1 << 7) | (1 << 9);
|
private const ushort DecimalMaskWhenLastIsSeven = (1 << 3) | (1 << 7) | (1 << 9);
|
||||||
private const ushort DecimalMaskWhenLastIsThree = (1 << 1) | (1 << 3) | (1 << 7) | (1 << 9);
|
private const ushort DecimalMaskWhenLastIsThree = (1 << 1) | (1 << 3) | (1 << 7) | (1 << 9);
|
||||||
@@ -34,7 +40,13 @@ internal static class DivisorGenerator
|
|||||||
out SmallPrimesLastThree,
|
out SmallPrimesLastThree,
|
||||||
out SmallPrimesPow2LastThree,
|
out SmallPrimesPow2LastThree,
|
||||||
out SmallPrimesLastNine,
|
out SmallPrimesLastNine,
|
||||||
out SmallPrimesPow2LastNine);
|
out SmallPrimesPow2LastNine,
|
||||||
|
out SmallPrimesLastOneWithoutLastThree,
|
||||||
|
out SmallPrimesPow2LastOneWithoutLastThree,
|
||||||
|
out SmallPrimesLastSevenWithoutLastThree,
|
||||||
|
out SmallPrimesPow2LastSevenWithoutLastThree,
|
||||||
|
out SmallPrimesLastNineWithoutLastThree,
|
||||||
|
out SmallPrimesPow2LastNineWithoutLastThree);
|
||||||
}
|
}
|
||||||
|
|
||||||
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||||
@@ -59,7 +71,13 @@ internal static class DivisorGenerator
|
|||||||
out uint[] lastThree,
|
out uint[] lastThree,
|
||||||
out ulong[] lastThreePow2,
|
out ulong[] lastThreePow2,
|
||||||
out uint[] lastNine,
|
out uint[] lastNine,
|
||||||
out ulong[] lastNinePow2)
|
out ulong[] lastNinePow2,
|
||||||
|
out uint[] lastOneWithoutThree,
|
||||||
|
out ulong[] lastOneWithoutThreePow2,
|
||||||
|
out uint[] lastSevenWithoutThree,
|
||||||
|
out ulong[] lastSevenWithoutThreePow2,
|
||||||
|
out uint[] lastNineWithoutThree,
|
||||||
|
out ulong[] lastNineWithoutThreePow2)
|
||||||
{
|
{
|
||||||
uint target = PerfectNumberConstants.PrimesLimit;
|
uint target = PerfectNumberConstants.PrimesLimit;
|
||||||
int targetInt = checked((int)target);
|
int targetInt = checked((int)target);
|
||||||
@@ -75,6 +93,13 @@ internal static class DivisorGenerator
|
|||||||
lastNinePow2 = new ulong[targetInt];
|
lastNinePow2 = new ulong[targetInt];
|
||||||
|
|
||||||
List<uint> primes = new(targetInt * 4 / 3 + 16);
|
List<uint> primes = new(targetInt * 4 / 3 + 16);
|
||||||
|
List<uint> lastOneWithoutThreeList = new(targetInt);
|
||||||
|
List<ulong> lastOneWithoutThreePow2List = new(targetInt);
|
||||||
|
List<uint> lastSevenWithoutThreeList = new(targetInt);
|
||||||
|
List<ulong> lastSevenWithoutThreePow2List = new(targetInt);
|
||||||
|
List<uint> lastNineWithoutThreeList = new(targetInt);
|
||||||
|
List<ulong> lastNineWithoutThreePow2List = new(targetInt);
|
||||||
|
|
||||||
uint candidate = 2U;
|
uint candidate = 2U;
|
||||||
int allCount = 0;
|
int allCount = 0;
|
||||||
int lastOneCount = 0;
|
int lastOneCount = 0;
|
||||||
@@ -82,24 +107,22 @@ internal static class DivisorGenerator
|
|||||||
int lastThreeCount = 0;
|
int lastThreeCount = 0;
|
||||||
int lastNineCount = 0;
|
int lastNineCount = 0;
|
||||||
|
|
||||||
while (allCount < targetInt || lastOneCount < targetInt || lastSevenCount < targetInt || lastThreeCount < targetInt || lastNineCount < targetInt)
|
while (allCount < targetInt ||
|
||||||
|
lastOneCount < targetInt ||
|
||||||
|
lastSevenCount < targetInt ||
|
||||||
|
lastThreeCount < targetInt ||
|
||||||
|
lastNineCount < targetInt)
|
||||||
{
|
{
|
||||||
bool isPrime = true;
|
bool isPrime = true;
|
||||||
int primesCount = primes.Count;
|
int primesCount = primes.Count;
|
||||||
for (int i = 0; i < primesCount; i++)
|
for (int i = 0; i < primesCount; i++)
|
||||||
{
|
{
|
||||||
uint p = primes[i];
|
uint p = primes[i];
|
||||||
// TODO: Reuse a single squared local for `p` within this loop so we avoid recalculating
|
|
||||||
// `p * p` repeatedly without storing every squared prime globally (the cache must remain
|
|
||||||
// limited to the small-divisor tables to honor the memory constraints).
|
|
||||||
if (p * p > candidate)
|
if (p * p > candidate)
|
||||||
{
|
{
|
||||||
break;
|
break;
|
||||||
}
|
}
|
||||||
|
|
||||||
// TODO: Replace this trial-division `%` with the sieve-based generator that avoids
|
|
||||||
// per-candidate modulo work so building the small-prime tables stops dominating
|
|
||||||
// startup time for large scans.
|
|
||||||
if (candidate % p == 0U)
|
if (candidate % p == 0U)
|
||||||
{
|
{
|
||||||
isPrime = false;
|
isPrime = false;
|
||||||
@@ -110,54 +133,80 @@ internal static class DivisorGenerator
|
|||||||
if (isPrime)
|
if (isPrime)
|
||||||
{
|
{
|
||||||
primes.Add(candidate);
|
primes.Add(candidate);
|
||||||
|
ulong candidateSquared = checked((ulong)candidate * candidate);
|
||||||
if (allCount < targetInt)
|
if (allCount < targetInt)
|
||||||
{
|
{
|
||||||
// TODO: Reuse a single cached `candidateSquared` per iteration so we avoid issuing
|
|
||||||
// three separate 64-bit multiplications when populating the pow2 arrays here and
|
|
||||||
// in the last-one/last-seven branches below.
|
|
||||||
all[allCount] = candidate;
|
all[allCount] = candidate;
|
||||||
allPow2[allCount] = checked((ulong)candidate * candidate);
|
allPow2[allCount] = candidateSquared;
|
||||||
allCount++;
|
allCount++;
|
||||||
}
|
}
|
||||||
|
|
||||||
if (candidate != 2U)
|
if (candidate != 2U)
|
||||||
{
|
{
|
||||||
LastDigit lastDigit = GetLastDigit(candidate);
|
LastDigit lastDigit = GetLastDigit(candidate);
|
||||||
|
bool allowedForLastOne = IsAllowedForLastOne(candidate, lastDigit);
|
||||||
|
bool allowedForLastSeven = IsAllowedForLastSeven(candidate, lastDigit);
|
||||||
|
bool allowedForLastThree = IsAllowedForLastThree(candidate, lastDigit);
|
||||||
|
bool allowedForLastNine = IsAllowedForLastNine(candidate, lastDigit);
|
||||||
|
bool lastDigitIsThree = lastDigit == LastDigit.Three;
|
||||||
|
|
||||||
if (lastOneCount < targetInt && IsAllowedForLastOne(candidate, lastDigit))
|
if (allowedForLastOne && lastOneCount < targetInt)
|
||||||
{
|
{
|
||||||
lastOne[lastOneCount] = candidate;
|
lastOne[lastOneCount] = candidate;
|
||||||
lastOnePow2[lastOneCount] = checked((ulong)candidate * candidate);
|
lastOnePow2[lastOneCount] = candidateSquared;
|
||||||
lastOneCount++;
|
lastOneCount++;
|
||||||
}
|
}
|
||||||
|
|
||||||
if (lastSevenCount < targetInt && IsAllowedForLastSeven(candidate, lastDigit))
|
if (allowedForLastSeven && lastSevenCount < targetInt)
|
||||||
{
|
{
|
||||||
lastSeven[lastSevenCount] = candidate;
|
lastSeven[lastSevenCount] = candidate;
|
||||||
lastSevenPow2[lastSevenCount] = checked((ulong)candidate * candidate);
|
lastSevenPow2[lastSevenCount] = candidateSquared;
|
||||||
lastSevenCount++;
|
lastSevenCount++;
|
||||||
}
|
}
|
||||||
|
|
||||||
if (lastThreeCount < targetInt && IsAllowedForLastThree(candidate, lastDigit))
|
if (allowedForLastThree && lastThreeCount < targetInt)
|
||||||
{
|
{
|
||||||
lastThree[lastThreeCount] = candidate;
|
lastThree[lastThreeCount] = candidate;
|
||||||
lastThreePow2[lastThreeCount] = checked((ulong)candidate * candidate);
|
lastThreePow2[lastThreeCount] = candidateSquared;
|
||||||
lastThreeCount++;
|
lastThreeCount++;
|
||||||
}
|
}
|
||||||
|
|
||||||
if (lastNineCount < targetInt && IsAllowedForLastNine(candidate, lastDigit))
|
if (allowedForLastNine && lastNineCount < targetInt)
|
||||||
{
|
{
|
||||||
lastNine[lastNineCount] = candidate;
|
lastNine[lastNineCount] = candidate;
|
||||||
lastNinePow2[lastNineCount] = checked((ulong)candidate * candidate);
|
lastNinePow2[lastNineCount] = candidateSquared;
|
||||||
lastNineCount++;
|
lastNineCount++;
|
||||||
}
|
}
|
||||||
|
|
||||||
|
if (!lastDigitIsThree && allowedForLastOne)
|
||||||
|
{
|
||||||
|
lastOneWithoutThreeList.Add(candidate);
|
||||||
|
lastOneWithoutThreePow2List.Add(candidateSquared);
|
||||||
|
}
|
||||||
|
|
||||||
|
if (!lastDigitIsThree && allowedForLastSeven)
|
||||||
|
{
|
||||||
|
lastSevenWithoutThreeList.Add(candidate);
|
||||||
|
lastSevenWithoutThreePow2List.Add(candidateSquared);
|
||||||
|
}
|
||||||
|
|
||||||
|
if (!lastDigitIsThree && allowedForLastNine)
|
||||||
|
{
|
||||||
|
lastNineWithoutThreeList.Add(candidate);
|
||||||
|
lastNineWithoutThreePow2List.Add(candidateSquared);
|
||||||
|
}
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
// TODO: Switch this increment to the Mod6 stepping table identified in the Mod6ComparisonBenchmarks
|
|
||||||
// so we skip the composite residues without relying on per-iteration `%` checks.
|
|
||||||
candidate = candidate == 2U ? 3U : candidate + 2U;
|
candidate = candidate == 2U ? 3U : candidate + 2U;
|
||||||
}
|
}
|
||||||
|
|
||||||
|
lastOneWithoutThree = lastOneWithoutThreeList.ToArray();
|
||||||
|
lastOneWithoutThreePow2 = lastOneWithoutThreePow2List.ToArray();
|
||||||
|
lastSevenWithoutThree = lastSevenWithoutThreeList.ToArray();
|
||||||
|
lastSevenWithoutThreePow2 = lastSevenWithoutThreePow2List.ToArray();
|
||||||
|
lastNineWithoutThree = lastNineWithoutThreeList.ToArray();
|
||||||
|
lastNineWithoutThreePow2 = lastNineWithoutThreePow2List.ToArray();
|
||||||
}
|
}
|
||||||
|
|
||||||
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||||
|
|||||||
@@ -3,7 +3,6 @@ namespace PerfectNumbers.Core.Gpu;
|
|||||||
public static class GpuPrimeWorkLimiter
|
public static class GpuPrimeWorkLimiter
|
||||||
{
|
{
|
||||||
private static SemaphoreSlim _semaphore = new(1, int.MaxValue);
|
private static SemaphoreSlim _semaphore = new(1, int.MaxValue);
|
||||||
private static readonly object _lock = new();
|
|
||||||
private static int _currentLimit = 1;
|
private static int _currentLimit = 1;
|
||||||
|
|
||||||
public static Lease Acquire()
|
public static Lease Acquire()
|
||||||
@@ -11,18 +10,11 @@ public static class GpuPrimeWorkLimiter
|
|||||||
// TODO: Inline the pooled limiter guard from the GpuLimiterThroughputBenchmarks so prime-sieve
|
// TODO: Inline the pooled limiter guard from the GpuLimiterThroughputBenchmarks so prime-sieve
|
||||||
// GPU jobs stop allocating Lease instances on every acquisition.
|
// GPU jobs stop allocating Lease instances on every acquisition.
|
||||||
var sem = _semaphore;
|
var sem = _semaphore;
|
||||||
sem.Wait();
|
// sem.Wait();
|
||||||
return new Lease(sem);
|
return new Lease(sem);
|
||||||
}
|
}
|
||||||
|
|
||||||
public static void SetLimit(int value)
|
public static void SetLimit(int value)
|
||||||
{
|
|
||||||
if (value < 1)
|
|
||||||
{
|
|
||||||
value = 1;
|
|
||||||
}
|
|
||||||
|
|
||||||
lock (_lock)
|
|
||||||
{
|
{
|
||||||
if (value == _currentLimit)
|
if (value == _currentLimit)
|
||||||
{
|
{
|
||||||
@@ -34,7 +26,6 @@ public static class GpuPrimeWorkLimiter
|
|||||||
_semaphore = new SemaphoreSlim(value, value);
|
_semaphore = new SemaphoreSlim(value, value);
|
||||||
_currentLimit = value;
|
_currentLimit = value;
|
||||||
}
|
}
|
||||||
}
|
|
||||||
|
|
||||||
public sealed class Lease
|
public sealed class Lease
|
||||||
{
|
{
|
||||||
@@ -47,7 +38,7 @@ public static class GpuPrimeWorkLimiter
|
|||||||
|
|
||||||
public void Dispose()
|
public void Dispose()
|
||||||
{
|
{
|
||||||
_sem.Release();
|
// _sem.Release();
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|||||||
@@ -94,7 +94,7 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
|
|||||||
var gpuLease = RentAccelerator();
|
var gpuLease = RentAccelerator();
|
||||||
var accelerator = gpuLease.Accelerator;
|
var accelerator = gpuLease.Accelerator;
|
||||||
|
|
||||||
Monitor.Enter(gpuLease.ExecutionLock);
|
// Monitor.Enter(gpuLease.ExecutionLock);
|
||||||
|
|
||||||
var kernel = GetKernel(accelerator);
|
var kernel = GetKernel(accelerator);
|
||||||
var resources = RentBatchResources(accelerator, batchCapacity);
|
var resources = RentBatchResources(accelerator, batchCapacity);
|
||||||
@@ -121,7 +121,7 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
|
|||||||
out processedCount);
|
out processedCount);
|
||||||
|
|
||||||
ReturnBatchResources(resources);
|
ReturnBatchResources(resources);
|
||||||
Monitor.Exit(gpuLease.ExecutionLock);
|
// Monitor.Exit(gpuLease.ExecutionLock);
|
||||||
ReturnAccelerator(gpuLease);
|
ReturnAccelerator(gpuLease);
|
||||||
|
|
||||||
if (composite)
|
if (composite)
|
||||||
@@ -560,7 +560,7 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
Monitor.Enter(_lease.ExecutionLock);
|
// Monitor.Enter(_lease.ExecutionLock);
|
||||||
|
|
||||||
EnsureExecutionResourcesLocked(length);
|
EnsureExecutionResourcesLocked(length);
|
||||||
|
|
||||||
@@ -601,7 +601,7 @@ public sealed class MersenneNumberDivisorByDivisorGpuTester : IMersenneNumberDiv
|
|||||||
ref byte hitRef = ref MemoryMarshal.GetReference(hitSlice);
|
ref byte hitRef = ref MemoryMarshal.GetReference(hitSlice);
|
||||||
hitView.CopyToCPU(ref hitRef, length);
|
hitView.CopyToCPU(ref hitRef, length);
|
||||||
|
|
||||||
Monitor.Exit(_lease.ExecutionLock);
|
// Monitor.Exit(_lease.ExecutionLock);
|
||||||
}
|
}
|
||||||
|
|
||||||
public void Dispose()
|
public void Dispose()
|
||||||
|
|||||||
@@ -30,7 +30,7 @@ public class MersenneNumberResidueGpuTester(bool useGpuOrder)
|
|||||||
var gpuLease = RentAccelerator();
|
var gpuLease = RentAccelerator();
|
||||||
var accelerator = gpuLease.Accelerator;
|
var accelerator = gpuLease.Accelerator;
|
||||||
|
|
||||||
Monitor.Enter(gpuLease.ExecutionLock);
|
// Monitor.Enter(gpuLease.ExecutionLock);
|
||||||
|
|
||||||
var stream = accelerator.CreateStream();
|
var stream = accelerator.CreateStream();
|
||||||
var kernel = GetPow2ModKernel(accelerator);
|
var kernel = GetPow2ModKernel(accelerator);
|
||||||
@@ -108,7 +108,7 @@ public class MersenneNumberResidueGpuTester(bool useGpuOrder)
|
|||||||
|
|
||||||
ReturnResources(accelerator, resources);
|
ReturnResources(accelerator, resources);
|
||||||
stream.Dispose();
|
stream.Dispose();
|
||||||
Monitor.Exit(gpuLease.ExecutionLock);
|
// Monitor.Exit(gpuLease.ExecutionLock);
|
||||||
ReturnAccelerator(gpuLease);
|
ReturnAccelerator(gpuLease);
|
||||||
limiter.Dispose();
|
limiter.Dispose();
|
||||||
}
|
}
|
||||||
|
|||||||
@@ -0,0 +1,667 @@
|
|||||||
|
using System;
|
||||||
|
using Open.Numeric.Primes;
|
||||||
|
using System.Collections.Generic;
|
||||||
|
using System.Numerics;
|
||||||
|
using System.Runtime.CompilerServices;
|
||||||
|
using System.Threading;
|
||||||
|
using ILGPU;
|
||||||
|
using ILGPU.Algorithms;
|
||||||
|
using ILGPU.Runtime;
|
||||||
|
using PerfectNumbers.Core.Gpu;
|
||||||
|
|
||||||
|
namespace PerfectNumbers.Core;
|
||||||
|
|
||||||
|
public sealed class HeuristicCombinedPrimeTester
|
||||||
|
{
|
||||||
|
static HeuristicCombinedPrimeTester()
|
||||||
|
{
|
||||||
|
HeuristicPrimeSieves.EnsureInitialized();
|
||||||
|
|
||||||
|
CombinedDivisorsEnding1TwoAOneB = BuildCombinedDivisors(1, CombinedDivisorPattern.TwoAOneB);
|
||||||
|
CombinedDivisorsEnding3TwoAOneB = BuildCombinedDivisors(3, CombinedDivisorPattern.TwoAOneB);
|
||||||
|
CombinedDivisorsEnding7TwoAOneB = BuildCombinedDivisors(7, CombinedDivisorPattern.TwoAOneB);
|
||||||
|
CombinedDivisorsEnding9TwoAOneB = BuildCombinedDivisors(9, CombinedDivisorPattern.TwoAOneB);
|
||||||
|
|
||||||
|
CombinedDivisorsEnding1OneAOneB = BuildCombinedDivisors(1, CombinedDivisorPattern.OneAOneB);
|
||||||
|
CombinedDivisorsEnding3OneAOneB = BuildCombinedDivisors(3, CombinedDivisorPattern.OneAOneB);
|
||||||
|
CombinedDivisorsEnding7OneAOneB = BuildCombinedDivisors(7, CombinedDivisorPattern.OneAOneB);
|
||||||
|
CombinedDivisorsEnding9OneAOneB = BuildCombinedDivisors(9, CombinedDivisorPattern.OneAOneB);
|
||||||
|
|
||||||
|
CombinedDivisorsEnding1ThreeAOneB = BuildCombinedDivisors(1, CombinedDivisorPattern.ThreeAOneB);
|
||||||
|
CombinedDivisorsEnding3ThreeAOneB = BuildCombinedDivisors(3, CombinedDivisorPattern.ThreeAOneB);
|
||||||
|
CombinedDivisorsEnding7ThreeAOneB = BuildCombinedDivisors(7, CombinedDivisorPattern.ThreeAOneB);
|
||||||
|
CombinedDivisorsEnding9ThreeAOneB = BuildCombinedDivisors(9, CombinedDivisorPattern.ThreeAOneB);
|
||||||
|
|
||||||
|
CombinedDivisorsEnding1ThreeATwoB = BuildCombinedDivisors(1, CombinedDivisorPattern.ThreeATwoB);
|
||||||
|
CombinedDivisorsEnding3ThreeATwoB = BuildCombinedDivisors(3, CombinedDivisorPattern.ThreeATwoB);
|
||||||
|
CombinedDivisorsEnding7ThreeATwoB = BuildCombinedDivisors(7, CombinedDivisorPattern.ThreeATwoB);
|
||||||
|
CombinedDivisorsEnding9ThreeATwoB = BuildCombinedDivisors(9, CombinedDivisorPattern.ThreeATwoB);
|
||||||
|
CombinedDivisorsEnding1 = CombinedDivisorsEnding1OneAOneB;
|
||||||
|
CombinedDivisorsEnding3 = CombinedDivisorsEnding3OneAOneB;
|
||||||
|
CombinedDivisorsEnding7 = CombinedDivisorsEnding7OneAOneB;
|
||||||
|
CombinedDivisorsEnding9 = CombinedDivisorsEnding9OneAOneB;
|
||||||
|
}
|
||||||
|
|
||||||
|
private const ulong Wheel210 = 210UL;
|
||||||
|
|
||||||
|
private static readonly bool UseHeuristicGroupBTrialDivision = false; // Temporary fallback gate for Group B.
|
||||||
|
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding1TwoAOneB;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding3TwoAOneB;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding7TwoAOneB;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding9TwoAOneB;
|
||||||
|
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding1OneAOneB;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding3OneAOneB;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding7OneAOneB;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding9OneAOneB;
|
||||||
|
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding1ThreeAOneB;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding3ThreeAOneB;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding7ThreeAOneB;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding9ThreeAOneB;
|
||||||
|
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding1ThreeATwoB;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding3ThreeATwoB;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding7ThreeATwoB;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding9ThreeATwoB;
|
||||||
|
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding1;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding3;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding7;
|
||||||
|
private static readonly uint[] CombinedDivisorsEnding9;
|
||||||
|
|
||||||
|
private static readonly ulong[] HeuristicSmallCycleSnapshot = MersenneDivisorCycles.Shared.ExportSmallCyclesSnapshot();
|
||||||
|
|
||||||
|
internal enum CombinedDivisorPattern : byte
|
||||||
|
{
|
||||||
|
TwoAOneB = 0,
|
||||||
|
OneAOneB = 1,
|
||||||
|
ThreeAOneB = 2,
|
||||||
|
ThreeATwoB = 3,
|
||||||
|
}
|
||||||
|
|
||||||
|
internal enum HeuristicDivisorGroup : byte
|
||||||
|
{
|
||||||
|
None = 0,
|
||||||
|
GroupAConstant = 1,
|
||||||
|
GroupAWheel = 2,
|
||||||
|
GroupB = 3,
|
||||||
|
}
|
||||||
|
|
||||||
|
internal readonly struct HeuristicDivisorCandidate(ulong value, HeuristicDivisorGroup group, byte ending, byte priorityIndex, ushort wheelResidue)
|
||||||
|
{
|
||||||
|
public readonly ulong Value = value;
|
||||||
|
public readonly HeuristicDivisorGroup Group = group;
|
||||||
|
public readonly byte Ending = ending;
|
||||||
|
public readonly byte PriorityIndex = priorityIndex;
|
||||||
|
public readonly ushort WheelResidue = wheelResidue;
|
||||||
|
}
|
||||||
|
|
||||||
|
internal readonly struct HeuristicDivisorPreparation(
|
||||||
|
in HeuristicDivisorCandidate candidate,
|
||||||
|
in MontgomeryDivisorData divisorData,
|
||||||
|
ulong cycleLengthHint,
|
||||||
|
bool hasCycleLengthHint)
|
||||||
|
{
|
||||||
|
public readonly HeuristicDivisorCandidate Candidate = candidate;
|
||||||
|
public readonly MontgomeryDivisorData DivisorData = divisorData;
|
||||||
|
public readonly ulong CycleLengthHint = cycleLengthHint;
|
||||||
|
public readonly bool HasCycleLengthHint = hasCycleLengthHint;
|
||||||
|
public readonly bool RequiresCycleComputation = !hasCycleLengthHint;
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
internal static int HeuristicGpuDivisorBatchSize = 4_096;
|
||||||
|
internal static int HeuristicDivisorInterleaveBatchSize = 64;
|
||||||
|
|
||||||
|
[ThreadStatic]
|
||||||
|
private static HeuristicCombinedPrimeTester? _tester;
|
||||||
|
|
||||||
|
public static HeuristicCombinedPrimeTester Exclusive => _tester ??= new();
|
||||||
|
|
||||||
|
public static bool IsPrimeCpu(ulong n)
|
||||||
|
{
|
||||||
|
byte nMod10 = (byte)n.Mod10();
|
||||||
|
ulong sqrtLimit = ComputeHeuristicSqrt(n);
|
||||||
|
|
||||||
|
if (sqrtLimit < 3UL)
|
||||||
|
{
|
||||||
|
return EvaluateWithOpenNumericFallback(n);
|
||||||
|
}
|
||||||
|
|
||||||
|
return HeuristicTrialDivisionCpu(n, sqrtLimit, nMod10);
|
||||||
|
}
|
||||||
|
|
||||||
|
public static bool IsPrimeGpu(ulong n)
|
||||||
|
{
|
||||||
|
byte nMod10 = (byte)n.Mod10();
|
||||||
|
ulong sqrtLimit = ComputeHeuristicSqrt(n);
|
||||||
|
return IsPrimeGpu(n, sqrtLimit, nMod10);
|
||||||
|
}
|
||||||
|
|
||||||
|
public static bool IsPrimeGpu(ulong n, ulong sqrtLimit, byte nMod10)
|
||||||
|
{
|
||||||
|
// TODO: Is this condition ever met on the execution path in EvenPerfectBitScanner?
|
||||||
|
if (sqrtLimit < 3UL)
|
||||||
|
{
|
||||||
|
return EvaluateWithOpenNumericFallback(n);
|
||||||
|
}
|
||||||
|
|
||||||
|
bool compositeDetected = HeuristicTrialDivisionGpuDetectsDivisor(n, sqrtLimit, nMod10);
|
||||||
|
return !compositeDetected;
|
||||||
|
}
|
||||||
|
|
||||||
|
private static bool HeuristicTrialDivisionCpu(ulong n, ulong sqrtLimit, byte nMod10)
|
||||||
|
{
|
||||||
|
ReadOnlySpan<uint> combinedDivisors = GetCombinedDivisors(nMod10);
|
||||||
|
if (combinedDivisors.IsEmpty)
|
||||||
|
{
|
||||||
|
return EvaluateWithOpenNumericFallback(n);
|
||||||
|
}
|
||||||
|
|
||||||
|
int length = combinedDivisors.Length;
|
||||||
|
for (int i = 0; i < length; i++)
|
||||||
|
{
|
||||||
|
uint entry = combinedDivisors[i];
|
||||||
|
if (entry > sqrtLimit)
|
||||||
|
{
|
||||||
|
break;
|
||||||
|
}
|
||||||
|
|
||||||
|
if (n % entry == 0UL)
|
||||||
|
{
|
||||||
|
return false;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
return EvaluateWithOpenNumericFallback(n);
|
||||||
|
}
|
||||||
|
|
||||||
|
private static bool HeuristicTrialDivisionGpuDetectsDivisor(ulong n, ulong sqrtLimit, byte nMod10)
|
||||||
|
{
|
||||||
|
int batchCapacity = HeuristicGpuDivisorBatchSize;
|
||||||
|
var divisorPool = ThreadStaticPools.UlongPool;
|
||||||
|
var hitPool = ThreadStaticPools.BytePool;
|
||||||
|
|
||||||
|
ulong[] divisorArray = divisorPool.Rent(batchCapacity);
|
||||||
|
byte[] hitFlags = hitPool.Rent(batchCapacity);
|
||||||
|
|
||||||
|
var limiter = GpuPrimeWorkLimiter.Acquire();
|
||||||
|
var gpu = PrimeTester.PrimeTesterGpuContextPool.Rent();
|
||||||
|
var accelerator = gpu.Accelerator;
|
||||||
|
var state = PrimeTester.GpuKernelState.GetOrCreate(accelerator);
|
||||||
|
var scratch = state.RentScratch(batchCapacity, accelerator);
|
||||||
|
|
||||||
|
bool compositeDetected = false;
|
||||||
|
int count = 0;
|
||||||
|
|
||||||
|
var input = scratch.Input;
|
||||||
|
var output = scratch.Output;
|
||||||
|
|
||||||
|
bool ProcessBatch(int length)
|
||||||
|
{
|
||||||
|
input.View.CopyFromCPU(ref divisorArray[0], length);
|
||||||
|
state.HeuristicTrialDivisionKernel(length, input.View, n, output.View);
|
||||||
|
accelerator.Synchronize();
|
||||||
|
output.View.CopyToCPU(ref hitFlags[0], length);
|
||||||
|
|
||||||
|
for (int i = 0; i < length; i++)
|
||||||
|
{
|
||||||
|
if (hitFlags[i] == 0)
|
||||||
|
{
|
||||||
|
continue;
|
||||||
|
}
|
||||||
|
|
||||||
|
return true;
|
||||||
|
// ulong divisor = divisorArray[i];
|
||||||
|
// if (divisor > 1UL && n % divisor == 0UL)
|
||||||
|
// {
|
||||||
|
// return true;
|
||||||
|
// }
|
||||||
|
}
|
||||||
|
|
||||||
|
return false;
|
||||||
|
}
|
||||||
|
|
||||||
|
ReadOnlySpan<uint> combinedDivisors = GetCombinedDivisors(nMod10);
|
||||||
|
if (combinedDivisors.IsEmpty)
|
||||||
|
{
|
||||||
|
goto Cleanup;
|
||||||
|
}
|
||||||
|
|
||||||
|
int length = combinedDivisors.Length;
|
||||||
|
for (int i = 0; i < length; i++)
|
||||||
|
{
|
||||||
|
uint entry = combinedDivisors[i];
|
||||||
|
if (entry > sqrtLimit)
|
||||||
|
{
|
||||||
|
break;
|
||||||
|
}
|
||||||
|
|
||||||
|
divisorArray[count] = entry;
|
||||||
|
count++;
|
||||||
|
|
||||||
|
if (count == batchCapacity)
|
||||||
|
{
|
||||||
|
if (ProcessBatch(count))
|
||||||
|
{
|
||||||
|
compositeDetected = true;
|
||||||
|
goto Cleanup;
|
||||||
|
}
|
||||||
|
|
||||||
|
count = 0;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
if (!compositeDetected && count > 0)
|
||||||
|
{
|
||||||
|
if (ProcessBatch(count))
|
||||||
|
{
|
||||||
|
compositeDetected = true;
|
||||||
|
goto Cleanup;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
Cleanup:
|
||||||
|
state.ReturnScratch(scratch);
|
||||||
|
divisorPool.Return(divisorArray);
|
||||||
|
hitPool.Return(hitFlags);
|
||||||
|
gpu.Dispose();
|
||||||
|
limiter.Dispose();
|
||||||
|
return compositeDetected;
|
||||||
|
}
|
||||||
|
|
||||||
|
private static bool EvaluateWithOpenNumericFallback(ulong n)
|
||||||
|
{
|
||||||
|
return Prime.Numbers.IsPrime(n);
|
||||||
|
}
|
||||||
|
|
||||||
|
internal static ulong ComputeHeuristicSqrt(ulong n)
|
||||||
|
{
|
||||||
|
ulong sqrt = (ulong)Math.Sqrt(n);
|
||||||
|
UInt128 square = (UInt128)sqrt * sqrt;
|
||||||
|
|
||||||
|
while (square > n)
|
||||||
|
{
|
||||||
|
sqrt--;
|
||||||
|
square = (UInt128)sqrt * sqrt;
|
||||||
|
}
|
||||||
|
|
||||||
|
ulong next = sqrt + 1UL;
|
||||||
|
UInt128 nextSquare = (UInt128)next * next;
|
||||||
|
while (nextSquare <= n)
|
||||||
|
{
|
||||||
|
sqrt = next;
|
||||||
|
next++;
|
||||||
|
nextSquare = (UInt128)next * next;
|
||||||
|
}
|
||||||
|
|
||||||
|
return sqrt;
|
||||||
|
}
|
||||||
|
|
||||||
|
private static ReadOnlySpan<uint> GetCombinedDivisors(byte nMod10)
|
||||||
|
{
|
||||||
|
return nMod10 switch
|
||||||
|
{
|
||||||
|
1 => CombinedDivisorsEnding1,
|
||||||
|
3 => CombinedDivisorsEnding3,
|
||||||
|
7 => CombinedDivisorsEnding7,
|
||||||
|
9 => CombinedDivisorsEnding9,
|
||||||
|
_ => ReadOnlySpan<uint>.Empty,
|
||||||
|
};
|
||||||
|
}
|
||||||
|
|
||||||
|
internal static ReadOnlySpan<uint> GetCombinedDivisors(byte nMod10, CombinedDivisorPattern pattern)
|
||||||
|
{
|
||||||
|
return pattern switch
|
||||||
|
{
|
||||||
|
CombinedDivisorPattern.TwoAOneB => nMod10 switch
|
||||||
|
{
|
||||||
|
1 => CombinedDivisorsEnding1TwoAOneB,
|
||||||
|
3 => CombinedDivisorsEnding3TwoAOneB,
|
||||||
|
7 => CombinedDivisorsEnding7TwoAOneB,
|
||||||
|
9 => CombinedDivisorsEnding9TwoAOneB,
|
||||||
|
_ => ReadOnlySpan<uint>.Empty,
|
||||||
|
},
|
||||||
|
CombinedDivisorPattern.OneAOneB => nMod10 switch
|
||||||
|
{
|
||||||
|
1 => CombinedDivisorsEnding1OneAOneB,
|
||||||
|
3 => CombinedDivisorsEnding3OneAOneB,
|
||||||
|
7 => CombinedDivisorsEnding7OneAOneB,
|
||||||
|
9 => CombinedDivisorsEnding9OneAOneB,
|
||||||
|
_ => ReadOnlySpan<uint>.Empty,
|
||||||
|
},
|
||||||
|
CombinedDivisorPattern.ThreeAOneB => nMod10 switch
|
||||||
|
{
|
||||||
|
1 => CombinedDivisorsEnding1ThreeAOneB,
|
||||||
|
3 => CombinedDivisorsEnding3ThreeAOneB,
|
||||||
|
7 => CombinedDivisorsEnding7ThreeAOneB,
|
||||||
|
9 => CombinedDivisorsEnding9ThreeAOneB,
|
||||||
|
_ => ReadOnlySpan<uint>.Empty,
|
||||||
|
},
|
||||||
|
CombinedDivisorPattern.ThreeATwoB => nMod10 switch
|
||||||
|
{
|
||||||
|
1 => CombinedDivisorsEnding1ThreeATwoB,
|
||||||
|
3 => CombinedDivisorsEnding3ThreeATwoB,
|
||||||
|
7 => CombinedDivisorsEnding7ThreeATwoB,
|
||||||
|
9 => CombinedDivisorsEnding9ThreeATwoB,
|
||||||
|
_ => ReadOnlySpan<uint>.Empty,
|
||||||
|
},
|
||||||
|
_ => ReadOnlySpan<uint>.Empty,
|
||||||
|
};
|
||||||
|
}
|
||||||
|
|
||||||
|
private static uint[] BuildCombinedDivisors(byte nMod10, CombinedDivisorPattern pattern)
|
||||||
|
{
|
||||||
|
ReadOnlySpan<int> groupA = HeuristicPrimeSieves.GroupADivisors;
|
||||||
|
ReadOnlySpan<uint> groupB = GetGroupBDivisors(nMod10);
|
||||||
|
|
||||||
|
var combined = new List<uint>(groupA.Length + groupB.Length);
|
||||||
|
|
||||||
|
int groupAIndex = 0;
|
||||||
|
while (groupAIndex < groupA.Length && groupA[groupAIndex] <= 19)
|
||||||
|
{
|
||||||
|
combined.Add((uint)groupA[groupAIndex]);
|
||||||
|
groupAIndex++;
|
||||||
|
}
|
||||||
|
|
||||||
|
int groupBIndex = GetGroupBStartIndex(groupB);
|
||||||
|
int groupAInterleaveCount;
|
||||||
|
int groupBInterleaveCount;
|
||||||
|
switch (pattern)
|
||||||
|
{
|
||||||
|
case CombinedDivisorPattern.OneAOneB:
|
||||||
|
groupAInterleaveCount = 1;
|
||||||
|
groupBInterleaveCount = 1;
|
||||||
|
break;
|
||||||
|
case CombinedDivisorPattern.TwoAOneB:
|
||||||
|
groupAInterleaveCount = 2;
|
||||||
|
groupBInterleaveCount = 1;
|
||||||
|
break;
|
||||||
|
case CombinedDivisorPattern.ThreeAOneB:
|
||||||
|
groupAInterleaveCount = 3;
|
||||||
|
groupBInterleaveCount = 1;
|
||||||
|
break;
|
||||||
|
case CombinedDivisorPattern.ThreeATwoB:
|
||||||
|
groupAInterleaveCount = 3;
|
||||||
|
groupBInterleaveCount = 2;
|
||||||
|
break;
|
||||||
|
default:
|
||||||
|
groupAInterleaveCount = 2;
|
||||||
|
groupBInterleaveCount = 1;
|
||||||
|
break;
|
||||||
|
}
|
||||||
|
|
||||||
|
while (groupAIndex < groupA.Length || groupBIndex < groupB.Length)
|
||||||
|
{
|
||||||
|
int addedA = 0;
|
||||||
|
while (addedA < groupAInterleaveCount && groupAIndex < groupA.Length)
|
||||||
|
{
|
||||||
|
combined.Add((uint)groupA[groupAIndex]);
|
||||||
|
groupAIndex++;
|
||||||
|
addedA++;
|
||||||
|
}
|
||||||
|
|
||||||
|
int addedB = 0;
|
||||||
|
while (addedB < groupBInterleaveCount && groupBIndex < groupB.Length)
|
||||||
|
{
|
||||||
|
combined.Add(groupB[groupBIndex]);
|
||||||
|
groupBIndex++;
|
||||||
|
addedB++;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
return combined.ToArray();
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
private static int GetGroupBStartIndex(ReadOnlySpan<uint> divisors)
|
||||||
|
{
|
||||||
|
int index = 0;
|
||||||
|
while (index < divisors.Length && IsGroupAPrefixValue(divisors[index]))
|
||||||
|
{
|
||||||
|
index++;
|
||||||
|
}
|
||||||
|
|
||||||
|
return index;
|
||||||
|
}
|
||||||
|
|
||||||
|
private static bool IsGroupAPrefixValue(uint value) => value is 3U or 7U or 11U or 13U;
|
||||||
|
|
||||||
|
private static HeuristicDivisorGroup ResolveGroup(uint value) => value switch
|
||||||
|
{
|
||||||
|
3U or 7U or 11U or 13U => HeuristicDivisorGroup.GroupAConstant,
|
||||||
|
_ when value % 10U == 3U => HeuristicDivisorGroup.GroupAWheel,
|
||||||
|
_ => HeuristicDivisorGroup.GroupB,
|
||||||
|
};
|
||||||
|
|
||||||
|
private static ReadOnlySpan<uint> GetGroupBDivisors(byte ending) => ending switch
|
||||||
|
{
|
||||||
|
1 => DivisorGenerator.SmallPrimesLastOneWithoutLastThree,
|
||||||
|
3 => ReadOnlySpan<uint>.Empty,
|
||||||
|
7 => DivisorGenerator.SmallPrimesLastSevenWithoutLastThree,
|
||||||
|
9 => DivisorGenerator.SmallPrimesLastNineWithoutLastThree,
|
||||||
|
_ => ReadOnlySpan<uint>.Empty,
|
||||||
|
};
|
||||||
|
|
||||||
|
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||||
|
internal static HeuristicDivisorEnumerator CreateHeuristicDivisorEnumerator(ulong sqrtLimit, byte nMod10, Span<HeuristicGroupBSequenceState> groupBBuffer)
|
||||||
|
{
|
||||||
|
return new HeuristicDivisorEnumerator(sqrtLimit, nMod10, groupBBuffer);
|
||||||
|
}
|
||||||
|
|
||||||
|
internal struct HeuristicGroupBSequenceState
|
||||||
|
{
|
||||||
|
}
|
||||||
|
|
||||||
|
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||||
|
internal static MersenneHeuristicDivisorEnumerator CreateMersenneDivisorEnumerator(ulong exponent, ulong maxDivisor)
|
||||||
|
{
|
||||||
|
return new MersenneHeuristicDivisorEnumerator(exponent, maxDivisor);
|
||||||
|
}
|
||||||
|
|
||||||
|
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||||
|
internal static HeuristicDivisorPreparation PrepareHeuristicDivisor(in HeuristicDivisorCandidate candidate)
|
||||||
|
{
|
||||||
|
ulong divisor = candidate.Value;
|
||||||
|
MontgomeryDivisorData divisorData = MontgomeryDivisorData.FromModulus(divisor);
|
||||||
|
bool hasCycleHint = TryGetCycleLengthHint(divisor, out ulong cycleLength);
|
||||||
|
return new HeuristicDivisorPreparation(candidate, divisorData, cycleLength, hasCycleHint);
|
||||||
|
}
|
||||||
|
|
||||||
|
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||||
|
internal static ulong ResolveHeuristicCycleLength(
|
||||||
|
ulong exponent,
|
||||||
|
in HeuristicDivisorPreparation preparation,
|
||||||
|
out bool cycleFromHint,
|
||||||
|
out bool cycleComputed,
|
||||||
|
out bool primeOrderFailed)
|
||||||
|
{
|
||||||
|
if (preparation.HasCycleLengthHint && preparation.CycleLengthHint != 0UL)
|
||||||
|
{
|
||||||
|
cycleFromHint = true;
|
||||||
|
cycleComputed = true;
|
||||||
|
primeOrderFailed = false;
|
||||||
|
return preparation.CycleLengthHint;
|
||||||
|
}
|
||||||
|
|
||||||
|
ulong divisor = preparation.Candidate.Value;
|
||||||
|
MontgomeryDivisorData divisorData = preparation.DivisorData;
|
||||||
|
|
||||||
|
bool trySuccess = MersenneDivisorCycles.TryCalculateCycleLengthForExponentCpu(
|
||||||
|
divisor,
|
||||||
|
exponent,
|
||||||
|
divisorData,
|
||||||
|
out ulong computedCycle,
|
||||||
|
out bool primeOrderFailedLocal);
|
||||||
|
|
||||||
|
if (trySuccess && computedCycle != 0UL)
|
||||||
|
{
|
||||||
|
cycleFromHint = false;
|
||||||
|
cycleComputed = true;
|
||||||
|
primeOrderFailed = primeOrderFailedLocal;
|
||||||
|
return computedCycle;
|
||||||
|
}
|
||||||
|
|
||||||
|
primeOrderFailed = primeOrderFailedLocal || !trySuccess || computedCycle == 0UL;
|
||||||
|
ulong resolvedCycle = MersenneDivisorCycles.CalculateCycleLength(divisor, divisorData, skipPrimeOrderHeuristic: primeOrderFailed);
|
||||||
|
cycleFromHint = false;
|
||||||
|
cycleComputed = resolvedCycle != 0UL;
|
||||||
|
return resolvedCycle;
|
||||||
|
}
|
||||||
|
|
||||||
|
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||||
|
internal static bool TryGetCycleLengthHint(ulong divisor, out ulong cycleLength)
|
||||||
|
{
|
||||||
|
ulong[] snapshot = HeuristicSmallCycleSnapshot;
|
||||||
|
if (divisor < (ulong)snapshot.Length)
|
||||||
|
{
|
||||||
|
cycleLength = snapshot[(int)divisor];
|
||||||
|
if (cycleLength != 0UL)
|
||||||
|
{
|
||||||
|
return true;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
cycleLength = 0UL;
|
||||||
|
return false;
|
||||||
|
}
|
||||||
|
|
||||||
|
internal ref struct HeuristicDivisorEnumerator
|
||||||
|
{
|
||||||
|
private readonly ulong sqrtLimit;
|
||||||
|
private readonly ReadOnlySpan<uint> combinedDivisors;
|
||||||
|
private int index;
|
||||||
|
|
||||||
|
public HeuristicDivisorEnumerator(ulong sqrtLimit, byte nMod10, Span<HeuristicGroupBSequenceState> groupBBuffer)
|
||||||
|
{
|
||||||
|
this.sqrtLimit = sqrtLimit;
|
||||||
|
_ = groupBBuffer;
|
||||||
|
combinedDivisors = GetCombinedDivisors(nMod10);
|
||||||
|
index = 0;
|
||||||
|
}
|
||||||
|
|
||||||
|
public bool TryGetNext(out HeuristicDivisorCandidate candidate)
|
||||||
|
{
|
||||||
|
while (index < combinedDivisors.Length)
|
||||||
|
{
|
||||||
|
uint entry = combinedDivisors[index++];
|
||||||
|
if (entry > sqrtLimit)
|
||||||
|
{
|
||||||
|
index = combinedDivisors.Length;
|
||||||
|
break;
|
||||||
|
}
|
||||||
|
|
||||||
|
HeuristicDivisorGroup group = ResolveGroup(entry);
|
||||||
|
ulong value = entry;
|
||||||
|
byte ending = (byte)(value % 10UL);
|
||||||
|
ushort residue = (ushort)(value % Wheel210);
|
||||||
|
candidate = new HeuristicDivisorCandidate(value, group, ending, 0, residue);
|
||||||
|
return true;
|
||||||
|
}
|
||||||
|
|
||||||
|
candidate = default;
|
||||||
|
return false;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
internal ref struct MersenneHeuristicDivisorEnumerator
|
||||||
|
{
|
||||||
|
private readonly GpuUInt128 step;
|
||||||
|
private readonly GpuUInt128 limit;
|
||||||
|
private GpuUInt128 current;
|
||||||
|
private MersenneDivisorResidueStepper residueStepper;
|
||||||
|
private bool active;
|
||||||
|
private ulong processedCount;
|
||||||
|
private ulong lastDivisor;
|
||||||
|
|
||||||
|
public MersenneHeuristicDivisorEnumerator(ulong exponent, ulong maxDivisor)
|
||||||
|
{
|
||||||
|
var stepLocal = new GpuUInt128(exponent);
|
||||||
|
stepLocal.ShiftLeft(1);
|
||||||
|
step = stepLocal;
|
||||||
|
|
||||||
|
var limitLocal = new GpuUInt128(maxDivisor);
|
||||||
|
limit = limitLocal;
|
||||||
|
|
||||||
|
var firstDivisor = stepLocal;
|
||||||
|
firstDivisor.Add(1UL);
|
||||||
|
|
||||||
|
bool hasCandidates = !stepLocal.IsZero && firstDivisor.CompareTo(limitLocal) <= 0;
|
||||||
|
|
||||||
|
current = hasCandidates ? firstDivisor : GpuUInt128.Zero;
|
||||||
|
residueStepper = hasCandidates ? new MersenneDivisorResidueStepper(exponent, stepLocal, firstDivisor) : default;
|
||||||
|
active = hasCandidates;
|
||||||
|
processedCount = 0UL;
|
||||||
|
lastDivisor = 0UL;
|
||||||
|
}
|
||||||
|
|
||||||
|
public bool TryGetNext(out HeuristicDivisorCandidate candidate)
|
||||||
|
{
|
||||||
|
while (active)
|
||||||
|
{
|
||||||
|
ulong value = current.Low;
|
||||||
|
processedCount++;
|
||||||
|
lastDivisor = value;
|
||||||
|
|
||||||
|
bool admissible = residueStepper.IsAdmissible();
|
||||||
|
|
||||||
|
Advance();
|
||||||
|
|
||||||
|
if (admissible)
|
||||||
|
{
|
||||||
|
candidate = CreateCandidate(value);
|
||||||
|
return true;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
candidate = default;
|
||||||
|
return false;
|
||||||
|
}
|
||||||
|
|
||||||
|
private void Advance()
|
||||||
|
{
|
||||||
|
GpuUInt128 next = current;
|
||||||
|
next.Add(step);
|
||||||
|
|
||||||
|
if (next.CompareTo(limit) > 0 || next.CompareTo(current) <= 0)
|
||||||
|
{
|
||||||
|
active = false;
|
||||||
|
return;
|
||||||
|
}
|
||||||
|
|
||||||
|
current = next;
|
||||||
|
residueStepper.Advance();
|
||||||
|
}
|
||||||
|
|
||||||
|
public ulong ProcessedCount => processedCount;
|
||||||
|
|
||||||
|
public ulong LastDivisor => lastDivisor;
|
||||||
|
|
||||||
|
public bool Exhausted => !active;
|
||||||
|
|
||||||
|
private static HeuristicDivisorCandidate CreateCandidate(ulong value)
|
||||||
|
{
|
||||||
|
byte ending = (byte)(value % 10UL);
|
||||||
|
HeuristicDivisorGroup group;
|
||||||
|
|
||||||
|
if (value == 3UL || value == 7UL || value == 11UL || value == 13UL)
|
||||||
|
{
|
||||||
|
group = HeuristicDivisorGroup.GroupAConstant;
|
||||||
|
}
|
||||||
|
else if (ending == 3)
|
||||||
|
{
|
||||||
|
group = HeuristicDivisorGroup.GroupAWheel;
|
||||||
|
}
|
||||||
|
else
|
||||||
|
{
|
||||||
|
group = HeuristicDivisorGroup.GroupB;
|
||||||
|
}
|
||||||
|
|
||||||
|
ushort residue = (ushort)(value % Wheel210);
|
||||||
|
return new HeuristicDivisorCandidate(value, group, ending, 0, residue);
|
||||||
|
}
|
||||||
|
}
|
||||||
|
}
|
||||||
@@ -122,10 +122,10 @@ public sealed class HeuristicPrimeTester
|
|||||||
return EvaluateWithOpenNumericFallback(n);
|
return EvaluateWithOpenNumericFallback(n);
|
||||||
}
|
}
|
||||||
|
|
||||||
if (!UseHeuristicGroupBTrialDivision)
|
// if (!UseHeuristicGroupBTrialDivision)
|
||||||
{
|
// {
|
||||||
return HeuristicTrialDivisionCpu(n, sqrtLimit, nMod10, includeGroupB: false);
|
// return HeuristicTrialDivisionCpu(n, sqrtLimit, nMod10, includeGroupB: false);
|
||||||
}
|
// }
|
||||||
|
|
||||||
bool compositeDetected = HeuristicTrialDivisionGpuDetectsDivisor(n, sqrtLimit, nMod10);
|
bool compositeDetected = HeuristicTrialDivisionGpuDetectsDivisor(n, sqrtLimit, nMod10);
|
||||||
return !compositeDetected;
|
return !compositeDetected;
|
||||||
@@ -260,8 +260,8 @@ public sealed class HeuristicPrimeTester
|
|||||||
|
|
||||||
private bool HeuristicTrialDivisionGpuDetectsDivisor(ulong n, ulong sqrtLimit, byte nMod10)
|
private bool HeuristicTrialDivisionGpuDetectsDivisor(ulong n, ulong sqrtLimit, byte nMod10)
|
||||||
{
|
{
|
||||||
int batchCapacity = Math.Max(1, HeuristicGpuDivisorBatchSize);
|
int batchCapacity = HeuristicGpuDivisorBatchSize;
|
||||||
int interleaveBatchSize = Math.Max(1, HeuristicDivisorInterleaveBatchSize);
|
int interleaveBatchSize = HeuristicDivisorInterleaveBatchSize;
|
||||||
var divisorPool = ThreadStaticPools.UlongPool;
|
var divisorPool = ThreadStaticPools.UlongPool;
|
||||||
var hitPool = ThreadStaticPools.BytePool;
|
var hitPool = ThreadStaticPools.BytePool;
|
||||||
|
|
||||||
@@ -277,7 +277,7 @@ public sealed class HeuristicPrimeTester
|
|||||||
bool compositeDetected = false;
|
bool compositeDetected = false;
|
||||||
int count = 0;
|
int count = 0;
|
||||||
|
|
||||||
lock (gpu.ExecutionLock)
|
// lock (gpu.ExecutionLock)
|
||||||
{
|
{
|
||||||
var input = scratch.Input;
|
var input = scratch.Input;
|
||||||
var output = scratch.Output;
|
var output = scratch.Output;
|
||||||
@@ -470,7 +470,7 @@ Cleanup:
|
|||||||
|
|
||||||
for (int i = 0; i < endingOrder.Length; i++)
|
for (int i = 0; i < endingOrder.Length; i++)
|
||||||
{
|
{
|
||||||
ReadOnlySpan<int> divisors = GetGroupBDivisors(endingOrder[i]);
|
ReadOnlySpan<uint> divisors = GetGroupBDivisors(endingOrder[i]);
|
||||||
int index = indices[i];
|
int index = indices[i];
|
||||||
if ((uint)index >= (uint)divisors.Length)
|
if ((uint)index >= (uint)divisors.Length)
|
||||||
{
|
{
|
||||||
@@ -507,7 +507,7 @@ Cleanup:
|
|||||||
{
|
{
|
||||||
for (int i = 0; i < endingOrder.Length; i++)
|
for (int i = 0; i < endingOrder.Length; i++)
|
||||||
{
|
{
|
||||||
ReadOnlySpan<int> divisors = GetGroupBDivisors(endingOrder[i]);
|
ReadOnlySpan<uint> divisors = GetGroupBDivisors(endingOrder[i]);
|
||||||
int index = indices[i];
|
int index = indices[i];
|
||||||
if ((uint)index >= (uint)divisors.Length)
|
if ((uint)index >= (uint)divisors.Length)
|
||||||
{
|
{
|
||||||
@@ -527,13 +527,13 @@ Cleanup:
|
|||||||
return false;
|
return false;
|
||||||
}
|
}
|
||||||
|
|
||||||
private static ReadOnlySpan<int> GetGroupBDivisors(byte ending) => ending switch
|
private static ReadOnlySpan<uint> GetGroupBDivisors(byte ending) => ending switch
|
||||||
{
|
{
|
||||||
1 => HeuristicPrimeSieves.GroupBDivisorsEnding1,
|
1 => DivisorGenerator.SmallPrimesLastOneWithoutLastThree,
|
||||||
3 => HeuristicPrimeSieves.GroupBDivisorsEnding3,
|
3 => ReadOnlySpan<uint>.Empty,
|
||||||
7 => HeuristicPrimeSieves.GroupBDivisorsEnding7,
|
7 => DivisorGenerator.SmallPrimesLastSevenWithoutLastThree,
|
||||||
9 => HeuristicPrimeSieves.GroupBDivisorsEnding9,
|
9 => DivisorGenerator.SmallPrimesLastNineWithoutLastThree,
|
||||||
_ => ReadOnlySpan<int>.Empty,
|
_ => ReadOnlySpan<uint>.Empty,
|
||||||
};
|
};
|
||||||
|
|
||||||
private static ReadOnlySpan<ushort> GetWheel210ResiduesForEnding(byte ending) => ending switch
|
private static ReadOnlySpan<ushort> GetWheel210ResiduesForEnding(byte ending) => ending switch
|
||||||
|
|||||||
@@ -1,11 +1,8 @@
|
|||||||
using System;
|
|
||||||
using System.Buffers;
|
using System.Buffers;
|
||||||
using System.Collections.Concurrent;
|
using System.Collections.Concurrent;
|
||||||
using System.Numerics;
|
using System.Numerics;
|
||||||
using System.Runtime.CompilerServices;
|
using System.Runtime.CompilerServices;
|
||||||
using System.Threading;
|
|
||||||
using ILGPU;
|
using ILGPU;
|
||||||
using ILGPU.Algorithms;
|
|
||||||
using ILGPU.Runtime;
|
using ILGPU.Runtime;
|
||||||
using PerfectNumbers.Core.Gpu;
|
using PerfectNumbers.Core.Gpu;
|
||||||
|
|
||||||
@@ -22,7 +19,7 @@ public sealed class PrimeTester
|
|||||||
|
|
||||||
public static PrimeTester Exclusive => _tester ??= new();
|
public static PrimeTester Exclusive => _tester ??= new();
|
||||||
|
|
||||||
public bool IsPrime(ulong n, CancellationToken ct)
|
public static bool IsPrime(ulong n, CancellationToken ct)
|
||||||
{
|
{
|
||||||
if (n <= 1UL)
|
if (n <= 1UL)
|
||||||
{
|
{
|
||||||
@@ -121,7 +118,7 @@ public sealed class PrimeTester
|
|||||||
public static bool IsPrimeCpu(ulong n, CancellationToken ct)
|
public static bool IsPrimeCpu(ulong n, CancellationToken ct)
|
||||||
{
|
{
|
||||||
// Preserve the legacy entry point for callers that bypass the thread-local caches.
|
// Preserve the legacy entry point for callers that bypass the thread-local caches.
|
||||||
return Exclusive.IsPrime(n, ct);
|
return IsPrime(n, ct);
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
@@ -137,8 +134,6 @@ public sealed class PrimeTester
|
|||||||
int totalLength = values.Length;
|
int totalLength = values.Length;
|
||||||
int batchSize = Math.Max(1, GpuBatchSize);
|
int batchSize = Math.Max(1, GpuBatchSize);
|
||||||
|
|
||||||
lock (gpu.ExecutionLock)
|
|
||||||
{
|
|
||||||
var scratch = state.RentScratch(batchSize, accelerator);
|
var scratch = state.RentScratch(batchSize, accelerator);
|
||||||
var input = scratch.Input;
|
var input = scratch.Input;
|
||||||
var output = scratch.Output;
|
var output = scratch.Output;
|
||||||
@@ -176,7 +171,6 @@ public sealed class PrimeTester
|
|||||||
|
|
||||||
pool.Return(temp, clearArray: false);
|
pool.Return(temp, clearArray: false);
|
||||||
state.ReturnScratch(scratch);
|
state.ReturnScratch(scratch);
|
||||||
}
|
|
||||||
|
|
||||||
gpu.Dispose();
|
gpu.Dispose();
|
||||||
limiter.Dispose();
|
limiter.Dispose();
|
||||||
@@ -190,8 +184,6 @@ public sealed class PrimeTester
|
|||||||
|
|
||||||
public Accelerator Accelerator { get; }
|
public Accelerator Accelerator { get; }
|
||||||
|
|
||||||
public object ExecutionLock { get; } = new();
|
|
||||||
|
|
||||||
private KernelState? _kernelState;
|
private KernelState? _kernelState;
|
||||||
|
|
||||||
public PooledContext()
|
public PooledContext()
|
||||||
@@ -226,7 +218,6 @@ public sealed class PrimeTester
|
|||||||
}
|
}
|
||||||
|
|
||||||
private static readonly ConcurrentQueue<PooledContext> Pool = new();
|
private static readonly ConcurrentQueue<PooledContext> Pool = new();
|
||||||
private static readonly object CreationLock = new();
|
|
||||||
|
|
||||||
internal static PrimeTesterGpuContextLease Rent()
|
internal static PrimeTesterGpuContextLease Rent()
|
||||||
{
|
{
|
||||||
@@ -235,11 +226,8 @@ public sealed class PrimeTester
|
|||||||
return new PrimeTesterGpuContextLease(ctx);
|
return new PrimeTesterGpuContextLease(ctx);
|
||||||
}
|
}
|
||||||
|
|
||||||
lock (CreationLock)
|
|
||||||
{
|
|
||||||
return new PrimeTesterGpuContextLease(new PooledContext());
|
return new PrimeTesterGpuContextLease(new PooledContext());
|
||||||
}
|
}
|
||||||
}
|
|
||||||
|
|
||||||
internal static void DisposeAll()
|
internal static void DisposeAll()
|
||||||
{
|
{
|
||||||
@@ -255,22 +243,16 @@ public sealed class PrimeTester
|
|||||||
Pool.Enqueue(ctx);
|
Pool.Enqueue(ctx);
|
||||||
}
|
}
|
||||||
|
|
||||||
internal struct PrimeTesterGpuContextLease
|
internal readonly struct PrimeTesterGpuContextLease(PrimeTesterGpuContextPool.PooledContext ctx)
|
||||||
{
|
{
|
||||||
private PooledContext? _ctx;
|
private readonly PooledContext? _ctx = ctx;
|
||||||
|
|
||||||
internal PrimeTesterGpuContextLease(PooledContext ctx)
|
public readonly Accelerator Accelerator => _ctx!.Accelerator;
|
||||||
{
|
|
||||||
_ctx = ctx;
|
|
||||||
}
|
|
||||||
|
|
||||||
public Accelerator Accelerator => _ctx!.Accelerator;
|
|
||||||
|
|
||||||
public object ExecutionLock => _ctx!.ExecutionLock;
|
public readonly KernelState State => _ctx!.KernelState;
|
||||||
|
|
||||||
public KernelState State => _ctx!.KernelState;
|
public readonly void Dispose()
|
||||||
|
|
||||||
public void Dispose()
|
|
||||||
{
|
{
|
||||||
Return(_ctx ?? throw new InvalidOperationException("GPU context lease is not initialized."));
|
Return(_ctx ?? throw new InvalidOperationException("GPU context lease is not initialized."));
|
||||||
}
|
}
|
||||||
@@ -280,20 +262,20 @@ public sealed class PrimeTester
|
|||||||
// Per-accelerator GPU state for prime sieve (kernel + uploaded primes).
|
// Per-accelerator GPU state for prime sieve (kernel + uploaded primes).
|
||||||
internal sealed class KernelState
|
internal sealed class KernelState
|
||||||
{
|
{
|
||||||
public Action<Index1D, ArrayView<ulong>, ArrayView<uint>, ArrayView<uint>, ArrayView<uint>, ArrayView<uint>, ArrayView<uint>, ArrayView<ulong>, ArrayView<ulong>, ArrayView<ulong>, ArrayView<ulong>, ArrayView<ulong>, ArrayView<byte>> Kernel { get; }
|
public Action<Index1D, ArrayView<ulong>, ArrayView<uint>, ArrayView<uint>, ArrayView<uint>, ArrayView<uint>, ArrayView<uint>, ArrayView<ulong>, ArrayView<ulong>, ArrayView<ulong>, ArrayView<ulong>, ArrayView<ulong>, ArrayView<byte>> Kernel;
|
||||||
public Action<Index1D, ArrayView<ulong>, ulong, ArrayView<byte>> HeuristicTrialDivisionKernel { get; }
|
public Action<Index1D, ArrayView<ulong>, ulong, ArrayView<byte>> HeuristicTrialDivisionKernel;
|
||||||
public MemoryBuffer1D<uint, Stride1D.Dense> DevicePrimesDefault { get; }
|
public MemoryBuffer1D<uint, Stride1D.Dense> DevicePrimesDefault;
|
||||||
public MemoryBuffer1D<uint, Stride1D.Dense> DevicePrimesLastOne { get; }
|
public MemoryBuffer1D<uint, Stride1D.Dense> DevicePrimesLastOne;
|
||||||
public MemoryBuffer1D<uint, Stride1D.Dense> DevicePrimesLastSeven { get; }
|
public MemoryBuffer1D<uint, Stride1D.Dense> DevicePrimesLastSeven;
|
||||||
public MemoryBuffer1D<uint, Stride1D.Dense> DevicePrimesLastThree { get; }
|
public MemoryBuffer1D<uint, Stride1D.Dense> DevicePrimesLastThree;
|
||||||
public MemoryBuffer1D<uint, Stride1D.Dense> DevicePrimesLastNine { get; }
|
public MemoryBuffer1D<uint, Stride1D.Dense> DevicePrimesLastNine;
|
||||||
public MemoryBuffer1D<ulong, Stride1D.Dense> DevicePrimesPow2Default { get; }
|
public MemoryBuffer1D<ulong, Stride1D.Dense> DevicePrimesPow2Default;
|
||||||
public MemoryBuffer1D<ulong, Stride1D.Dense> DevicePrimesPow2LastOne { get; }
|
public MemoryBuffer1D<ulong, Stride1D.Dense> DevicePrimesPow2LastOne;
|
||||||
public MemoryBuffer1D<ulong, Stride1D.Dense> DevicePrimesPow2LastSeven { get; }
|
public MemoryBuffer1D<ulong, Stride1D.Dense> DevicePrimesPow2LastSeven;
|
||||||
public MemoryBuffer1D<ulong, Stride1D.Dense> DevicePrimesPow2LastThree { get; }
|
public MemoryBuffer1D<ulong, Stride1D.Dense> DevicePrimesPow2LastThree;
|
||||||
public MemoryBuffer1D<ulong, Stride1D.Dense> DevicePrimesPow2LastNine { get; }
|
public MemoryBuffer1D<ulong, Stride1D.Dense> DevicePrimesPow2LastNine;
|
||||||
public bool IsDisposed { get; private set; }
|
public bool IsDisposed;
|
||||||
private readonly System.Collections.Concurrent.ConcurrentBag<ScratchBuffers> _scratchPool = [];
|
private readonly ConcurrentBag<ScratchBuffers> _scratchPool = [];
|
||||||
// TODO: Replace this ConcurrentBag with the lock-free ring buffer variant validated in
|
// TODO: Replace this ConcurrentBag with the lock-free ring buffer variant validated in
|
||||||
// GpuModularArithmeticBenchmarks so renting scratch buffers stops contending on the bag's internal locks when
|
// GpuModularArithmeticBenchmarks so renting scratch buffers stops contending on the bag's internal locks when
|
||||||
// thousands of GPU batches execute per second.
|
// thousands of GPU batches execute per second.
|
||||||
@@ -369,7 +351,6 @@ public sealed class PrimeTester
|
|||||||
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
[MethodImpl(MethodImplOptions.AggressiveInlining)]
|
||||||
public ScratchBuffers RentScratch(int minCapacity, Accelerator accel)
|
public ScratchBuffers RentScratch(int minCapacity, Accelerator accel)
|
||||||
{
|
{
|
||||||
minCapacity = Math.Max(1, minCapacity);
|
|
||||||
while (_scratchPool.TryTake(out var sb))
|
while (_scratchPool.TryTake(out var sb))
|
||||||
{
|
{
|
||||||
if (sb.Capacity >= minCapacity)
|
if (sb.Capacity >= minCapacity)
|
||||||
@@ -428,7 +409,7 @@ public sealed class PrimeTester
|
|||||||
{
|
{
|
||||||
var lazy = States.GetOrAdd(
|
var lazy = States.GetOrAdd(
|
||||||
accelerator,
|
accelerator,
|
||||||
acc => new Lazy<KernelState>(() => new KernelState(acc), System.Threading.LazyThreadSafetyMode.ExecutionAndPublication));
|
acc => new Lazy<KernelState>(() => new KernelState(acc), LazyThreadSafetyMode.ExecutionAndPublication));
|
||||||
return lazy.Value;
|
return lazy.Value;
|
||||||
}
|
}
|
||||||
|
|
||||||
@@ -525,4 +506,4 @@ public sealed class PrimeTester
|
|||||||
|
|
||||||
return u << shift;
|
return u << shift;
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
Reference in new issue
Block a user