Exposure of backward error

[NAThompson]

I notice that the solution object given by the solve call has a number of useful utilities:

julia> sol.
alg        dense       destats     errors      interp      k           prob        retcode     t           tslocation  u           u_analytic

However, since the solver comes packaged with an interpolant, would it be sensible to also compute a backward error? This would allow easy comparison of various solution methods.

(Was going to tag this with "feature request" but looks like I need more permissions for that.)

Discussion topic; also am willing to submit a PR. The code is pretty trivial, at least in the scalar case:

julia> using DifferentialEquations
julia> using PlotThemes
julia> using Plots
julia> theme(:solarized)
julia> f(u,p,t) = 1.01*u;
julia> u0 = 1/2;
julia> tspan = (0.0, 1.0);
julia> prob = ODEProblem(f,u0,tspan);
julia> sol1 = solve(prob, Tsit5(), reltol=1e-8, abstol=1e-8)
julia> sol2 = solve(prob, DP5(), reltol=1e-8, abstol=1e-8)
julia> residual1(t) = sol1(t,Val{1}) - sol1.prob.f(sol1(t), sol1.prob.p, t)
julia> residual2(t) = sol2(t,Val{1}) - sol2.prob.f(sol2(t), sol2.prob.p, t)

Obviously if this was added to the solution object, you'd want to do something like expose the residual in some norm.

[ChrisRackauckas]

Is that a standard usage? It is not the way backward error analysis is normally considered since it's the pullback of a perturbation of the solution, and that pullback is going to be a differential equation itself that would need to be solved.

https://publish.uwo.ca/~rmoir2/docs/Reconsidering%20Backward%20Error%20Analysis%20for%20Ordinary%20Differential%20Equations.pdf

This is interesting though, and would be useful to have as a getproperty overload I think. Or maybe just a function that could be used over the solution object.

[NAThompson]

I'm using the language from Corless's book A Graduate Introduction to Numerical Analysis: From a viewpoint of backward error; also it's Definition 1.4 of your link, no? If the usage is nonstandard, I think "residual" would be language that everyone finds acceptable. ("Defect" might also work?)

[ChrisRackauckas]

I thought in that work it's the pullback of errors at the final value? @YingboMa

[NAThompson]

Here's Corless on the topic:

As we do throughout this book, from our a posteriori backward error analysis point of view, we then interpret the residual as a backward error, and henceforth say that our numerical solution gives us the exact solution of a nearby problem. We will look at more sophisticated backward errors later.

The "more sophisticated backward errors" are precisely what you had identified earlier.

[YingboMa]

You are computing the interpolant's residual, not the solution residual. They are not the same.

[NAThompson]

The interpolant is the solution, in Corless's language.