(* Stability of ODE's and methods *) 

(* An ODE is said to be stable over a certain interval if two solutions starting out close to each other remain closer over the entire interval *) 

(* An ODE can be stable in one region, and unstable in another *)

(* For example, the one below is unstable if t < 1, and unstable if t > 1 *)

(* Even though a particular ODE {x' == f(t, x)} may be stable, a method applied to it may be unstable. For the 1-st order Euler method, stability condition is | 1 + h*df/dx | > 1 *)

solutionA = NDSolve[{x'[t] == -10(t - 1)x[t], x[1] == 1.0}, x, {t, 0, 2}]

plotA = Plot[Evaluate[x[t] /. solutionA ], {t, 0, 2}, PlotStyle->{Hue[1], Thickness[0.01]}, PlotLabel -> FontForm["x(0)= 1.0", {"Courier", 15}], AxesLabel -> {t , x[t]}]


solutionB = NDSolve[{x'[t] == -10(t - 1)x[t], x[1] == 0.05}, x, {t, 0, 2}]

plotB = Plot[Evaluate[x[t] /. solutionB ], {t, 0, 2}, PlotStyle->{Hue[0.5], Thickness[0.01]}, PlotLabel -> FontForm["x(0) = 0.05", {"Courier", 15}], AxesLabel -> {t, x[t]}]

Show[plotA, plotB, PlotLabel -> FontForm["both solutions", {"Courier", 15}], AxesLabel -> {t, x[t]} ] 


Plot[Evaluate[x[t] /. solutionB ] - Evaluate[x[t] /. solutionA ], {t, 0, 2}, PlotLabel -> FontForm["Difference between them", {"Courier", 15}], AxesLabel -> {t, x[t]}]