1. 程式人生 > >[數學建模(六)]使用MATLAB實現插值

[數學建模(六)]使用MATLAB實現插值

數值 font 多項式 not new pie 插值法 包含 data

常用的插值:拉格朗日多項式插值、牛頓插值、分段線性插值、Hermite 插值和三次樣條插值。

1.拉格朗日插值法

function y=lagrange(x0,y0,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;

  end

  y(i)=s;

end

測試:

>> x0=1:5;

>> y0=x0.^2;

>> x=2.5;

>> y=lagrange(x0,y0,x)

y =

6.2500

>> plot(x0,y0)

>> hold on

>> plot(x,y,‘*‘)

2.牛頓插值

function yi=newton(x,y,xi)

n=length(x);

m=length(y);

if n~=m

    error(The lengths of X ang Y must be equal!);

    return;

end

Y=zeros(n);

Y(:,
1)=y; for k=1:n-1 for i=1:n-k if abs(x(i+k)-x(i))<eps error(the DATA is error!); return; end Y(i,k+1)=(Y(i+1,k)-Y(i,k))/(x(i+k)-x(i)); end end yi=0; for i=1:n z=1; for k=1:i-1 z=z*(xi-x(k)); end yi=yi+Y(1,i)*z; end

測試:

x0=1:5;

y0=x0.^2;

x=2.5;

y=newton(x0,y0,x)

plot(x0,y0)

hold on

plot(x,y,‘*‘)

y =

6.2500

3.分段線性插值

y=interp1(x0,y0,x,‘method‘)

method 指定插值的方法,默認為線性插值。其值可為:

‘nearest‘ 最近項插值

‘linear‘ 線性插值

‘spline‘ 逐段3 次樣條插值

‘cubic‘ 保凹凸性3 次插值。

所有的插值方法要求 x0 是單調的。

測試:

x0=1:5;

y0=x0.^2;

x=2.5;

y=interp1(x0,y0,x,‘linear‘)

plot(x0,y0)

hold on

plot(x,y,‘*‘)

y =

6.5000

4. 埃爾米特(Hermite)插值

function y=hermite(x0,y0,y1,x); %y1是導數值

n=length(x0);m=length(x);

for k=1:m

  yy=0.0;

  for i=1:n

    h=1.0;

    a=0.0;

    for j=1:n

      if j~=i

        h=h*((x(k)-x0(j))/(x0(i)-x0(j)))^2;

        a=1/(x0(i)-x0(j))+a;

      end

    end

    yy=yy+h*((x0(i)-x(k))*(2*a*y0(i)-y1(i))+y0(i));

  end

  y(k)=yy;

end

測試:

>> x0=1:5;

y0=x0.^2;

x=2.5;

y1=[1 1 1 1 1];

y=hermite(x0,y0,y1,x)

plot(x0,y0)

hold on

plot(x,y,‘*‘)

y =

7.4228

5.樣條插值

Matlab 中三次樣條插值也有現成的函數:

(1)y=interp1(x0,y0,x,‘spline‘);

(2)y=spline(x0,y0,x);

(3)pp=csape(x0,y0,conds);

y=ppval(pp,x)

說明:

csape 的返回值是pp 形式,要求插

值點的函數值,必須調用函數ppval。

pp=csape(x0,y0):使用默認的邊界條件,即Lagrange 邊界條件。

pp=csape(x0,y0,conds)中的conds 指定插值的邊界條件,其值可為:

‘complete‘ 邊界為一階導數,即默認的邊界條件

‘not-a-knot‘ 非扭結條件

‘periodic‘ 周期條件

‘second‘ 邊界為二階導數,二階導數的值[0, 0]。

‘variational‘ 設置邊界的二階導數值為[0,0]。

測試:

x0=1:5;

y0=x0.^2;

x=2.5;

y1=interp1(x0,y0,x,‘spline‘)

y2=spline(x0,y0,x)

pp=csape(x0,y0,‘second‘)

y3=ppval(pp,x)

plot(x0,y0)

hold on

plot(x,y1,‘*‘)

plot(x,y2,‘*‘)

plot(x,y3,‘*‘)

y1 =

6.2500

y2 =

6.2500

pp =

包含以下字段的 struct:

form: ‘pp‘

breaks: [1 2 3 4 5]

coefs: [4×4 double]

pieces: 4

order: 4

dim: 1

y3 =

6.2321

6.二維插值

一維插值:節點為一維變量,插值函數是一元函數(曲線)。

二維插值:節點是二維的,插值函數就是二元函數,即曲面。

6.1插值節點為網格節點

(1)z=interp2(x0,y0,z0,x,y,‘method‘)

(2)pp=csape({x0,y0},z0,conds,valconds); %三次樣條插值

z=fnval(pp,{x,y})

6.2 插值節點為散亂節點

技術分享

[數學建模(六)]使用MATLAB實現插值