3.1: Euler's Method
- Page ID
- 103481
\( \newcommand{\vecs}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)
\( \newcommand{\vecd}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash {#1}}} \)
\( \newcommand{\id}{\mathrm{id}}\) \( \newcommand{\Span}{\mathrm{span}}\)
( \newcommand{\kernel}{\mathrm{null}\,}\) \( \newcommand{\range}{\mathrm{range}\,}\)
\( \newcommand{\RealPart}{\mathrm{Re}}\) \( \newcommand{\ImaginaryPart}{\mathrm{Im}}\)
\( \newcommand{\Argument}{\mathrm{Arg}}\) \( \newcommand{\norm}[1]{\| #1 \|}\)
\( \newcommand{\inner}[2]{\langle #1, #2 \rangle}\)
\( \newcommand{\Span}{\mathrm{span}}\)
\( \newcommand{\id}{\mathrm{id}}\)
\( \newcommand{\Span}{\mathrm{span}}\)
\( \newcommand{\kernel}{\mathrm{null}\,}\)
\( \newcommand{\range}{\mathrm{range}\,}\)
\( \newcommand{\RealPart}{\mathrm{Re}}\)
\( \newcommand{\ImaginaryPart}{\mathrm{Im}}\)
\( \newcommand{\Argument}{\mathrm{Arg}}\)
\( \newcommand{\norm}[1]{\| #1 \|}\)
\( \newcommand{\inner}[2]{\langle #1, #2 \rangle}\)
\( \newcommand{\Span}{\mathrm{span}}\) \( \newcommand{\AA}{\unicode[.8,0]{x212B}}\)
\( \newcommand{\vectorA}[1]{\vec{#1}} % arrow\)
\( \newcommand{\vectorAt}[1]{\vec{\text{#1}}} % arrow\)
\( \newcommand{\vectorB}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)
\( \newcommand{\vectorC}[1]{\textbf{#1}} \)
\( \newcommand{\vectorD}[1]{\overrightarrow{#1}} \)
\( \newcommand{\vectorDt}[1]{\overrightarrow{\text{#1}}} \)
\( \newcommand{\vectE}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash{\mathbf {#1}}}} \)
\( \newcommand{\vecs}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)
\( \newcommand{\vecd}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash {#1}}} \)
\(\newcommand{\avec}{\mathbf a}\) \(\newcommand{\bvec}{\mathbf b}\) \(\newcommand{\cvec}{\mathbf c}\) \(\newcommand{\dvec}{\mathbf d}\) \(\newcommand{\dtil}{\widetilde{\mathbf d}}\) \(\newcommand{\evec}{\mathbf e}\) \(\newcommand{\fvec}{\mathbf f}\) \(\newcommand{\nvec}{\mathbf n}\) \(\newcommand{\pvec}{\mathbf p}\) \(\newcommand{\qvec}{\mathbf q}\) \(\newcommand{\svec}{\mathbf s}\) \(\newcommand{\tvec}{\mathbf t}\) \(\newcommand{\uvec}{\mathbf u}\) \(\newcommand{\vvec}{\mathbf v}\) \(\newcommand{\wvec}{\mathbf w}\) \(\newcommand{\xvec}{\mathbf x}\) \(\newcommand{\yvec}{\mathbf y}\) \(\newcommand{\zvec}{\mathbf z}\) \(\newcommand{\rvec}{\mathbf r}\) \(\newcommand{\mvec}{\mathbf m}\) \(\newcommand{\zerovec}{\mathbf 0}\) \(\newcommand{\onevec}{\mathbf 1}\) \(\newcommand{\real}{\mathbb R}\) \(\newcommand{\twovec}[2]{\left[\begin{array}{r}#1 \\ #2 \end{array}\right]}\) \(\newcommand{\ctwovec}[2]{\left[\begin{array}{c}#1 \\ #2 \end{array}\right]}\) \(\newcommand{\threevec}[3]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \end{array}\right]}\) \(\newcommand{\cthreevec}[3]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \end{array}\right]}\) \(\newcommand{\fourvec}[4]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \\ #4 \end{array}\right]}\) \(\newcommand{\cfourvec}[4]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \\ #4 \end{array}\right]}\) \(\newcommand{\fivevec}[5]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \\ #4 \\ #5 \\ \end{array}\right]}\) \(\newcommand{\cfivevec}[5]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \\ #4 \\ #5 \\ \end{array}\right]}\) \(\newcommand{\mattwo}[4]{\left[\begin{array}{rr}#1 \amp #2 \\ #3 \amp #4 \\ \end{array}\right]}\) \(\newcommand{\laspan}[1]{\text{Span}\{#1\}}\) \(\newcommand{\bcal}{\cal B}\) \(\newcommand{\ccal}{\cal C}\) \(\newcommand{\scal}{\cal S}\) \(\newcommand{\wcal}{\cal W}\) \(\newcommand{\ecal}{\cal E}\) \(\newcommand{\coords}[2]{\left\{#1\right\}_{#2}}\) \(\newcommand{\gray}[1]{\color{gray}{#1}}\) \(\newcommand{\lgray}[1]{\color{lightgray}{#1}}\) \(\newcommand{\rank}{\operatorname{rank}}\) \(\newcommand{\row}{\text{Row}}\) \(\newcommand{\col}{\text{Col}}\) \(\renewcommand{\row}{\text{Row}}\) \(\newcommand{\nul}{\text{Nul}}\) \(\newcommand{\var}{\text{Var}}\) \(\newcommand{\corr}{\text{corr}}\) \(\newcommand{\len}[1]{\left|#1\right|}\) \(\newcommand{\bbar}{\overline{\bvec}}\) \(\newcommand{\bhat}{\widehat{\bvec}}\) \(\newcommand{\bperp}{\bvec^\perp}\) \(\newcommand{\xhat}{\widehat{\xvec}}\) \(\newcommand{\vhat}{\widehat{\vvec}}\) \(\newcommand{\uhat}{\widehat{\uvec}}\) \(\newcommand{\what}{\widehat{\wvec}}\) \(\newcommand{\Sighat}{\widehat{\Sigma}}\) \(\newcommand{\lt}{<}\) \(\newcommand{\gt}{>}\) \(\newcommand{\amp}{&}\) \(\definecolor{fillinmathshade}{gray}{0.9}\)If an initial value problem
\[\label{eq:3.1.1} y'=f(x,y),\quad y(x_0)=y_0 \]
cannot be solved analytically, it is necessary to resort to numerical methods to obtain useful approximations to a solution of Equation \ref{eq:3.1.1}. We will consider such methods in this chapter.
We’re interested in computing approximate values of the solution of Equation \ref{eq:3.1.1} at equally spaced points \(x_0\), \(x_1\), …, \(x_n=b\) in an interval \([x_0,b]\). Thus,
\[x_i=x_0+ih,\quad i=0,1, \dots,n, \nonumber \]
where
\[h={b-x_0\over n}.\nonumber\]
we will denote the approximate values of the solution at these points by \(y_0\), \(y_1\), …, \(y_n\); thus, \(y_i\) is an approximation to \(y(x_i)\). we will call
\[e_i=y(x_i)-y_i \nonumber\]
the error at the \(i\)th step. Because of the initial condition \(y(x_0)=y_0\), we will always have \(e_0=0\). However, in general \(e_i\ne0\) if \(i>0\).
We encounter two sources of error in applying a numerical method to solve an initial value problem:
- The formulas defining the method are based on some sort of approximation. Errors due to the inaccuracy of the approximation are called truncation errors.
- Computers do arithmetic with a fixed number of digits, and therefore make errors in evaluating the formulas defining the numerical methods. Errors due to the computer’s inability to do exact arithmetic are called roundoff errors.
That said, a careful analysis of error is beyond the scope of this book.
Euler’s Method
The simplest numerical method for solving Equation \ref{eq:3.1.1} is Euler’s method. Euler’s method is based on the assumption that the tangent line to the integral curve of Equation \ref{eq:3.1.1} at \((x_i,y(x_i))\) approximates the integral curve over the interval \([x_i,x_{i+1}]\). Since the slope of the integral curve of Equation \ref{eq:3.1.1} at \((x_i,y(x_i))\) is \(y'(x_i)=f(x_i,y(x_i))\), the equation of the tangent line to the integral curve at \((x_i,y(x_i))\) is
\[\label{eq:3.1.2} y=y(x_i)+f(x_i,y(x_i))(x-x_i).\]
Setting \(x=x_{i+1}=x_i+h\) in Equation \ref{eq:3.1.2} yields
\[\label{eq:3.1.3} y_{i+1}=y(x_i)+hf(x_i,y(x_i))\]
as an approximation to \(y(x_{i+1})\). Since \(y(x_0)=y_0\) is known, we can use Equation \ref{eq:3.1.3} with \(i=0\) to compute
\[y_1=y_0+hf(x_0,y_0).\nonumber\]
However, setting \(i=1\) in Equation \ref{eq:3.1.3} yields
\[y_2=y(x_1)+hf(x_1,y(x_1)),\nonumber\]
which isn’t useful, since we don’t know \(y(x_1)\). Therefore we replace \(y(x_1)\) by its approximate value \(y_1\) and redefine
\[y_2=y_1+hf(x_1,y_1).\nonumber\]
Having computed \(y_2\), we can compute
\[y_3=y_2+hf(x_2,y_2).\nonumber\]
In general, Euler’s method starts with the known value \(y(x_0)=y_0\) and computes \(y_1\), \(y_2\), …, \(y_n\) successively by with the formula
\[\label{eq:3.1.4} y_{i+1}=y_i+hf(x_i,y_i),\quad 0\le i\le n-1.\]
The next example illustrates the computational procedure indicated in Euler’s method.
Use Euler’s method with \(h=0.1\) to find approximate values for the solution of the initial value problem
\[\label{eq:3.1.5} y'+2y=x^3e^{-2x},\quad y(0)=1\]
at \(x=0.1,0.2,0.3\).
Solution
We rewrite Equation \ref{eq:3.1.5} as
\[y'=-2y+x^3e^{-2x},\quad y(0)=1, \nonumber\]
which is of the form Equation \ref{eq:3.1.1}, with
\[f(x,y)=-2y+x^3e^{-2x}, \, x_0=0, \, \text{and} \, y_0=1. \nonumber\]
Euler’s method yields
\[\begin{align*} y_1 &= y_0+hf(x_0,y_0) \\ &= 1+(0.1)f(0,1)=1+(0.1)(-2)=0.8,\\[4pt] y_2 & = y_1+hf(x_1,y_1)\\ & = 0.8+(0.1)f(0.1,0.8)=0.8+(0.1)\left(-2(0.8)+(0.1)^3e^{-0.2}\right)= 0.640081873,\\[4pt] y_3 & = y_2+hf(x_2,y_2)\\ & = 0.640081873+(0.1)\left(-2(0.640081873)+(0.2)^3e^{-0.4}\right)= 0.512601754. \end{align*}\]
We’ve written the details of these computations to ensure that you understand the procedure. However, in the rest of the examples as well as the exercises in this chapter, we will assume that you can use a calculator to carry out the necessary computations.
Examples Illustrating The Error in Euler’s Method
Use Euler’s method with step sizes \(h=0.1\), \(h=0.05\), and \(h=0.025\) to find approximate values of the solution of the initial value problem
\[y'+2y=x^3e^{-2x},\quad y(0)=1\nonumber \]
at \(x=0\), \(0.1\), \(0.2\), \(0.3\), …, \(1.0\). Compare these approximate values with the values of the exact solution
\[\label{eq:3.1.6} y={e^{-2x}\over4}(x^4+4),\]
which can be obtained by the method of Section 2.3 (verify this).
Table 3.1.1 shows the values of the exact solution Equation \ref{eq:3.1.6} at the specified points, and the approximate values of the solution at these points obtained by Euler’s method with step sizes \(h=0.1\), \(h=0.05\), and \(h=0.025\). In examining this table, keep in mind that the approximate values in the column corresponding to \(h=0.05\) are actually the results of 20 steps with Euler’s method. We haven’t listed the estimates of the solution obtained for \(x=0.05\), \(0.15\), …, since there’s nothing to compare them with in the column corresponding to \(h=0.1\). Similarly, the approximate values in the column corresponding to \(h=0.025\) are actually the results of 40 steps with Euler’s method.
| \(x\) | \(h=0.1\) | \(h=0.05\) | \(h=0.025\) | Exact |
|---|---|---|---|---|
| 0.0 | 1.000000000 | 1.000000000 | 1.000000000 | 1.000000000 |
| 0.1 | 0.800000000 | 0.810005655 | 0.814518349 | 0.818751221 |
| 0.2 | 0.640081873 | 0.656266437 | 0.663635953 | 0.670588174 |
| 0.3 | 0.512601754 | 0.532290981 | 0.541339495 | 0.549922980 |
| 0.4 | 0.411563195 | 0.432887056 | 0.442774766 | 0.452204669 |
| 0.5 | 0.332126261 | 0.353785015 | 0.363915597 | 0.373627557 |
| 0.6 | 0.270299502 | 0.291404256 | 0.301359885 | 0.310952904 |
| 0.7 | 0.222745397 | 0.242707257 | 0.252202935 | 0.261398947 |
| 0.8 | 0.186654593 | 0.205105754 | 0.213956311 | 0.222570721 |
| 0.9 | 0.159660776 | 0.176396883 | 0.184492463 | 0.192412038 |
| 1.0 | 0.139778910 | 0.154715925 | 0.162003293 | 0.169169104 |
You can see from Table 3.1.1 that decreasing the step size improves the accuracy of Euler’s method. For example,
\[y_{exact}(1)-y_{approx}(1)\approx \left\{\begin{array}{l} 0.0293 \text{with} h=0.1,\\ 0.0144\mbox{ with }h=0.05,\\ 0.0071\mbox{ with }h=0.025. \end{array}\right.\nonumber \]
Based on this scanty evidence, you might guess that the error in approximating the exact solution at a fixed value of \(x\) by Euler’s method is roughly halved when the step size is halved. You can find more evidence to support this conjecture by examining Table 3.1.2 , which lists the approximate values of \(y_{exact}-y_{approx}\) at \(x=0.1\), \(0.2\), …, \(1.0\).
| \(x\) | \(h=0.1\) | \(h=0.05\) | \(h=0.025\) |
|---|---|---|---|
| 0.1 | 0.0187 | 0.0087 | 0.0042 |
| 0.2 | 0.0305 | 0.0143 | 0.0069 |
| 0.3 | 0.0373 | 0.0176 | 0.0085 |
| 0.4 | 0.0406 | 0.0193 | 0.0094 |
| 0.5 | 0.0415 | 0.0198 | 0.0097 |
| 0.6 | 0.0406 | 0.0195 | 0.0095 |
| 0.7 | 0.0386 | 0.0186 | 0.0091 |
| 0.8 | 0.0359 | 0.0174 | 0.0086 |
| 0.9 | 0.0327 | 0.0160 | 0.0079 |
| 1.0 | 0.0293 | 0.0144 | 0.0071 |
Tables 3.1.3 and 3.1.4 show analogous results for the initial value problem
\[\label{eq:3.1.7} y'=-2y^2+xy+x^2,\ y(0)=1,\]
except in this case we cannot solve Equation \ref{eq:3.1.7} exactly. The results in the “Exact” column were obtained by using a more accurate numerical method known as the Runge-Kutta method (discussed later) with a small step size. They are exact to eight decimal places.
| \(x\) | \(h=0.1\) | \(h=0.05\) | \(h=0.025\) | "Exact" |
|---|---|---|---|---|
| 0.0 | 1.000000000 | 1.000000000 | 1.000000000 | 1.000000000 |
| 0.1 | 0.800000000 | 0.821375000 | 0.829977007 | 0.837584494 |
| 0.2 | 0.681000000 | 0.707795377 | 0.719226253 | 0.729641890 |
| 0.3 | 0.605867800 | 0.633776590 | 0.646115227 | 0.657580377 |
| 0.4 | 0.559628676 | 0.587454526 | 0.600045701 | 0.611901791 |
| 0.5 | 0.535376972 | 0.562906169 | 0.575556391 | 0.587575491 |
| 0.6 | 0.529820120 | 0.557143535 | 0.569824171 | 0.581942225 |
| 0.7 | 0.541467455 | 0.568716935 | 0.581435423 | 0.593629526 |
| 0.8 | 0.569732776 | 0.596951988 | 0.609684903 | 0.621907458 |
| 0.9 | 0.614392311 | 0.641457729 | 0.654110862 | 0.666250842 |
| 1.0 | 0.675192037 | 0.701764495 | 0.714151626 | 0.726015790 |
Since we think it is important in evaluating the accuracy of the numerical methods that we will be studying in this chapter, we often include a column listing values of the exact solution of the initial value problem, even if the directions in the example or exercise don’t specifically call for it. If quotation marks are included in the heading, the values were obtained by applying the Runge-Kutta method in a way that’s explained in Section 3.3. If quotation marks are not included, the values were obtained from a known formula for the solution. In either case, the values are exact to eight places to the right of the decimal point.
| \(x\) | \(h=0.1\) | \(h=0.05\) | \(h=0.025\) |
|---|---|---|---|
| 0.1 | 0.0376 | 0.0162 | 0.0076 |
| 0.2 | 0.0486 | 0.0218 | 0.0104 |
| 0.3 | 0.0517 | 0.0238 | 0.0115 |
| 0.4 | 0.0523 | 0.0244 | 0.0119 |
| 0.5 | 0.0522 | 0.0247 | 0.0121 |
| 0.6 | 0.0521 | 0.0248 | 0.0121 |
| 0.7 | 0.0522 | 0.0249 | 0.0122 |
| 0.8 | 0.0522 | 0.0250 | 0.0122 |
| 0.9 | 0.0519 | 0.0248 | 0.0121 |
| 1.0 | 0.0508 | 0.0243 | 0.0119 |
In the next two sections we will study two other numerical methods for solving initial value problems, the improved Euler method and the Runge- Kutta method.