Skip to content

Commit

Permalink
feat: allow remade initialization systems to be incomplete
Browse files Browse the repository at this point in the history
  • Loading branch information
AayushSabharwal committed Jan 9, 2025
1 parent a976298 commit a932605
Show file tree
Hide file tree
Showing 4 changed files with 44 additions and 3 deletions.
9 changes: 8 additions & 1 deletion src/systems/diffeqs/abstractodesystem.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1292,6 +1292,7 @@ function InitializationProblem{iip, specialize}(sys::AbstractSystem,
fully_determined = nothing,
check_units = true,
use_scc = true,
allow_incomplete = false,
kwargs...) where {iip, specialize}
if !iscomplete(sys)
error("A completed system is required. Call `complete` or `structural_simplify` on the system before creating an `ODEProblem`")
Expand Down Expand Up @@ -1333,7 +1334,13 @@ function InitializationProblem{iip, specialize}(sys::AbstractSystem,
# TODO: throw on uninitialized arrays
filter!(x -> !(x isa Symbolics.Arr), uninit)
if !isempty(uninit)
throw(IncompleteInitializationError(uninit))
allow_incomplete || throw(IncompleteInitializationError(uninit))
# for incomplete initialization, we will add the missing variables as parameters.
# they will be updated by `update_initializeprob!` and `initializeprobmap` will
# use them to construct the new `u0`.
newparams = map(toparam, uninit)
append!(get_ps(isys), newparams)
isys = complete(isys)
end

neqs = length(equations(isys))
Expand Down
9 changes: 8 additions & 1 deletion src/systems/nonlinear/initializesystem.jl
Original file line number Diff line number Diff line change
Expand Up @@ -317,6 +317,12 @@ function SciMLBase.remake_initialization_data(
u0map[dvs[i]] = newu0[i]
end
end
# ensure all unknowns have guesses in case they weren't given one
# and become solvable
for i in eachindex(dvs)
haskey(guesses, dvs[i]) && continue
guesses[dvs[i]] = newu0[i]
end
if p === missing
# the user didn't pass `p` to `remake`, so they want to retain
# existing values. Fill all parameters in `pmap` so that none of
Expand All @@ -341,7 +347,8 @@ function SciMLBase.remake_initialization_data(
op, missing_unknowns, missing_pars = build_operating_point(
u0map, pmap, defs, cmap, dvs, ps)
kws = maybe_build_initialization_problem(
sys, op, u0map, pmap, t0, defs, guesses, missing_unknowns; use_scc, initialization_eqs)
sys, op, u0map, pmap, t0, defs, guesses, missing_unknowns;
use_scc, initialization_eqs, allow_incomplete = true)
return get(kws, :initialization_data, nothing)
end

Expand Down
5 changes: 4 additions & 1 deletion src/systems/problem_utils.jl
Original file line number Diff line number Diff line change
Expand Up @@ -545,7 +545,10 @@ function maybe_build_initialization_problem(
(!is_time_dependent(sys) || t !== nothing)
initializeprob = ModelingToolkit.InitializationProblem(
sys, t, u0map, pmap; guesses, kwargs...)
initializeprobmap = getu(initializeprob, unknowns(sys))

all_init_syms = Set(all_symbols(initializeprob))
solved_unknowns = filter(var -> var in all_init_syms, unknowns(sys))
initializeprobmap = getu(initializeprob, solved_unknowns)

punknowns = [p
for p in all_variable_symbols(initializeprob)
Expand Down
24 changes: 24 additions & 0 deletions test/initializationsystem.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1282,3 +1282,27 @@ end

@test SciMLBase.successful_retcode(solve(newprob))
end

@testset "Issue#3295: Incomplete initialization of pure-ODE systems" begin
@variables X(t) Y(t)
@parameters p d
eqs = [
D(X) ~ p - d * X,
D(Y) ~ p - d * Y
]
@mtkbuild osys = ODESystem(eqs, t)

# Make problem.
u0_vals = [X => 4, Y => 5.0]
tspan = (0.0, 10.0)
p_vals = [p => 1.0, d => 0.1]
oprob = ODEProblem(osys, u0_vals, tspan, p_vals)
integ = init(oprob)
@test integ[X] 4.0
@test integ[Y] 5.0
# Attempt to `remake`.
rp = remake(oprob; u0 = [Y => 7])
integ = init(rp)
@test integ[X] 4.0
@test integ[Y] 7.0
end

0 comments on commit a932605

Please sign in to comment.