Convergence analysis

Compute errors and estimate convergence

example_gni_project.convergence_analysis(ivp, scheme, hs, duration, function=<function <lambda>>)

Assess the convergence of the given function/observable of the discrete trajectory generated by the given scheme applied to the given initial value problem.

Parameters:

  • ivp (IVP) – initial value problem

  • scheme (function) – function defining the integration scheme

  • hs (list of floats) – time step sizes The smallest time step is used to compute the reference.

  • duration (float) – duration of simulation The duration needs to be divisible by all time step sizes in hs.

  • function (function) – function of the state As a default, the identity function yields the global truncation error.

Returns:

  • hs (list of floats) – time steps (excluding the reference)

  • errors (list of floats) – corresponding to every time step in hs: norm of the deviation of the function value for the final state of the simulation from the corresponding value for the final state of the reference simulation

  • p (float) – mean of the slopes of all segments in the log-log plot of error over time step

Example

Study the convergence of the pendulum energy computed with the symplectic Euler method.

>>> import autograd.numpy as np
>>> from example_gni_project import *
>>> def pendulum_vector_field(x):
...     q, p = x
...     return np.array([p, -np.sin(q)])
>>> pendulum_ivp = IVP(
...     ("q", "p"),
...     pendulum_vector_field,
...     np.array([1., 0.])
... )
>>> def pendulum_energy(x):
...     q, p = x
...     return 1/2 * p**2 + (1 - np.cos(q))
>>> hs, errors, ps, p = convergence_analysis(
...     pendulum_ivp,
...     symplectic_euler,
...     hs=[0.01, 0.05, 0.125, 0.25],
...     duration=15.25,
...     function=pendulum_energy
... )
>>> plot_errors(hs, errors)
>>> print(f"estimated order of convergence: {round(p, 1)}")
estimated order of convergence: 1.4

Log-log plot

example_gni_project.plot_errors(hs, errors, filename=None)

Create a log-log plot of errors over timesteps. The plot can be used to empirically analyze the convergence of a numerical scheme.

Parameters:

  • hs (list of floats) – time steps (excluding the reference)

  • errors (list of floats) – corresponding to every time step in hs: norm of the deviation of the function value for the final state of the simulation from the corresponding value for the final state of the reference simulation

  • filename (str) – If a filename (path) is provided, the plot is saved.