Volterra Integral Equations

First Kind

This package provides the function inteq.SolveVolterra() which approximates the solution to the Volterra Integral Equation of the first kind using the method in Betto and Thomas (2021).

inteq.SolveVolterra(k: Callable[[numpy.ndarray, numpy.ndarray], numpy.ndarray], f: Callable[[numpy.ndarray], numpy.ndarray] = <function <lambda>>, a: float = 0.0, b: float = 1.0, num: int = 1000, method: str = 'midpoint') → numpy.ndarray

Approximate the solution, g(x), to the Volterra Integral Equation of the first kind:

\[f(s) = \int_a^s K(s,y) g(y) dy\]

using the method in Betto and Thomas (2021).

Parameters
kfunction

The kernel function that takes two arguments.

ffunction

The left hand side (free) function with f(a) = 0.

afloat

Lower bound of the integral, defaults to 0.

bfloat

Upper bound of the estimate, defaults to 1.

numint

Number of estimation points between zero and b.

methodstring

Use either the ‘midpoint’ (default) or ‘trapezoid’ rule.

Returns
grid2-D array

Input values are in the first row and output values are in the second row.

Example

example of first kind vie

from inteq import SolveVolterra
import numpy as np
import matplotlib.pyplot as plt

# define kernel
def k(s, t):
    return np.cos(s - t)


# define true solution
def trueg(s):
    return (2 + s ** 2) / 2


# look decent with just 8 grid points but we'll use 25
s, g = SolveVolterra(k, num=25, method="trapezoid")

# plot functions
figure = plt.figure()

plt.plot(s, g)
plt.plot(s, trueg(s))

plt.legend(["Estimate", "True Value"])

plt.xlabel("s")
plt.ylabel("g(s)")

figure.set_dpi(100)
figure.set_size_inches(8, 5)

figure.savefig("..\\docs\\volterra\\volterra-example.svg", bbox_inches="tight")

Trapezoid vs Midpoint Rule

Volterra integral equations are typically solved using the midpoint rule. However, the trapezoid rule often converges faster. See below an example of the trapezoid rule performing well with just six grid points.

example of trapezoid rule converging faster

Thus, the trapezoid rule typically performs better. However, the trapezoid rule is less stable than the midpoint rule. An example where this this instability is an issue is provided below.

example of trapezoid rule having issues

This can be remedied by smoothing the function. For example, with inteq.helpers.smooth().

Second Kind

Implementation forthcoming.