新冠插值法matlab_多项式函数插值：全域多项式插值（一）单项式基插值、拉格朗日插值、牛顿插值 [MATLAB]...

全域多项式插值指的是在整个插值区域内形成一个多项式函数作为插值函数。关于多项式插值的基本知识，见“计算基本理论”。

在单项式基插值和牛顿插值形成的表达式中，求该表达式在某一点处的值使用的Horner嵌套算法啊，见"Horner嵌套算法"。

1. 单项式(Monomial)基插值

1)插值函数基　　单项式基插值采用的函数基是最简单的单项式：$$\phi_j(t)=t^{j-1}, j=1,2,...n;\quad f(t)=p_{n-1}(t)=x_1+x_2t+x_3t^2+...x_nt^{n-1}=\sum\limits_{j=1}^nx_jt^{j-1}$$　　所要求解的系数即为单项式系数 $x_1,x_2,...x_n$ ，在这里仍然采用1,2,...n的下标记号而不采用和单项式指数对应的0,1,2,...,n-1的下标仅仅是出于和前后讨论一致的需要。

2)叠加系数

单项式基插值采用单项式函数基，若有m个离散数据点需要插值，设使用n项单项式基底：

$$x_1+t_1x_2+t_1^2x_3+...+t_1^{n-1}x_n=y_1\\ x_1+t_2x_2+t_2^2x_3+...+t_2^{n-1}x_n=y_2\\ ......   ......   ......   ......   ......   ......\\ x_1+t_mx_2+t_m^2x_3+...+t_m^{n-1}x_n=y_m$$　　系数矩阵为一 $m\times n$ 的矩阵($m\leq n$)，范德蒙(Vandermonde)矩阵：

$$\begin{bmatrix}1&t_1&t_1^2&...&t_1^{n-1}\\1&t_2&t_2^2&...&t_2^{n-1}\\...&...&...&...&...\\1&t_n&t_n^2&...&t_n^{n-1}\end{bmatrix} \begin{bmatrix}x_1\\x_2\\...\\x_n\end{bmatrix}=\begin{bmatrix}y_1\\y_2\\...\\y_n\end{bmatrix}$$　　根据计算基本理论中的讨论，多项式插值的函数基一定线性无关，且只要离散数据点两两不同，所构成的矩阵行也一定线性无关，这保证了矩阵一定行满秩。此时当且仅当m=n时，叠加系数有且仅有一组解。因此，n=插值基底的个数=离散数据点的个数=多项式次数+1。

3)问题条件和算法复杂度

生成范德蒙矩阵的复杂度大致在 $O(n^2)$ 量级；由于范德蒙矩阵并没有什么好的性质，既没有稀疏性，也没有对称性，因此只能使用标准的稠密矩阵分解(如LU)来解决，复杂度在 $O(n^3)$ 量级。因此，问题的复杂度在 $O(n^3)$ 量级。

范德蒙矩阵存在的问题是，当矩阵的维数越来越高的时候，解线性方程组时问题将越来越病态，条件数越来越大；从另一个角度来说，单项式基底本身趋势相近，次数增大时将越来越趋于平行(见下图)。这将造成随着数据点的增加，确定的叠加系数的不确定度越来越大，因此虽然单项式基很简单，进行插值时却往往不用这一方法。如果仍然采用单项式基底，有时也会进行两种可以改善以上问题的操作：平移(shifting)和缩放(scaling)，即将 $((t-c)/d)^{j-1}$ 作为基底。常见的平移和缩放将所有数据点通过线性变换转移到区间[-1,1]中，即：$c=(t_1+t_n)/2,d=(t_n-t_1)/2$ 。

4)算法实现

使用MATLAB实现单项式插值代码如下：

function [ polyfunc, vecx, condition ] = MonoPI( vect, vecy, shift, scale )

% 计算单项式型插值多项式系数

% 输入四个参数：插值点行向量vect，插值点函数值行向量vecy，平移shift，压限scale；

% 输出两个参数：插值多项式各项系数行向量vecx，矩阵条件数condition；

% 设置缺省值：若只输入两个参数，则不平移不缩放

if nargin==2

shift = 0; scale = 1;

end

% 求解系数

vecsize = size(vect, 2);

basis = (vect - shift * ones(1, vecsize))/scale; % 确定基底在各个数据点的取值向量basis

Mat = vander(basis); condition = cond(Mat); % 用vander命令生成basis的范德蒙矩阵并求条件数

[L, U] = lu(Mat); vecx = (U\(L\vecy.')).'; vecx = fliplr(vecx); % 标准lu分解解矩阵方程

% 生成句柄函数polyfunc

syms t;

monomial = (t - shift)/scale; vecsize = size(vecx, 2); funeval = vecx(vecsize);

for i = vecsize:-1:2 % 生成函数的算法采用Horner算法提高效率

funeval = vecx(i - 1) + monomial*funeval;

end

polyfunc = matlabFunction(funeval, 'Vars', t);

end

比如对