GameCult/CultLib
Watch
9
0
Fork
You've already forked CultLib
 0
Code Issues Pull requests 1 Projects Releases Packages Wiki Activity Actions

CultMath BoundedLeastSquares: the free-set stationarity test is a traced decision; the stall boundary is pinned to the ulp

Browse source 
FreeSetViolated is the one owner of the free-set test and records each violating column with its gradient magnitude and the threshold it was judged against (SolveEvent.FreeViolated). Tests pin the column and the threshold by exact equality, one ulp each side of the tolerance at four scales, and the stall comparison f >= mark - allowance at mark == allowance and one ulp each way.

Builds are unverified until the batch's single verify job; this commit has not been compiled.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
This commit is contained in:
Meta Vens 2026-10-11 02:37:05 +02:00
commit f8253837ed
2 changed files with 158 additions and 7 deletions
Download patch file Download diff file Expand all files Collapse all files

25
packages/cultmath/src/CultMath/BoundedLeastSquares.cs
View file

@ -46,6 +46,9 @@ internal enum SolveEvent
/// <summary>A bound released: column, the stall mark after the release.</summary>
Release,
/// <summary>A free column whose gradient is above the tolerance, so the free set is not stationary: column, the gradient's magnitude, the tolerance it was judged against. One per such column, before the pass's other decisions.</summary>
FreeViolated,
/// <summary>KKT holds and the solve is over: the tolerance.</summary>
Kkt,
@ -348,11 +351,7 @@ public static class BoundedLeastSquares
Gradient(ata, atb, xd, g, n);
// Stationarity of the free set is part of KKT: a Newton step on a near-singular free set can leave a residual gradient.
var freeViolated = false;
for (var i = 0; i < freeCount; i++)
if (Math.Abs(g[freeIndex[i]]) > kktTol)
freeViolated = true;
var freeViolated = FreeSetViolated(g, freeIndex, freeCount, kktTol, trace);
var release = SelectRelease(ata, atb, xd, g, state, lo, hi, kktTol, n, trace);
if (JudgePass(release, freeViolated, kktTol, ref fMark, xd, g, atb, n, trace)) { status = BoundedLeastSquaresStatus.Converged; break; }
if (iterations >= maxIterations) { trace?.Record(SolveEvent.Cap, -1, 0.0, 0.0, 0.0); status = BoundedLeastSquaresStatus.IterationLimit; break; }
@ -366,6 +365,22 @@ public static class BoundedLeastSquares
return status;
}
// Whether any free column's gradient is above the tolerance; every such column is recorded with the magnitude and the threshold it was judged against.
internal static bool FreeSetViolated(ReadOnlySpan<double> g, ReadOnlySpan<int> freeIndex, int freeCount, double kktTol, ISolveTrace? trace)
{
var violated = false;
for (var i = 0; i < freeCount; i++)
{
var magnitude = Math.Abs(g[freeIndex[i]]);
if (magnitude > kktTol)
{
trace?.Record(SolveEvent.FreeViolated, freeIndex[i], magnitude, kktTol, 0.0);
violated = true;
}
}
return violated;
}
// The bound to release: the one with the most violated gradient sign, the first of equal violations, and only if
// the violation exceeds both the tolerance (the starting `worst`) and the rounding bound of its own gradient.
// -1 when none qualifies.

140
packages/cultmath/tests/CultMath.Tests/BoundedLeastSquaresTraceTests.cs
View file

@ -11,7 +11,7 @@ namespace CultMath.Tests;
// equality, never a status or a tolerance band. Two sources of truth:
// - BlsBaseTraces: the same decisions as the base commit dbcb2a8e took them on the same problems, so Solve's defaults are
// unchanged decision by decision, not only in the end;
// - hand-built inputs to the two judgement functions, for states a solve may never reach (a release after a stationarity
// - hand-built inputs to the three judgement functions, for states a solve may never reach (a release after a stationarity
// miss with a finite stall mark) where only the decision itself can be asked.
public sealed class BoundedLeastSquaresTraceTests
{
@ -23,11 +23,58 @@ public sealed class BoundedLeastSquaresTraceTests
foreach (var c in BlsTrace.Cases())
{
Assert.True(BlsBaseTraces.Expected.TryGetValue(c.Name, out var expected), $"{c.Name}: no base trace");
Assert.Equal(Normal(expected!), BlsTrace.Trace(c));
Assert.Equal(Normal(expected!), WithoutFreeViolated(BlsTrace.Trace(c)));
}
Assert.Equal(BlsTrace.Cases().Count(), BlsBaseTraces.Expected.Count);
}
// The base commit took the free-set test inline and recorded nothing for it; its record is the one thing the base traces lack.
private static string WithoutFreeViolated(string trace)
=> string.Join("
", trace.Split('
').Where(l => !l.StartsWith("FreeViolated ", StringComparison.Ordinal)));
// In a real solve the free-set test is a decision of its own, taken before the pass's release and outcome: a miss or a
// stall always follows the violating columns that caused it, KKT never does, and each recorded magnitude is above the
// threshold recorded beside it.
[Fact]
public void EveryMissAndStallIsPrecededByTheFreeColumnsThatBrokeStationarityAndKktByNone()
{
var misses = 0;
foreach (var c in BlsTrace.Cases())
{
var violated = 0;
foreach (var line in BlsTrace.Trace(c).Split('
'))
{
var parts = line.Split(' ');
switch (parts[0])
{
case "FreeViolated":
violated++;
Assert.InRange(int.Parse(parts[1], System.Globalization.CultureInfo.InvariantCulture), 0, c.N - 1);
Assert.True(double.Parse(parts[2], System.Globalization.CultureInfo.InvariantCulture) > double.Parse(parts[3], System.Globalization.CultureInfo.InvariantCulture), c.Name);
break;
case "Miss":
case "Stall":
misses++;
Assert.True(violated > 0, $"{c.Name}: {parts[0]} without a violating free column");
violated = 0;
break;
case "Kkt":
Assert.True(violated == 0, $"{c.Name}: KKT after a violating free column");
break;
case "Release":
case "Step":
case "Cap":
violated = 0;
break;
}
}
}
Assert.True(misses > 0);
}
// Each kind of decision, and the ones the stop rule's variants differ on, occurs in the recorded problems: a case table
// that stopped exercising a decision would pin nothing about it.
[Fact]
@ -130,6 +177,95 @@ public sealed class BoundedLeastSquaresTraceTests
Assert.Equal($"Stall -1 2 2 {BlsTrace.F(allowance)}", trace);
}
// The stall comparison is f >= mark - allowance, pinned at its boundary to one ulp each side like the release gate. The point has
// cost exactly 0 (so the mark is the allowance itself and the subtraction is exact) and scale 4, allowance 4e-15: a stall at mark ==
// allowance, progress one ulp above it, whatever factor, offset or clamp the comparison puts on the allowance.
[Fact]
public void TheStallComparisonIsExactAtTheAllowanceToTheUlp()
{
var x = new[] { 1.0, 1.0 };
var g = new[] { 3.0, 1.0 };
var atb = new[] { 1.0, 3.0 };
var allowance = 1e-15 * 4.0;
var mark = allowance; // 0 >= allowance - allowance
var trace = Judge(-1, true, ref mark, x, g, atb, out var over);
Assert.True(over);
Assert.Equal($"Stall -1 0 {BlsTrace.F(allowance)} {BlsTrace.F(allowance)}", trace);
mark = Math.BitIncrement(allowance); // 0 < the smallest positive difference
trace = Judge(-1, true, ref mark, x, g, atb, out over);
Assert.False(over);
Assert.Equal(0.0, mark);
Assert.Equal($"Miss -1 0 {BlsTrace.F(Math.BitIncrement(allowance))} {BlsTrace.F(allowance)}", trace);
mark = Math.BitDecrement(allowance); // 0 >= a negative difference
trace = Judge(-1, true, ref mark, x, g, atb, out over);
Assert.True(over);
Assert.Equal($"Stall -1 0 {BlsTrace.F(Math.BitDecrement(allowance))} {BlsTrace.F(allowance)}", trace);
}
private static string FreeViolated(double[] g, int[] freeIndex, int freeCount, double kktTol, out bool violated)
{
var rec = new BlsTrace.Recorder();
violated = BoundedLeastSquares.FreeSetViolated(g, freeIndex, freeCount, kktTol, rec);
return string.Join("
", rec.Lines);
}
// Every free column is judged, and each one that breaks stationarity is named with its gradient magnitude and the threshold
// it was judged against: not the first or last free column only, a column that is not free never, and the columns past
// freeCount never.
[Fact]
public void EveryFreeColumnIsJudgedAndTheViolatingOnesAreNamedWithTheirThreshold()
{
var g = new[] { 0.1, -0.3, 9.0, 0.2, -0.5, 7.0 };
// Column 1 only, in the middle of the free list.
var trace = FreeViolated(g, new[] { 0, 1, 3, 4, 5 }, 3, 0.25, out var violated);
Assert.True(violated);
Assert.Equal("FreeViolated 1 0.3 0.25 0", trace);
// The free list is not the column order: the recorded column is the column, in list order, and column 2 (a bound) is not listed.
trace = FreeViolated(g, new[] { 4, 3, 1, 0, 5 }, 4, 0.25, out violated);
Assert.True(violated);
Assert.Equal("FreeViolated 4 0.5 0.25 0
FreeViolated 1 0.3 0.25 0", trace);
// Only the last listed column violates; only the first does.
trace = FreeViolated(g, new[] { 0, 3, 4 }, 3, 0.25, out violated);
Assert.Equal("FreeViolated 4 0.5 0.25 0", trace);
trace = FreeViolated(g, new[] { 4, 0, 3 }, 3, 0.25, out violated);
Assert.Equal("FreeViolated 4 0.5 0.25 0", trace);
trace = FreeViolated(g, new[] { 0, 3 }, 2, 0.25, out violated);
Assert.False(violated);
Assert.Equal("", trace);
trace = FreeViolated(g, new int[0], 0, 0.0, out violated);
Assert.False(violated);
}
// The threshold is the tolerance itself, at any scale: equal is stationary, one ulp above it (either sign) is not, one ulp
// below it is, and the threshold recorded is the one compared against.
[Theory]
[InlineData(1e-6)]
[InlineData(0.25)]
[InlineData(3.0)]
[InlineData(1e6)]
public void TheFreeStationarityThresholdIsTheToleranceToTheUlp(double kktTol)
{
var freeIndex = new[] { 0, 1, 2 };
var above = Math.BitIncrement(kktTol);
var g = new[] { kktTol, -above, Math.BitDecrement(kktTol) };
var trace = FreeViolated(g, freeIndex, 3, kktTol, out var violated);
Assert.True(violated);
Assert.Equal($"FreeViolated 1 {BlsTrace.F(above)} {BlsTrace.F(kktTol)} 0", trace);
trace = FreeViolated(new[] { kktTol, -kktTol, Math.BitDecrement(kktTol) }, freeIndex, 3, kktTol, out violated);
Assert.False(violated);
Assert.Equal("", trace);
}
// The cost and its scale are 0.5 sum x_j (g_j - atb_j) and 0.5 sum |x_j| (|g_j| + |atb_j|) over every column.
[Fact]
public void TheCostAndItsScaleCountEveryColumnAndTheirOwnSigns()