一直想学习数值分析,无奈每次看到复杂的公式,各种证明就头疼,对于一个不是数学专业的人来说,其实能够掌握基本原理会实现代码就够了。最近看到王能超老师编写的一本《计算方法——算法设计及其Matlab》,感觉挺好的。我喜欢这种薄薄的书,能讲到要点,思路也很清楚,所以打算写下一点读书笔记。
另外,重要能够在博客中写latex公式了,为了这再怎么也要写点东西。
首先要讲的是插值算法,我们再日常工作和研究中常常需要做差值运算,例如图像插值,虽然是二维情况,但是跟一维原理是一样的。所以,插值问题可归结为以下情景,即已知一些列离散的点,$({x}_i,{y}_i),i=1,2,...,n$,目的是利用这些已知点求解任意给定的x所对应的y值,这些已知的${x}_i$称为插值节点。
那么,如何充分利用这些已知点的信息呢?一种自然的想法是对插值点x,取与它邻近的已知点${x}_i$的${y}_i$值作为结果,只是精度不高而已。那么进一步向能够利用已知点的线性组合来得到插值结果呢?于是问题就变为如何选择系数${lambda}_i$使得一下组合具有较高的精度,
$y = sum_{i=0}^{n}{{lambda}_iy_i}$
在这里,${lambda}_i$是关于x的函数。
由于$y_i = f(x_i)$,所以,上式求得的y是f(x)的近似,那么插值所要解决的问题就是如何选取${lambda}_i$所得y能够较准确地近似f(x)。另外,一般函数都能够用多项式近似代替,因此,为保证y具有较高的精度,只要令他对于尽可能多的代数多项式f(x)尽可能准确成立。
定义1 如果近似关系$yapprox f(x)$对于幂函数$f = 1,x,x^2,...,x^m$均能准确成立,而对于$f = x^{m+1}$不准确,则可以说该近似关系具有m阶精度。
下面以两点插值公式开始探讨Lagrange插值公式。
一、Lagrange插值公式
1、两点插值
给定两个已知点$(x_0,y_0),(x_1,y_1)$,求形如$y = {lambda}_0y_0+{lambda}_1y_1$的插值公式。
为此,令近似关系$yapprox f(x)$对y = 1,x均准确成立,即近似关系具有一阶精度。
egin{cases}
f(1) = {lambda}_0+{lambda}_1 = 1& \
f(x) = {lambda}_0x_0+{lambda}_1y_0 = x&
end{cases}
由该方程组可解得${lambda}_0,{lambda}_1$,则插值结果为
$y = frac{x-x_1}{x_{0} - x_1}y_0+frac{x-x_0}{x1-x_0}y_1$
2、多点插值
给定n个离散点$({x}_i,{y}_i),i=1,2,...,n$,求插值公式$y = sum_{i=0}^{n}{{lambda}_iy_i}$,由两点插值公式可以看到,n点插值具有n阶精度,即对于$f = 1,x,x^2,...,x^n$均准确成立,可得方程组
其系数是范德蒙特行列式,运用Gramer法则,可求得根为
${lambda}_i = prod_{j=0 \ j eq i}^{n}{frac{x-x_j}{x_i-x_j}},i=0,1,2,...,n$
拉格朗日插值Matlab实现:
%输入插值节点x(数组),及对应的y(数组),xh为待插值点 function lagrange( x,y,xh ) %lagrange 拉格朗日多项式插值 %lagrange的公式为yh = ∑(y(i)∏((xh-xj)/(xi-xj))) (i != j ) n = length(x); m = length(xh); yh = zeros(1,m); c1 = ones(n-1,1); c2 = ones(1,m); for i=1:n xtemp = x([1:i-1 i+1:n]); yh = yh + y(i)*prod((c1*xh - xtemp'*c2)./(x(i)-xtemp'*c2)); end plot(xh,yh,'b'); hold on; plot(x,y,'r'); end