From f1563817370335b6462487fbafa2e60aa69fbaa9 Mon Sep 17 00:00:00 2001 From: Herman Sletmoen Date: Sun, 8 Sep 2024 12:03:19 +0200 Subject: [PATCH 1/5] Prioritize filtered_u0 deferral when creating initialization system --- src/systems/nonlinear/initializesystem.jl | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/src/systems/nonlinear/initializesystem.jl b/src/systems/nonlinear/initializesystem.jl index 737beda9c1..ad10520dc7 100644 --- a/src/systems/nonlinear/initializesystem.jl +++ b/src/systems/nonlinear/initializesystem.jl @@ -44,7 +44,10 @@ function generate_initializesystem(sys::ODESystem; y = get(schedule.dummy_sub, x[1], x[1]) y = ModelingToolkit.fixpoint_sub(y, full_diffmap) - if y isa Symbolics.Arr + if y ∈ set_full_states + # defer initialization until defaults are merged below + push!(filtered_u0, y => x[2]) + elseif y isa Symbolics.Arr _y = collect(y) # TODO: Don't scalarize arrays @@ -55,8 +58,6 @@ function generate_initializesystem(sys::ODESystem; # y is a derivative expression expanded # add to the initialization equations push!(eqs_ics, y ~ x[2]) - elseif y ∈ set_full_states - push!(filtered_u0, y => x[2]) else error("Initialization expression $y is currently not supported. If its a higher order derivative expression, then only the dummy derivative expressions are supported.") end From d648bc1da6441ddf7d419cef4114413b629370dd Mon Sep 17 00:00:00 2001 From: Herman Sletmoen Date: Sun, 8 Sep 2024 12:07:57 +0200 Subject: [PATCH 2/5] Test that defaults are overridden when sys is initialized with map from unknowns(sys) --- test/initializationsystem.jl | 19 +++++++++++++++++++ 1 file changed, 19 insertions(+) diff --git a/test/initializationsystem.jl b/test/initializationsystem.jl index ca41bcd072..9bfb50a470 100644 --- a/test/initializationsystem.jl +++ b/test/initializationsystem.jl @@ -486,3 +486,22 @@ sys = extend(sysx, sysy) ssys = structural_simplify(sys) @test_throws ArgumentError ODEProblem(ssys, [x => 3], (0, 1), []) # y should have a guess end + +# https://github.com/SciML/ModelingToolkit.jl/issues/3025 +@testset "Override defaults when setting initial conditions with unknowns(sys) or similar" begin + @variables x(t) y(t) + + # system 1 should solve to x = 1 + ics1 = [x => 1] + sys1 = ODESystem([D(x) ~ 0], t; defaults = ics1, name = :sys1) |> structural_simplify + prob1 = ODEProblem(sys1, [], (0.0, 1.0), []) + sol1 = solve(prob1, Tsit5()) + @test all(sol1[x] .== 1) + + # system 2 should solve to x = y = 2 + sys2 = extend(sys1, ODESystem([D(y) ~ 0], t; initialization_eqs = [y ~ 2], name = :sys2)) |> structural_simplify + ics2 = unknowns(sys1) .=> 2 # should be equivalent to "ics2 = [x => 2]" + prob2 = ODEProblem(sys2, ics2, (0.0, 1.0), []; fully_determined = true) + sol2 = solve(prob2, Tsit5()) + @test all(sol2[x] .== 2) && all(sol2[y] .== 2) +end From aa3f485e293356aa78f26775efd617fd7d3c4649 Mon Sep 17 00:00:00 2001 From: Herman Sletmoen Date: Sun, 8 Sep 2024 12:37:00 +0200 Subject: [PATCH 3/5] Fix small inconsistencies --- src/systems/nonlinear/initializesystem.jl | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/src/systems/nonlinear/initializesystem.jl b/src/systems/nonlinear/initializesystem.jl index ad10520dc7..47119b7cbd 100644 --- a/src/systems/nonlinear/initializesystem.jl +++ b/src/systems/nonlinear/initializesystem.jl @@ -48,13 +48,12 @@ function generate_initializesystem(sys::ODESystem; # defer initialization until defaults are merged below push!(filtered_u0, y => x[2]) elseif y isa Symbolics.Arr + # scalarize array # TODO: don't scalarize arrays _y = collect(y) - - # TODO: Don't scalarize arrays - for i in 1:length(_y) + for i in eachindex(_y) push!(filtered_u0, _y[i] => x[2][i]) end - elseif y isa ModelingToolkit.BasicSymbolic + elseif y isa Symbolics.BasicSymbolic # y is a derivative expression expanded # add to the initialization equations push!(eqs_ics, y ~ x[2]) From a4a5ccfc4abc34b8f3e2d1641d856c376ca3882b Mon Sep 17 00:00:00 2001 From: Herman Sletmoen Date: Sun, 8 Sep 2024 12:59:20 +0200 Subject: [PATCH 4/5] filtered_u0 cannot be a Pair here --- src/systems/nonlinear/initializesystem.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/systems/nonlinear/initializesystem.jl b/src/systems/nonlinear/initializesystem.jl index 47119b7cbd..cc0b8098a1 100644 --- a/src/systems/nonlinear/initializesystem.jl +++ b/src/systems/nonlinear/initializesystem.jl @@ -61,7 +61,7 @@ function generate_initializesystem(sys::ODESystem; error("Initialization expression $y is currently not supported. If its a higher order derivative expression, then only the dummy derivative expressions are supported.") end end - filtered_u0 = filtered_u0 isa Pair ? todict([filtered_u0]) : todict(filtered_u0) + filtered_u0 = todict(filtered_u0) end else dd_guess = Dict() From 937ef537944239f0b203e7d7f41a91df3f49e99a Mon Sep 17 00:00:00 2001 From: Herman Sletmoen Date: Sun, 8 Sep 2024 13:12:08 +0200 Subject: [PATCH 5/5] Format --- test/initializationsystem.jl | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/test/initializationsystem.jl b/test/initializationsystem.jl index 9bfb50a470..482147beba 100644 --- a/test/initializationsystem.jl +++ b/test/initializationsystem.jl @@ -499,7 +499,10 @@ end @test all(sol1[x] .== 1) # system 2 should solve to x = y = 2 - sys2 = extend(sys1, ODESystem([D(y) ~ 0], t; initialization_eqs = [y ~ 2], name = :sys2)) |> structural_simplify + sys2 = extend( + sys1, + ODESystem([D(y) ~ 0], t; initialization_eqs = [y ~ 2], name = :sys2) + ) |> structural_simplify ics2 = unknowns(sys1) .=> 2 # should be equivalent to "ics2 = [x => 2]" prob2 = ODEProblem(sys2, ics2, (0.0, 1.0), []; fully_determined = true) sol2 = solve(prob2, Tsit5())