Example 11-5: Global error

Example 11-5: Global error#

Given the many tools in our toolbox for solving ODEs, it is important to understand their various advantages and disadvantages. One way to quantify the accuracy of different algorithms is using the global error of the approximation. Here we will show how to calculate the global error numerically for Euler’s method.

Summary of commands#

In this exercise, we will demonstrate the following:

Global error#

The global error can be calculated using the following formula:

\[ \text{global error} = \sqrt{\frac{\sum_{n=1}^{N} |y^{\mathrm{exact}} (t_n) - y^{\mathrm{numerical}} (t_n) |^2}{N} } \]

While this expression looks complex, it is nothing more than the Euclidean norm (2-norm or magnitude) of the vector that is the difference between \(y^{\mathrm{exact}}\) and \(y^{\mathrm{numerical}}\). Therefore, we can calculate it very quickly by using routine linear algebra methods.

Now, consider the same IVP that we saw previously:

\[ y' = -y + 10 \sin(3t), \quad y(0) = 0, \quad 0 \le t \le 2 \]

The exact solution is \(y = -3 \cos (3t) + \sin (3t) + 3e^{-t}\).

We’ll calculate the global error of forward Euler using a step size of \(h = 0.2\).

# import libraries
import numpy as np
import matplotlib.pyplot as plt

# user-defined function
def myODE(t, y):
    return -y + 10 * np.sin(3 * t)

# Euler's method function
def my_euler(f, t0, tf, y0, h):
    t = np.arange(t0, tf+h, h)
    y = [y0]
    for n in range(len(t) - 1):
        y.append(y[n] + h * f(t[n], y[n]))
    return t, np.array(y)

# initialize
h = 0.2
t0 = 0
tf = 2
y0 = 0

# forward Euler
t, y_euler = my_euler(myODE, t0, tf, y0, h)
y_exact = -3 * np.cos(3 * t) + np.sin(3 * t) + 3 * np.exp(-t)

# compute global error using norm - very fast!
err_global = np.linalg.norm(y_euler - y_exact) / np.sqrt(len(t))
print(f"The global error is εG = {err_global:.5f}")

The global error is εG = 0.59359

Tip

What’s going on in the last line?? That is special syntax for f-strings, or formatted string literals, one of the coolest features of Python. It is like a regular string, but everything inside the curly braces {} is evaluated on-the-fly, allowing you to insert variables and even functions dynamically into strings. Here we insert the global error and ask it to round the decimal (float) to 5 places.

Contents