Can we ignore the error function in numerical solution of ODEs?

Question

In 'Numerical methods for engineers' by Chapra and Canale there was an example in order to prove that Runge-Kutta-Fehlberg method is not always accurate, where given the equation $$\frac{dy}{dx}+0.6y=10e^{-(x-2)^2/[2(0.075)^2]}$$ with initial condition $y(0)=0.5$ the general solution is said to be $$y=0.5e^{-0.6x}$$

That way it is true (as easily seen by the graphs) that the approximate and the true solution are not in close agreement, but why is the error function term ignored by the general solution? If we take in consideration this term then the two solutions are not so different.

Answer 1

Look at a graph of the right side. Or consider that its exponent is $$ -\frac{(x-2)^2}{2\left(\frac3{40}\right)^2}=-\frac{1600(x-2)^2}{18} $$ and $\ln(10^{-17})=-39.14..$ which level the expression will reach for $|x-2|>0.66..$.

Thus from a numerical perspective, the right side does not exist before $x=1.3$ and after $x=2.7$, and on the segment $[0,1.3]$ the solution to the now homogeneous equation is indeed $0.5·e^{-0.6x}$.

---

Using classical RK4 estimating the error by repeating steps with half the step size (or more precisely, for every two steps one step for the same time interval using the double step size) I get the following picture for a relative error threshold of $10^{-4},10^{-5},10^{-6}$.

with a computation log for the $10^{-4}$ integration of

new step size    h=  0.06666667  at t=  0.00000000       with fac= 10.0000
internal step to t=  0.06666667, y[0]=  0.48039472,      abs. error =  1.99859e-14      rel. error =  2.99788e-13       curve length =     0.066667 
new step size    h=  0.66666667  at t=  0.06666667       with fac= 10.0000
internal step to t=  0.73333333, y[0]=  0.32202024,      abs. error =  1.91098e-09      rel. error =  2.60588e-09       curve length =     0.733333 
internal step to t=  1.40000000, y[0]=  0.21585798,      abs. error =  1.21963e-05      rel. error =  8.71164e-06       curve length =     1.400000 
new step size    h=  0.22417033  at t=  1.40000000       with fac=  0.3363
internal step to t=  1.62417033, y[0]=  0.18869247,      abs. error =   0.00176841      rel. error =   0.00108881       curve length =     1.624170 
new step size    h=  1.15353919  at t=  1.62417033       with fac=  5.1458
internal step to t=  2.77770952, y[0]=  0.09701207,      abs. error =    0.0038085      rel. error =   0.00137109       curve length =     2.777710 
internal step to t=  3.93124872, y[0]=  0.04856108,      abs. error =   0.00825655      rel. error =   0.00210023       curve length =     3.931249 
internal step to t=  5.08478791, y[0]=  0.02430809,      abs. error =    0.0178082      rel. error =   0.00350225       curve length =     5.084788 
internal step to t=  6.23832710, y[0]=  0.01216784,      abs. error =    0.0383644      rel. error =   0.00614979       curve length =     6.238327 
internal step to t=  7.39186630, y[0]=  0.00609082,      abs. error =    0.0826262      rel. error =     0.011178       curve length =     7.391866 
internal step to t=  8.54540549, y[0]=  0.00304887,      abs. error =     0.177942      rel. error =    0.0208231       curve length =     8.545405 
internal step to t=  9.69894468, y[0]=  0.00152616,      abs. error =     0.383207      rel. error =    0.0395102       curve length =     9.698945 
new step size    h=  3.92357535  at t=  9.69894468       with fac=  3.4013
internal step to t= 13.62252003, y[0]=  0.00016007,      abs. error =      1.88675      rel. error =     0.138502       curve length =    13.622520 

where one sees that at $x=1.62417033$, directly before the jump, the step size is increased to $dx=1.15353919$, thus neatly stepping over the bump even in the half step evaluations.

But this is more artefact of the implementation and a rather low error threshold. Already tightening the relative error margin to $10^{-6}$ corrects the solution to be as good as the best numerical solution.