Numerical approximation when exact solutions are impossible
After completing this section, you should be able to:
Many differential equations lack closed-form analytical solutions. When we cannot solve $y' = f(x, y)$ explicitly, we turn to numerical methods that approximate the solution step by step, following the direction field at each point. Euler's method is the simplest and most intuitive approach to this problem.
Given the initial value problem
we want to find approximate values of $y$ at points $x_1, x_2, \ldots$ ahead of $x_0$.
At each current point $(x_n, y_n)$, the slope of the solution curve is $f(x_n, y_n)$. Instead of following the actual curve, Euler's method simply follows the tangent line for a small step of size $h$:
where $x_{n+1} = x_n + h$ and $h$ is the step size (also called mesh width).
Starting from $(x_0, y_0)$, compute successive approximations via: $$x_{n+1} = x_n + h$$ $$y_{n+1} = y_n + h \cdot f(x_n, y_n)$$ for $n = 0, 1, 2, \ldots$ until reaching the target $x$ value.
At each point $(x_n, y_n)$, we draw the tangent line to the solution curve (which has slope $f(x_n, y_n)$). We then move horizontally a distance $h$, and vertically an amount $h \cdot f(x_n, y_n)$ to find the next approximation point. This is like walking along the direction field with a fixed step size.
Local Truncation Error (LTE): The error introduced in a single step, assuming all previous steps are exact. For Euler's method, LTE $= O(h^2)$.
Global Truncation Error (GTE): The total accumulated error over all steps from $x_0$ to $x_f$. For Euler's method, GTE $= O(h)$.
This makes Euler a first-order method. Doubling the number of steps (halving $h$) roughly halves the global error, but at the cost of doing twice the work.
Halving the step size $h$ improves accuracy significantly, but it doubles the computational cost. In practice, we often need to balance accuracy against computational feasibility.
A better strategy is to use a two-stage corrector method. First, we predict using Euler's method, then we correct using the average slope between the current point and the predicted next point.
Predictor step: $$\tilde{y}_{n+1} = y_n + h \cdot f(x_n, y_n)$$
Corrector step: $$y_{n+1} = y_n + \frac{h}{2} \left[ f(x_n, y_n) + f(x_{n+1}, \tilde{y}_{n+1}) \right]$$
The improved Euler method (also called Heun's method) has:
For certain differential equations (especially stiff equations), if the step size $h$ is too large, the numerical solution can become unstable or diverge. Checking that your computed values remain reasonable is an important sanity check. If in doubt, use a smaller step size.
To test stability, compute the solution with step sizes $h$ and $h/2$. If the two results agree to your desired precision, you can trust the smaller-$h$ result. This is the basis of adaptive step-size control in modern codes.
This is the exponential growth equation. The exact solution is $y(x) = e^x$. We'll take 5 steps with $h = 0.1$ to reach $x = 0.5$.
Setup: $f(x, y) = y$, $(x_0, y_0) = (0, 1)$, $h = 0.1$
Computations:
| $n$ | $x_n$ | $y_n$ (Euler) | $y_{\text{exact}}$ | Error |
|---|---|---|---|---|
| 0 | 0.0 | 1.0000 | 1.0000 | 0.0000 |
| 1 | 0.1 | 1.1000 | 1.1052 | 0.0052 |
| 2 | 0.2 | 1.2100 | 1.2214 | 0.0114 |
| 3 | 0.3 | 1.3310 | 1.3499 | 0.0189 |
| 4 | 0.4 | 1.4641 | 1.4918 | 0.0277 |
| 5 | 0.5 | 1.6105 | 1.6487 | 0.0382 |
Step-by-step:
Step 0 to 1: $m_0 = f(0, 1) = 1 \Rightarrow y_1 = 1 + 0.1(1) = 1.1$
Step 1 to 2: $m_1 = f(0.1, 1.1) = 1.1 \Rightarrow y_2 = 1.1 + 0.1(1.1) = 1.21$
Continue similarly for steps 3, 4, and 5.
Result: $y(0.5) \approx 1.6105$ (Euler) vs. $e^{0.5} \approx 1.6487$ (exact). The error grows because we're accumulating local errors.
This is a linear first-order ODE. The exact solution is $y(x) = 2e^x - x - 1$. We'll take 5 steps with $h = 0.2$ from $x = 0$ to $x = 1$.
Setup: $f(x, y) = x + y$, $(x_0, y_0) = (0, 1)$, $h = 0.2$
| $n$ | $x_n$ | $f(x_n, y_n)$ | $y_n$ (Euler) | $y_{\text{exact}}$ | Error |
|---|---|---|---|---|---|
| 0 | 0.0 | 1.0 | 1.0000 | 1.0000 | 0.0000 |
| 1 | 0.2 | 1.2 | 1.2000 | 1.2428 | 0.0428 |
| 2 | 0.4 | 1.6 | 1.4400 | 1.5595 | 0.1195 |
| 3 | 0.6 | 2.0400 | 1.7680 | 1.9739 | 0.2059 |
| 4 | 0.8 | 2.5680 | 2.1816 | 2.5526 | 0.3710 |
| 5 | 1.0 | 3.1816 | 2.7979 | 3.4366 | 0.6387 |
Result: $y(1) \approx 2.7979$ (Euler) vs. $2e - 1 - 1 \approx 3.4366$ (exact). The error is substantial with this step size, showing why finer grids or better methods are needed.
Same problem as Example 1, but now using the improved Euler (Heun's method) to show the dramatic improvement in accuracy.
Setup: $f(x, y) = y$, $(x_0, y_0) = (0, 1)$, $h = 0.1$
| $n$ | $x_n$ | $\tilde{y}_{n+1}$ (Pred.) | $y_n$ (Improved) | $y_{\text{exact}}$ | Error |
|---|---|---|---|---|---|
| 0 | 0.0 | — | 1.0000 | 1.0000 | 0.0000 |
| 1 | 0.1 | 1.1000 | 1.1050 | 1.1052 | 0.0002 |
| 2 | 0.2 | 1.2155 | 1.2214 | 1.2214 | 0.0000 |
| 3 | 0.3 | 1.3436 | 1.3499 | 1.3499 | 0.0000 |
| 4 | 0.4 | 1.4849 | 1.4918 | 1.4918 | 0.0000 |
| 5 | 0.5 | 1.6410 | 1.6487 | 1.6487 | 0.0000 |
Detail for Step 1:
Predictor: $\tilde{y}_1 = 1 + 0.1 \cdot 1 = 1.1$
Corrector: $y_1 = 1 + \frac{0.1}{2}[1 + 1.1] = 1 + 0.05(2.1) = 1.105$
Observation: The improved Euler method is dramatically more accurate! The errors are on the order of $10^{-4}$ or better, compared to $10^{-2}$ for basic Euler. This is the power of second-order methods.
This equation has the exact solution $y(x) = e^{-x^2}$ (a Gaussian curve). We'll use Euler's method to see how it tracks this smooth, decaying function.
Setup: $f(x, y) = -2xy$, $(x_0, y_0) = (0, 1)$, $h = 0.1$
| $n$ | $x_n$ | $m_n = -2x_n y_n$ | $y_n$ (Euler) | $y_{\text{exact}}$ | Error |
|---|---|---|---|---|---|
| 0 | 0.0 | 0.0 | 1.0000 | 1.0000 | 0.0000 |
| 1 | 0.1 | 0.0 | 1.0000 | 0.9900 | 0.0100 |
| 2 | 0.2 | -0.0200 | 0.9980 | 0.9608 | 0.0372 |
| 3 | 0.3 | -0.0599 | 0.9920 | 0.9139 | 0.0781 |
| 4 | 0.4 | -0.0794 | 0.9801 | 0.8521 | 0.1280 |
| 5 | 0.5 | -0.0980 | 0.9603 | 0.7788 | 0.1815 |
Result: The method captures the general decay, but as $x$ increases, the slope becomes steeper (in absolute value), and Euler's method underestimates the decay. Using a smaller step size would help.
To illustrate the convergence of Euler's method as we decrease the step size, we solve the same problem with different values of $h$. The exact answer is $e \approx 2.71828$.
| Step Size $h$ | Number of Steps | $y_{\text{approx}}(1)$ | Error | Error Ratio |
|---|---|---|---|---|
| 0.5 | 2 | 2.2500 | 0.4683 | — |
| 0.25 | 4 | 2.4414 | 0.2769 | 1.69 |
| 0.1 | 10 | 2.5937 | 0.1246 | 2.22 |
| 0.05 | 20 | 2.6533 | 0.0650 | 1.92 |
Observation: As $h$ decreases, the error roughly decreases proportionally to $h$ (the error ratio is close to $h_{\text{old}} / h_{\text{new}}$). This confirms that Euler's method is first-order: error $\propto h$.
For example, from $h = 0.5$ to $h = 0.25$ (halving $h$), the error ratio is approximately $2$. From $h = 0.25$ to $h = 0.1$, the ratio is about $2.2$. This linear convergence is characteristic of first-order methods.
Apply Euler's method to $y' = 2x$, $y(0) = 0$, with $h = 0.5$ for one step. What is $y_1$ (the approximation at $x = 0.5$)?
For the equation $y' = xy$ with $y(0) = 1$ and step size $h = 0.1$, what is the Euler approximation $y_2$ (after two steps)?
Which of the following best describes the global truncation error of Euler's method?
In the improved Euler (Heun's) method, what is the purpose of the corrector step?
If you halve the step size $h$ in Euler's method, approximately how much does the global error improve?
For improved Euler, what order method is it (in terms of convergence)?
$$y_{n+1} = y_n + h \cdot f(x_n, y_n)$$ $$x_{n+1} = x_n + h$$
Order: 1st (error $\propto h$)
$$\tilde{y}_{n+1} = y_n + h f(x_n, y_n)$$ $$y_{n+1} = y_n + \frac{h}{2}[f(x_n, y_n) + f(x_{n+1}, \tilde{y}_{n+1})]$$
Order: 2nd (error $\propto h^2$)
Euler: LTE $O(h^2)$, GTE $O(h)$
Improved: LTE $O(h^3)$, GTE $O(h^2)$
Better accuracy requires smaller $h$ or higher-order methods.
• Smaller $h$ = more accurate but slower
• Test with $h$ and $h/2$ to verify stability
• Improved Euler gives $\approx 4\times$ better accuracy for same work