diff --git a/Project.toml b/Project.toml index 1c18f31..b7b39d9 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "OrdinaryDiffEqOperatorSplitting" uuid = "760fc936-9fa0-4281-9c1a-468957eedc89" -version = "0.4.4" +version = "0.4.5" authors = ["Dennis Ogiermann and contributors"] [deps] diff --git a/src/integrator.jl b/src/integrator.jl index d4b20e4..8c09afc 100644 --- a/src/integrator.jl +++ b/src/integrator.jl @@ -1146,7 +1146,11 @@ end return current_retcode end -_child_retcode(child::DEIntegrator) = SciMLBase.check_error(child) +# The stored retcodes, for leaves as for nodes: every child was examined through +# `child_failed` right after it stepped, and a failure found there was either retried +# from a clean rollback or escalated into this node's own retcode. Running the leaves' +# `check_error` again here would scan the same state a second time. +_child_retcode(child::DEIntegrator) = child.sol.retcode _child_retcode(child::SplitSubIntegrator) = child.status.retcode function setup_u(prob::OperatorSplittingProblem, solver, alias_u0) diff --git a/src/utils.jl b/src/utils.jl index c547863..c57918f 100644 --- a/src/utils.jl +++ b/src/utils.jl @@ -134,7 +134,7 @@ function forward_sync_subintegrator!( end _reinit_child_controller!(child) end - forward_sync_internal!(parent.u, parent.uprev, child, solution_indices) + forward_sync_internal!(parent.u, child, solution_indices) @timeit_debug "external sync" forward_sync_external!(parent, child, sync) return nothing end @@ -184,20 +184,14 @@ function _suppress_duplicate_fsal_reeval!(child) return nothing end -# Shared internal helper: copy the parent u slice → child DEIntegrator u/uprev. -# `uprev_parent` is the parent's rollback buffer: the refresh of the child's own -# rollback anchor must never write through a child `uprev` that aliases it, because -# that buffer has to survive the whole step untouched (rejection and the palindromic -# mid-step rewind restore from it). When they alias, the slice already holds the -# interval start state and no copy is needed. -function forward_sync_internal!(u_source, uprev_parent, child::DEIntegrator, solution_indices) +# Shared internal helper: copy the parent u slice → child DEIntegrator u. The child's +# rollback anchor `uprev` is the child's own business: `mark_state_modified!` tells it +# that `u` changed from the outside, and a leaf refreshes `uprev` from `u` in the header +# of its next step anyway (OrdinaryDiffEqCore in `on_derivative_discontinuity_at_init!` +# and `apply_step!`). Copying it here as well would write the same state twice. +function forward_sync_internal!(u_source, child::DEIntegrator, solution_indices) @views usrc = u_source[solution_indices] - @timeit_debug "sync vectors" begin - sync_vectors!(child.u, usrc) - if need_sync(child.uprev, uprev_parent) - sync_vectors!(child.uprev, child.u) - end - end + @timeit_debug "sync vectors" sync_vectors!(child.u, usrc) mark_state_modified!(child) return nothing end diff --git a/test/failure_escalation.jl b/test/failure_escalation.jl index 4c50d9e..bac737b 100644 --- a/test/failure_escalation.jl +++ b/test/failure_escalation.jl @@ -110,6 +110,41 @@ PPLTG = PalindromicPairLieTrotterGodunov @test integ.t < 1.0 # aborted mid-solve, nobody could retry end + @testset "A failed step leaves the last accepted state behind" begin + fail_from(du, u, p, t) = t < 0.25 ? odeB(du, u, p, t) : (du .= NaN; nothing) + f = GenericSplitFunction((fA, ODEFunction(fail_from)), (dofs, dofs)) + prob = OperatorSplittingProblem(f, copy(u0), (0.0, 1.0)) + integ = DiffEqBase.init( + prob, LieTrotterGodunov((Euler(), Euler())); dt = 0.1, verbose = false + ) + for _ in 1:3 + DiffEqBase.step!(integ) + end + @test integ.sol.retcode == ReturnCode.Default + ulast, tlast = copy(integ.u), integ.t + DiffEqBase.step!(integ) + @test integ.sol.retcode == ReturnCode.Unstable + @test integ.t == tlast + @test integ.u == ulast + end + + @testset "check_error at the root reads the leaves' verdicts" begin + # The leaves were examined by `child_failed` right after they stepped; the + # root's own check must not scan their state again. + scans = Ref(0) + counting = (dt, u, p, t) -> (scans[] += 1; any(isnan, u)) + f = GenericSplitFunction((fA, ODEFunction(odeB)), (dofs, dofs)) + prob = OperatorSplittingProblem(f, copy(u0), (0.0, 1.0)) + integ = DiffEqBase.init( + prob, LieTrotterGodunov((Euler(), Euler())); dt = 0.1, unstable_check = counting + ) + DiffEqBase.step!(integ) + n = scans[] + @test n > 0 + @test SciMLBase.check_error(integ) == ReturnCode.Default + @test scans[] == n + end + @testset "Rollback clears sticky Newton-style leaf failure flags" begin # SciMLBase's generic check_error re-derives ConvergenceFailure from # `last_stepfail` on non-adaptive leaves even after the retcode is reset, diff --git a/test/sync.jl b/test/sync.jl index 2a8fd4a..b92a95e 100644 --- a/test/sync.jl +++ b/test/sync.jl @@ -82,29 +82,28 @@ end end @testset "forward_sync_internal!" begin - @testset "independent child buffers are synced" begin + @testset "the child state is synced, its rollback anchor is left to the child" begin u_parent = [1.0, 2.0, 3.0] - uprev_parent = [10.0, 20.0, 30.0] idxs = [1, 2] child = MockLeaf([0.0, 0.0], [0.0, 0.0], false) - OS.forward_sync_internal!(u_parent, uprev_parent, child, idxs) + OS.forward_sync_internal!(u_parent, child, idxs) @test child.u == u_parent[idxs] - @test child.uprev == u_parent[idxs] + @test child.uprev == [0.0, 0.0] @test child.u_modified end - @testset "a child uprev aliasing the parent rollback buffer survives" begin - # A child whose uprev is a view into the parent's uprev: the sync must not - # scribble on the parent's rollback anchor mid-step. - u_parent = [1.0, 2.0, 3.0] - uprev_parent = [10.0, 20.0, 30.0] - idxs = [1, 2] - child = MockLeaf([0.0, 0.0], view(uprev_parent, idxs), false) - - OS.forward_sync_internal!(u_parent, uprev_parent, child, idxs) - @test child.u == u_parent[idxs] - @test uprev_parent == [10.0, 20.0, 30.0] # untouched - @test child.u_modified + @testset "a leaf steps from the synced state" begin + # The leaf refreshes its own `uprev` from the synced `u` before it steps: with + # overlapping slices the second operator has to start from what the first left. + ode1(du, u, p, t) = (du .= 1.0; nothing) + ode2(du, u, p, t) = (du .= u; nothing) + u0 = [1.0, 2.0] + f = GenericSplitFunction((ODEFunction(ode1), ODEFunction(ode2)), ([1, 2], [1, 2])) + prob = OperatorSplittingProblem(f, u0, (0.0, 1.0)) + integ = DiffEqBase.init(prob, LieTrotterGodunov((Euler(), Euler())); dt = 0.5) + DiffEqBase.step!(integ) + u1 = u0 .+ 0.5 + @test integ.u ≈ u1 .+ 0.5 .* u1 end end