diff --git a/src/systems/problem_utils.jl b/src/systems/problem_utils.jl index 4cba787852..4525b0e46b 100644 --- a/src/systems/problem_utils.jl +++ b/src/systems/problem_utils.jl @@ -938,8 +938,7 @@ function maybe_build_initialization_problem( is_parameter_solvable(p, pmap, defs, guesses) || continue get(op, p, missing) === missing || continue p = unwrap(p) - stype = symtype(p) - op[p] = get_temporary_value(p, floatT) + op[p] = getu(initializeprob, p)(initializeprob) if iscall(p) && operation(p) === getindex arrp = arguments(p)[1] op[arrp] = collect(arrp) @@ -948,7 +947,7 @@ function maybe_build_initialization_problem( if is_time_dependent(sys) for v in missing_unknowns - op[v] = get_temporary_value(v, floatT) + op[v] = getu(initializeprob, v)(initializeprob) end empty!(missing_unknowns) end diff --git a/test/initial_values.jl b/test/initial_values.jl index 79b6b8e067..b3614de0f4 100644 --- a/test/initial_values.jl +++ b/test/initial_values.jl @@ -281,3 +281,31 @@ end @test prob.p isa Vector{Float64} @test length(prob.p) == 5 end + +@testset "Temporary values for solved variables are guesses" begin + @parameters σ ρ β=missing [guess = 8 / 3] + @variables x(t) y(t) z(t) w(t) w2(t) + + eqs = [D(D(x)) ~ σ * (y - x), + D(y) ~ x * (ρ - z) - y, + D(z) ~ x * y - β * z, + w ~ x + y + z + 2 * β, + 0 ~ x^2 + y^2 - w2^2 + ] + + @mtkbuild sys = ODESystem(eqs, t) + + u0 = [D(x) => 2.0, + x => 1.0, + y => 0.0, + z => 0.0] + + p = [σ => 28.0, + ρ => 10.0] + + tspan = (0.0, 100.0) + prob = ODEProblem(sys, u0, tspan, p, jac = true, guesses = [w2 => -1.0], + warn_initialize_determined = false) + @test prob[w2] ≈ -1.0 + @test prob.ps[β] ≈ 8 / 3 +end diff --git a/test/initializationsystem.jl b/test/initializationsystem.jl index 4274f31cc7..1ccc70607a 100644 --- a/test/initializationsystem.jl +++ b/test/initializationsystem.jl @@ -1186,10 +1186,10 @@ end @mtkbuild sys = ODESystem([D(x) ~ x * p + q, x^3 + y^3 ~ 3], t) prob = ODEProblem( sys, [], (0.0, 1.0), [p => 1.0]; guesses = [x => 1.0, y => 1.0, q => 1.0]) - @test prob[x] == 0.0 - @test prob[y] == 0.0 + @test prob[x] == 1.0 + @test prob[y] == 2.0 @test prob.ps[p] == 1.0 - @test prob.ps[q] == 0.0 + @test prob.ps[q] == 3.0 integ = init(prob) @test integ[x] ≈ 1 / cbrt(3) @test integ[y] ≈ 2 / cbrt(3)