EULERS METHOD

Consider the problem of calculating the shape of an unknown curve which starts at a given point and satisfies a given differential equation. Here, a differential equation can be thought of as a formula by which the slope of the tangent line to the curve can be computed at any point on the curve, once the position of that point has been calculated.

The idea is that while the curve is initially unknown, its starting point, which we denote by  A_{0}, is known (see the picture on top right). Then, from the differential equation, the slope to the curve at  A_{0} can be computed, and so, the tangent line.

Take a small step along that tangent line up to a point  A_{1}. Along this small step, the slope does not change too much, so  A_{1} will be close to the curve. If we pretend that  A_{1} is still on the curve, the same reasoning as for the point  A_{0} above can be used. After several steps, a polygonal curve  A_{0}A_{1}A_{2}A_{3}\dots  is computed. In general, this curve does not diverge too far from the original unknown curve, and the error between the two curves can be made small if the step size is small enough and the interval of computation is finite:

 {\displaystyle y'(t)=f(t,y(t)),\qquad y(t_{0})=y_{0}.}

Choose a value h for the size of every step and set    t_{n}=t_{0}+nh. Now, one step of the Euler method from

t_{n} to  t_{n+1}=t_{n}+h

y_{n+1}=y_{n}+hf(t_{n},y_{n}).

The value of  y_{n} is an approximation of the solution to the ODE at time  t_{n} y_{n}\approx y(t_{n}). The Euler method is  explicit i.e. the solution  y_{n+1} is an explicit function of y_{i} for  i\leq n.

While the Euler method integrates a first-order ODE, any ODE of order N can be represented as a first-order ODE: to treat the equation

 y^{(N)}(t)=f(t,y(t),y'(t),\ldots ,y^{(N-1)}(t)),

we introduce auxiliary variables  z_{1}(t)=y(t),z_{2}(t)=y'(t),\ldots ,z_{N}(t)=y^{(N-1)}(t) and obtain the equivalent equation

 \mathbf {z} '(t)={\begin{pmatrix}z_{1}'(t)\\\vdots \\z_{N-1}'(t)\\z_{N}'(t)\end{pmatrix}}={\begin{pmatrix}y'(t)\\\vdots \\y^{(N-1)}(t)\\y^{(N)}(t)\end{pmatrix}}={\begin{pmatrix}z_{2}(t)\\\vdots \\z_{N}(t)\\f(t,z_{1}(t),\ldots ,z_{N}(t))\end{pmatrix}}

This is a first-order system in the variable  \mathbf {z} (t) and can be handled by Euler’s method or, in fact, by any other scheme for first-order systems

Example 

Given the initial value problem

 y'=y,\quad y(0)=1,

we would like to use the Euler method to approximate  y(4) 

Using step size equal to 1 (h = 1) 

Illustration of numerical integration for the equation  y'=y,y(0)=1.Blue is the Euler method; green, the  midpoint method; red, the exact solution, y=e^{t}. The step size is h = 1.0.

The Euler method is

 y_{n+1}=y_{n}+hf(t_{n},y_{n}).\qquad \qquad

so first we must compute  f(t_{0},y_{0}). In this simple differential equation, the function f is defined by  f(t,y)=y. We have

 f(t_{0},y_{0})=f(0,1)=1.\qquad \qquad

By doing the above step, we have found the slope of the line that is tangent to the solution curve at the point  (0,1). Recall that the slope is defined as the change in  y divided by the change in  t, or  \Delta y/\Delta t.

The next step is to multiply the above value by the step size  h, which we take equal to one here:

 h\cdot f(y_{0})=1\cdot 1=1.\qquad \qquad

Since the step size is the change in t, when we multiply the step size and the slope of the tangent, we get a change in  y value. This value is then added to the initial  y value to obtain the next value to be used for computations.

 y_{0}+hf(y_{0})=y_{1}=1+1\cdot 1=2.\qquad \qquad

The above steps should be repeated to find  y_{2},  y_{3} and  y_{4}.

 {\begin{aligned}y_{2}&=y_{1}+hf(y_{1})=2+1\cdot 2=4,\\y_{3}&=y_{2}+hf(y_{2})=4+1\cdot 4=8,\\y_{4}&=y_{3}+hf(y_{3})=8+1\cdot 8=16.\end{aligned}}

Due to the repetitive nature of this algorithm, it can be helpful to organize computations in a chart form, as seen below, to avoid making errors.

 n  y_{n}  t_{n}  f(t_{n},y_{n})  h  \Delta y  y_{n+1}
0 1 0 1 1 1 2
1 2 1 2 1 2 4
2 4 2 4 1 4 8
3 8 3 8 1 8 16

The conclusion of this computation is that  y_{4}=16. The exact solution of the differential equation is y(t)=e^{t}, so y(4)=e^{4}\approx 54.598. Thus, the approximation of the Euler method is not very good in this case. However, as the figure shows, its behaviour is qualitatively right.

Using other step sizes 

The same illustration for h = 0.25.

As suggested in the introduction, the Euler method is more accurate if the step size h is smaller. The table below shows the result with different step sizes. The top row corresponds to the example in the previous section, and the second row is illustrated in the figure.

step size result of Euler’s method error
1 16 38.598
0.25 35.53 19.07
0.1 45.26 9.34
0.05 49.56 5.04
0.025 51.98 2.62
0.0125 53.26 1.34

The error recorded in the last column of the table is the difference between the exact solution at  t=4 and the Euler approximation. In the bottom of the table, the step size is half the step size in the previous row, and the error is also approximately half the error in the previous row.

Other methods, such as the  mid point method also illustrated in the figures, behave more favourably : the error of the midpoint method is roughly proportional to the square of the step size. For this reason, the Euler method is said to be a first-order method, while the midpoint method is second order.

We can extrapolate from the above table that the step size needed to get an answer that is correct to three decimal places is approximately 0.00001, meaning that we need 400,000 steps. This large number of steps entails a high computational cost. For this reason, people usually employ alternative, higher-order methods such as  runge kutta method or linear multi step methods , especially if a high accuracy is desired

Derivation 

The Euler method can be derived in a number of ways. Firstly, there is the geometrical description mentioned above.

Another possibility is to consider the taylor expansion of the function y around  t_{0}:

 y(t_{0}+h)=y(t_{0})+hy'(t_{0})+{\frac {1}{2}}h^{2}y''(t_{0})+O(h^{3}).

The differential equation states that  y'=f(t,y). If this is substituted in the Taylor expansion and the quadratic and higher-order terms are ignored, the Euler method arises  The Taylor expansion is used below to analyze the error committed by the Euler method, and it can be extended to produce runge-kutta method

A closely related derivation is to substitute the forward finite differences formula for the derivative,

 y'(t_{0})\approx {\frac {y(t_{0}+h)-y(t_{0})}{h}}

in the differential equation  y'=f(t,y). Again, this yields the Euler method.  A similar computation leads to the midpoint rule and the backward euler method

Finally, one can integrate the differential equation from  t_{0} to  t_{0}+h and apply the  fundamental theorem of calcus to get:

 y(t_{0}+h)-y(t_{0})=\int _{t_{0}}^{t_{0}+h}f(t,y(t))\,\mathrm {d} t.

Now approximate the integral by the left-hand  rectangle method (with only one rectangle):

 \int _{t_{0}}^{t_{0}+h}f(t,y(t))\,\mathrm {d} t\approx hf(t_{0},y(t_{0})).