Lagrange插值MATLAB源程序(14页).doc

-Lagrange插值MATLAB源程序 function y=lagrange(x0,y0,x); % x0自变量取值向量已知 y0为已知对应x0的函数取值，x为要求插值点坐标n=length(x0); m=length(x); for i=1:m z=x(i); s=0.0; for k=1:n p=1.0; for j=1:n if j=k p=p*(z-x0(j)/(x0(k)-x0(j);  %插值基函数 end end s=p*y0(k)+s; %lagrange插值多项式 end y(i)=s;end 测试：x0=0:2 y0=2 3 5 x=0.5Lagrange(x0,y0,x)x0=0:2 y0=2 3 5 x=0:0.01:2Lagrange(x0,y0,x)  Newton插值MATLAB源程序function f=Newton(x,y,x0) syms t; n = length(x); c(1:n) = 0.0; f = y(1); y1 = 0; l = 1; for(i=1:n-1) for(j=i+1:n) y1(j) = (y(j)-y(i)/(x(j)-x(i); end c(i) = y1(i+1); l = l*(t-x(i); f = f + c(i)*l; simplify(f); y = y1; if(i=n-1) if(nargin = 3) f = subs(f,'t',x0); else f = collect(f); end end endtest: x=1 -1 2 y=0 -3 4 x0=-1:0.1:2Newton(x,y,x0)分段插值MATLAB源程序function y = div_linear(x0,y0,x,n) for i = 1:n-1 if (x >= x0(i) && (x <= x0(i+1) y = (x - x0(i+1)/(x0(i) - x0(i+1)*y0(i) + ( x - x0(i)/(x0(i+1) -x0(i)*y0(i+1); else continue; end endtest：（例题3）x0=-1:0.2:1;y0= 1./(25*x0.2+1); y=interp1(x0,y0,x0,'linear') plot(x0,y0,x0,y,'p');Newton插值MATLAB源程序(2)function = newton(x,y,v)x=input(“X数组=”); y=input(“Y数组=”);v=input(“插值点数值=”); n=length(x); t=zeros(n,n);u=0; for i=1:n t(i,1)=y(i); end for j=2:n for i=2:n if i>=j t(i,j)=(t(i,j-1)-t(i-1,j-1)/(x(i)-x(i-j+1); end end end for k=1:n s=1; m=1; for j=1:k if j<k s=s*(v-x(j); end end m=s*t(k,k); u=u+m; end disp(“插值结果=”); disp(u); end>>newton Newton插值MATLAB源程序  function yi=newton(x,y,xi) m=length(x);n=length(y); if m=n error('x and y must same');end f=zeros(n+1,1); k=2; f(1)=y(1) while k=n+1 f(1)=y(k);k,x(k) for i=1:k-1 if i=k-1 f(i+1)=(f(i)-y(i)/(x(k)-x(i); end end cs(i)=f(i+1); y(k)=f(k); k=k+1;end cfwh=0; for i=1:n-2 w=1; for j=1:i w=w*(xi-x(j); end cfwh=cfwh+cs(i)*w; end yi=y(1)+cfwh;x=0:2 y=2 3 5 xi=0.5newton(x,y,xi)Newton2)function f = Newton(x,y,x0) syms t;  if(length(x) = length(y)       n = length(x);  c(1:n) = 0.0;  else       disp('x和y的维数不相等！');      return;  end  f = y(1);  y1 = 0;  l  = 1;  for(i=1:n-1)        for(j=i+1:n)           y1(j) = (y(j)-y(i)/(x(j)-x(i);      end       c(i) = y1(i+1);           l = l*(t-x(i);        f = f + c(i)*l;  simplify(f);      y = y1;      if(i=n-1)           if(nargin = 3)               f = subs(f,'t',x0);          else             f = collect(f);                %将插值多项式展开              f = vpa(f, 6);          end     end endx=0:2 y=2 3 5 x0=0.5Newton(x,y,x0) function  = Newton(x,y,v) x=input('X数组=：'); y=input('Y数组=：'); v=input('插值点数值=：'); n=length(x); t=zeros(n,n); u=0;   for i=1:n       t(i,1)=y(i);   end   for j=2:n     for i=2:n         if i>=j       t(i,j)=(t(i,j-1)-t(i-1,j-1)/(x(i)-x(i-j+1);         end end   end     for k=1:n     s=1;     m=1;      for j=1:k        if j<k        s=s*(v-x(j);        end      end     m=s*t(k,k);     u=u+m;   end disp('插值结果=');disp(u); endx=0:2 y=2 3 5 v=0.5newton(x,y,v) x0=1 -1 2 y0=0 -3 4 x=-1:0.1:2Lagrange(x0,y0,x)第 15 页-