函数逼近——两点三次Hermite(埃尔米特)、分段Hermite插值法 | 北太天元 or Matlab

一、两点三次Hermite插值

已知: 两个插值节点机器对应的导数值 $(x_0,y_0),(x_1,y_1),y_0^{\prime },y_1^{\prime }.$

$$H(x)={\alpha }_0(x)y_0+{\alpha }_1(x)y_1+{\beta }_0(x)y_0^{\prime }+{\beta }_1(x)y_1^{\prime }$$

其中
$$\{\begin{array}{l}{\alpha }_0(x)=\left(1+2\frac{x-x_0}{x_1-x_0}\right){\left(\frac{x-x_1}{x_0-x_1}\right)}^2\\ {\alpha }_1(x)=\left(1+2\frac{x-x_1}{x_0-x_1}\right){\left(\frac{x-x_0}{x_1-x_0}\right)}^2\end{array}\{\begin{array}{l}{\beta }_0(x)=(x-x_0){\left(\frac{x-x_1}{x_0-x_1}\right)}^2\\ {\beta }_1(x)=(x-x_1){\left(\frac{x-x_0}{x_1-x_0}\right)}^2\end{array}$$

二、分段三次Hermite插值

$$H(x_i)=f_i,H^{\prime }(x_i)=f_i^{\prime },i=0,1,\cdots ,n.$$

在每个$[x_i,x_{i-1}]$上利用 $(x_i,y_i),(x_{i-1},y_{i-1}),y_i^{\prime },y_{i-1}^{\prime }$ 进行三次Hermite插值

三、北太天元 or Matlab 实现

两点三次Hermite插值

function [output]= Hermite_interp(x0,y0,dy0,x)
% 两点三次 Hermite 插值
% x0 : [x1 x2]
% y0 : [y1 y2]
% dy0 : [y1' y2']
% length(x)=1
%
%   Version:            1.0
%   last modified:      05/16/2023
    h = x0(2)-x0(1);
    h1 = x - x0(1);
    h2 = x - x0(2);
    X1 = (h2/h)^2;
    X2 = (h1/h)^2;
    alpha1 = (1+2*h1/h)*X1;
    alpha2 = (1+2*h2/(-h))*X2;
    beta1 = h1*X1;
    beta2 = h2*X2;
    output = alpha1*y0(1)+alpha2*y0(2)+beta1*dy0(1)+beta2*dy0(2);
end

将上述代码保存为 Hermite_interp.m 文件.

分段三次Hermite插值

function [out1] = piecewise_Hermite_interp(x0,y0,dy0,x)
    % 分段Hermite插值
    % x0 : [x1 x2 ...]
    % y0 : [y1 y2...]
    % dy0 : [y1' y2'...]
    %
    %   Version:            1.0
    %   last modified:      09/13/2023
    %   file need: Hermite_interp.m
    m = length(x);
    n = length(x0);
    y = zeros(1,m);
    for k =1:1:m
        for i =1:1:n-1
        	   if x0(i)<=x(k) && x(k)<= x0(i+1)
        	   	 y(k) = Hermite_interp(x0([i,i+1]),y0([i,i+1]),dy0([i,i+1]),x(k));
        	   end
        end
    end
    out1 = y;
end

将上述代码保存为 piecewise_Hermite_interp 文件.

这里命名比较长,不过也没什么关系,我比较喜欢把意思表示清楚,另外也有代码的自动补全,也没有怎么不方便

四、简单实现一下

例1

% Hermite_v0 例子1
%   last modified:      04/17/2024
%   file need: Hermite_interp.m   piecewise_Hermite_interp.m
clc;clear all;
%样本点
x0 = [0 pi/2];%x0 = [0 pi];
y0 = sin(x0);
dy0 = cos(x0);

x = 0:pi/20:pi
y = zeros(1,length(x));
for i= 1:1:length(x)
    y(i) = Hermite_interp(x0,y0,dy0,x(i));
end

plot(x,y);
hold on 
    plot(x,sin(x),'r');
    legend('H(x)','sin(x)');
    title('H(x)与sin(x)的图像')
hold off

运行后得到
上面这张图中设置的插值节点是在 $0,\frac{\pi }{2}$，所以在前半部分比较贴合。

例2

%% 分段Hermite_v0 例子2
clc;clear all;
Runge = @(x)1./(x.^2+1);  % 定义Runge函数
dRunge = @(x) -(2*x)./((x.^2+1).^2); % 定义导函数

x0 = linspace(-5,5,11);
y0 = Runge(x0);
dy0 = dRunge(x0);
x = linspace(-5,5,100);
yy = piecewise_Hermite_interp(x0,y0,dy0,x);
plot(x,yy);
hold on
plot(x,Runge(x),'r');
legend('分段Hermite','Runge函数')
hold off

下面这张图为分段三次Hermite插值在Runge函数上的表现，也算是分段上的一种优势。

这一小节内容比较简单，就介绍到这里。至于详细的理论证明，书本上已经写的很好了，我为什么要再抄一遍呢。