Fix Float32/static-array mass-matrix OOP Rosenbrock JacReuseState W typing

[ChrisRackauckas-Claude]

Note

This PR should be ignored until reviewed by @ChrisRackauckas.

Fixes #3721.

Problem

A Float32 static-array mass-matrix problem (e.g. DiffEqGPU's CPU-reference ROBER solve) solved with Rosenbrock23() and dt = 0.1 (a Float64, which promotes the integrator's time type to Float64) fails with

MethodError: Cannot `convert` an object of type
  StaticWOperator{true, SMatrix{3,3,Float64,9}, ...} to
  StaticWOperator{true, SMatrix{3,3,Float32,9}, ...}
 [1] setproperty!(x::JacReuseState{Float64, SMatrix{3,3,Float32,9}, ...}, ...)
 [2] calc_rosenbrock_differentiation @ derivative_utils.jl:922

Two distinct failure modes in the OOP Jacobian-reuse path:

1. Seed/solve type mismatch (the reported error). JacReuseState.cached_W is seeded in alg_cache with build_J_W's W, whose eltype follows the state (Float32). At solve time calc_W builds StaticWOperator(J - mass_matrix * inv(dtgamma)), whose eltype is promoted by typeof(dtgamma) (Float64). The fresh-Jacobian branch's jac_reuse.cached_W = W then throws, because StaticWOperator has no eltype-converting convert. (The dense path silently survives the same mismatch only because LinearAlgebra defines convert between LU eltypes.)

2. Reuse branch built the wrong wrapper for static caches. The cached-J reuse branch in calc_rosenbrock_differentiation rebuilt W via DiffEqBase.default_factorize, producing a StaticArrays.LU — not the StaticWOperator the cached_W slot holds. Any OOP static-array mass-matrix solve with at least one rejected step (the mass-matrix decision returns (false, true) on rejection retries) hits this and throws, even in homogeneous Float64. Existing tests passed only because their runs happened to have zero rejections.

Fix

• Add build_jac_reuse_W_seed(W, J, u, mass_matrix, dtgamma) in OrdinaryDiffEqDifferentiation, which builds the cached_W seed with the same concrete type calc_W produces at solve time (mirroring its StaticWOperator / WOperator / AbstractSciMLOperator / factorized-dense branches), given a dtgamma prototype carrying the solve-time dt * gamma type.
• Seed JacReuseState in the three OOP alg_caches (Rosenbrock23, Rosenbrock32, RodasTableauAlgorithms) with that seed and with zero(dt * tab.d) / zero(dt * tab.gamma) for the dtgamma fields (previously zero(dt)), so the field types match solve-time dtgamma exactly. Nested-Dual behavior (Second order differentiation of solve with ForwardDiff errors due to incorrect types in JacReuseState #3486) is preserved: the product keeps dt's full Dual type.
• Make the reuse branch in calc_rosenbrock_differentiation mirror calc_W's branches so it builds a StaticWOperator for static caches.

Verification

Ran locally (all previously failing, now passing):

• The issue's repro: Float32 SVector ROBER + Rosenbrock23() + dt = 0.1 → Success, endpoint matches Float64 reference to <1e-4, with 10 step rejections also exercising the static reuse branch.
• Homogeneous Float64 static ROBER at tolerances with rejections (the previously-latent failure mode 2) → Success, matches reference.
• Rodas4/Rodas23W/ROS34PW2 on the Float32 static problem (previously the same MethodError through the RodasTableau OOP cache) → Success.
• Float64 static result is bit-identical to pre-fix master on the no-rejection run (no behavior change where the old code worked).
• Full Pkg.test Core groups of OrdinaryDiffEqRosenbrock and OrdinaryDiffEqDifferentiation.

Pre-existing issues found while testing, NOT addressed here (all reproduce identically on unmodified master):

• Default CheckInit for OOP static-array mass-matrix problems broadcasts in-place into an SVector (SciMLBase initialization.jl:216); the tests use BrownFullBasicInit() like the DiffEqGPU test does.
• Float32 SVector state with Float64 tspan on a non-mass-matrix problem errors with "Detected non-constant types" (first step promotes u to Float64).
• Rodas5P on Float32 ROBER hits DtNaN (precision limitation; identical on the dense path).

Versions: OrdinaryDiffEqDifferentiation 3.2.0 → 3.2.1 (new function), OrdinaryDiffEqRosenbrock 2.3.0 → 2.3.1 with compat OrdinaryDiffEqDifferentiation = "3.2.1" (it now requires the new function — matters for downgrade CI).

🤖 Generated with Claude Code

[ChrisRackauckas-Claude]

CI triage (none of the failures trace to this PR's changes)

Jobs covering the changed libraries:

Job	Result	Notes
Differentiation Core (julia 1, pre), Sparse, QA, ModelingToolkit	✅ pass
downgrade-sublibraries (lib/OrdinaryDiffEqDifferentiation)	✅ pass	the raised compat lower bound resolves correctly
Rosenbrock Core (lts, pre), QA (1, lts)	❌	tests passed in all four (log shows Testing OrdinaryDiffEqRosenbrock tests passed); the jobs then fail at the coverage step with directory "lib/OrdinaryDiffEqRosenbrock/ext" not found! — the SciML/.github reusable workflow unconditionally passes lib/X/ext to julia-processcoverage and this lib has no ext/. Fix: SciML/.github#84
Rosenbrock Core (julia 1)	❌	runner received a shutdown signal mid-run ("The operation was canceled") — infra flake, needs a re-run (I lack rerun permissions)
Rosenbrock GPU (julia 1)	❌	job cancelled at runner level — same

Failures unrelated to this PR (pre-existing on master):

• Spell Check: fails on unmodified master (oscilations in docs, Lamba/Strang false positives, plus typos from the StochasticDiffEq merge) — fixed in Fix Spell Check (typos) failures on master #3733.
• Downgrade / InterfaceI: Unsatisfiable requirements ... ADTypes driven by OrdinaryDiffEqExtrapolation's ADTypes = 1.22.0 bound — downgrade-CI re-enablement fallout, none of those compat entries are touched here.
• 17 downgrade-sublibraries jobs and ~30 sublibrary-ci jobs in libraries this PR doesn't touch (BDF, Default, SDIRK, IMEXMultistep, Extrapolation, StochasticDiffEq*, …): master's own latest Sublibrary CI run (27214025857) shows 39 failures / 2 successes, so these reproduce without this PR.

Local verification (full output in my run logs): Pkg.test Core groups of both modified libs pass — Rosenbrock 16+56+169 (includes the new 11-assertion regression testset), Differentiation all green.

🤖 Generated with Claude Code

[ChrisRackauckas-Claude]

CI complete (4 checks still queued, all in untouched libraries/downstream). Final tally for the two modified libraries — every job that ran tests passed them, including both downgrade jobs:

• ✅ downgrade-sublibraries (lib/OrdinaryDiffEqRosenbrock) and (lib/OrdinaryDiffEqDifferentiation) pass — this is the job class that motivated Rosenbrock23 on Float32 mass-matrix problem: JacReuseState Float32/Float64 StaticWOperator convert MethodError #3721 (DiffEqGPU downgrade re-enablement).
• ✅ Differentiation Core (julia 1, pre), Sparse, QA (julia 1), ModelingToolkit pass.
• The remaining red on the changed libs is all post-test infra, with the test logs showing tests passed where logs exist:
  • 4× Rosenbrock (lts/pre Core, QA×2): coverage ext bug → Skip missing directories in coverage processing instead of erroring .github#84
  • Rosenbrock Core (julia 1): runner shutdown signal; Rosenbrock GPU + Differentiation QA (lts): runner-level cancelled with no failed steps
  • Differentiation Core (lts): tests passed, then could not lock config file ~/.gitconfig: File exists in a post step — concurrent self-hosted runners (demeter3-15) racing on the shared global gitconfig

A re-run of the failed jobs (needs admin) should clear everything except the known pre-existing failures triaged above.

🤖 Generated with Claude Code

[ChrisRackauckas-Claude]

Reworked per review (98bf3cb): the W typing is now fixed at the source in build_J_W rather than by re-deriving a seed afterwards.

build_J_W computes a dtgamma_prototype = promote(t, dt)[2] * one(constvalue(uEltypeNoUnits)) (Dual-safe; constvalue strips eltypes the same way the tableau constructors do, and non-Number eltypes fall back to dt's type) and builds every non-DAE W with the type calc_W produces at solve time:

• OOP static: StaticWOperator(J - f.mass_matrix * inv(oneunit(dtgamma_prototype)), false)
• OOP dense: lu_instance(J - f.mass_matrix * inv(oneunit(dtgamma_prototype)))
• All WOperator constructions get the dtgamma prototype as the gamma slot (previously promote(t, dt)[2] or raw dt)
• DAE branches untouched — their W isn't built from J − M/(dt·γ)

The build_jac_reuse_W_seed helper is deleted and the OOP Rosenbrock alg_caches pass build_J_W's W straight through. A side benefit: Newton-based OOP caches that store build_J_W's W and later assign calc_W results had the same latent mismatch under mixed precision, now also fixed. For homogeneous-precision problems the constructed types are identical to before, so this is type-neutral where things already worked.

Re-verified locally:

• Issue repro and full battery: bit-identical results to the previous revision (Float32 static Success with rejections, Float64 static endpoint identical to pre-fix master, Rodas4/Rodas23W/ROS34PW2 Float32 static all Success).
• The 11-assertion regression testset passes.
• Full Pkg.test Core suites, all EXIT 0 with tests passed: OrdinaryDiffEqRosenbrock, OrdinaryDiffEqDifferentiation, and — since build_J_W is shared infrastructure — OrdinaryDiffEqSDIRK and OrdinaryDiffEqBDF (the latter covering the DAE branches and the heaviest Newton-cache usage).

🤖 Generated with Claude Code