Skip to content

Use stiff solver in DEIM test#164

Merged
ChrisRackauckas merged 2 commits into
SciML:mainfrom
ChrisRackauckas-Claude:agent/fix-deim-stiff-solver
Jul 12, 2026
Merged

Use stiff solver in DEIM test#164
ChrisRackauckas merged 2 commits into
SciML:mainfrom
ChrisRackauckas-Claude:agent/fix-deim-stiff-solver

Conversation

@ChrisRackauckas-Claude

Copy link
Copy Markdown
Member

This PR should be ignored until reviewed by @ChrisRackauckas.

Summary

  • Switch the DEIM PDE test from Tsit5() to Rodas5P() for both the full and reduced solves.
  • Keeps the test on the current dependency stack while avoiding the one-timepoint MaxIters solution that makes POD fail before DEIM is exercised.

Root Cause

test/deim.jl discretizes a stiff PDE-derived system but used the explicit Tsit5() solver. With current compatible dependencies, Tsit5() exhausts maxiters, MethodOfLines returns a one-timepoint solution, and POD(snapshot, 3) errors because the snapshot has only one time column.

Dependency boundary evidence from local probes:

  • Current resolved stack fails: MethodOfLines 0.11.14, ModelingToolkit 11.31.1, PDEBase 0.1.27, OrdinaryDiffEq 7.1.2.
  • Compatible older stack passes: MethodOfLines 0.11.12, ModelingToolkit 11.23.0, ModelingToolkitBase 1.41.0, ModelingToolkitTearing 1.15.0, PDEBase 0.1.24, with the current ODE/SciMLBase stack.
  • First compatible post-floor stack fails: MethodOfLines 0.11.13, ModelingToolkit 11.27.0, ModelingToolkitBase 1.51.0, ModelingToolkitTearing 1.18.2, PDEBase 0.1.27.

Validation

  • timeout 3600 env GROUP=Core /home/crackauc/.juliaup/bin/julia +release --project=. --startup-file=no -e 'using Pkg; Pkg.test(; coverage=false)' passes.
  • Runic.main(["--check", "."]) with Runic v1.7.0 passes.

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

Co-Authored-By: Chris Rackauckas <accounts@chrisrackauckas.com>
The full-order solve returned InitialFailure on the current MTK v11 /
MethodOfLines stack, so the snapshot had a single time column and
`POD(snapshot, 3)` errored with "Number of modes should be in {1,...,1}".

Root cause is initialization, not the solver. Passing
`missing_guess_value = Constant(0.0)` left every discretized state as a
missing guess, so MethodOfLines' initialization system was structurally
singular and underdetermined (1 equation for 11 unknowns). MTK fell back
to a least-squares initialization that diverged from the all-zero guess
(residual norm ~4.9e4, MaxIters), producing a garbage u0 (values ~1e19)
and an InitialFailure return code. Switching the solver to Rodas5P did
not help because the failure happens during initialization, before any
step is taken.

The boundary conditions fix every state to 0 at t=0 (v(x,0)=w(x,0)=0 and
i₀(0)=0), so the correct initial condition is simply all zeros. Providing
it explicitly makes the initialization well-posed and the solve succeeds,
yielding a 15-column snapshot that exercises POD/DEIM as intended.

Verified on the CI-failing stack (MethodOfLines 0.11.14, MTK 11.31.1,
ModelingToolkitBase 1.51.0, PDEBase 0.1.27): GROUP=Core Pkg.test() passes
(Core/deim.jl 4/4) and Runic --check is clean.

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

Copy link
Copy Markdown
Member Author

Last CI failure diagnosed and fixed (commit 0e4059c)

The Rodas5P() switch alone did not fix CI — Core (julia 1) still failed with:

AssertionError: Number of modes should be in {1,2,...,1}.

Root cause: initialization, not the solver

Reproduced locally on the exact failing stack (MethodOfLines 0.11.14, ModelingToolkit 11.31.1, ModelingToolkitBase 1.51.0, PDEBase 0.1.27). The full-order solve returns InitialFailure before any step is taken, so sol.original_sol has a single time column and POD(snapshot, 3) errors (only 1 mode available).

init(ode_prob, Rodas5P()) returns ReturnCode.InitialFailure directly. Instrumenting the initialization shows why:

  • missing_guess_value = Constant(0.0) leaves every discretized state as a missing guess, so MethodOfLines' initialization system is structurally singular / underdetermined — 1 equation for 11 unknowns.
  • MTK falls back to a least-squares initialization, which diverges from the all-zero guess: MaxIters, residual norm ≈ 4.9e4, giving a garbage u0 (values ≈ 1e19).
  • Because the failure is in initialization, changing the stepping solver (Tsit5 → Rodas5P) can't help.

This is also why LTS passed but julia 1 failed: the older resolved stack produced a well-posed init; the post-floor stack does not.

Fix

The BCs fix every state to 0 at t=0 (v(x,0)=w(x,0)=0, and i₀(0)=0), so the correct initial condition is simply all zeros. Providing it explicitly makes initialization well-posed:

u0 = [u => 0.0 for u in ModelingToolkit.get_unknowns(simp_sys)]
ode_prob = ODEProblem(simp_sys, u0, tspan)

The solve then succeeds with a 15-column snapshot that exercises POD/DEIM as intended.

Validation (on the CI-failing stack)

  • GROUP=Core Pkg.test() passes: Core/deim.jl 4/4, DataReduction 15/15, utils 10/10, "ModelOrderReduction tests passed".
  • Runic.main(["--check", "test/deim.jl"]) (v1.7.0) clean.

Note: docs/src/tutorials/deim_FitzHugh-Nagumo.md still uses the missing_guess_value hack but with N=15, where the docs build currently passes — left unchanged to keep this PR focused. It carries the same latent fragility and could get the same treatment in a follow-up.

Should be ignored until reviewed by @ChrisRackauckas.

@ChrisRackauckas
ChrisRackauckas marked this pull request as ready for review July 12, 2026 05:42
@ChrisRackauckas
ChrisRackauckas merged commit aff82c5 into SciML:main Jul 12, 2026
10 checks passed
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