2.145.0 causes MaxIters in MTK optimization problem

[cstjean]

Updating my packages suddenly had my optimization code take 70X longer than before, and fail with MaxIters. Bisecting the registry designates

  [0bca4576] ↑ SciMLBase v2.144.2 ⇒ v2.145.0

as the culprit.

Before that update (on SciMLBase v2.144.2):

julia> @time res, _, sol = run_optimization(sim);
  0.730308 seconds (10.20 M allocations: 291.045 MiB, 27.10% gc time, 9.77% compilation time)

with a successful return code, while after that update (on SciMLBase v2.145.0),

julia> @time res, _, sol = run_optimization(sim);
 53.272994 seconds (1.42 G allocations: 34.099 GiB, 9.24% gc time)

julia> res.retcode
ReturnCode.MaxIters = 4

julia> res.stats
SciMLBase.OptimizationStats
Number of iterations:                              3000
Time in seconds:                                   53.215000
Number of function evaluations:                    44051
Number of gradient evaluations:                    3002
Number of hessian evaluations:                     0

All timings are for warm JIT.

(worryingly, I also saw return codes of SciMLBase.ReturnCode.Infeasible pop up during the bisection)

run_optimization performs constrained optimization of a fairly vanilla MTK model using Ipopt and AutoForwardDiff, following https://docs.sciml.ai/ModelingToolkit/stable/examples/remake/. I am quite confident that the problem has a single optimum.

Registry bisect log

git bisect start 'eb858dc5b9d' '5f08671057b'
git bisect good 347f4cc60990dfc9dff686ad2f6b7f81e2bea24d
git bisect bad 6600297d1d03fdc3b619056f04547332a24491e5
git bisect good cc83a43779c53abe1172dca7c78ee4e61593767f
git bisect good 8f8b07927f42c56ff47582e0a967b6267b0c4077
git bisect bad ffdb90c3bcefe4137426d90955ee5653163e711d
git bisect good d02eebb97ff81debd95b1a95324c94ee7e759991
git bisect good 64774d9d0862ec28d402fe0e0c574c344525a3e6
git bisect good ac02530260f4005ffdf630b08e00ec4b31b7c13f
git bisect good 812bf925db3446d7d225116f556dc2a5913550c9
git bisect bad dda9eb05c7ced14e680903dbf964b1b3efa2ba31
git bisect good 2a0a9e18e42a524e07fbd6b4c83aa770c210adaa
dda9eb05c7ced14e680903dbf964b1b3efa2ba31 is the first bad commit

[cstjean]

The diff logs look completely unrelated (I don't have Zygote in my Manifest), so I tried again and got the same results as above. Repeated runs (in the same session) consistently yield <1s on v2.144.2, and >50s on v2.145.0. There are no other changes to the Manifest.

This looks like a weird Julia / precompilation bug...

The ODE's sol.alg is the same in both cases, and so is Optimization's solution's alg:

julia> res.alg
OptimizationIpopt.IpoptOptimizer(1.0e-6, 15, 1.0, 0.0001, 0.0001, "no", "no", "mumps", "none", "", "", "yes", "gradient-based", 100.0, "no", "no", "monotone", "quality-function", 0.1, "obj-constr-filter", "no", "exact", 6, "bfgs", "no", "filter", "no", Dict{String, Any}())

Comparing res.stats (the output from solve(::OptimizationProblem)):

julia> res.stats
SciMLBase.OptimizationStats
Number of iterations:                              35
Time in seconds:                                   0.424000
Number of function evaluations:                    253
Number of gradient evaluations:                    38
Number of hessian evaluations:                     0

vs.

julia> res.stats
SciMLBase.OptimizationStats
Number of iterations:                              3000
Time in seconds:                                   52.662000
Number of function evaluations:                    44051
Number of gradient evaluations:                    3002
Number of hessian evaluations:                     0

The optimization objective is 886.0337593381878 on 2.144.2, but slightly higher on 2.145.0, at 886.033767262743. They both find the same solution, to 10 digits (it's simple unidimensional optimization). This is odd; perhaps the termination conditions end up different?

My constraint on the ODE solution has lcons = [10] and ucons=[Inf]. I inserted a print statement in my constraint evaluation function, and on 2.145.0, it has no problem finding a solution:

res = [9.99999991038167]
res = [9.99999991813885]
res = [10.000000001025608]
res = [9.999999894788086]
res = [9.999999897886447]
res = [9.999999899102514]

My optimization code is

    optfn = OptimizationFunction(optim_loss, Optimization.AutoForwardDiff(); cons=cons_fn)

    p = create_optim_p(odeprob, optimized_params, target, constraints, d_discontinuities)
    lcons = Float64.(get_lb.(last.(constraints)))
    ucons = Float64.(get_ub.(last.(constraints)))

    optprob = OptimizationProblem(optfn, initial_values, p; lcons, ucons)
    return solve(optprob, IpoptOptimizer(); maxtime=60)

[cstjean]

Precompiling everything (I'd been running on --compiled=existing) made the problem disappear on 2.45.0, so clearly 2.45.0 is not the true issue, based on how inoffensive of an update it is. I'm quite confident that my initial bisection was correct.

I don't understand what's going on. Something is very wrong with my system (/ Julia 1.12?), and I can't explain it.

For the past few weeks, our Github CI has been consistently failing this optimization test on our package, with ReturnCode.MaxIters (as described above), which prompted me to update my local environment to reproduce the issue.

Today I made a PR upper-bounding SciMLBase to 2.44.2, and it made the precompilation of OrdinaryDiffEq fail with supports_solve_rng not defined in SciMLBase, because its compat bounds for SciMLBase have not been updated.

I don't see how I can reasonably fix my CI without bounding all the SciML dependencies on our package. Any clue is appreciated. Otherwise I'll just pray that future updates somehow magically resolve the CI situation.

[ChrisRackauckas]

Why do you think it's SciMLBase, not something in LinearSolve, Optimization, or OrdinaryDiffEq?

[ChrisRackauckas]

I think it's the LinearSolve residual safety check change?

[cstjean]

The original bisección pointed at scimlbase, but subsequent experiments (and lack of meaningful change in that release) suggest that it's something else. Given that it doesn't look cleanly reproducible, and that my last problem was solved by upgrading Julia, I'm going to wait it out I think. Bisecting takes forever.

[ChrisRackauckas]

Okay because I assumed that SciML/LinearSolve.jl#912 would fix it, so when you said you had it on most recent versions I was confused, and then the fact that it was LinearSolve related 😅 I was like uuhhhh are you sure? Because it sounds like this.

[cstjean]

Locally, I can run the same script repeatedly and get ReturnCode.Infeasible or ReturnCode.Success randomly (or maybe alternatingly?) Below, I got Infeasible, then Success, then Infeasible again.

PS C:\Users\a1058035\Documents\SimDenkai\debug_env> julia --project=. --compiled=existing
               _
   _       _ _(_)_     |  Documentation: https://docs.julialang.org
  (_)     | (_) (_)    |
   _ _   _| |_  __ _   |  Type "?" for help, "]?" for Pkg help.
  | | | | | | |/ _` |  |
  | | |_| | | | (_| |  |  Version 1.12.5 (2026-02-09)
 _/ |\__'_|_|_|\__'_|  |  Official https://julialang.org release
|__/                   |

julia> include("bisect.jl")
┌ Warning: Using arrays or dicts to store parameters of different types can hurt performance.
│ Consider using tuples instead.
└ @ SciMLBase C:\Users\a1058035\.julia\packages\SciMLBase\eTBDr\src\performance_warnings.jl:32

******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
 Ipopt is released as open source code under the Eclipse Public License (EPL).
         For more information visit https://github.com/coin-or/Ipopt
******************************************************************************

  1.085537 seconds (6.53 M allocations: 196.564 MiB, 66.10% gc time, 6.60% compilation time)
Test Failed at C:\Users\a1058035\Documents\SimDenkai\debug_env\bisect.jl:19
  Expression: res.retcode == ReturnCode.Success
   Evaluated: SciMLBase.ReturnCode.Infeasible == SciMLBase.ReturnCode.Success

ERROR: LoadError: There was an error during testing
in expression starting at C:\Users\a1058035\Documents\SimDenkai\debug_env\bisect.jl:19

julia> include("bisect.jl")
  0.400397 seconds (6.53 M allocations: 196.484 MiB, 19.77% gc time)
Test Failed at C:\Users\a1058035\Documents\SimDenkai\debug_env\bisect.jl:19
  Expression: res.retcode == ReturnCode.Success
   Evaluated: SciMLBase.ReturnCode.Infeasible == SciMLBase.ReturnCode.Success

ERROR: LoadError: There was an error during testing
in expression starting at C:\Users\a1058035\Documents\SimDenkai\debug_env\bisect.jl:19

julia>
PS C:\Users\a1058035\Documents\SimDenkai\debug_env> julia --project=. --compiled=existing
               _
   _       _ _(_)_     |  Documentation: https://docs.julialang.org
  (_)     | (_) (_)    |
   _ _   _| |_  __ _   |  Type "?" for help, "]?" for Pkg help.
  | | | | | | |/ _` |  |
  | | |_| | | | (_| |  |  Version 1.12.5 (2026-02-09)
 _/ |\__'_|_|_|\__'_|  |  Official https://julialang.org release
|__/                   |

julia> include("bisect.jl")
┌ Warning: Using arrays or dicts to store parameters of different types can hurt performance.
│ Consider using tuples instead.
└ @ SciMLBase C:\Users\a1058035\.julia\packages\SciMLBase\eTBDr\src\performance_warnings.jl:32

******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
 Ipopt is released as open source code under the Eclipse Public License (EPL).
         For more information visit https://github.com/coin-or/Ipopt
******************************************************************************

  0.568899 seconds (10.08 M allocations: 287.425 MiB, 11.99% gc time, 11.64% compilation time)
Test Passed

julia> include("bisect.jl")
  0.452374 seconds (10.07 M allocations: 287.392 MiB, 6.28% gc time)
Test Passed

julia> include("bisect.jl")
  0.441829 seconds (10.07 M allocations: 287.398 MiB, 6.30% gc time)
Test Passed

julia> include("bisect.jl")
  0.440739 seconds (10.07 M allocations: 287.346 MiB, 6.78% gc time)
Test Passed

julia> include("bisect.jl")
  0.457158 seconds (10.07 M allocations: 287.431 MiB, 6.10% gc time)
Test Passed

PS C:\Users\a1058035\Documents\SimDenkai\debug_env> julia --project=. --compiled=existing
               _
   _       _ _(_)_     |  Documentation: https://docs.julialang.org
  (_)     | (_) (_)    |
   _ _   _| |_  __ _   |  Type "?" for help, "]?" for Pkg help.
  | | | | | | |/ _` |  |
  | | |_| | | | (_| |  |  Version 1.12.5 (2026-02-09)
 _/ |\__'_|_|_|\__'_|  |  Official https://julialang.org release
|__/                   |

julia> include("bisect.jl")
┌ Warning: Using arrays or dicts to store parameters of different types can hurt performance.
│ Consider using tuples instead.
└ @ SciMLBase C:\Users\a1058035\.julia\packages\SciMLBase\eTBDr\src\performance_warnings.jl:32

******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
 Ipopt is released as open source code under the Eclipse Public License (EPL).
         For more information visit https://github.com/coin-or/Ipopt
******************************************************************************

  1.090115 seconds (6.53 M allocations: 196.587 MiB, 66.32% gc time, 6.03% compilation time)
Test Failed at C:\Users\a1058035\Documents\SimDenkai\debug_env\bisect.jl:19
  Expression: res.retcode == ReturnCode.Success
   Evaluated: SciMLBase.ReturnCode.Infeasible == SciMLBase.ReturnCode.Success

ERROR: LoadError: There was an error during testing
in expression starting at C:\Users\a1058035\Documents\SimDenkai\debug_env\bisect.jl:19

julia> include("bisect.jl")
  0.408373 seconds (6.53 M allocations: 196.595 MiB, 17.25% gc time)
Test Failed at C:\Users\a1058035\Documents\SimDenkai\debug_env\bisect.jl:19
  Expression: res.retcode == ReturnCode.Success
   Evaluated: SciMLBase.ReturnCode.Infeasible == SciMLBase.ReturnCode.Success

ERROR: LoadError: There was an error during testing
in expression starting at C:\Users\a1058035\Documents\SimDenkai\debug_env\bisect.jl:19

This is unedited, and there was no code change or anything. My startup script is empty (ergo no Revise usage). I double-checked that I do not use random numbers anywhere (or anything that might be non-deterministic), and I'm careful to use OrderedDict too. No serialization / deserialization, no HTTP connection or anything.

This explains why my bisection from Friday blamed an innocent SciMLBase release.

Is anything at all initialized randomly in the SciML universe?

bisect.jl is very far from a MWE but FWIW it looks like

using Pkg, SimDenkai, Test, IntervalSets, Optimization

sim = define_simulation(...)

run_optimization(sim);  # JIT warm-up to avoid time-out
@time res, _, sol = run_optimization(sim);

@test res.retcode == ReturnCode.Success
@test sol(20, idxs=[sim._electro_room.production_tank.mass_NaOH])[1] ≈ 10 rtol=0.0001

[cstjean]

Further tests in the same terminal for successive Julia sessions yielded Success, Success, Infeasible.

Github CI for some PR on my project passed this morning, for the first time in a few weeks, suggesting that the LinearSolve update from last week really solved at least part of the problem. However, restarting CI on that PR yielded a failure with SciMLBase.ReturnCode.Infeasible. It exhibits the same non-deterministism as my local environment. CI is on Ubuntu.

[cstjean]

I went back to my environment from January 23rd, and reran bisect.jl. There repeated Julia sessions yielded Success, Success, Success, Success, Success. Furthermore, runtime is about 4X faster and allocates 4X less

julia> include("..\\debug_env\\bisect.jl")
  0.093769 seconds (1.08 M allocations: 42.747 MiB)
Test Passed

compared to an environment updated this morning. Looks like there's some commit on the SciML side, between January 23rd and March 9th, that causes this issue.

[oscardssmith]

@ChrisRackauckas I'm pretty sure this is one of yours.

[ChrisRackauckas]

Is there code I can run?

[cstjean]

Sorry, I won't have bandwidth for an MWE in the near future, especially with the non-determinism making it that much more painful. If you have a hunch about which package might be problematic I could test it on my side.

My bisection log had bad commits on February 27th, so presumably the bad change is between January 23rd and February 27th. Though I suppose it's possible that the non-determinism was always there, and only surfaced because of SciML/LinearSolve.jl#911

[oscardssmith]

Would a non-minimal example be easy?

[ChrisRackauckas]

The non-determinism should only be deterministic now?