用Matlab实现了3次样条曲线插值的算法边界条件取为自然.docVIP

用Matlab实现了3次样条曲线插值的算法。边界条件取为自然边界条件，即：两个端点处的2阶导数等于0； 共包含3各个函数文件，主函数所在文件（即使用的时候直接调用的函数）为spline3.m,另外两个函数文件是在splin3函数文件中被调用的自定义函数。一个是GetParam.m,一个是GetM.m。 %GetParam.m文件的内容： %根据给定的离散点的横坐标所构成的向量，计算各个区间段的h值； function GetParam(Vx,Vy) global gh; global gf; global gu; global gr; global gd; global gff; global gM; %global gn; %n=length(Vx);%length()为向量Vx所含元素的个数； %n=legth(Vx); %gn=n; %n=gn; n=length(Vx); gh(1)=Vx(2)-Vx(1); gf(1)=(Vy(2)-Vy(1))/gh(1); for i=2:1:n-1%从区间0到区间n-1; gh(i)=Vx(i+1)-Vx(i); gf(i)=(Vy(i+1)-Vy(i))/gh(i); gu(i)=gh(i-1)/(gh(i-1)+gh(i)); gr(i)=1-gu(i); gff(i)=(gf(i-1)-gf(i))/(Vx(i-1)-Vx(i+1)); gd(i)=6*gff(i); end %设置与边界条件有关的参数； gM(1)=0;%起点的2阶导数； gM(n)=0;%终点的2阶导数； end %GetM.m文件的内容： function GetM(Vx) global gh; global gf; global gu; global gr; global gd; global gff; global gM; %global gn; nn=length(Vx); %nn=gn; n=nn-2; b=zeros(n,1); A=zeros(n,n); A(1,1)=2;A(1,2)=gr(2); b(1)=gd(2)-gu(2)*gM(1); for i=2:1:n-1 A(i,i)=2; A(i,i-1)=gu(i+1); A(i,i+1)=gr(i+1); b(i)=gd(i+1); end A(n,n-1)=gu(n);A(n,n)=2; b(n)=gd(nn-1)-gr(nn-1)*gM(nn); X=(inv(A))*b; for i=2:1:nn-1 gM(i)=X(i-1); end end %主函数文件spline3.m的内容： function result=spline3(x,Vx,Vy) global gh; global gf; global gu; global gr; global gd; global gff; global gM; %global gn; GetParam(Vx,Vy); GetM(Vx); %n=length(Vx); %n=gn; n=length(Vx); nn=length(x); y=zeros(1,nn); for j=1:1:nn i=1; while(x(j)Vx(i+1)) i=i+1; end sn=i; t1=(Vx(sn+1)-x(j))^3/(6*gh(sn)); t1=t1*gM(sn); t2=(x(j)-Vx(sn))^3/(6*gh(sn)); t2=t2*gM(sn+1); t3=Vy(sn)-gM(i)*((gh(i))^2)/6; t3=t3*(Vx(sn+1)-x(j))/gh(sn); t4=Vy(sn+1)-gM(sn+1)*((gh(sn))^2)/6; t4=t4*(x(j)-Vx(sn))/gh(sn); y(j)=t1+t2+t3+t4; end result=y; end 函数调用的时候，result=spline3(x,Vx,Vy)，x为代求点的横坐标向量， (Vx,Vy)为已知的点的坐标。