Skip to content

Guard Hegedüs warm starts with residual predictiveness - #1123

Draft
ChrisRackauckas-Claude wants to merge 1 commit into
SciML:mainfrom
ChrisRackauckas-Claude:agent/guard-hegedus-predictiveness
Draft

Guard Hegedüs warm starts with residual predictiveness#1123
ChrisRackauckas-Claude wants to merge 1 commit into
SciML:mainfrom
ChrisRackauckas-Claude:agent/guard-hegedus-predictiveness

Conversation

@ChrisRackauckas-Claude

Copy link
Copy Markdown
Member

This PR should be ignored until reviewed by @ChrisRackauckas.

Companion to SciML/OrdinaryDiffEq.jl#4040 and SciML/OrdinaryDiffEq.jl#4034.

What changed

WarmStart.Hegedus now uses its projected previous solution only when that
projection reduces the raw initial residual by at least a factor of two. If the
previous direction is not predictive, the solve starts cold.

The public documentation now describes that fallback and no longer says that
OrdinaryDiffEq resolves WarmStart.Auto to Hegedüs. It also says explicitly that
warm starting is workload-dependent and should be benchmarked against a cold
start. The threshold remains an implementation detail; no public API was added.
This is a bug-fix release bump to 5.3.1.

Why

The Hegedüs scalar

$$\xi = \frac{\langle Au,b\rangle}{\lVert Au\rVert^2}$$

minimizes the residual along the previous solution direction. When Au is
nearly orthogonal to b, however, that direction predicts essentially none of
the new solution. In the OrdinaryDiffEq index-1 DAE reproducer, the previous
Newton increment had converged to round-off, |ξ| spanned roughly 63 orders of
magnitude, and GMRES sometimes accepted the rescaled noise at iteration zero.
That made each Newton update depend on the entire prior solve history and
amplified a one-ulp input perturbation into a 3.73e-3 solution change.

For the least-squares scalar, the projected residual ratio is

$$\frac{\lVert b-\xi Au\rVert}{\lVert b\rVert} = \sqrt{1-\cos^2\angle(Au,b)}.$$

The implementation therefore checks the equivalent cosine threshold using the
already-computed Au and inner product. It adds no operator application or
workspace allocation. For an unpreconditioned solve it also reuses the RHS norm
that is already needed for the stopping criterion.

Robustness and performance evidence

On the small index-1 DAE that exposed the bug:

configuration Hairer4 one-ulp drift Hairer4 nf
cold 6.64e-8 842
unguarded Hegedüs 3.73e-3 670
guarded Hegedüs 1.33e-15 617

The guard both removes the round-off amplification and improves the iteration
count on that correlated tiny-DAE workload.

It is deliberately not used to re-enable OrdinaryDiffEq's Auto default.
On a matrix-free 1D Brusselator (n=800, default diagonal preconditioner),
guarded Hegedüs still regressed against cold starts:

algorithm nf change wall-time change
KenCarp4 +12.7% +4.4%
TRBDF2 +30.0% +39.3%
FBDF +7.5% +11.3%

The guard makes explicit Hegedüs safer; it does not make warm starting a
universal speedup.

Local verification

  • Focused test/Core/warm_start.jl: 41/41 assertions passed on Julia 1.12.6.
  • GROUP=Core: passed.
  • GROUP=QA: passed, including Allocation QA (48/48) and SupernodalLU Allocation QA (8/8).
  • Repository-wide Runic 1.7 --check: passed.
  • OrdinaryDiffEq end-to-end guarded-Hegedüs reproduction: Hairer4 drift
    1.33e-15 (nf=617), Hairer42 drift 1.33e-15 (nf=625), and Cash4
    drift 2.77e-9 (nf=740).

Co-Authored-By: Chris Rackauckas <accounts@chrisrackauckas.com>

Copy link
Copy Markdown
Member Author

Cross-repository CI note: any immediate OrdinaryDiffEq downstream resolver failure is the clean OrdinaryDiffEq master NonlinearSolveBase compatibility conflict already fixed independently in SciML/OrdinaryDiffEq.jl#4051. It is unrelated to this three-file Hegedüs guard; local LinearSolve Core and QA (including allocation checks) are green.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants