Skip to main content

Section 19.2 Euler's method

To solve an ODE numerically, we “discretize” it: instead of variable \(t\) varying continuously over some interval, we choose step size \(h\) and work with \(t_0\text{,}\) \(t_1=t_0+h\text{,}\) \(t_2=t_1+h\text{,}\) and so on. (An adaptive method will change the step size \(h\) based on the outcome of computation; more on this later.) We try to find approximate values of \(y\) at the points \(t_k\text{,}\) starting with \(y(t_0)= y_0\) which is known (initial condition). The idea of Euler's method is geometric: draw the line segment starting at \((t_0,y_0)\) with the slope \(f(t_0,y_0)\text{.}\) When it reaches the abscissa \(t_1\text{,}\) we have new point \((t_1,y_1)\text{.}\) Then process is repeated. In a formula, it is expressed as

\begin{equation} y_{k+1} = y_k + h f(t_k, y_k)\label{eq-euler-method}\tag{19.2.1} \end{equation}

Implement Euler's method to solve the differential equation \(y'=-ty\) with initial condition \(y(0)=1\text{.}\) Plot the result on the interval \([0, 3]\) and compare it with the exact solution, which is \(y(t) = \exp(-t^2/2)\text{.}\) Try different step sizes: 0.3, 0.1, 0.01.

Repeat the above for \(y'= ty\text{,}\) where the exact solution is \(y(t) = \exp(t^2)\text{.}\)

Answer
f = @(t, y) -t*y;   % or without -
h = 0.3;   %  or 0.1, 0.01 
t = 0:h:3;
y0 = 1;
y = y0*ones(size(t));

for k = 1:numel(t)-1
    y(k+1) = y(k) + h*f(t(k), y(k));
end

exact = exp(-t.^2 /2);   % or without -
plot(t, y, t, exact)
legend('Euler', 'exact')

As Example 19.2.1 shows, a weakness of Euler's method is its tendency to fall behind the solution when the solution is increasing and concave up, or decreasing and concave down. The reason is that it uses the rate of change from point \(t_k\) for the entire interval \([t_k, t_{k+1}]\text{.}\) When applied to an ODE of the form \(y'=f(t)\text{,}\) Euler's method becomes the left endpoint method of integration. So it is not very accurate.