西安交大 计算方法b上机作业

《西安交大 计算方法b上机作业》由会员上传分享，免费在线阅读，更多相关内容在行业资料-天天文库。

1、实用标准计算方法（B）上机作业一、三次样条拟合某通信公司在一次施工中，需要在水面宽度为20米的河沟底部沿直线走向铺设一条沟底光缆。在铺设光缆之前需要对沟底的地形进行初步探测，从而估计所需光缆的长度，为工程预算提供依据。已探测到一组等分点位置的深度数据(单位:米)如下表所示：分点0123456深度9.018.967.967.978.029.0510.13分点78910111213深度11.1812.2613.2813.3212.6111.2910.22分点14151617181920深度9.157.907.958.869.8110.8010

2、.93(1)请用合适的曲线拟合所测数据点；(2)估算所需光缆长度的近似值，并作出铺设河底光缆的曲线图；解：1、算法实现的思想及依据题目（1）为曲线拟合问题多项式插值、分段插值和最小二乘法。多项式插值，随着插值数据点的数目增多，误差也会随之增大，因此不选用。最小二乘法适于数据点较多的场合，在此也不适用。故选用分段插值。分段插值又分为分段线性插值、分段二次插值、三次样条插值及更高阶的多项式插值。由本题的物理背景知，光缆正常工作时各点应该是平滑过渡，因此至少选用三次样条插值法。对于更高阶的多项式插值，由于“龙格现象”而不选用。题目（2）求光缆长

3、度，即求拟合曲线在0到20的长度，对弧长进行积分即可。光缆长度的第一型线积分表达式为。文案大全实用标准2、算法实现的结构参考教材给出的SPLINEM算法和TTS算法，在选定边界条件和选定插值点等距分布后，可以先将数据点的二阶差商求出来并赋值给右端向量d，再根据TSS解法求解M。光缆长度的第一型线积分表达式为。3、程序运行结果及分析图1.1三种边界条件下三次样条插值图1.2光缆长度4、MATLAB代码：1）自己编程实现代码clear;clc;I=input('你想使用第几种边界条件？请输入1、2、3之一:');文案大全实用标准x=0:20;

4、y=[9.018.967.967.978.029.0510.1311.1812.2613.2813.3212.6111.2910.229.157.907.958.869.8110.810.93];plot(x,-y,'k.','markersize',15)%y为深度，取负号holdon%%计算一阶差商y1=ones(1,21);fori=2:1:21y1(i)=(y(i)-y(i-1))/(x(i)-x(i-1));end%%计算二阶差商y2=ones(1,21);fori=3:1:21y2(i)=(y1(i)-y1(i-1))/(x(

5、i)-x(i-2));end%%计算三阶差商y3=ones(1,21);fori=4:1:21y3(i)=(y2(i)-y2(i-1))/(x(i)-x(i-3));end%%选择边界条件（I）ifI==1d(1)=0;d(21)=0;a(21)=0;c(1)=0;%第一个点和最后一个点的二阶差商为0endifI==2d(1)=6*y1(1);d(21)=-6*y1(21);a(1)=1;c(1)=1;endifI==3d(1)=-12*y3(1);d(21)=-12*y3(21);a(21)=-2;c(1)=-2;%endfori=2:

6、20d(i)=6*y2(i+1);end文案大全实用标准%%构造带状矩阵求解（追赶法）b=2*ones(1,21);a=0.5*ones(1,21);%a(21)=-2;c=0.5*ones(1,21);%c(1)=-2;u(1)=b(1);r(1)=c(1);%%追yz(1)=d(1);fori=2:21l(i)=a(i)/u(i-1);u(i)=b(i)-l(i)*r(i-1);r(i)=c(i);yz(i)=d(i)-l(i)*yz(i-1);end%%赶xg(21)=yz(21)/u(21);fori=20:-1:1xg(i)=(

7、yz(i)-r(i)*xg(i+1))/u(i);endM=xg;%%所有点的二阶导数值%%求函数表达式并积分t=1;h=1;N=1000x1=0:20/(N-1):20;length=0;fori=1:Nforj=2:20ifx1(i)<=x(j)t=j;break;elset=j+1;endendf1=x(t)-x1(i);f2=x1(i)-x(t-1);S(i)=(M(t-1)*f1^3/6/h+M(t)*f2^3/6/h+(y(t-1)-M(t-1)*h^2/6)*f1+(y(t)-M(t)*h^2/6)*f2)/h;Sp(i)=

8、-M(t-1)*f1^2/2/h+M(t)*f2^2/2/h+(y(t)-y(t-1))/h-(M(t)-M(t-1))*h/6;文案大全实用标准length(i+1)=sqrt(1+Sp(i