计算方法B上机报告材料
《计算方法B上机报告材料》由会员分享,可在线阅读,更多相关《计算方法B上机报告材料(15页珍藏版)》请在装配图网上搜索。
1、word计算方法B上机报告第1题某通信公司在一次施工中,需要在水面宽度为20米的河沟底部沿直线走向铺设一条沟底光缆。在铺设光缆之前需要对沟底的地形进展初步探测,从而估计所需光缆的长度,为工程预算提供依据。已探测到一组等分点位置的深度数据(单位:米)如下表所示:分点0123456深度分点78910111213深度分点14151617181920深度 (1)请用适宜的曲线拟合所测数据点;(2)估算所需光缆长度的近似值,并作出铺设河底光缆的曲线图;问题分析和算法思想:此题的主要目的是对21个测量数据进展拟合,同时对拟合曲线进展线积分即可得到河底光缆长度的近似值,可以用的插值方法很多:多项式插值、La
2、grange插值、Newton插值、三次样条插值等。由于数值点较多时,采用高次多项式插值将产生很大的误差,用拉格朗日插值多项式会出现龙格现象。故为了将所有的数据点都用上,且题中光缆为柔性,可光滑铺设于水底,鉴于此特性,采用三次样条插值的方法较为适宜。计算光缆长度近似值,只需将每两点之间的距离算出,然后依次相加,所得的折线长度,即为光缆长度的近似值。光缆长度计算公式:算法结构:三次样条算法结构见计算方法教程P110。源程序:clear;clc;x=0:20;y=9.01 8.96 7.96 7.97 8.02 9.05 10.13 11.18 12.26 13.28 13.32 12.61 11
3、.29 10.22 9.15 7.90 7.95 8.86 9.81 10.80 10.93;d=y;plot(x,y,k.,markersize,15)hold on%计算二阶差商for k=1:2 for i=21:-1:(k+1) d(i)=(d(i)-d(i-1)/(x(i)-x(i-k); endend%假定d的边界条件,采用自然三次样条for i=2:20 d(i)=6*d(i+1);endd(1)=0;d(21)=0;%追赶法求解带状矩阵的m值a=0.5*ones(1,21);b=2*ones(1,21);c=0.5*ones(1,21);a(1)=0;c(21)=0;u=one
4、s(1,21);u(1)=b(1);r=c;yy(1)=d(1);%追的过程for k=2:21 l(k)=a(k)/u(k-1); u(k)=b(k)-l(k)*r(k-1); yy(k)=d(k)-l(k)*yy(k-1);end%赶的过程m(21)=yy(21)/u(21);for k=20:-1:1 m(k)=(yy(k)-r(k)*m(k+1)/u(k);end%利用插值点画出拟合曲线k=1;nn=100;xx=linspace(0,20,nn);l=0;for j=1:nn for i=2:20 if xx(j)0 %判断在区间内是否有error(No answer,because
5、 f0*f10);end for i=1:201 x=(x0+x1)/2; if abs(x1-x)e %判断解是否符合误差 disp(方程的解,num2str(x); disp(迭代步数,num2str(i); disp(误差,num2str(abs(x1-x); break; end f=T3Sub(x,n); if (f*f0)0 x1=x; f1=f; else x0=x; f0=f; end endfunction f=T3Sub(x,n)%复化梯形求积公式h=pi/n; f=0;for i=1:n u=i*h; f=f+cos(x*sin(u-h)+cos(x*sin(u);end
6、 f=(h*f)/(2*pi); end 运行结果:第4题线性方程组求解1编写程序实现大规模方程组的高斯消去法程序,并对所附的方程组进展求解。所附方程组的类型为对角占优的带状方程组。2针对本专业中所碰到的实际问题,提炼一个使用方程组进展求解的例子,并对求解过程进展分析、求解。算法思想:高斯消去法是利用现行方程组初等变换中的一种变换,将一个不为零的数乘到一个方程后加到另一个方程,使方程组变成同解的上三角方程组,然后再自下而上对上三角方程组求解。算法结构:1.读取二进制文件,存入计算矩阵2.对矩阵进展初等变换,然后求解(计算方法教程高斯消去法与列主元高斯消去法算法)源代码:clear;clc;%
7、读取系数矩阵f,p=uigetfile(*.dat,选择数据文件); %读取数据文件num=5; %输入系数矩阵文件头的个数name=strcat(p,f);file=fopen(name,r);head=fread(file,num,uint); %读取二进制头文件id=dec2hex(head(1); %读取标识符fprintf(文件标识符为:);idver=dec2hex(head(2); %读取版本号fprintf(文件版本号为:);vern=head(3); %读取阶数fprintf(矩阵A的阶数:);nq=head(4); %上带宽fprintf(矩阵A的上带宽:);q p=hea
8、d(5); %下带宽fprintf(矩阵A的下带宽:);p dist=4*num;fseek(file,dist,bof); %把句柄值转向第六个元素开头处A,count=fread(file,inf,float); %读取二进制文件,获取系数矩阵fclose(file); %关闭二进制头文件% 对非压缩带状矩阵进展求解if ver=102, a=zeros(n,n); for i=1:n, for j=1:n, a(i,j)=A(i-1)*n+j); %求系数矩阵a(i,j) end end b=zeros(n,1); for i=1:n, b(i)=A(n*n+i); end for k=
9、1:n-1, %列主元高斯消去法 m=k; for i=k+1:n, %寻找主元 if abs(a(m,k)abs(a(i,k) m=i; end end if a(m,k)=0 %遇到条件终止 disp(错误!) return end for j=1:n, %交换元素位置得主元 t=a(k,j); a(k,j)=a(m,j); a(m,j)=t; t=b(k); b(k)=b(m); b(m)=t; end for i=k+1:n, %计算l(i,k)并将其放到a(i,k)中 a(i,k)=a(i,k)/a(k,k); for j=k+1:n a(i,j)=a(i,j)-a(i,k)*a(k
10、,j); end b(i)=b(i)-a(i,k)*b(k); end end x=zeros(n,1); %回代过程 x(n)=b(n)/a(n,n); for k=n-1:-1:1, x(k)=(b(k)-sum(a(k,k+1:n)*x(k+1:n)/a(k,k); endend% 对压缩带状矩阵进展求解if ver=202, %高斯消去法 m=p+q+1; a=zeros(n,m); for i=1:1:n for j=1:1:m a(i,j)=A(i-1)*m+j); end end b=zeros(n,1); for i=1:1:n b(i)=A(n*m+i); %求b(i) en
11、d for k=1:1:(n-1) %开始消去 if a(k,(p+1)=0 disp(错误!); break; end st1=n; if (k+p)n st1=k+p; end for i=(k+1):1:st1 a(i,(k+p-i+1)=a(i,(k+p-i+1)/a(k,(p+1); for j=(k+1):1:(k+q) a(i,j+p-i+1)=a(i,j+p-i+1)-a(i,k+p-i+1)*a(k,j+p-k+1); end b(i)=b(i)-a(i,k+p-i+1)*b(k); end end x=zeros(n,1); %回代 x(n)=b(n)/a(n,p+1);
12、sum=0; for k=(n-1):-1:1 sum=b(k); st2=n; if (k+q)abs(A(max,i) max=k; end end temp=A(i,i:m);A(i,i:m)=A(max,i:m);A(max,i:m)=temp; for j=(i+1):m A(j,:)=A(j,:)-A(i,:)*A(j,i)/A(i,i); endenddisp(回代求解) x(m)=A(m,n)/A(m,m); for k=(m-1):-1:1 x(k)=(A(k,n)-A(k,k+1:m)*x(k+1:m)/A(k,k); end x运算结果:T1=28.7368,T2=24.2632,T3=20.6842,T4=18.3158。运算结果和正确结果吻合。精彩文档
- 温馨提示:
1: 本站所有资源如无特殊说明,都需要本地电脑安装OFFICE2007和PDF阅读器。图纸软件为CAD,CAXA,PROE,UG,SolidWorks等.压缩文件请下载最新的WinRAR软件解压。
2: 本站的文档不包含任何第三方提供的附件图纸等,如果需要附件,请联系上传者。文件的所有权益归上传用户所有。
3.本站RAR压缩包中若带图纸,网页内容里面会有图纸预览,若没有图纸预览就没有图纸。
4. 未经权益所有人同意不得将文件中的内容挪作商业或盈利用途。
5. 装配图网仅提供信息存储空间,仅对用户上传内容的表现方式做保护处理,对用户上传分享的文档内容本身不做任何修改或编辑,并不能对任何下载内容负责。
6. 下载文件中如有侵权或不适当内容,请与我们联系,我们立即纠正。
7. 本站不保证下载资源的准确性、安全性和完整性, 同时也不承担用户因使用这些下载资源对自己和他人造成任何形式的伤害或损失。