计算方法第一次上机报告.docx
- 文档编号:10820578
- 上传时间:2023-05-27
- 格式:DOCX
- 页数:20
- 大小:490.76KB
计算方法第一次上机报告.docx
《计算方法第一次上机报告.docx》由会员分享,可在线阅读,更多相关《计算方法第一次上机报告.docx(20页珍藏版)》请在冰点文库上搜索。
计算方法第一次上机报告
《计算方法》实验报告一
题目:
班级:
计1402
学号:
1403140126
姓名:
应成龙
日期:
2016.3.17
程序名:
一、问题和要求:
1.土木工程和环境工程师在设计一条排水渠道时必须考虑渠道的各种参数(如宽度,深度,渠道内壁光滑度)及水流速度、流量、水深等物理量之间的关系。
假设修一条横断面为矩形的水渠,其宽度为B,假定水流是定常的,也就是说水流速度不随时间而变化。
为了不同的工业目的(比如说要把污染物稀释到一定的浓度以下,或者为某工厂输入一定量的水),需要指定流量Q和B,求出水的深度。
这样,就需要求解
其中Q是水的流量(
),U是流速(
),H是水的深度(
)。
一个具体的案例是
.求出渠道中水的深度H。
2.分析用迭代法求解方程组
的收敛性,并求出使||x(k+1)–x(k)||1≤0.0001的近似解及相应的迭代次数。
考虑用
(1)雅可比迭代法;
(2)赛德尔迭代法;(3)超松驰迭代法(取1.2,1.3,1.9,0.9)。
3.对下面的9阶对称矩阵分别做LU分解和LDLT分解
二、设计总体方案及分析:
2.1问题分析
第1题考查的是非线性方程组的迭代法(我用的是牛顿迭代法)
第2题考查的是线性代数方程组的迭代法,分别是
(1)雅可比迭代法,
(2)赛德尔迭代法,(3)超松驰迭代法
第3题考查的是线性代数方程组的矩阵三角分解法,分别是LU分解和LDLT分解
2.2概要设计
2.3详细设计(二合一)
第1题用牛顿迭代法,先判断f’(x)是否收敛,从x=0开始
,利用公式不断迭代,直至Xk+1-Xk的值小于0.000001,才停止。
并输出所有的迭代值
第2-1题用雅可比迭代法时,由于
公式中,Xj迭代的次数是k次,所以需要用一个额外的数组来保存之前第k次迭代结果,我用了x[6][N]={0},来存放那个数组,通过控制每个Xk+1次迭代与Xk次迭代差<=0.001控制循环次数,当每个差都同时<=0.001,循环停止,迭代结束
第2-2题用的是赛德尔迭代法,
由于已经进行过的迭代,直接覆盖原值,所以赛德尔迭代法程序设计比雅可比迭代要简单,不用在设计数组存放X的k次迭代值,而且只需在在雅可比迭代的基础上稍作修改就可以,而且循环控制方法与雅可比迭代相同
第2-3题用的是超松弛迭代法,
超松弛法在赛德尔迭代法基础上,加入w,只需稍作修改就好,循环控制与上面几题相同,这3题互相关联,只要写出一个程序,其余2个就一定会写了。
在难度上雅可比迭代稍难一些。
在程序中只要修改宏定义w就行了。
第3-1题用的是LU分解
这题我用了三个嵌套循环,在最内层循环有两个并列循环,先给U赋值,再给L赋值,只有U赋值了,L才能赋值,最外层两个循环实现a[][],l[][],u[][]的遍历,最终实现全部赋值。
第3-2题用的是LDLT分解
实际用到的是这3个公式
这一题我通过修改我自己写的LU分解,而且只用L和D两个公式行不通,一直出不来结果,借鉴了同学的代码,也是通过3个循环嵌套,实际用的是这3个公式,方法同LU分解一样。
2.4调试分析
第1题结果我就算出来两个0,和0.206213,估计是0这个数太接近实际结果了,所以迭代较快。
第2-1这题我是这里面用时比较长的一题,一开始是不知道怎么写,然后从公式入手,还得保存X的k次迭代,通过数组的方法来保存,最后终于写出来了。
我还把迭代次数及到该迭代完成后X的值在同一行表示出来
第2-2,2-3,都是在2-1的基础上稍作修改就很快算出了结果。
第3-1也是从公式入手,通过观察公式,发现规律,有了前面的经验之后很快就能写出来。
第3-2是由于公式选择的问题,导致一直编译出错,最后在借鉴了同学的代码之后完成。
2.5测试结果
第1题
第2-1
第2-2题
第2-3题
第3-1
第3-2
三、源程序及注释
第1题
#include
#include
#include
#definen0.03
#defineB20
#defineS0.0002
#defineQ5
usingnamespacestd;
boolshoulian(doubleH)
{
if(0.8/B*pow(n*Q/pow(S,0.5)/(B+2*H),0.6)<1)
returntrue;
elsereturnfalse;
}
doublediedai(doubleH)
{
return(pow(n*Q/pow(S,0.5)*pow(B+2*H,2/3),0.6)/B);
}
intmain()
{
doubleH;
for(doubleh=0;h<=100;h=h+0.5)
{
if(shoulian(h))
{
H=h;
break;
}
}
while(fabs(H-diedai(H))>=0.00000001)
{
cout< H=diedai(H); } cout< } 第2-1题 #include #include #include usingnamespacestd; #defineN1000 intmain() { intcounter=0; inta[][6]={4,-1,0,-1,0,0, -1,4,-1,0,-1,0, 0,-1,4,-1,0,-1, -1,0,-1,4,-1,0, 0,-1,0,-1,4,-1, 0,0,-1,0,-1,4}; doublex[6][N]={0}; intb[6]={0,5,0,6,-2,6}; for(intn=0;n { for(inti=0;i<6;i++) { doubles=0; for(intj=0;j<6;j++) { if(j! =i) { s-=a[i][j]*x[j][n]; } } s+=b[i]; s/=a[i][i]; x[i][n+1]=s; } if(fabs(x[0][n+1]-x[0][n])<=0.0001&&fabs(x[1][n+1]-x[1][n])<=0.0001&&fabs(x[2][n+1]-x[2][n])<=0.0001&&fabs(x[3][n+1]-x[3][n])<=0.0001&&fabs(x[4][n+1]-x[4][n])<=0.0001&&fabs(x[5][n+1]-x[5][n])<=0.0001) { counter=n+1; break; } } for(inti=0;i<=counter;i++) { cout< } } 第2-2题 #include #include #include usingnamespacestd; #defineN1000 intmain() { intcounter=0; inta[][6]={4,-1,0,-1,0,0, -1,4,-1,0,-1,0, 0,-1,4,-1,0,-1, -1,0,-1,4,-1,0, 0,-1,0,-1,4,-1, 0,0,-1,0,-1,4}; doublex[6]={0}; doublex1[6]={0}; intb[6]={0,5,0,6,-2,6}; for(intn=0;n { cout< if(n! =0) { if(fabs(x[0]-x1[0])<=0.0001&&fabs(x[1]-x1[1])<=0.0001&&fabs(x[2]-x1[2])<=0.0001&&fabs(x[3]-x1[3])<=0.0001&&fabs(x[4]-x1[4])<=0.0001&&fabs(x[5]-x1[5])<=0.0001) { counter=n+1; break; } } for(inti=0;i<6;i++) { x1[i]=x[i]; doubles=0; for(intj=0;j<6;j++) { if(j! =i) { s-=a[i][j]*x[j]; } } s+=b[i]; s/=a[i][i]; x[i]=s; } } } 第2-3题 #include #include #include usingnamespacestd; #defineN1000 #definew1.2//w取1.2,1.3,1.9,0.9 intmain() { intcounter=0; inta[][6]={4,-1,0,-1,0,0, -1,4,-1,0,-1,0, 0,-1,4,-1,0,-1, -1,0,-1,4,-1,0, 0,-1,0,-1,4,-1, 0,0,-1,0,-1,4}; doublex[6]={0}; doublex1[6]={0}; intb[6]={0,5,0,6,-2,6}; for(intn=0;n { cout< if(n! =0) { if(fabs(x[0]-x1[0])<=0.0001&&fabs(x[1]-x1[1])<=0.0001&&fabs(x[2]-x1[2])<=0.0001&&fabs(x[3]-x1[3])<=0.0001&&fabs(x[4]-x1[4])<=0.0001&&fabs(x[5]-x1[5])<=0.0001) { counter=n+1; break; } } for(inti=0;i<6;i++) { x1[i]=x[i]; doubles=0; for(intj=0;j<6;j++) { if(j! =i) { s-=a[i][j]*x[j]; } } s+=b[i]; s/=a[i][i]; s*=w; s+=(1-w)*x[i]; x[i]=s; } } } 第3-1题 #include #include #include usingnamespacestd; intmain() { doublea[9][9]={4,-1,0,-1,0,0,0,0,0, -1,4,-1,0,-1,0,0,0,0, 0,-1,4,0,0,-1,0,0,0, -1,0,0,4,-1,0,-1,0,0, 0,-1,0,-1,4,-1,0,-1,0, 0,0,-1,0,-1,4,-1,0,-1, 0,0,0,-1,0-1,4,-1,0, 0,0,0,0,-1,0-1,4,-1, 0,0,0,0,0,-1,0,-1,4}; doublel[9][9]={0}; doubleu[9][9]={0}; for(inti=0;i<9;i++) { u[0][i]=a[0][i]; l[i][i]=1; if(i! =0) { l[i][0]=a[i][0]/a[0][0]; } } for(intr=1;r<9;r++) { inti=r; doubles2=0; for(;i<9;i++) { doubles1=0; for(intk=0;k { s1-=l[r][k]*u[k][i]; } s1+=a[r][i]; u[r][i]=s1; } for(intk=0;k { s2-=l[r+1][k]*u[k][r]; } s2+=a[r+1][r]; s2/=u[r][r]; l[r+1][r]=s2; } cout<<"l[9][9]"< for(inti=0;i<9;i++) { for(intj=0;j<9;j++) { cout< } cout< } cout<<"u[9][9]"< for(inti=0;i<9;i++) { for(intj=0;j<9;j++) { cout< } cout< } } 第3-2题 #include #include #include #definen9 usingnamespacestd; intmain() { doublea[n][n]={4,-1,0,-1,0,0,0,0,0, -1,4,-1,0,-1,0,0,0,0, 0,-1,4,0,0,-1,0,0,0, -1,0,0,4,-1,0,-1,0,0, 0,-1,0,-1,4,-1,0,-1,0, 0,0,-1,0,-1,4,-1,0,-1, 0,0,0,-1,0,-1,4,-1,0, 0,0,0,0,-1,0,-1,4,-1, 0,0,0,0,0,-1,0,-1,4}; doublet[n][n]={0}; doublel[n][n]={0}; doubled[n][n]={0}; inti,j,k; for(i=0;i d[i][i]=a[i][i]; for(k=0;k d[i][i]-=t[i][k]*l[i][k]; for(j=0;j t[i][j]=a[i][j]; for(k=0;k t[i][j]-=t[i][k]*l[j][k]; l[i][j]=t[i][j]/d[i][i]; } } } cout<<"l[9][9]"< for(inti=0;i<9;i++) { for(intj=0;j<9;j++) { cout< } cout< } cout<<"d[9][9]"< for(inti=0;i<9;i++) { for(intj=0;j<9;j++) { cout< } cout< } } 四、设计总结(程序设计过程中的收获、遇到问题、解决问题过程的思考、程序调试的思考、对数据结构这门课程的思考认识等.对算法的程序的讨论、分析,改进设想,其它经验教训等) 通过这次上机,我发觉要正确编写计算方法代码,就要把公式理解透,有助于编写。 我发现了计算方法很多编程的相似之处,只要通过修改其中一些代码,就能得出其他题目的代码,前提是要深入理解公式。 这题来讲,这个上机,难度不是很大,只是纯粹的公式转换成代码的过程。
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 计算方法 第一次 上机 报告