SciML / StochasticDiffEq.jl

Solvers for stochastic differential equations which connect with the scientific machine learning (SciML) ecosystem
Other
248 stars 66 forks source link

SKenCarp failing for simple non-scalar SDEs #463

Closed spcornelius closed 2 years ago

spcornelius commented 2 years ago

Before I complain, let me first say thanks for your fantastic work! The DifferentialEquations suite is really a killer app for the Julia ecosystem.

After a recent upgrade to StochasticDiffEq@6.44,I'm noticing that SDE systems that previously integrated fine with SKenCarp are now failing with:

Warning: dt <= dtmin. Aborting. There is either an error in your model specification or the true solution is unstable

Two minimal examples:

1) Uncoupled OU processes with small additive noise:

using StochasticDiffEq

σ = 0.05

function f(du, u, p, t)
    @. du = -u
end

function g(du, u, p, t)
    @. du = σ
end

u0 = [1.0, 1.0]
prob = SDEProblem(f, g, u0, (0.0, 10.0))

sol = solve(prob, SKenCarp())

2) This example from ModelingToolkit (but with additive noise):

using ModelingToolkit, StochasticDiffEq

# Define some variables
@parameters σ ρ β
@variables t x(t) y(t) z(t)
D = Differential(t)

eqs = [D(x) ~ σ*(y-x),
       D(y) ~ x*(ρ-z)-y,
       D(z) ~ x*y - β*z]

noiseeqs = [0.1,
            0.1,
            0.1]

@named de = SDESystem(eqs,noiseeqs,t,[x,y,z],[σ,ρ,β])

u0map = [
    x => 1.0,
    y => 0.0,
    z => 0.0
]

parammap = [
    σ => 10.0,
    β => 26.0,
    ρ => 2.33
]

prob = SDEProblem(de,u0map,(0.0,100.0),parammap)
sol = solve(prob,SKenCarp())

Both of these examples work fine as of StochasticDiffEq@6.42 but fail as of v6.44.

I couldn't say what the exact issue is, but it seems related to #450 and the broader change to the linear solver interface across DifferentialEquations.

Note: Scalar (out-of-place) systems like Example 1 here still seem to work.

ChrisRackauckas commented 2 years ago

Yes, it's an issue from that linsolve change, now fixed in #466. Thanks! Good test for it.