Add Kantorovich homotopy continuation

[ChrisRackauckas-Claude]

Please ignore this draft until it has been reviewed by @ChrisRackauckas.

Summary

• Add KantorovichHomotopy, a documented HomotopyProblem solver whose parameter-step controller is adapted from ImplicitDiscreteSolve.KantorovichTypeController.
• Reuse the HomotopySweep path driver, including cache reuse, endpoint polishing, predictor trust, tracking limits, Jacobian forwarding, and warm handoff to ArcLengthContinuation.
• Measure inner-corrector residual contraction when the selected solver exposes an iterative cache; use the controller's Θmin convention for polyalgorithms and no-init solvers.
• Add strict contraction rejection, the Deuflhard step-size formula, constant/secant predictor orders, validation, tests, docs, and the required minor version bumps.

Why

ImplicitDiscreteSolve uses a Newton--Kantorovich continuation-style controller across accepted discrete states, but it is not itself one fixed HomotopyProblem: its residual mutates the previous accepted state. The reusable part is the contraction-based path-step controller. This PR implements that controller for a genuine fixed homotopy while retaining NonlinearSolve's existing continuation machinery.

OrdinaryDiffEq already has a generic HomotopyNonlinearSolveAlg adapter that builds the normalized per-stage homotopy from the known zero-step anchor. No adapter-specific source change is needed: callers can pass HomotopyNonlinearSolveAlg(KantorovichHomotopy(...)) once this solver is available.

This feature branch is rebased directly onto current master; the former #1075 and #1079 prerequisites are now included upstream.

Comparison

After correcting the controller to use the actual first residual contraction (predictor residual to first corrected residual), explicit finite-difference Newton comparisons are mixed. Values below are homotopy residual-function evaluations.

Problem	Sweep (secant)	Kantorovich (constant)	Kantorovich (secant)	ArcLength
Smooth cubic	66	51	57	101
Curved exponential	96	150	99	230
Narrow transition	78	216	285	1497
Coupled cubic (n=50)	538	792	483	707
Two-fold S curve	478 success via path jump	3163 strict failure	2843 strict failure	320 success

The predictor and controller effects are distinct. Secant prediction makes the Kantorovich controller competitive on the curved exponential and best on the coupled path, but it remains much worse through the narrow transition. The constant predictor matches ImplicitDiscreteSolve's method, yet is worse than HomotopySweep on three of these five standalone paths. Natural-parameter Kantorovich continuation cannot traverse a fold; strict rejection detects the bad contraction, while ArcLengthContinuation follows the branch.

The existing OrdinaryDiffEq adapter shows real wins on some stiff workloads but the same lack of uniform robustness. Values are sol.stats.nf from ImplicitEuler runs.

Workload	Sweep (secant)	Kantorovich (constant, strict)	Result
Stiff van der Pol	139749	99531	both succeed; 28.8% fewer calls
Robertson	3225477	2366489	both succeed; 26.6% fewer calls
Fixed-step ignition, dt=5	2138 success	1564 then failure at step 36/40	strict controller does not finish
Fixed-step ignition, dt=10	1160 success	882 then failure at step 6/20	strict controller does not finish
Fixed-step ignition, dt=20	3686 success	5738 success	55.7% more calls

At ignition dt=10, strict=false finishes in 751 calls, confirming that the apparent speedup comes from accepting a corrector whose initial contraction violates the strict path-safety criterion. That permissive setting can jump branches and is not a safe general default.

The conclusion is therefore conservative: this controller is useful as an opt-in experiment and can reduce work substantially on some smooth stiff-stage homotopies, but it is not a replacement for HomotopySweep, and ArcLengthContinuation remains the fold solver. The default homotopy polyalgorithm is intentionally unchanged.

Validation

• Rebased head (63f359999): 50 passed / 50 total in 46.5s, including a retained shrunken-span Kantorovich → arclength warm handoff.
• On the patch-equivalent pre-rebase head, all 40 continuation regression files on Julia 1.10: 527 passed / 527 total in 6m29.8s.
• The same pre-rebase aggregate on Julia 1.12: 527 passed, 2 failed / 529 total in 8m54.8s; the only failures are the two pre-existing allocation assertions reproduced on clean master and tracked in #1076.
• GROUP=QA: 15 passed / 15 total.
• End-to-end OrdinaryDiffEq adapter smoke with stiff van der Pol: ReturnCode.Success, endpoint error 1.450765923494135e-4 against a tight TRBDF2 reference.
• Repository-wide Runic check: exit 0.
• git diff --check: exit 0.

GROUP=Core reaches an unrelated failure in the pre-existing Broyden test before the continuation section. The failure was reproduced on clean master and bisected separately; #1057 is the generalized-Rosenbrock fix and #1056 tracks the platform-dependent Brown behavior. This PR does not alter or suppress that test.

Process notes

The work audited the ImplicitDiscreteSolve controller and residual mutation semantics, reviewed the existing OrdinaryDiffEq homotopy adapter, extracted the shared continuation driver hooks, tested exact controller factors and rejection behavior, benchmarked smooth, curved, stiff-ODE, and fold cases, then verified documentation and the full continuation regression surface. Reproduction details and clean-master CI investigations are recorded in the PR comments.

[ChrisRackauckas-Claude]

Commit affc2ecd8 implements the planned controller extraction and comparison.

Execution notes:

1. Audited ImplicitDiscreteSolve and confirmed it is continuation-like across mutable accepted states, not one fixed HomotopyProblem.
2. Reused the existing HomotopySweep driver and isolated controller-specific corrector measurement, accepted-step sizing, rejected-step sizing, and strict rejection hooks.
3. Kept the current default polyalgorithm unchanged after benchmarking showed mixed results and fold sensitivity.
4. Verified the existing OrdinaryDiffEq HomotopyNonlinearSolveAlg adapter accepts this solver without adapter source changes.

Local evidence on the committed tree:

• Kantorovich homotopy | 43 pass / 43 total.
• Continuation regression files | 520 pass / 520 total (all homotopy and arclength Core files).
• GROUP=QA ... Pkg.test() → 15 pass / 15 total.
• Repository-wide Runic --check . → exit 0, no output.
• End-to-end stiff van der Pol adapter solve → retcode = Success, endpoint error 0.0001450765923494135.

The unrelated Core Broyden stop was reproduced and bisected separately; details and links are in the PR body.

[ChrisRackauckas-Claude]

CI note: the red Downgrade / Downgrade Tests - Core job is a current-master baseline failure, not introduced by this PR. The same SciMLOperators resolver conflict appears in master run 29393716642, job 87282441811: the downgrade pins SciMLOperators 1.22.0 while the in-tree NonlinearSolveQuasiNewton requires 1.24+. The focused fix is already under review in #1070; I am independently reproducing/bisecting it per repository instructions and will rerun this PR after the baseline fix lands.

[ChrisRackauckas-Claude]

Additional CI note: downgrade-sublibraries / test (lib/SimpleNonlinearSolve) is the same current-master dependency-floor failure. Master run 29393716591, job 87282459268 pins SciMLOperators 1.22 while the in-tree NonlinearSolveBase requires 1.24; this PR’s job differs only in the expected Base patch version (2.36 instead of 2.35). The #1070 baseline investigation is being extended to cover both the root and SimpleNonlinearSolve downgrade environments without skipping lower-bound tests.

[ChrisRackauckas-Claude]

The third red check, PolyAlgorithms (Julia 1, macos-latest), is runner networking after successful tests: the job reports NonlinearSolve tests passed, then coverage processing cannot resolve pkg.julialang.org or github.com while installing CoverageTools. No source/test failure occurred. The current token cannot rerun upstream Actions jobs (GitHub returns HTTP 403); the next branch update for the separately verified downgrade prerequisite will naturally retrigger it.

[ChrisRackauckas-Claude]

The NonlinearSolveBase [QA] red is another exact current-master failure: both this PR job and master job 87282502519 report that NonlinearSolveBaseForwardDiffExt imports non-public FunctionWrappersWrappers.SingleCacheStorage. The focused local QA run used Julia 1.10, where Julia public checks are inactive; CI on Julia 1.12 exposed it. A separate clean-master dependency/API investigation is queued per repository instructions; this PR will not add an ignore or reach further into the dependency internal.

[ChrisRackauckas-Claude]

Supplemental current-stable verification: after instantiating the repository-local Julia 1.12 environment, timeout 3600 julia +1.12 --project=. test/Core/kantorovich_homotopy_tests__item1.jl exited 0. git status --short and git diff --check remained empty afterward. This is in addition to the previously reported 43/43 focused suite and 520/520 continuation regression run.

[ChrisRackauckas-Claude]

Rebased onto current master after #1070 merged, after fetching both remotes and confirming the remote feature branch still contained only the prior agent commit. Post-rebase local verification: Julia 1.12 focused suite 43 passed / 43 total; repository-wide Runic check exited 0; git diff --check and worktree status are clean. The branch was updated with --force-with-lease (fad6ad6ab). CI has restarted against the corrected downgrade floors.

[ChrisRackauckas-Claude]

Post-rebase CI note: Downgrade / Downgrade Tests - Core reached the pre-existing Broyden block and stopped with the exact clean-master baseline: 94 passed, the generalized Rosenbrock case failed for algorithm 4, the Brown almost-linear case produced one platform-dependent unexpected pass, and 19 remained broken. This is the separately reproduced/bisected issue described in the PR body; #1057 contains the generalized-Rosenbrock fix and #1056 tracks Brown. No homotopy test failed, and this PR does not alter, skip, or loosen the baseline tests.

[ChrisRackauckas-Claude]

Follow-up commit a3bd3b1 corrects the contraction observation to measure the first Newton ratio from the predictor residual, rather than starting at the second ratio. The new regression deliberately has first contraction 0.99 and later contraction 0, so the prior implementation would select the wrong value. Exact post-change local run: Kantorovich homotopy 48 passed / 48 total in 2m47.6s on Julia 1.12; repository-wide Runic and git diff --check both exited 0. The branch was current with upstream/master and was pushed normally without force.

[ChrisRackauckas-Claude]

Corrected benchmark note: measuring the true first contraction materially changed the comparison, so I reran rather than reusing the earlier numbers and replaced the PR table. Default strict Kantorovich is 28.8% lower in residual calls on stiff van der Pol and 26.6% lower on Robertson, but it is not uniformly good: it fails to finish fixed-step ignition at dt=5 and dt=10, and costs 55.7% more at dt=20. On standalone paths it wins the smooth cubic, loses badly on a narrow transition, and correctly fails at a fold where ArcLength succeeds. At ignition dt=10, strict=false recovers a 751-vs-1160 call win, but only by accepting the poor initial contraction, which reintroduces path-jump risk. The conclusion is now explicitly that this remains opt-in and must not replace the existing default. Temporary benchmark environments/scripts were removed and both worktrees are clean.

[ChrisRackauckas-Claude]

Benchmark reproduction details: standalone runs used Julia 1.12, explicit NewtonRaphson(autodiff=AutoFiniteDiff()), abstol=reltol=1e-10, maxiters=50, and counted every H evaluation. Paths were H=u^3-(1+lambda); H=u^3-exp(2lambda)^3; H=u^3-(1+0.5tanh(25(lambda-0.5)))^3; a 50-variable tridiagonally coupled cubic with known path u(lambda)=1+lambda*v; and H=u^3-3u-(-3+6lambda) for the two-fold curve. ODE runs used the existing HomotopyNonlinearSolveAlg adapter with ImplicitEuler and sol.stats.nf: the repository stiff-VDP problem (mu=1000, tspan 0:0.3, tolerances 1e-6), repository Robertson problem (tspan 0:1e5, tolerances 1e-8), and u prime = u^2(1-u), u0=0.01, tspan 0:200 at fixed dt 5/10/20. Each algorithm was reconstructed for every run.

[ChrisRackauckas-Claude]

Post-fix aggregate continuation verification is complete. Each test/Core/*homotopy* / *arclength* regression file was included in an isolated module to avoid cross-file constant collisions.

Julia 1.10 exited 0:

Test Summary:                 | Pass  Total      Time
Continuation regression files |  525    525  10m17.0s

Julia 1.12 produced:

Test Summary:                 | Pass  Fail  Total      Time
Continuation regression files |  525     2    527  15m02.7s

The two Julia 1.12 failures are exactly the pre-existing homotopy_alloc_tests__item1.jl allocation assertions reproduced on clean master (Sweep 154.4 B/step < 48 and ArcLength 80 B/step < 48). They are tracked separately in #1076; no test was disabled, loosened, or changed. The disposable test environments have been removed and the feature worktree is clean.

[ChrisRackauckas-Claude]

The newly reported downgrade-sublibraries / test (lib/NonlinearSolveHomotopyContinuation) resolver failure is an exact clean-master baseline, not introduced by this branch.

• Feature job: https://github.com/SciML/NonlinearSolve.jl/actions/runs/29413260099/job/87349445725 (NonlinearSolveBase 2.36.0)
• Clean master 6264446f4 job: https://github.com/SciML/NonlinearSolve.jl/actions/runs/29408197365/job/87330819019 (NonlinearSolveBase 2.35.0)

Both fail with the same constraint chain: #1070 requires LinearSolve 5 in the source-tree Base package, while the downgraded SciMLBase 3.33 and SciMLLogging 2.0.3 constrain the test-only registered root NonlinearSolve to <=4.20.3, which requires LinearSolve <=4.3. The mandated clean-master investigator is reproducing and repairing this separately; this PR does not suppress or alter the downgrade test.

[ChrisRackauckas-Claude]

Rebased the feature commits onto the focused clean-master downgrade fix in #1075. The remote branch was audited first (only the two known feature commits), and the rewritten push used an explicit lease against old head a3bd3b1026c384ac86e2072922d7ccdcf0016aa7.

New head: 57e36d7a1cce59ddb8057a095aa2a1e0b1c7fca1.

Post-rebase local verification:

Test Summary:        | Pass  Total   Time
Kantorovich homotopy |   48     48  51.1s

• Julia Runic repository check: exit 0.
• git diff --check: exit 0.
• Worktree: clean.

This stack includes #1075 only. The independent root dependency-floor PR #1079 is not a prerequisite for this feature.

[ChrisRackauckas-Claude]

Final documentation audit found that HomotopyPolyAlgorithm's implementation already included KantorovichHomotopy in the warm natural-parameter → arclength handoff, while its public docstring still described the handoff as Sweep-only. Commit 254829887 corrects that documentation and the corresponding source comments; behavior is unchanged.

Local validation after the documentation edit:

Test Summary:        | Pass  Total     Time
Kantorovich homotopy |   48     48  2m30.5s

• Whole-repository Julia Runic check: exit 0.
• git diff --check: exit 0.
• Worktree: clean.

The exact #1071 + #1075 downgrade sandbox was also reproduced locally: NLS 4.21.0 / local Base 2.36.0 / SciMLBase 3.34.0 / SciMLLogging 2.0.3 / LinearSolve 5.0.0, with AllRoots 169/169 and Single Root 50/50, exit 0.

[ChrisRackauckas-Claude]

Added direct regression coverage for the optimized warm handoff in c8541601f: the fold polyalgorithm retains its warm-stage solution, and the test asserts that the fallback's λspan starts strictly inside the original span rather than restarting cold.

Test Summary:        | Pass  Total   Time
Kantorovich homotopy |   50     50  50.3s

Whole-repository Runic and git diff --check both exited 0; worktree clean.

[ChrisRackauckas-Claude]

Latest-head aggregate validation on Julia 1.10.11 (c8541601f): all 40 Core continuation files passed, including sweep, arclength, polyalgorithm, allocation, Jacobian, tolerance, retention, and the new Kantorovich warm-handoff checks.

Test Summary:                 | Pass  Total     Time
Continuation regression files |  527    527  6m29.8s

Exit code 0. No test was skipped, disabled, loosened, or silenced. The isolated temporary environment used to resolve the exact local source tree on Julia 1.10 was removed after the run; the feature worktree is clean.

[ChrisRackauckas-Claude]

Clean-master Adjoint CI failure investigation is complete and confirms an unrelated dependency-floor regression.

• Exact current master 6264446f4, Julia 1.10.11, deployed source-tree helper, GROUP=Adjoint, and CI Pkg.test flags: resolver exits 1 because local NonlinearSolveQuasiNewton 1.14.1 fixes LinearSolve 5.0.0, while released SciMLSensitivity <= 7.116.0 allows only LinearSolve through 4.3.
• Exact immediate parent 247a8bfcb: resolves SciMLSensitivity 7.116.0 + LinearSolve 4.3.0; Adjoint sets pass 3/3, 3/3, and 1/1; package tests pass, exit 0.
• Formal bisect identifies 6264446f4ce3ad82fea62c71590c64ac04f69bac (Align source-tree downgrade dependency floors #1070) as first bad.
• The proper existing fix is SciMLSensitivity draft #1555, which publishes the LinearSolve 5 compatibility already merged in SciMLSensitivity #1549. Injecting that candidate into a disposable clean-master environment resolves SciMLSensitivity 7.116.1 + LinearSolve 5.0.0; the same three Adjoint sets pass and the package exits 0.

No test or dependency constraint in this homotopy PR was changed to mask the baseline. All disposable investigation worktrees, clones, and depots were removed.

[ChrisRackauckas-Claude]

Latest-head Julia 1.12.6 aggregate (c8541601f) completed all 40 continuation files:

Test Summary:                 | Pass  Fail  Total     Time
Continuation regression files |  527     2    529  8m54.8s

kantorovich_homotopy_tests__item1.jl passed 50/50. Every file except homotopy_alloc_tests__item1.jl passed completely. Its two failures are exactly the clean-master allocation baselines tracked in #1076:

• sweep_per_step(prob_big): 154.4 < 48
• arclength_per_step(prob_big): 80.0 < 48

No test was skipped, disabled, loosened, or silenced. The disposable Julia 1.12 environment was removed and the feature worktree is clean.

[ChrisRackauckas-Claude]

Rebased directly onto current clean upstream/master (74405d303) after #1075 and #1079 merged.

Remote safety audit before rewriting:

• fetched upstream and origin;
• remote feature head was exactly the known c8541601f;
• local and remote feature branches had no one-sided commits;
• rebased only the four feature commits, dropping the obsolete stacked prerequisite hash;
• pushed with an explicit force-with-lease bound to c8541601f.

New head: 63f3599995de5865de2feff192ec4f767e7eef01.

git range-diff shows the four feature commits are patch-equivalent; the only first-commit differences are the expected new master baselines (4.21.1 → 4.22.0 and NonlinearSolveBase 2.35 → 2.36). Post-rebase local validation:

Test Summary:        | Pass  Total   Time
Kantorovich homotopy |   50     50  46.5s

• Repository-wide Runic: exit 0.
• git diff --check: exit 0.
• Worktree clean and remote head verified exactly.
• GitHub now reports the PR mergeable; CI has restarted.

[ChrisRackauckas-Claude]

Scheduled one-hour CI audit for rebased head 63f3599995de5865de2feff192ec4f767e7eef01 at 2026-07-20 07:04 EDT:

• Re-fetched both upstream and the fork. Upstream master remains 74405d3035f759a5631ea5f5776d27d4091063bc; this branch is 0 behind / 4 ahead, its worktree is clean, and GitHub reports the draft MERGEABLE.
• Current CI totals: 68 successful, 3 skipped, 7 failed, 2 in progress, 13 queued. The unfinished jobs are macOS matrix lanes.
• Every completed failure was independently reproduced or matched on the exact current clean-master tip:
  • root downgrade resolver conflict, first bad Require operator-capable NonlinearSolveBase #1079: focused fix Restore root LinearSolve 5 floor #1087;
  • duplicate yuan2015recent bibliography key: focused fix Remove duplicate BibTeX entry #1088 (docs CI green);
  • Base QA non-public SingleCacheStorage access, first bad Make AutoSpecialize infer through functor residuals; zero-alloc homotopy stepping #1046: focused public-API fix Use public cache construction for ForwardDiff wrappers #1073, reapplied to current master and QA 18/18;
  • CUDA Broyden DimensionMismatch, first bad Quasi-Newton J⁻¹ init as a linsolve workspace: factorization + identity solve through the LinearSolve default #1039: focused fix Restore GPU Broyden inverse initialization #1066, whose V100 job passes 22/22;
  • Julia 1.12/1.13 continuation allocation signature 154.4/80.0 B per step: focused fix Eliminate dense continuation step allocations on Julia 1.12 #1061; its current-master merge ref passes 15/15 at 0.0/0.0 on Julia 1.12 and 15/15 on Julia prerelease;
  • macOS prerelease Brown/Broyden residual 5.5: exact current-master push failure tracked by Brown bad-Broyden broken test is environment-sensitive #1056; the lane stops before homotopy tests.
• No unrelated compatibility, documentation, public-API, allocation, CUDA, or Broyden fix was stacked into Add Kantorovich homotopy continuation #1071; doing so would defeat the requested direct clean-master rebase.
• Rebased-head local validation remains: Kantorovich homotopy 50/50 passed in 46.5s, repository-wide Runic exit 0, and git diff --check exit 0.

No completed CI failure implicates the new Kantorovich controller.