复旦大学：《Matlab Math》（双语版）CHAPTER 3 插值多项式

Matlab math 多项式和样条 Cleve morler著 一陈文斌(wbchen@fudan.edu.cn) 复旦大学2002

Matlab Math Cleve Morler著 陈文斌（wbchen@fudan.edu.cn) 复旦大学2002 多项式和样条

插值多项式 在平面上给定n个点(x)可以唯一确定一个最多n 1次的多项式通过这些点,这个多项式叫插值多项式 P(xk)=yk,k=1,2, Lagrange插值形式 X-X P(x)=∑I k(j≠Ek X-X

在平面上给定n个点(xk ,yk ),可以唯一确定一个最多n- 1次的多项式通过这些点，这个多项式叫插值多项式 插值多项式 P(xk ) = yk , k = 1,2,…,n Lagrange插值形式           − − = k  k j k k j j y x x x x P(x)

插值多项式例子 MATLAB X=0:3; y=[-5-6-116]; disp(lx; y D P(x) (x-1)(x-2)(x-3 (-5)+x(x-2)x-3) (-6)+ (-6) x(x-1)(x-3 (-1)+ (x-1)(x (16) (6) P(x)=x3-2x-5

插值多项式例子 x =0:3; y = [-5 –6 –1 16]; disp([x; y]) (16) (6) ( 1)( 2) ( 1) ( 2) ( 1)( 3) ( 6) (2) ( 2)( 3) ( 5) ( 6) ( 1)( 2)( 3) ( ) − − − + − − − − + − − − + − − − − = x x x x x x x x x x x x P x ( ) 2 5 3 P x = x − x −

Monomial基 P(x)=Cx"+C2x+.+Cnx+c VI 2 y2 2 Vandermonde矩阵 h-J MATLAB x=03;y=[-56-116,P(x)=x3-2x-5 V= vander(x) Vandermonde矩阵是非 奇异的,但条件数是非 常坏的

Monomial基 n n n n P x = c x + c x + + c x + c − − − 1 2 2 1 1 ( ) ...               =                             − − − − − − n n n n n n n n n n n y y y c c c x x x x x x x x x           2 1 2 1 1 2 2 2 2 1 2 1 2 1 1 1 1 1 1 Vandermonde矩阵 n j k j k v x − , = x =0:3; y = [-5 –6 –1 16]; V = vander(x) c = V\y’ ( ) 2 5 3 P x = x − x − Vandermonde矩阵是非 奇异的，但条件数是非 常坏的

Polyinterp (Lagrange插值形式) MATLAB function v= polyinterp(x, y, u) n=length(x for k= 1: n w=ones(size(u)) forj=[1: k-1 k+l: n w=(u-O)./(x(k-x).*w er +w*y(k); end

Polyinterp(Lagrange插值形式) function v = polyinterp(x,y,u) n = length(x); for k = 1:n w = ones(size(u)); for j = [1:k-1 k+1:n] w = (u – x(j))./(x(k) – x(j)).*w; end v = v + w*y(k); end

matlab u=- 25:01: 3.25,v=polyinterp(x,y, u) plot(, y,o', u,v,-) 0.5

-0.5 0 0.5 1 1.5 2 2.5 3 3.5 -10 -5 0 5 10 15 20 25 u = -.25:.01:3.25; v = polyinterp(x,y,u); plot(x,y,’o’,u,v,’-’)

Polyinterp(符号运算) MATLAB smx =sym(X) P= polyinterp(x, y, symx) pretty(P) P= simplify(p) X^3-2*x-5

Polyinterp(符号运算) symx =sym(‘x’) P = polyinterp(x,y,symx) pretty(P) P = simplify(P) P = x^3-2*x-5

Polyinterpe(另外的例子) MATLAB x=1:6;y=[161821171512]; u=.75: 05: 6.25; v=polyinterp(x, y, u) plot(x,y,o', u,v,-) 086

0 1 2 3 4 5 6 7 6 8 10 12 14 16 18 20 22 Polyinterp(另外的例子) x = 1:6; y = [16 18 21 17 15 12]; u = .75:.05:6.25; v = polyinterp(x,y,u); plot(x,y,’o’,u,v,’-’);

分片线性插值 MATLAB x=1:6;y=[161821171512]; plot(xy,o, u,v, - k+1 n2 L(x)=Dk+(x-xxi+l- XK

0 1 2 3 4 5 6 7 10 12 14 16 18 20 22 分片线性插值 x = 1:6; y = [16 18 21 17 15 12]; plot(x,y,’o’,u,v,’-’); k k k k k k x x y y L x y x x − − = + − + + 1 1 ( ) ( )

MATLAB function v, Sigma= piecelin(x y, u) d=diff(y). /diff(x); First divided difference %o Find subinterval indices, x(k)=x()=j end %o Evaluate interpolant s=u-x(k) y(k)+S. d(k)

function [v,sigma] = piecelin(x,y,u) d = diff(y)./diff(x); % First divided difference % Find subinterval indices, x(k) = x(j)) = j; end % Evaluate interpolant s = u - x(k); v = y(k) + s.*d(k);