From 8c0b98e491a123ed21a7aeb307f5e8b0eb0d9632 Mon Sep 17 00:00:00 2001 From: Keno Fischer Date: Wed, 25 Oct 2023 04:40:43 +0000 Subject: [PATCH 1/2] Correct initial saved value for `dures` We were saving `u` twice rather than `du` at the start. While we're at it, also defer the choice of the interpolator to the SciMLBase code (which does the same thing currently, but let's have this in one place and give the SciMLBase constructor access to the `du` vector if we have it). --- src/common_interface/solve.jl | 8 ++------ 1 file changed, 2 insertions(+), 6 deletions(-) diff --git a/src/common_interface/solve.jl b/src/common_interface/solve.jl index 33434b4..c178640 100644 --- a/src/common_interface/solve.jl +++ b/src/common_interface/solve.jl @@ -1254,12 +1254,8 @@ function DiffEqBase.__init(prob::DiffEqBase.AbstractDAEProblem{uType, duType, tu sol = DiffEqBase.build_solution(prob, alg, ts, - ures; + ures, dense ? dures : nothing; dense = dense, - interp = dense ? - DiffEqBase.HermiteInterpolation(ts, ures, - dures) : - DiffEqBase.LinearInterpolation(ts, ures), calculate_error = false, timeseries_errors = timeseries_errors, retcode = retcode, @@ -1323,7 +1319,7 @@ function DiffEqBase.__init(prob::DiffEqBase.AbstractDAEProblem{uType, duType, tu if save_start save_value!(ures, integrator.u, uType, save_idxs) - save_value!(dures, integrator.u, uType, save_idxs) + save_value!(dures, integrator.du, duType, save_idxs) end initialize_callbacks!(integrator) From 248081fefe8ef99e99cc22903ce0fe04abf46db8 Mon Sep 17 00:00:00 2001 From: Keno Fischer Date: Thu, 26 Oct 2023 18:01:08 -0400 Subject: [PATCH 2/2] Add testcase for initial du save --- test/common_interface/ida.jl | 10 ++++++++++ 1 file changed, 10 insertions(+) diff --git a/test/common_interface/ida.jl b/test/common_interface/ida.jl index 91057fd..3a8dcb0 100644 --- a/test/common_interface/ida.jl +++ b/test/common_interface/ida.jl @@ -121,3 +121,13 @@ f_noconverge(out, du, u, p, t) = out .= [du[1] + u[1] / (t - 1)] prob = DAEProblem(f_noconverge, [1.0], [1.0], (0, 2); differential_vars = [true]) sol = solve(prob, IDA()) @test !(sol.retcode in (ReturnCode.Success, ReturnCode.MaxIters)) + +# Test that we're saving the correct initial data for du +function f_inital_data(du, u, p, t) + return [du[1] - (u[1] + 10.0)] +end +prob = DAEProblem(f_inital_data, [0.0], [1.0], (0.0, 1.0); differential_vars = [true]) +sol = solve(prob, IDA()) +# If this is one, it incorrectly saved u, if it is 0., it incorrectly solved +# the pre-init value rather than the post-init one. +@test sol(0.0, Val{1})[1] ≈ 11.0