# Showing this iterative ODE solver converges quadratically

Given the ODE:

$$y'(t) = f(y(t)), y(0) = y_0,$$

And the following method to solve the ODE:

$$y_{n+1} = y_n +\frac{h}{2}(f(y_n)+f(y_n+hf(y_n))),$$

I am trying to show the method converges quadratically.

I looked at the error at a time $$t_n=hn$$ which I denoted as $$e_n = |y_n-y(t_n)|$$ and then tried to find $$p>0$$ s.t. :

$$\lim_{n\rightarrow{\infty}}\frac{|e_{n+1}|}{|e_n|^p}=\lambda$$

for some $$\lambda>0,$$

but I ran into the problem in the calculation that is: how do I use $$y(t_n)$$, or at least simplify it or manipulate it to make it useful for finding such a $$\lambda$$?

Look for Heun's (2nd order) method, you will find more details. What quadratic or second order convergence means for ODE solutions is that the error of the numerical solution against the exact solution is $$O(h^2)$$, which is different from quadratic convergence for the Newton method.
Or more precisely, let $$y(t)$$ be the exact solution for the IVP with $$y(t_0)=y_0$$ and $$y_n\approx y(t_n)$$ the numerical solution for step size $$h$$. Further $$L$$ the Lipschitz constant of $$f$$. Then there exists some constant $$M\sim \sup_t |y'''(t)|$$ so that the individual local truncation errors are bounded by $$Mh^3$$.
Propagation of previous errors happens at a rate bounded by the Lipschitz constant $$L$$, giving an error propagation inequality of $$e_{n+1}\le(1+Lh)e_n+Mh^3,$$ which solves to $$e_n=|y_n-y(t_n)|\le \frac{e^{L|t_n-t_0|}-1}L\cdot Mh^2,$$ using $$(1+Lh)^n\le e^{Lhn}$$.