Skip to content

Commit 2de171c

Browse files
fix: copy initials to u0 if u0 not provided to remake
1 parent 30bf372 commit 2de171c

File tree

2 files changed

+33
-0
lines changed

2 files changed

+33
-0
lines changed

src/systems/nonlinear/initializesystem.jl

Lines changed: 12 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -593,6 +593,18 @@ function SciMLBase.late_binding_update_u0_p(
593593
prob, sys::AbstractSystem, u0, p, t0, newu0, newp)
594594
supports_initialization(sys) || return newu0, newp
595595
u0 === missing && return newu0, (p === missing ? copy(newp) : newp)
596+
# If the user passes `p` to `remake` but not `u0` and `u0` isn't empty,
597+
# and if the system supports initialization (so it has initial parameters),
598+
# and if the initialization solves for `u0`,
599+
# THEN copy the values of `Initial`s to `newu0`.
600+
if u0 === missing && newu0 !== nothing && p !== missing && supports_initialization(sys) && prob.f.initialization_data !== nothing && prob.f.initialization_data.initializeprobmap !== nothing
601+
if ArrayInterface.ismutable(newu0)
602+
copyto!(newu0, getu(sys, Initial.(unknowns(sys)))(newp))
603+
else
604+
T = StaticArrays.similar_type(newu0)
605+
newu0 = T(getu(sys, Initial.(unknowns(sys)))(newp))
606+
end
607+
end
596608
# non-symbolic u0 updates initials...
597609
if !(eltype(u0) <: Pair)
598610
# if `p` is not provided or is symbolic

test/initializationsystem.jl

Lines changed: 21 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1512,3 +1512,24 @@ end
15121512
@inferred remake(prob; u0 = 2 .* prob.u0, p = prob.p)
15131513
@inferred solve(prob)
15141514
end
1515+
1516+
@testset "Issue#3570: `Initial`s are copied to `u0` if `u0` not provided to `remake`" begin
1517+
@parameters g
1518+
@variables x(t) [state_priority = 10] y(t) λ(t)
1519+
eqs = [D(D(x)) ~ λ * x
1520+
D(D(y)) ~ λ * y - g
1521+
x^2 + y^2 ~ 1]
1522+
@mtkbuild pend = ODESystem(eqs, t)
1523+
1524+
prob = ODEProblem(
1525+
pend, [x => (2 / 2)], (0.0, 1.5), [g => 1], guesses ==> 1, y => 2 / 2])
1526+
sol = solve(prob)
1527+
1528+
setter = setsym_oop(prob, [Initial(x)])
1529+
(u0, p) = setter(prob, [0.8])
1530+
1531+
new_prob = remake(prob; p, initializealg = BrownFullBasicInit())
1532+
@test new_prob[x] 0.8
1533+
new_sol = solve(new_prob)
1534+
@test new_sol[x, 1] 0.8
1535+
end

0 commit comments

Comments
 (0)