[ChrisRackauckas-Claude]

🚧 UNREVIEWED — awaiting review by @ChrisRackauckas. Filed by an AI agent running as @ChrisRackauckas; Chris has not seen this.
Harness: Claude Code 2.1.301 · Model: claude-opus-5-5
Conversation: https://claude.ai/code/session_01LPHREnnonfLg1VcE1EJovv

Please ignore until reviewed by @ChrisRackauckas.

Stacked on PR 1365. This branch includes that PR's Jacobian residual-buffer changes, first-iteration reuse, and the public API and compatibility corrections described below. The independent review below assesses the full stacked diff.

What changed and why

construct_jacobian_cache evaluates the Jacobian at u0 during init (it needs the matrix type to build the linear-solver cache). The first step! then evaluated it again at the same point, because make_new_jacobian starts out true. This PR reuses the init-time Jacobian on the first step when it is valid:

• NonlinearSolveBase. JacobianCache records whether construction actually evaluated J. That is true only on the AD path with no analytic jac, no jac_prototype and no operator Jacobian, and never for scalar problems. Solver packages read the flag through the public NonlinearSolveBase.jacobian_computed_at_construction(cache), which defaults to false for every other cache type, has a docstring and is rendered in docs/src/devdocs/jacobian.md. Both solvers import the helper directly from NonlinearSolveBase; their QA ignore lists have no additions for it. Base is versioned at 2.56.0 and both consumers require that minor.
• NonlinearSolveFirstOrder (NewtonRaphson, TrustRegion, LevenbergMarquardt, PseudoTransient, GaussNewton, …). The cache flag initial_jacobian_valid is set at init and cleared by reinit! and after the first step. Step 1 reuses J only while the flag is set and u == u_cache. The second check matters because a polyalgorithm with alias_u0 = true shares u between its subcaches, so an earlier subalgorithm may have moved it. The descent still receives new_jacobian = true, so refactorization, damping and the Jacobian-reuse policy behave as before.
• NonlinearSolveQuasiNewton. TrueJacobianInitialization (Broyden/Klement with init_jacobian = Val(:true_jacobian) or Val(:true_jacobian_diagonal)) is marked initialized when construction evaluated J. The same aliasing guard applies.
• DampedNewtonDescent (:simple mode, used by PseudoTransient). Its init damped the aliased Jacobian-cache storage in place. A reused first-step J would then have been damped twice: before this fix, PT took 90 steps instead of the base's 89. init now damps a copy. This adds one n×n copy at PseudoTransient init.

Not reused, by design: lazy/Krylov JacobianOperator (J is never formed), an analytic jac and a jac_prototype (J is only allocated at init), and scalar problems.

Bitwise identical to #1365

Comparison on the same 98 configurations as #1365:

• Problems: dense 100 IIP and OOP, dense IIP with a user jac, sparse 1000 with a jac_prototype and with TracerSparsityDetector, SVector{4}, scalar.
• Algorithms: NR, NR + FiniteDiff, NR + GMRES (lazy and concrete), TR, TR(Hei), LM, PT, Broyden, Broyden/Klement with true J, Klement with the true diagonal, and PolyAlg(PT, NR) and PolyAlg(PT, Broyden-true-J) with alias_u0 = true.
• Each configuration covers a one-shot solve, every iterate of a step! loop, and solve! after reinit! with a new u0 and p, with a new p, and with the original values.

98 / 98 IDENTICAL: every iterate, residual, retcode and stat other than njacs (nf, nfactors, nsolve, nsteps) is the same. Change in njacs:

njacs change (solve, cache solve, three reinit solves)	configurations
−1, −1, 0, 0, 0	40 (every direct solve where J is computed at init)
0 everywhere	50 (user jac, jac_prototype, Krylov operator, scalar, identity-init Broyden)
−2, −1, 0, 0, 0	6 (polyalgorithms: both subcaches save one; with alias_u0 the second subcache's state has moved on the cache path, so the guard recomputes)
−2, −2, 0, 0, 0	2 (SVector polyalgorithms: immutable u, so no aliasing)

After reinit!, the Jacobian is always re-evaluated (delta 0).

Before → after

Pinned to one core with taskset -c 100, @benchmark solve(prob, alg) minimum. Base is #1365 (7dd45902), and the two runs were done back to back on the same core:

problem	alg	njacs	min µs	change
dense 100 IIP	NewtonRaphson	3 → 2	3820 → 2623	−31%
dense 100 IIP	TrustRegion	3 → 2	5719 → 4533	−21%
dense 100 IIP	LevenbergMarquardt	4 → 3	12491 → 11393	−9%
sparse 1000 IIP (sparsity detector)	NewtonRaphson	3 → 2	1119 → 1100	−2%
sparse 1000 IIP (sparsity detector)	TrustRegion	6 → 5	3624 → 3564	−2%
sparse 1000 IIP (sparsity detector)	LevenbergMarquardt	29 → 28	1.21e6 → 1.40e6	noise (see below)
sparse 1000 IIP (jac_prototype)	NR / TR / LM	2 / 5 / 28 unchanged	758 / 3292 / 1.30e6 → 753 / 3367 / 1.29e6	none expected: J is not computed at init
SVector{4} OOP	NewtonRaphson	6 → 5	1.97 → 1.44
SVector{4} OOP	TrustRegion	6 → 5	3.01 → 3.16	noise
SVector{4} OOP	LevenbergMarquardt	10 → 9	69.9 → 69.9	noise

Dense 100 is the only case where the effect is large and clearly above noise. There, one ForwardDiff Jacobian (100 columns at chunk size 1) dominates a Newton step. For the sparse Bratu problem, a colored Jacobian costs about three function evaluations, so dropping one is within run-to-run variation. The sparse LM row is identical in njacs and nsteps, and its 1.2 s vs 1.4 s gap is load noise: the machine had a load average of about 70 during these runs. An earlier pair of runs against master on the same core gave 3808 → 2597 / 5699 → 4513 / 12517 → 11284 µs for dense NR / TR / LM.

Verification

Regression tests that fail on the base (#1365, run in the same environment):

• lib/NonlinearSolveFirstOrder/test/misc_tests__item11.jl counts dual-number evaluations of f and the njacs delta in the first step!:
  • 4 algorithms × (IIP, OOP, SVector, sparsity-detector) plus chunked ForwardDiff on the dense problems.
  • Not-computed-at-init paths: jac_prototype, user jac, Krylov.
  • A moved u before the first step.
  • reinit! with a new u0, p, or both, for NR/TR/LM/PT: the first step after reinit evaluates the Jacobian, and the result equals a fresh solve (u, retcode, njacs).
• lib/NonlinearSolveQuasiNewton/test/core_tests__item14.jl is the same for Broyden(true J), Klement(true J) and Klement(true diagonal).

On the base:

first step reuses the Jacobian from init: Test Failed at lib/NonlinearSolveFirstOrder/test/misc_tests__item11.jl:51
  Expression: njacs == 0
   Evaluated: 1 == 0
a moved or reinitialized state invalidates the Jacobian from init: Test Failed at lib/NonlinearSolveFirstOrder/test/misc_tests__item11.jl:97
  Expression: reinit_sol.stats.njacs == fresh.stats.njacs
   Evaluated: 7 == 8
Test Summary:                                                       | Pass  Fail  Total     Time
first-step reuse (base)                                             |   61    50    111  1m22.2s

first step reuses the true Jacobian from init: Broyden: Test Failed at lib/NonlinearSolveQuasiNewton/test/core_tests__item14.jl:33
  Expression: ad_calls[] == calls
   Evaluated: 16 == 8
Test Summary:                                            | Pass  Fail  Total   Time
B-QN                                                     |   39    27     66  28.0s

On this branch, in the full Core groups (Julia 1.12):

GROUP=Core lib/NonlinearSolveFirstOrder
  First step reuses the Jacobian from init |  111    111  1m07.5s
  TrustRegion reinit! resets trust_region  |   11     11  1m18.4s
  Adaptive Jacobian reuse                  |  125    125  17.6s
  Native bounded methods                   | 6001   6001  17m38.9s
  SciMLOperator Jacobians                  |   14     14  17.8s
     Testing NonlinearSolveFirstOrder tests passed
GROUP=Core lib/NonlinearSolveQuasiNewton
  True Jacobian initialization is not re-evaluated on the first step |   66     66  18.7s
     Testing NonlinearSolveQuasiNewton tests passed
GROUP=Core lib/NonlinearSolveBase         Testing NonlinearSolveBase tests passed
GROUP=QA   lib/NonlinearSolveBase         QA |   20     20   Testing NonlinearSolveBase tests passed
GROUP=QA   lib/NonlinearSolveFirstOrder   QA |   20     20  53.2s   Testing NonlinearSolveFirstOrder tests passed
GROUP=QA   lib/NonlinearSolveQuasiNewton  QA |   20     20  44.9s   Testing NonlinearSolveQuasiNewton tests passed

The reinit assertions now pass for all four FirstOrder algorithms, including PseudoTransient after the double-damping fix. Runic (Runic --check) and typos are clean on every changed file.

Not verified

• Not run: the other sublibraries (SpectralMethods, SimpleNonlinearSolve, SCCNonlinearSolve, top-level NonlinearSolve), Julia LTS and pre, GPU, and the docs build. SimpleNonlinearSolve shares none of this code and is untouched.
• Not changed:
  • The native bounded algorithms (native_bounded.jl) also evaluate J twice around init. That is a separate path.
  • A user jac still adds 1 to njacs at init without evaluating anything (jacobian.jl).
• Not detected: mutating p in place between init and the first step without reinit!. Mutating u in place is caught by the u == u_cache guard. On the base, such a mutation also leaves fu stale, so it is not a supported path.

For reviewers to push back on

• PseudoTransient init now copies J once.

🤖 Generated with Claude Code (model: claude-opus-5-5) https://claude.ai/code/session_01LPHREnnonfLg1VcE1EJovv

Review findings addressed

Findings from the independent review (codex-cli, gpt-6-astra, verdict CHANGES):

1. P1, raise the NonlinearSolveBase compat floors of FirstOrder and QuasiNewton. Not changed in this PR. In this repo, versions and floors move in separate release PRs: Base is released first, then the dependent package is released with its floor raised (Release NonlinearSolveBase v2.54.0 #1312 then Release NonlinearSolveFirstOrder v2.10.0 #1313; Release NonlinearSolveBase v2.52.0 #1295 then Raise NonlinearSolveFirstOrder's NonlinearSolveBase floor to 2.52 #1297). A floor that includes this function has to name an unreleased Base version (2.56), so raising it here would also mean bumping Base here. The release requirement is spelled out under "Release boundary" below.
2. P2, make the helper public instead of exempting it from QA. Fixed. jacobian_computed_at_construction is now a top-level NonlinearSolveBase function, declared public with @compat(public, ...) next to reused_jacobian. It has a docstring and an entry in docs/src/devdocs/jacobian.md. It is gone from InternalAPI and from its module docstring. Both sublibraries import it the way they import reused_jacobian, and the two QA ignore-list additions are reverted, so both test/qa/qa.jl files match master again. New unit test: lib/NonlinearSolveBase/test/jacobian_computed_at_construction.jl. It checks true for in-place and out-of-place AD with reused_jacobian(cache, u0) equal to the true Jacobian. It checks false for an analytic jac, a jac_prototype, a Krylov JacobianOperator, a scalar problem and a non-cache argument.

This supersedes three statements in the description above. The helper is no longer reached as InternalAPI.jacobian_computed_at_construction. It is no longer in the QA ignore lists. The PR is no longer free of public API changes: it adds one public, documented NonlinearSolveBase name.

Evidence (Julia 1.13.0, 4 threads)

QA with the ignore entries removed but the helper still in InternalAPI (the source before this commit):

NonPublicQualifiedAccessException
  Module `NonlinearSolveFirstOrder` has explicit imports of names from modules in which they are not public (i.e. exported or declared public in Julia 1.11+):
  - `jacobian_computed_at_construction` is not public in `NonlinearSolveBase.InternalAPI` but it was imported from `NonlinearSolveBase.InternalAPI` at `lib/NonlinearSolveFirstOrder/src/solve.jl:352:25`
QA                                             |   19      1     20  42.7s
  - `jacobian_computed_at_construction` is not public in `NonlinearSolveBase.InternalAPI` but it was imported from `NonlinearSolveBase.InternalAPI` at `lib/NonlinearSolveQuasiNewton/src/initialization.jl:129:31`
QA                                             |   19      1     20  34.0s

With this commit:

GROUP=QA   lib/NonlinearSolveFirstOrder   QA |   20     20  40.5s   Testing NonlinearSolveFirstOrder tests passed
GROUP=QA   lib/NonlinearSolveQuasiNewton  QA |   20     20  31.8s   Testing NonlinearSolveQuasiNewton tests passed
GROUP=QA   lib/NonlinearSolveBase         QA |   20     20  36.8s   Testing NonlinearSolveBase tests passed
GROUP=Core lib/NonlinearSolveBase         jacobian_computed_at_construction | 9  9   Testing NonlinearSolveBase tests passed
GROUP=Core lib/NonlinearSolveQuasiNewton  True Jacobian initialization is not re-evaluated on the first step | 66  66   Testing NonlinearSolveQuasiNewton tests passed
GROUP=Core lib/NonlinearSolveFirstOrder  First step reuses the Jacobian from init | 111  111; Native bounded methods | 6001  6001; SciMLOperator Jacobians | 14  14   Testing NonlinearSolveFirstOrder tests passed
Julia 1.10.12, in-repo Base developed as CI does:
GROUP=QA   lib/NonlinearSolveFirstOrder   QA | 18  18   Testing NonlinearSolveFirstOrder tests passed

A strict Documenter build (doctest = false) of a page with the devdocs/jacobian.md @docs block, plus the docstrings those entries cross-reference, renders the new docstring and resolves its @ref links to reused_jacobian and construct_jacobian_cache.
Runic (--check) and typos are clean on every changed file.

Not run: the full docs/make.jl build (it needs Sundials, HomotopyContinuation, SciPy through Python, and linkcheck), Julia pre, and the other sublibraries. On Julia 1.10, only FirstOrder QA was run, not the Core groups or the other QA groups.

Release boundary

No breaking change. The new public function means NonlinearSolveBase needs a minor bump (2.55.0 → 2.56.0) when it is released. A release of NonlinearSolveFirstOrder or NonlinearSolveQuasiNewton cut from a tree that contains this PR must raise its floor to NonlinearSolveBase = "2.56". The current floors are "2.54" and "2.43". Both admit Base releases without jacobian_computed_at_construction or the new JacobianCache field, and with those the packages fail to load. The reviewer reproduced that failure with Base 2.55.0, and so did I. Julia 1.10 ignores [sources], so a plain GROUP=QA Pkg.test() of lib/NonlinearSolveFirstOrder on 1.10 resolves the registered NonlinearSolveBase v2.55.0, and precompilation then fails with UndefVarError: jacobian_computed_at_construction not defined. CI does not hit this, because the sublibrary workflow develops the in-repo [sources] on 1.10. The repo already runs DowngradeSublibraries.yml. I did not check whether it resolves NonlinearSolveBase from the registry or from the in-repo path, so it is not established that it would catch this.

---

🤖 Posted by an AI agent — harness: Claude Code 2.1.301 · model: claude-opus-5-5
Conversation: local session ID ff2c686d-08da-4517-80f8-8868d37c73f8 (headless, no URL)

Breaking change

These version bumps are required, not optional. This PR deliberately changes no version or compat entry: when and how to release is Chris's decision. This section restates the "Release boundary" section above as a hard release condition.

• NonlinearSolveBase needs a minor bump, 2.55.0 → 2.56.0. It gains the public function jacobian_computed_at_construction and a new JacobianCache field. No existing public NonlinearSolveBase API is removed or changed.
• NonlinearSolveFirstOrder and NonlinearSolveQuasiNewton must require NonlinearSolveBase = "2.56" in the first release that contains this code. Their current floors are "2.54" and "2.43". Both packages import the new function at load time. With any NonlinearSolveBase older than 2.56, they fail to precompile with UndefVarError: jacobian_computed_at_construction not defined; both independent reviewers and I reproduced this against Base 2.55.0. Neither package may be registered from a tree containing this PR until its floor is raised. Their own public API is unchanged.
• The floors cannot be raised in this PR without the Base version bump. With NonlinearSolveBase = "2.56" in FirstOrder and the in-repo Base still at 2.55.0, Pkg.resolve() fails: ERROR: empty intersection between NonlinearSolveBase@2.55.0 and project compatibility 2.56.0 - 2.
• No existing CI job enforces this boundary. Sublibrary CI develops the in-repo [sources], so it always tests against the in-repo Base. The shared downgrade workflow (sublibrary-downgrade.yml) adds every in-repo lib/* sublibrary to its skip list, so it never pins NonlinearSolveBase to a floor. An enforceable check would be a job that loads each sublibrary against the registered versions of its in-repo siblings. That would be a change to the shared SciML workflow, outside this PR.

Review findings addressed (re-review)

1. P1, the compat floors still allow a combination that fails to load. Not changed in code. Fixing it means bumping NonlinearSolveBase's version, which is left to Chris; the required bumps are stated unconditionally under "Breaking change" above.
2. Non-blocking: the early sections of the description still contain superseded claims. Not rewritten, because edits to this description are append-only. "Review findings addressed" above lists the corrections: the helper is a public top-level NonlinearSolveBase function, both QA ignore lists match master, and the PR adds one public name.

No code changed in this round. The head is still 5f24a5421a11799c1efc0eed9fff186cddb5f7bf.

---

🤖 Posted by an AI agent — harness: Claude Code 2.1.301 · model: claude-opus-5-5
Conversation: local session ID ff2c686d-08da-4517-80f8-8868d37c73f8 (headless, no URL)

Compat floors resolved

NonlinearSolveBase is bumped from 2.55.0 to 2.56.0, the next minor above upstream master at 8d5e163e. NonlinearSolveFirstOrder and NonlinearSolveQuasiNewton now require NonlinearSolveBase = "2.56", so their dependency constraints exclude Base releases without the public jacobian_computed_at_construction helper. A repository-wide search found no other packages importing or calling the helper; the top-level package and other sublibraries need no floor change. No fallback shim is added.

Compatibility verification (Julia 1.12.7 and 1.10.12)

On Julia 1.12.7, timeout 3600 julia +1.12 --startup-file=no --project=lib/<package> -e 'using Pkg; Pkg.resolve(); Pkg.instantiate()' succeeds for Base, FirstOrder and QuasiNewton. Both solvers select the in-repo NonlinearSolveBase 2.56.0.

Julia 1.10.12 was checked with clean copies of each Project.toml, without manifests or developed siblings, so [sources] is ignored. Before this fix, both solver projects resolve/instantiate with registered Base 2.55.0. With the raised floors, both reject that combination at resolution:

ERROR: Unsatisfiable requirements detected for package NonlinearSolveBase [be0214bd]:
 ├─possible versions are: 1.0.0-2.55.0 or uninstalled
 └─restricted to versions 2.56.0-2 by NonlinearSolveFirstOrder [5959db7a] — no versions left

QuasiNewton reports the same restriction from NonlinearSolveQuasiNewton [9a2c21bd]. Base's own project resolves/instantiates on 1.10.12. Base 2.56.0 is not registered in the registry used for these checks; it must be registered before these solver projects can resolve on Julia 1.10 without developed siblings. No workaround was applied.

Local tests and formatting

All test invocations use GROUP=<group> NONLINEARSOLVE_TEST_GROUP=<group> timeout 3600 julia +1.12 --startup-file=no --project=lib/<package> -e 'using Pkg; Pkg.test()', four Julia threads, four precompile tasks and one OpenBLAS thread.

Package	Core assertions passed	QA checks passed
NonlinearSolveBase	558	20
NonlinearSolveFirstOrder	9089	20
NonlinearSolveQuasiNewton	1132	20

Each complete invocation ends with Testing <package> tests passed. The public predicate test passes 9/9, FirstOrder reuse 111/111, QuasiNewton reuse 66/66, native bounded methods 6001/6001, and operator Jacobians 14/14.

The first FirstOrder Core attempt reached the 3600-second timeout during BandedMatrices precompilation for a sparsity test, after 9059 passing assertions and no failed assertions. After Base Core completed that dependency's precompilation, the identical full FirstOrder Core command passed. No tests or assertions were changed.

Real Runic 1.10.0 --check passes on all 15 Julia files in the full PR. The requested $HOME/sandbox/goals/performance/runic-env project is absent; the command used a local copy of an existing Runic environment through JULIA_LOAD_PATH, without modifying a global environment. typos and git diff --check pass on all four follow-up files (three Project.toml entries and the repository compatibility convention in AGENTS.md).

Not rerun: full docs build, GPU, downstream tests, or Julia 1.10 runtime suites. The three existing-head CI failures (top-level LTS Core on Ubuntu/macOS and SCCNonlinearSolve Core) have not been compared with master, and this follow-up has not yet received independent review.

The independent review classified the full stacked diff as MH1 (medium-high); this follow-up updates three Project.toml entries and records the sibling-API compatibility convention in AGENTS.md. A fresh independent review and comparison of the outstanding CI failures with master are still pending.

---

🤖 Posted by an AI agent — harness: Claude Code 2.1.301 · model: claude-opus-5-5
Conversation: local session ID ff2c686d-08da-4517-80f8-8868d37c73f8 (headless, no URL)

Independent review (Codex CLI 0.162.1 gpt-6.1-sol, head 5f24a54)

🤖 Automated review section from an AI agent running as Chris Rackauckas's account — not written or reviewed by Chris.

VERDICT: CHANGES
Reviewed the full 16-file diff and supplied discussion; no linked issue was supplied. The Jacobian reuse regression tests have real failing-before/passing-after evidence, but the dependency compatibility blocker remains unresolved.

Local verification used Julia 1.12.7, JULIA_NUM_THREADS=4, JULIA_NUM_PRECOMPILE_TASKS=4 and timeout 3600. Against merge-base dd8a550, the new FirstOrder test has 61 passes/50 failures, QuasiNewton 39 passes/27 failures, and residual preservation 2 passes/2 failures. At head these pass 111/111, 66/66 and 4/4; the public predicate test passes 9/9. An independent analytic-root probe covering six algorithms, in-place/out-of-place residuals and parameter-changing reinit passes 48/48 assertions on each tree. All 12 first iterates are identical in that probe, with first-step njacs decreasing from 1 to 0.

G1: analytic Jacobian expectations and the independent root checks provide accuracy evidence. G2: the supplied rollup reports passing CUDA and Adjoint groups, and their test wiring covers the affected solvers; these were not rerun locally. G3: the new moved-state and reinit tests exercise cache invalidation and pass. The separate requirement to establish no CI regressions against master remains unresolved.

Matched criteria include L1 (documentation), L2 (stronger tests), ML7 (new public helper), and ML8 (cache validity state). The highest is MH1: the full diff includes the private-residual-buffer change at lib/NonlinearSolveBase/src/jacobian.jl:132, deliberately changing residual handling during Jacobian evaluation on the default Newton/trust-region lanes. The reuse optimization alone targets ML11, but this review includes its stacked prerequisite.

1. P1 — Enforce the new NonlinearSolveBase requirement. lib/NonlinearSolveFirstOrder/src/NonlinearSolveFirstOrder.jl:46 and lib/NonlinearSolveQuasiNewton/src/NonlinearSolveQuasiNewton.jl:41 import jacobian_computed_at_construction, while lib/NonlinearSolveFirstOrder/Project.toml:50 still permits Base 2.54 and lib/NonlinearSolveQuasiNewton/Project.toml:46 permits Base 2.43. Confirmed by running the head solver packages against the base-tree NonlinearSolveBase 2.55.0: FirstOrder precompilation fails with UndefVarError: jacobian_computed_at_construction not defined in NonlinearSolveFirstOrder; initializing Broyden with true-Jacobian initialization fails with the corresponding UndefVarError in NonlinearSolveQuasiNewton. These are declared-compatible combinations. The release notes acknowledge this but do not enforce it. Coordinate the Base version bump and both compat floors before shipping the dependent changes, and verify against permitted sibling versions. A future release reminder does not resolve this finding.

Not verified: full Core/QA suites, timing benchmarks, docs, GPU/adjoint execution locally, Julia LTS/pre, or tests on the latest master tip (the before tests used merge-base dd8a550). The three failing supplied checks—Core on LTS Ubuntu, Core on LTS macOS, and SCCNonlinearSolve Core on Julia 1—could not be compared with master: no failure logs or master run evidence were supplied, and credential-free GitHub CLI access was unavailable. None is classified as pre-existing or as a confirmed regression. The supplied green affected-library, downstream, docs and accelerator checks do not settle those failures.

🤖 Posted by an AI agent — harness: Codex CLI 0.162.1 · model: gpt-6.1-sol
Conversation: Codex exec session 01a14b26-dcfc-7b23-964d-a8740ab171ea (no public URL)

Codex CLI gpt-6.1-sol rated it medium-high (MH1) at head 5f24a54

Risk assessment

• Risk: high
• Blast radius: the first iteration of every GeneralizedFirstOrderAlgorithm solve and of true-Jacobian QuasiNewton solves with an AD Jacobian computed at init. Also polyalgorithm subcaches (aliasing guard), reinit! paths, and the DampedNewtonDescent :simple init. One public, documented NonlinearSolveBase helper is added, with the minor version bumped to 2.56.0.
• Evidence: all three libraries pass Core and QA locally on Julia 1.12.7. Clean resolution selects Base 2.56.0 on Julia 1.12 and rejects registered Base 2.55.0 under the new solver floors on Julia 1.10. The original regression evidence and independent review are retained above.
• Independent review: Codex CLI 0.162.1, gpt-6.1-sol rated the full stacked diff medium-high (MH1) at head 5f24a54; see the independent review below. Review of the compat-floor follow-up is pending.
• Merge: needs human review. It is stacked on Evaluate Jacobians into a private residual buffer #1365 and changes behavior on the hot path.

---

🤖 Updated by an AI agent — harness: Codex CLI 0.163.0 · model: gpt-6.1-sol
Session: local session, transcript on a fleet worker host (headless, no public URL)

🤖 Automated review section from an AI agent running as Chris Rackauckas's account — not written or reviewed by Chris.

Independent executed findings at head 850fbf9 — session 01a15e30-7596-7a03-9c72-51d911618e85

The first-step cache reuse ignores explicit recompute_jacobian=true, violating the documented override. Raising the required Base floor also causes a reproduced Julia1.10 benchmark bootstrap resolution regression.

Preserve dependency floors while supplying complete local sibling graph to LTS benchmark setup; do not hide the red regression or release packages here. Add discriminating explicit-true first-step test. Original current commit session attribution still names path/host; no identity/history repair is inferred. Full stack and higher existing human/API risk hold remain.

Risk assessment

• Risk: medium-high (MH1); higher existing parent-stack or author assessments remain controlling.
• Blast radius: The first-step cache reuse ignores explicit recompute_jacobian=true, violating the documented override. Raising the required Base floor also causes a reproduced Julia1.10 benchmark bootstrap resolution regression.
• Evidence: Exact changed override: head dual-evaluations0/njacs0, clean master and prerequisite1/1; explicitfalse0onall. Added tests181/181 and independent prescribed-root/first-step50/50 alltrees pass; Base/FirstOrder/QuasiNewton QA20each. Clean LTSconsumer headcannotresolveBase2.56, masterresolves, matching actual failed benchmark log. FirstOrder Core timedout after1hour with8933passing assertions; completeCore not claimed. Actual CI CUDA22+termination20,adjoint6+wrapper3+Enzyme1+structural3 logs read. SCCdefault and head12pass1fail separately investigated.
• Independent review: Codex CLI0.163.0, model gpt-6.1-sol, local session 01a15e30-7596-7a03-9c72-51d911618e85 (no session URL), reports CHANGES/medium-high (MH1) at head 850fbf9. Findings-only: supplied intake exposed historical verdict/risk labels, so this rating does not count toward merging.
• Merge: Held until current-head counted-review requirements, no-regression gates, open findings and applicable human/API/provenance decisions are satisfied. No old rating is repinned; no approval, merge, closure or release is performed by this review.

Links:

• PR: Reuse the init-time Jacobian on the first Newton step (stacked on #1365) #1367
• Reviewed source: 850fbf9

---

🤖 Posted by an AI agent — harness: Codex CLI0.163.0 · model:gpt-6.1-sol
Conversation: local session 01a154f0-8986-78e1-8af2-1119e8cdea4a (no session URL)

[ChrisRackauckas]

🤖 Automated comment from an AI agent running as @ChrisRackauckas — not written or reviewed by Chris.

Independent review (Devin CLI 3000.11.3, model fusion-claude-opus-5-5-high-sidekick-swe-2-medium): MERGE, risk medium.

Full review

VERDICT: MERGE
RISK: medium

Merge only after #1365 (#1365), which is still open. The risk rating is medium because the change sits on the first iteration of every GeneralizedFirstOrderAlgorithm solve and every true-Jacobian QuasiNewton solve. No public API changes.

Blocking findings

None.

Non-blocking findings

1. Test file names collide with sibling PRs (confirmed with gh pr list). These open PRs also add the same test files:

   • lib/NonlinearSolveFirstOrder/test/misc_tests__item11.jl: Opt-in bounds_handling = BoundsProjection() for Newton-type solvers #1354 and Fix EisenstatWalkerForcing2 two-iteration-stale residual ratio #1351.
   • lib/NonlinearSolveQuasiNewton/test/core_tests__item14.jl: Opt-in bounds_handling = BoundsProjection() for Newton-type solvers #1354, Jacobian reuse across reinit! and for matrix-free preconditioners #1353 and Add a gradient-stationarity termination criterion and make it the least-squares default #1181.

   Whichever PR merges second gets a file conflict and has to renumber its test file. Jacobian reuse across reinit! and for matrix-free preconditioners #1353 adds reinit!(cache; reuse_jacobian = true). Its semantics look compatible with this PR: here reinit! only clears initial_jacobian_valid, and Jacobian reuse across reinit! and for matrix-free preconditioners #1353's reuse decision stays with make_new_jacobian. Still, check the two together when rebasing.

2. Mutating p in place between init and the first step! now also leaves J stale (confirmed by running). The PR body acknowledges this. scratch/edge.jl does cache.p .= [20.0, 3.0] after init:

   • Base: step 1 uses a stale fu with a fresh J.
   • Head: step 1 uses a stale fu and a stale J. The first iterate changes ([0.590, 0.759] → [1.182, 1.062]).
   • The final u matches a fresh solve to 1e-16 (NR) and 7e-15 (Broyden) in both base and head.
   • TrustRegion stalls on both base and head, so that problem is pre-existing.

   fu was already stale on this path, so it was never supported. It would only matter if some caller (an SII setp on the cache, SCC-style reuse) mutates p without reinit!. I grepped SCCNonlinearSolve's src for reinit!, .p =, setp and set_parameter and found no hits. I did not trace how SCCNonlinearSolve actually hands p to its subproblem solvers (see "What I did not verify").

3. Missing regression test for the real polyalgorithm case. The u == u_cache guard is only tested by mutating cache.u by hand. No test in the repo runs the scenario that motivated the guard: NonlinearSolvePolyAlgorithm with alias_u0 = true, where the first subalgorithm moves u. My probe shows this case is bitwise identical to base:

   • dense_iip/Poly_PT_NR: cache-path njacs Δ = −1, not −2, so the guard fires.
   • svec/Poly_PT_NR: Δ −2 / −2, because an immutable u is not aliased.

   A test along the lines of Poly_PT_NR in scratch/probe.jl would pin this down.

4. Suspected from reading only, pre-existing, outside this PR: double damping in :simple mode on reuse steps. In damped_newton.jl (the :simple branch of solve!, around line 296), J is always re-damped, because J !== nothing holds in FirstOrder. cache.J aliases the Jacobian-cache storage when alias_J = true and D is not a matrix, so dampen_jacobian!!(cache.J, J, D) damps that storage in place. On a step where the reuse policy keeps J (PseudoTransient(; jacobian_reuse = true)), the storage would then be damped a second time. This PR fixes exactly that hazard for init → step 1 only. Base and head behave identically for PT_reuse in my probe, so this is not a regression. It is worth a separate issue if it is confirmed.

5. Attribution/process format doesn't follow CLAUDE.md (confirmed from gh pr view and git log):

   • The PR body's first line is "Please ignore until reviewed…", not the required > 🚧 **UNREVIEWED …** banner.
   • The footer is "Generated with Claude Code (model…)", not the required footer format.
   • There is no Links section.
   • The commit's Agent-Harness: Claude Code has no version.

6. NonlinearSolveBase compat floors aren't raised. FirstOrder (floor 2.54) and QuasiNewton (floor 2.43) now call InternalAPI.jacobian_computed_at_construction. With a NonlinearSolveBase release older than this PR, they would throw UndefVarError at init. I checked git log on master: repo convention raises the floor in a separate release PR (for example Release NonlinearSolveFirstOrder v2.10.0 #1313), and the PR body says so. Whoever does the release needs to remember this.

7. Minor cost. In the aliased :simple case, DampedNewtonDescent.init now makes an n×n copy(J). Factorization solvers with alias_A_for_refactorization copy A again inside construct_linear_solver, so PT init does two copies. This is acceptable, as the PR notes.

What I ran

All of this was in scratch/, with two Julia 1.12 envs:

• env_head: the PR tree, 3c0bb26.
• env_base: git archive 7dd45902, i.e. Evaluate Jacobians into a private residual buffer #1365 without this PR's commit.

pathof confirms each env loads its own tree. isdefined(InternalAPI, :jacobian_computed_at_construction) is true in head and false in base.

Rerun of the PR's new tests (confirmed):

Test file	head	base
misc_tests__item11.jl	111/111 pass	61 pass / 50 fail. First failure at :50, calls == 0 evaluated 8 == 0
core_tests__item14.jl	66/66 pass	39 pass / 27 fail. First failure at :33, 16 == 8

The tests fail before the change and pass after it.

The assertions have real force:

• The reinit test checks reinit_sol.u == fresh.u bitwise. The fresh solve uses the reused init-time J and the reinit solve recomputes J, so this cross-checks that the reused J equals the recomputed one.
• It also checks the damping fix for PT.
• No tolerances were fitted, and no tests were skipped or loosened.

Neighbouring tests on head (confirmed):

• misc_tests__item9.jl (adaptive Jacobian reuse): 125/125.
• operator_jacobian.jl: 14/14.
• NonlinearSolveBase/test/jacobian_private_fu.jl: 4/4.
• QuasiNewton core_tests__item12.jl: 63/63.

Independent bitwise A/B, base vs head (confirmed): scratch/probe.jl, compared with scratch/compare.jl; output is in scratch/compare.out. I wrote this harness myself rather than reusing the author's.

• Problems:
  • 30-dim dense with sin: IIP, OOP, and IIP with a user jac.
  • 40-dim sparse Bratu with TracerSparsityDetector.
  • SVector{4}.
  • A 3→7 NonlinearLeastSquaresProblem.
• Algorithms:
  • NewtonRaphson variants: plain, FiniteDiff, jacobian_reuse = true, GMRES with concrete J, lazy GMRES, ForwardDiff chunksize 2.
  • TrustRegion: default, Hei and Yuan radius schemes.
  • LevenbergMarquardt, GaussNewton.
  • PseudoTransient, with and without jacobian_reuse.
  • Broyden(true J), Klement(true J), Klement(true diagonal).
  • PolyAlg(PT, NR) and PolyAlg(PT, Broyden-true-J), both with alias_u0 = true, maxiters = 8.
• Each configuration compares, bitwise with isequal:
  • a one-shot solve;
  • every iterate u and fu of a step! loop;
  • init+solve!;
  • reinit! + solve! with a new u0 and p, a new p only, and the original values.

Results:

• 41/41 configurations are bitwise identical. That covers every iterate, residual, retcode, nf, nfactors, nsolve and nsteps.
• Klement(true J) on dense IIP throws the same DomainError on both base and head. That is pre-existing.
• njacs drops by 1 on the solve paths exactly where J is computed at init. It doesn't change for lazy GMRES, a user jac, or after reinit!.
• The polyalgorithm rows match the PR's table: Δ −2/−1 with an aliased Vector u, and −2/−2 for SVector.

This reproduces the PR's "bitwise identical, only njacs changes" claim on an independent configuration set.

Edge probes (scratch/edge.jl, confirmed):

• In-place p mutation: see non-blocking finding 2.
• PseudoTransient with a mass_matrix (matrix D, so the descent works on an unaliased copy): same steps and solution on base and head (387 steps; njacs 388 → 387).

Formatting, spelling and QA (confirmed):

• Runic --check --diff on the 12 changed files: clean.
• typos on the same files: clean.
• GROUP=QA Pkg.test() for FirstOrder, QuasiNewton and Base, run on an rsync copy of the PR tree: 20/20 pass for each.

Code reading:

• computed_at_construction is true exactly on the branch where construct_jacobian_cache calls DI.jacobian, and false for the scalar method.
• OverrideInit goes through update_initial_values! → reinit!, which clears both flags (FirstOrder reinit_self!; InitializedApproximateJacobianCache.reinit_self!).
• Nothing between construct_jacobian_cache and the first step mutates J's storage except the :simple damping that this PR fixes:
  • Newton descent and construct_linear_solver copy or don't touch A.
  • linsolve_workspace copies J.
  • TrustRegion and LineSearch init only read J.
• No other code calls the JacobianCache or GeneralizedFirstOrderAlgorithmCache positional constructors.
• About 12 of the 238 added lines are comments (~5%), and none of them narrate history.

What I did not verify

• I did not rerun the timing benchmarks (−31% / −21% / −9% on dense 100). I only confirmed the njacs reductions they rest on.
• I did not run the full GROUP=Core suites for FirstOrder, QuasiNewton or Base. Notably that leaves out the 6001-test native-bounded set, which the PR claims passes. I ran only the targeted files above.
• I did not run the top-level NonlinearSolve tests, SCCNonlinearSolve, SpectralMethods, Julia LTS or pre-release, GPU arrays (the u == u_cache guard on device arrays), the docs build, or Enzyme/Zygote backends.
• I did not test the MTK OverrideInit path end to end. I traced it by reading code only.
• Finding 4 (double damping on reuse steps) is from reading only.

---

🤖 Posted by an AI agent — harness: Devin CLI 3000.11.3 (review), Claude Code 2.1.285 (posting) · model: fusion-claude-opus-5-5-high-sidekick-swe-2-medium
Conversation: local Claude session 3c311569-dda6-4687-bc9f-65febb212c33 (fleet master)