Numerical analysis 1–Lagrangian interpolation method

Directory

  • Principle introduction
  • Code
  • Results display

Principle introduction

1. Definition

2. Is Ln(x) unique that satisfies the interpolation conditions?
Assuming that Ln(x) is a polynomial of degree n, then the following n + 1 equations can be obtained by bringing in the conditions of each interpolation point:

We want to prove the uniqueness of Ln(x), which is to prove the uniqueness of a0 to an. 1, x0, x0^2… can be regarded as a coefficient matrix, and its coefficient determinant is as follows:

Obviously this is a Vandermond determinant. According to its calculation method, the solution of this determinant is:

And because x is selected from different points on the coordinate system curve, if i != j, then x_i != x_j, then the determinant is not 0.
Finally, the only condition for the solution of the non-homogeneous equation system AX = B is R(A)=R(A,B)=n, which is obviously satisfied, therefore:

3. Find Ln(x) using the undetermined coefficient method
As shown in the figure below, the only thing to note is that the question gives three zero points and two derivative zero points, so there are five conditions in total. We know that the five conditions can solve the five-variable function, that is, a0-a4. At this time we The polynomial set should be a fourth-order polynomial, and the other content is primary school mathematics.

4. Lagrangian polynomials
Definition:
Description: Treat n + 1 y_i as each coefficient of the polynomial.

Example: I wrote the point-slope form of L1(x), then converted it into a two-point form, and successfully converted it into the above form.

When n >= 1:

example:

5. Lagrangian polynomial interpolation remainder

It is defined as the interpolation remainder, which is the truncation error.
Using Rolle’s theorem and its extension, we get the following definition.

6. The formula of Lagrangian general term is as follows:

7.quiz

Given six interpolation points, we can know from the general expression of l(x) that l(x) is a fifth-order polynomial curve, as shown in A.

Code implementation

The definition of Lagrangian class:

import numpy as np
import sympy as sym
import matplotlib.pyplot as plt
from Lagrangeinterpolation.utils import interp_utils


class LagrangeInterpolation:
    def __init__(self, x, y):
        self.x = np.asarray(x, dtype=float) # as_array is the view, array is the copy
        self.y = np.asarray(y, dtype=float)
        if len(self.x) > 1 and len(self.x) == len(self.y): # For robustness judgment, you must know what you want to judge, and there must be no errors in data transmission.
            self.n = len(self.x)
        else:
            raise ValueError("Interpolation data (x, y) dimensions do not match")
        self.polynomial = None # final polynomial
        self.poly_coefficient = None #Polynomial coefficient
        self.coefficient_order = None # Corresponding order of polynomial coefficients
        self.y0 = None #The value of the interpolation point sought

    def fit_interp(self):
        """
        Core algorithm: Generating Lagrange polynomials
        :return:
        """
        x = sym.Symbol('x')
        self.polynomial = 0.0 # Interpolation polynomial instantiation
        for i in range(self.n): # Pass in 10 interpolation points, 0-9, len(x) is 10, so traversing all interpolation points is range(self.n)
            basic_fuc = self.y[i]
            for j in range(i):
                basic_fuc *= (x - self.x[j]) / (self.x[i] - self.x[j])
            for j in range(i + 1, self.n):
                basic_fuc *= (x - self.x[j]) / (self.x[i] - self.x[j])
            self.polynomial + = basic_fuc
        self.polynomial = sym.expand(self.polynomial)
        polynomial = sym.Poly(self.polynomial, x)
        self.poly_coefficient = polynomial.coeffs()
        self.coefficient_order = polynomial.monoms()

    def cal_interp_x0(self, x0):
        self.y0 = interp_utils.cal_interp_x0(self.polynomial, x0)
        return self.y0

    def plt_interpolation(self, x0, y0):
        params = (self.polynomial, self.x, self.y, 'Lagrange', x0, y0)
        interp_utils.plt_interpolation(params)

Tool functions:

import numpy as np
import matplotlib.pyplot as plt


def cal_interp_x0(polynomial, x0):
    """
    Calculate the value of a given interpolation point, that is, interpolation
    :param polynomial:
    :param x0: interpolated abscissa
    :return:
    """
    x0 = np.asarray(x0, dtype=np.float64)
    n0 = len(x0) #The number of interpolation points required
    y_0 = np.zeros(n0) # Store the interpolation value corresponding to the interpolation point
    t = polynomial.free_symbols.pop()
    for i in range(n0):
        y_0[i] = polynomial.evalf(subs={<!-- -->t: x0[i]})
    return y_0


def plt_interpolation(params):
    polynomial, x, y, title, x0, y0 = params
    plt.figure(figsize=(8, 6))
    plt.plot(x, y, "yo", label="Interpolation base point")
    xi = np.linspace(min(x), max(x), 100)
    yi = cal_interp_x0(polynomial, xi)
    plt.plot(xi, yi, "r--", label="Interpolation polynomial")
    if x0 is not None and y0 is not None:
        plt.plot(x0, y0, "g*", label="Interpolation point values")
    plt.legend()
    plt.xlabel("x", fontdict={<!-- -->"fontsize": 12})
    plt.ylabel("y", fontdict={<!-- -->"fontsize": 12})
    plt.title(title + "interpolation polynomial and values")
    plt.grid(ls=":")
    plt.show()

main function

from Lagrangeinterpolation.lagrange import LagrangeInterpolation
import numpy as np

x = np.linspace(0, 2 * np.pi, 10, endpoint=True)
y = np.sin(x)
x0 = np.array([1.2, 1.3, 1.5, 2, 3, 4, 5, 6])
lag_interp = LagrangeInterpolation(x=x, y=y)
lag_interp.fit_interp()
print("The Lagrange polynomial is as follows:")
print(lag_interp.polynomial)
print("Lagrangian interpolation polynomial coefficient vector and corresponding order:")
print(lag_interp.poly_coefficient)
print(lag_interp.coefficient_order)
y0 = lag_interp.cal_interp_x0(x0)
lag_interp.plt_interpolation(x0, y0)

Result display



Excerpted from: (1) https://www.bilibili.com/video/BV1Lq4y1U7Hj/?spm_id_from=333.337.search-card.all.click
(2) https://www.bilibili.com/video/BV1AK4y1k7Px/?spm_id_from=333.337.search-card.all.click & amp;vd_source=35775f5151031f11ee2a799855c8e368