Skip to content

Commit

Permalink
fix: propagate initializeprobpmap and update_initializeprob! to `…
Browse files Browse the repository at this point in the history
…DAEFunction`
  • Loading branch information
AayushSabharwal committed Nov 27, 2024
1 parent f3caa10 commit 5b702ca
Show file tree
Hide file tree
Showing 2 changed files with 22 additions and 1 deletion.
6 changes: 5 additions & 1 deletion src/systems/diffeqs/abstractodesystem.jl
Original file line number Diff line number Diff line change
Expand Up @@ -498,6 +498,8 @@ function DiffEqBase.DAEFunction{iip}(sys::AbstractODESystem, dvs = unknowns(sys)
checkbounds = false,
initializeprob = nothing,
initializeprobmap = nothing,
initializeprobpmap = nothing,
update_initializeprob! = nothing,
kwargs...) where {iip}
if !iscomplete(sys)
error("A completed system is required. Call `complete` or `structural_simplify` on the system before creating a `DAEFunction`")
Expand Down Expand Up @@ -551,7 +553,9 @@ function DiffEqBase.DAEFunction{iip}(sys::AbstractODESystem, dvs = unknowns(sys)
jac_prototype = jac_prototype,
observed = observedfun,
initializeprob = initializeprob,
initializeprobmap = initializeprobmap)
initializeprobmap = initializeprobmap,
initializeprobpmap = initializeprobpmap,
update_initializeprob! = update_initializeprob!)
end

function DiffEqBase.DDEFunction(sys::AbstractODESystem, args...; kwargs...)
Expand Down
17 changes: 17 additions & 0 deletions test/initializationsystem.jl
Original file line number Diff line number Diff line change
Expand Up @@ -958,3 +958,20 @@ end
@test_warn ["structurally singular", "initialization", "Guess", "heuristic"] ODEProblem(
pend, [x => 1, y => 0], (0.0, 1.5), [g => 1], guesses ==> 1])
end

@testset "DAEProblem initialization" begin
@variables x(t) [guess = 1.0] y(t) [guess = 1.0]
@parameters p=missing [guess = 1.0] q=missing [guess = 1.0]
@mtkbuild sys = ODESystem(
[D(x) ~ p * y + q * t, x^3 + y^3 ~ 5], t; initialization_eqs = [p^2 + q^3 ~ 3])

# FIXME: solve for du0
prob = DAEProblem(
sys, [D(x) => cbrt(4), D(y) => -1 / cbrt(4)], [x => 1.0], (0.0, 1.0), [p => 1.0])

integ = init(prob, DImplicitEuler())
@test integ[x] 1.0
@test integ[y]cbrt(4) rtol=1e-6
@test integ.ps[p] 1.0
@test integ.ps[q]cbrt(2) rtol=1e-6
end

0 comments on commit 5b702ca

Please sign in to comment.