Gill's Method

A formula for numerical solution of differential equations,

$\begin{displaymath} y_{n+1}=y_n+{\textstyle{1\over 6}}[k_1+(2-\sqrt{2}\,)k_2+(2+\sqrt{2}\,)k_3+k_4]+{\mathcal O}(h^5), \end{displaymath}$

where

$\displaystyle k_1$	$\textstyle =$	$\displaystyle hf(x_n,y_n)$
$\displaystyle k_2$	$\textstyle =$	$\displaystyle hf(x_n+{\textstyle{1\over 2}}h, y_n+{\textstyle{1\over 2}}k_1)$
$\displaystyle k_3$	$\textstyle =$	$\displaystyle hf[x_n+{\textstyle{1\over 2}}h, y_n+{\textstyle{1\over 2}}(-1+\sqrt{2}\,)k_1+(1-{\textstyle{1\over 2}}\sqrt{2}\,)k_2]$
$\displaystyle k_4$	$\textstyle =$	$\displaystyle hf[x_n+h, y_n-{\textstyle{1\over 2}}\sqrt{2}\,k_2+(1+{\textstyle{1\over 2}}\sqrt{2}\,)k_3].$

References