资源描述:
《三次样条插值自然边界条件》由会员上传分享,免费在线阅读,更多相关内容在行业资料-天天文库。
1、例:已知一组数据点,编写一程序求解三次样条插值函数满足并针对下面一组具体实验数据0.250.30.390.450.530.50000.54770.62450.67080.7280求解,其中边界条件为.1)三次样条插值自然边界条件源程序:functions=spline3(x,y,dy1,dyn)%x为节点,y为节点函数值,dy1,dyn分别为x=0.25,0.53处的二阶导m=length(x);n=length(y);ifm~=nerror('xory输入有误')returnendh=zeros(1,n-1);h(n-1)=x(n)-x(n-1);fork=1:n-2h(k)=
2、x(k+1)-x(k);v(k)=h(k+1)/(h(k+1)+h(k));u(k)=1-v(k);endg(1)=3*(y(2)-y(1))/h(1)-h(1)/2*dy1;g(n)=3*(y(n)-y(n-1))/h(n-1)+h(n-1)/2*dyn;fori=2:n-1g(i)=3*(u(i-1)*(y(i+1)-y(i))/h(i)+v(i-1)*(y(i)-y(i-1))/h(i-1));endfori=2:n-1;A(i,i-1)=v(i-1);A(i,i+1)=u(i-1);endA(n,n-1)=1;A(1,2)=1;A=A+2*eye(n);M=zhuigf(
3、A,g);%调用函数,追赶法求Mfprintf('三次样条(三对角)插值的函数表达式');symsX;fork=1:n-1fprintf('S%d--%d:',k,k+1);s(k)=(h(k)+2*(X-x(k)))./h(k).^3.*(X-x(k+1)).^2.*y(k)...+(h(k)-2*(X-x(k+1)))./h(k).^3.*(X-x(k)).^2.*y(k+1)...+(X-x(k)).*(X-x(k+1)).^2./h(k).^2*M(k)+(X-x(k+1)).*...(X-x(k)).^2./h(k).^2*M(k+1);ends=s.';s=v
4、pa(s,4);%画三次样条插值函数图像fori=1:n-1X=x(i):0.01:x(i+1);st=(h(i)+2*(X-x(i)))./(h(i)^3).*(X-x(i+1)).^2.*y(i)...+(h(i)-2.*(X-x(i+1)))./(h(i)^3).*(X-x(i)).^2.*y(i+1)...+(X-x(i)).*(X-x(i+1)).^2./h(i)^2*M(i)+(X-x(i+1)).*...(X-x(i)).^2./h(i)^2*M(i+1);plot(x,y,'o',X,st);holdonEndplot(x,y);gridon%%%%%%%%%%%
5、%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%调用的函数:%追赶法functionM=zhuigf(A,g)n=length(A);L=eye(n);U=zeros(n);fori=1:n-1U(i,i+1)=A(i,i+1);endU(1,1)=A(1,1);fori=2:nL(i,i-1)=A(i,i-1)/U(i-1,i-1);U(i,i)=A(i,i)-L(i,i-1)*A(i-1,i);endY(1)=g(1);fori=2:nY(i)=g(i)-L(i,i-1)*Y(i-1);endM(n)=Y(n)/U(n,n);fori=n-1:-1:1M(i)
6、=(Y(i)-A(i,i+1)*M(i+1))/U(i,i);end2)在命令窗口输入x,y,dy1,dyn,得到三次样条函数:>>x=[0.25,0.3,0.39,0.45,0.53];>>y=[0.5,0.5477,0.6245,0.6708,0.7280];>>dy1=0;>>dyn=0;>>s=spline3(x,y,dy1,dyn)运行结果:三次样条(三对角)插值的函数表达式S1--2:S2--3:S3--4:S4--5:.5000*(-3600.+.1600e5*X)*(X-.3000)^2+.5477*(5200.-.1600e5*X)*(X-.2500)^2+39
7、4.8*(X-.2500)*(X-.3000)^2+355.2*(X-.3000)*(X-.2500)^2.5477*(-699.6+2743.*X)*(X-.3900)^2+.6245*(1193.-2743.*X)*(X-.3000)^2+109.6*(X-.3000)*(X-.3900)^2+96.77*(X-.3900)*(X-.3000)^2.6245*(-3333.+9259.*X)*(X-.4500)^2+.6708*(4444.-9259.*X)*(X-.3900