cmlabs.interpolate.remainder

cmlabs.interpolate.remainder(xvals, M, x=None, method='auto')[source]

Remainder in interpolation formulas.

Common notation for the remainder value is:

\[\begin{split}\begin{aligned} R_n(x) &= f(x) - L_n(x) = \frac{f^{(n+1)}(\xi)}{(n+1)!} \prod_{i=0}^{n} (x - x_i), \quad \xi \in (a, b) \\ \end{aligned}\end{split}\]

where \(L_n(x)\) is the interpolation polynomial of degree n, \(R_n(x)\) is the remainder term and n = len(xvals)-1.

Parameters:

  • xvals (array_like, 1-D) – The x-coordinates of the data points.

  • M (float) – The bound for the (n+1)-th derivative of the function.

  • x (float, optional, default: None) – The x-coordinate at which to evaluate the polynomial. Only used if method is ‘exact’.

  • mehtod ({'auto', 'exact', 'bound'}, optional, default: 'auto') –

    Defines the method to compute the remainder. The following options are available (default is ‘bound’):

    • ’auto’ : uses ‘exact’ if \(x\) is provided, otherwise uses ‘bound’.

    • ’exact’ : uses the exact formula for the remainder term.

    • ’bound’ : uses the bound formula for the remainder term.

Returns:

res – The value of the remainder term.

Return type:

float

Notes

The method parameter determines how the remainder term is computed. If method is ‘exact’, the function will compute the exact value of the remainder term using the formula:

\[|R_n(x)| \leq \frac{M_{n+1}}{(n+1)!} \prod_{i=0}^{n} |x - x_i|\]

where \(n\) is the length of the xvals array, \(M_{n+1}\) is the bound for the (n+1)-th derivative of the function.

If method is ‘bound’, the function will compute the bound for the remainder term using the formula:

\[|f(x) - L_n(x)| \leq \frac{M_{n+1}(x_{[n]}-x_{[0]})^{n+1}}{2^{2n+1}(n+1)!}\]

where \(x_{[i]}\) is the i-th element of the sorted xvals array

This formula follows from the properties of Chebyshev polynomials.

Examples

>>> import numpy as np
>>> from cmlabs.interpolate import lagrange, remainder
>>> # f(x) = sin(x)
>>> xvals = np.array([0, np.pi/2, np.pi/2])
>>> yvals = np.array([0, 1/2, 1])
>>> # M_3 = max |f'''(x)| = 1
>>> M_3 = 1.0
>>> x = np.pi / 8
>>> abs(np.sin(x) - lagrange(xvals, yvals, x))
0.007941567634910107
>>> # Exact remainder
>>> remainder(xvals, M, x)
0.010093189023535093
>>> # Bound remainder
>>> remainder(xvals, M)
0.020186378047070193