KROGH_INTERPOLATE

Overview

The KROGH_INTERPOLATE function performs polynomial interpolation using the Krogh algorithm, which constructs an interpolating polynomial that passes exactly through all provided data points. This method is particularly useful when both function values and derivatives need to be matched at interpolation nodes, making it a form of Hermite interpolation.

The algorithm is based on Fred Krogh’s 1970 paper “Efficient Algorithms for Polynomial Interpolation and Numerical Differentiation,” which describes an efficient method for computing interpolating polynomials and their derivatives. Unlike Lagrange or Newton interpolation, the Krogh method allows specifying derivative values by repeating x-coordinates in the input data. For example, providing the same x-value twice indicates that both the function value and first derivative are specified at that point.

Given n+1 data points (x_0, y_0), (x_1, y_1), \ldots, (x_n, y_n) with distinct x-coordinates, the function constructs the unique polynomial p(x) of degree at most n satisfying:

p(x_i) = y_i \quad \text{for } i = 0, 1, \ldots, n

This implementation uses the SciPy library’s krogh_interpolate function. For more details, see the SciPy interpolation documentation and the underlying KroghInterpolator class.

The function also supports computing derivatives of the interpolating polynomial at specified points through the der parameter. Setting der=0 returns interpolated values, while der=1 returns first derivatives, and so on.

Note on numerical stability: Polynomial interpolation can become numerically unstable for large numbers of points or poorly chosen node distributions. The Runge phenomenon demonstrates that equidistant nodes may produce oscillations near interval boundaries. For improved stability, consider using Chebyshev nodes or spline-based methods for high-degree interpolation.

This example function is provided as-is without any representation of accuracy.

Excel Usage

=KROGH_INTERPOLATE(xi, yi, x, der)
  • xi (list[list], required): X-coordinates of the data points (column vector)
  • yi (list[list], required): Y-coordinates of the data points (column vector)
  • x (list[list], required): X-coordinates at which to evaluate the interpolated values
  • der (int, optional, default: 0): Order of derivative to compute (0 for interpolated values)

Returns (list[list]): A 2D list (column vector) of interpolated values or derivatives, or an error message (str) if invalid.

Examples

Example 1: Quadratic polynomial interpolation

Inputs:

xi yi x
0 0 0.5
1 1 1.5
2 4

Excel formula:

=KROGH_INTERPOLATE({0;1;2}, {0;1;4}, {0.5;1.5})

Expected output:

Result
0.25
2.25

Example 2: Polynomial interpolation single point

Inputs:

xi yi x
0 1 0.5
1 2
2 5

Excel formula:

=KROGH_INTERPOLATE({0;1;2}, {1;2;5}, {0.5})

Expected output:

Result
1.25

Example 3: First derivative computation

Inputs:

xi yi x der
0 0 1 1
1 1
2 4

Excel formula:

=KROGH_INTERPOLATE({0;1;2}, {0;1;4}, {1}, 1)

Expected output:

Result
2

Example 4: Evaluate at multiple points

Inputs:

xi yi x der
1 1 1.5 0
2 4 2.5
3 9 3.5
4 16

Excel formula:

=KROGH_INTERPOLATE({1;2;3;4}, {1;4;9;16}, {1.5;2.5;3.5}, 0)

Expected output:

Result
2.25
6.25
12.25

Python Code

import math
from scipy.interpolate import krogh_interpolate as scipy_krogh_interpolate

def krogh_interpolate(xi, yi, x, der=0):
    """
    Krogh polynomial interpolation.

    See: https://docs.scipy.org/doc/scipy/reference/generated/scipy.interpolate.krogh_interpolate.html

    This example function is provided as-is without any representation of accuracy.

    Args:
        xi (list[list]): X-coordinates of the data points (column vector)
        yi (list[list]): Y-coordinates of the data points (column vector)
        x (list[list]): X-coordinates at which to evaluate the interpolated values
        der (int, optional): Order of derivative to compute (0 for interpolated values) Default is 0.

    Returns:
        list[list]: A 2D list (column vector) of interpolated values or derivatives, or an error message (str) if invalid.
    """
    def to2d(val):
        return [[val]] if not isinstance(val, list) else val

    def flatten(arr):
        return [item for sublist in arr for item in sublist]

    try:
        # Normalize inputs to 2D lists
        xi = to2d(xi)
        yi = to2d(yi)
        x = to2d(x)

        # Validate that inputs are 2D lists
        if not isinstance(xi, list) or not all(isinstance(row, list) for row in xi):
            return "Error: xi must be a 2D list."
        if not isinstance(yi, list) or not all(isinstance(row, list) for row in yi):
            return "Error: yi must be a 2D list."
        if not isinstance(x, list) or not all(isinstance(row, list) for row in x):
            return "Error: x must be a 2D list."

        # Flatten the 2D lists
        xi_flat = flatten(xi)
        yi_flat = flatten(yi)
        x_flat = flatten(x)

        # Validate that all elements are numeric
        try:
            xi_flat = [float(val) for val in xi_flat]
            yi_flat = [float(val) for val in yi_flat]
            x_flat = [float(val) for val in x_flat]
        except (ValueError, TypeError):
            return "Error: all elements must be numeric."

        # Check for non-finite values
        for val in xi_flat + yi_flat + x_flat:
            if math.isnan(val) or math.isinf(val):
                return "Error: all values must be finite."

        # Validate der parameter
        if not isinstance(der, (int, float)):
            return "Error: der must be an integer."
        der_int = int(der)
        if der_int < 0:
            return "Error: der must be non-negative."

        # Check that xi and yi have the same length
        if len(xi_flat) != len(yi_flat):
            return "Error: xi and yi must have the same length."

        # Check minimum number of points
        if len(xi_flat) < 1:
            return "Error: at least one data point is required."

        # Check for duplicate xi values
        if len(xi_flat) != len(set(xi_flat)):
            return "Error: xi values must be unique."

        # Perform interpolation
        result = scipy_krogh_interpolate(xi_flat, yi_flat, x_flat, der=der_int)

        # Convert result to 2D list (column vector)
        if hasattr(result, '__iter__'):
            result_list = [[float(val)] for val in result]
        else:
            result_list = [[float(result)]]

        # Check for non-finite values in result
        for row in result_list:
            for val in row:
                if math.isnan(val) or math.isinf(val):
                    return "Error: result contains non-finite values."

        return result_list

    except Exception as exc:
        return f"Error: {exc}"

Online Calculator