SciML / OrdinaryDiffEq.jl

High performance ordinary differential equation (ODE) and differential-algebraic equation (DAE) solvers, including neural ordinary differential equations (neural ODEs) and scientific machine learning (SciML)
https://diffeq.sciml.ai/latest/
Other
536 stars 205 forks source link

Switch RK methods to use an inf norm for `eigen_est` #2204

Closed oscardssmith closed 4 months ago

oscardssmith commented 4 months ago

fixes https://github.com/SciML/OrdinaryDiffEq.jl/issues/2198.

Originally, something seemed to be wrong since

function exrober(du, u, p, t)
    y₁, y₂, y₃ = u
    k₁, k₂, k₃ = p
    du .= vcat([-k₁ * y₁ + k₃ * y₂ * y₃,
      k₁ * y₁ - k₃ * y₂ * y₃ - k₂ * y₂^2,
      k₂ * y₂^2, ], u[4:end])
end
n = 100
jac_prototype = sparse(I(n+3))
jac_prototype[1:3, 1:3] .= 1.0
prob_ex_rober = ODEProblem(ODEFunction(exrober; jac_prototype), vcat([1.0,0.0,0.0], ones(n)),(0.0,100.0),(0.04,3e7,1e4))

julia> solve(prob_ex_rober, AutoTsit5(FBDF())).stats
SciMLBase.DEStats
Number of function 1 evaluations:                  4553
Number of function 2 evaluations:                  0
Number of W matrix evaluations:                    842
Number of linear solves:                           3278
Number of Jacobians created:                       3
Number of nonlinear solver iterations:             3258
Number of nonlinear solver convergence failures:   4
Number of fixed-point solver iterations:                     0
Number of fixed-point solver convergence failures:           0
Number of rootfind condition calls:                0
Number of accepted steps:                          859
Number of rejected steps:                          277
Maximum eigenvalue recorded:                       3623401289626744889992244603839849986785280

julia> solve(prob_ex_rober, AutoTsit5(FBDF(;autodiff=false))).stats
SciMLBase.DEStats
Number of function 1 evaluations:                  767728
Number of function 2 evaluations:                  0
Number of W matrix evaluations:                    179601
Number of linear solves:                           597222
Number of Jacobians created:                       33647
Number of nonlinear solver iterations:             498502
Number of nonlinear solver convergence failures:   111620
Number of fixed-point solver iterations:                     0
Number of fixed-point solver convergence failures:           0
Number of rootfind condition calls:                0
Number of accepted steps:                          34616
Number of rejected steps:                          240
Maximum eigenvalue recorded:                       3623401289626744889992244603839849986785280

but further investigation shows that this is expected given the magnitude of u.

ChrisRackauckas commented 4 months ago

Rebase?

oscardssmith commented 4 months ago

rebased. This isn't quite ready to merge yet due to the issues with autodiff=false on AutoSwitch solvers for this problem (I suspect that the jacobian isn't getting refreshed correctly for some reason, but haven't figure out why yet).

oscardssmith commented 4 months ago

This now is actually working. The problem I was seeing before is that du = u with a long timespan leads to very large u values that make finite differencing understandably mad at you.