数值分析上机作业最强版.docx
《数值分析上机作业最强版.docx》由会员分享,可在线阅读,更多相关《数值分析上机作业最强版.docx(28页珍藏版)》请在冰豆网上搜索。
![数值分析上机作业最强版.docx](https://file1.bdocx.com/fileroot1/2023-1/26/54a53dad-8c68-446c-95aa-01e4d20bcfa3/54a53dad-8c68-446c-95aa-01e4d20bcfa31.gif)
数值分析上机作业最强版
数值分析上机作业
姓名:
唐皓
学号:
142460
专业:
道路与铁道工程
院系:
交通学院
授课教师:
吴宏伟
日期:
2015年1月
习题一
1题目
17.(上机题)舍入误差与有效数
设
,其精确值为
。
(1)编制按从大到小的顺序
,计算
的通用程序;
(2)编制按从小到大的顺序
,计算
的通用程序;
(3)按两种顺序分别计算
,
,
,并指出有效位数。
(编制程序时用单精度);
(4)通过本上机题你明白了什么?
2通用程序代码
2.1按从小到大的顺序计算
voidAscendSum(unsignedlongintN)//计算从大到小的总和
{
for(unsignedlongintj=2;j<=N;j++)
ascendSum+=(float)1.0/(j*j-1);
cout<<"SumFrom1toN(Ascend)is:
"<Error(ascendSum);
Delimiter();
}
2.2按从大到小的顺序计算
voidDescendSum(unsignedlongintN)//计算从小到大的总和
{
for(unsignedlongintj=N;j>=2;j--)
descendSum+=(float)1.0/(j*j-1);
cout<<"SumFromNto1(Descend)is:
"<Error(descendSum);
Delimiter();
}
3计算结果展示
图1N=100时的计算结果
图2N=10000时的计算结果
图3N=1000000时的计算结果
表1-1计算结果汇总
精确值
按从小到大
按从大到小
值
有效位数
值
有效位数
0.7400494814
0.7400494814
10
0.740049541
6
0.7498999834
0.7498521209
4
0.7498999834
10
0.7499989867
0.751856029
2
0.7529925108
2
4计算结果分析
(1)如果采用单精度数据结构进行计算,则相较于双精度的数据结果,由于数据存储字长的限制导致计算机存在较大的舍入误差,因此本程序采用的是双精度数据存储方式。
(2)由计算结果可知,正序计算和逆序计算的精度是不稳定的。
由计算结果可以发现,当N=100时,正序计算(1-N)的精度较高;当N=10000时,逆序计算(N-1)的精度较高;当N=1000000时,正序计算和逆序计算的精度一样。
当然,和其他同学做出来的结果对比,结论并不一致。
我个人分析这是因为电脑的硬件、软件(位数)不同等原因导致的。
但总体而言,在N较小时,正序计算精度高于逆序计算的精度。
当N较大时,正序和逆序计算的精度接近。
(3)由于计算机的实际计算过程是一种舍入机制,故对于我们计算所采用的加法交换律是不成立的。
计算机中若干数相加时,先要进行对阶操作,即将两数的阶数统一为绝对值较大的数的阶数。
这样一来将导致绝对值较小的数的有效数字可能会大量损失,增大舍入误差,即所谓的“大数吃小数”现象。
为了避免这种现象的出现,在进行加减法的时候应该先将绝对值较小的数相加,再与绝对值较大的数相加这样按阶逐步递增的相加。
5完整代码
#include
#include
#include
usingnamespacestd;
floataccurateSum=0,ascendSum=0,descendSum=0;
voidDelimiter()//输出一系列星号以间隔
{
for(inti=1;i<=50;i++)cout<<"*";
cout<}
voidError(floatSum)//计算绝对误差
{
floaterror;
error=fabs(Sum-accurateSum);
intflag;
for(flag=0;flag<10;flag++)
{
error=error*10;
if(error>0.5)
break;
}
cout<<"Thereare"<}
voidAccurateSum(unsignedlongintN)//计算精确值
{
accurateSum=0.5*(1.5-(float)1/N-(float)1/(N+1));
cout<<"Accuratesumis:
"<Delimiter();
}
voidAscendSum(unsignedlongintN)//计算从大到小的总和
{
for(unsignedlongintj=2;j<=N;j++)
ascendSum+=(float)1.0/(j*j-1);
cout<<"SumFrom1toN(Ascend)is:
"<Error(ascendSum);
Delimiter();
}
voidDescendSum(unsignedlongintN)//计算从小到大的总和
{
for(unsignedlongintj=N;j>=2;j--)
descendSum+=(float)1.0/(j*j-1);
cout<<"SumFromNto1(Descend)is:
"<Error(descendSum);
Delimiter();
}
voidmain()//主程序
{
longintN;
while
(1)
{
cout<<"InputanintegerN(N>=2):
";
cin>>N;
if(N<2)
cout<<"Invalidinput,Pleasetryitagain!
\n";
else
break;
}
AccurateSum(N);
AscendSum(N);
DescendSum(N);
}
习题二
1题目
20.(上机题)Newton迭代法
(1)给定初值
及容许误差
,编制Newton法解方程
根的通用程序。
(2)给定方程
,易知其有三个根
,
,
。
①由Newton方法的局部收敛性可知存在
,当
时,Newton迭代序列收敛于根
。
试确定尽可能大的
。
②试取若干初始值,观察当
,
,
,
,
时Newton序列是否收敛以及收敛于哪一个根。
(3)通过本上机题,你明白了什么?
2通用程序代码
2.1Newton法解方程通用程序
doubleNewtonIteration(doublex0,doubleeps)//Newton迭代法求解子程序
{
doublex1,x2;
x1=x0;
x2=x1-f(x1)/df(x1);
while(fabs(x1-x2)>=eps)
{
x1=x2;
x2=x1-f(x1)/df(x1);
}
returnx1;
}
2.2求解尽可能大
doubleMaximalDeviateRange()//求解尽可能大的范围
{
doublestep=1e-5;//steplength
intcnt=1;//stepcount
doubledelta;
cout<<"**********************NewtonIteration(eps=1e*5)**********************"<while(fabs(NewtonIteration(step*cnt,eps))<=eps)
{
cnt++;
}
delta=step*cnt;
cout<<"Themaximaldeviaterangeforxconvergingtox2*=0is(-"<returndelta;
}
6计算结果展示
图2-1计算结果
在取步长为
的情况下,允许误差
时有
轴上的一个小区间
为迭代序列在
处的尽可能大的局部收敛区间。
当
,
,
,
,
时牛顿迭代序列分别收敛于
。
7计算结果分析
(1)通过本次上机编程并通过多次的调试,可以发现运行结果很好的验证了教材上牛顿迭代法具有局部收敛性这一重要性质。
(2)选择不同的初值区间,迭代序列会收敛于不同的根。
所以为了收敛到需要的计算结果,就需要选择合适的牛顿迭代法大范围收敛区间。
(3)产生上述结果的原因是所选取的部分区间不满足大范围收敛的条件,导致并没有收敛到理想的结果。
8完整代码
#include
#include
#include
usingnamespacestd;
doubleeps=1e-5;
doublef(doublex)//函数f(x)的录入
{
return(x*x*x)/3-x;
}
doubledf(doublex)//函数f(x)的导数
{
return(x*x)-1;
}
doubleNewtonIteration(doublex0,doubleeps)//Newton迭代法求解子程序
{
doublex1,x2;
x1=x0;
x2=x1-f(x1)/df(x1);
while(fabs(x1-x2)>=eps)
{
x1=x2;
x2=x1-f(x1)/df(x1);
}
returnx1;
}
doubleMaximalDeviateRange()//求解尽可能大的范围
{
doublestep=1e-5;//steplength
intcnt=1;//stepcount
doubledelta;
cout<<"**********************NewtonIteration(eps=1e*5)**********************"<while(fabs(NewtonIteration(step*cnt,eps))<=eps)
{
cnt++;
}
delta=step*cnt;
cout<<"Themaximaldeviaterangeforxconvergingtox2*=0is(-"<returndelta;
}
voidCalculate(doublex0)//计算Newton迭代法和相应的子程序
{
printf("Ifx0=%11.4f,xconvergestothevalueof:
%8.6f\n",x0,NewtonIteration(x0,eps));
}
voidmain()//主程序
{
doubledelta=MaximalDeviateRange();
Calculate(-10000);
Calculate((-1-delta)/2);
Calculate(-delta/2);
Calculate(delta/2);
Calculate((1+delta)/2);
Calculate(10000);
}
习题三
9题目
39.(上机题)列主元消去法
对于某电路的分析,归结为求解线性方程组。
其中
(1)编制解
阶线性方程组的列主元消去法的通用程序;
(2)用所编程序解线性方程组,并打印出解向量,保留5位有效数;
(3)本题编程之中,你提高了哪些编程能力?
10通用程序代码
voidGaussElimination(floatmatrix[ROW][COL])//高斯列主元消去法主程序
{
introw,col;
inttmpRow,maxRow,tmpCol;
floatratio;
floatroot[ROW];
//Calculatetheechelonmatrix
for(row=0;row{
col=row;
maxRow=row;
for(tmpRow=row;tmpRow{
if(fabs(matrix[tmpRow][col])>fabs(matrix[maxRow][col]))
maxRow=tmpRow;
}
SwapRow(matrix[row],matrix[maxRow]);
for(tmpRow=row+1;tmpRow{
ratio=matrix[tmpRow][col]/matrix[row][col];
if(ratio==0)
continue;
for(tmpCol=0;tmpCol
{
matrix[tmpRow][tmpCol]-=matrix[row][tmpCol]*ratio;
}
}
}
printf("TheEchelonMatrixis:
\n");
PrintMatrix(AugmentedMatrix);
PrintDelimiter();
//Substitutebackforthefinalroots.
root[ROW-1]=AugmentedMatrix[ROW-1][COL-1]/AugmentedMatrix[ROW-1][ROW-1];
floatsum;
intcnt;
for(tmpRow=ROW-2;tmpRow>=0;tmpRow--)
{
sum=0;
for(cnt=tmpRow+1;cntsum+=AugmentedMatrix[tmpRow][cnt]*root[cnt];
root[tmpRow]=(AugmentedMatrix[tmpRow][COL-1]-sum)/AugmentedMatrix[tmpRow][tmpRow];
}
printf("TheSolutionVectoris:
\n[");
for(cnt=0;cnt{
printf("%8.5f",root[cnt]);
}
printf("]\n");
}
11计算结果展示
图3-1计算结果
解向量为:
12结果分析和心得体会
(1)在本题编程过程中,由于考虑到列主元Gauss消去法的算法特点:
程序的过程重复性高,故本程序采用模块化程序设计方法。
这样就可以将每步都要重复进行的消去、交换列主元等过程的功能以模块函数来实现。
这样可以使程序结构清晰,从而调试变得相对容易。
(2)在输入/输出的控制方面,本程序采用格式化输出,使输出结果自动对齐,提高了程序的结果的可阅读性。
(3)为了能更方便快捷地进行矩阵数据的输入。
将需要求解的线性方程增广矩阵存储在文件中,程序运行时自动调用并输出。
因此,对于任意线性矩阵,只需将其增广矩阵存储于其中,并修改ROW和COL为实际行数与列数即可。
(4)通过该程序的编写,提高了我对循环语句的运用、子程序调用以及找bug(错位)的能力。
此外,我也对列主元Gauss消去法的步骤有了更加深刻的理解。
13完整代码
#include
#include
#include
#include
constintROW=9,COL=10;
floatAugmentedMatrix[ROW][COL];
voidPrintMatrix(floatmatrix[ROW][COL])//输出解子程序
{
introw,col;
for(row=0;row{
for(col=0;col
{
printf("%9.5f",matrix[row][col]);
}
printf("\n");
}
}
voidPrintDelimiter()
{
for(inti=0;i<100;i++)
printf("*");
printf("\n");
}
voidReadData(charfilename[],floatmatrix[ROW][COL])//读入数据程序
{
ifstreaminfile(filename);
if(!
infile)
cout<<"CannotopentheAugmentedMatrixdatafile!
\n";
cout<<"TheAugmentedMatrixis:
"<for(inti=0;i{
for(intj=0;j<(COL);j++)
{
infile>>matrix[i][j];
cout<}
cout<}
PrintDelimiter();
infile.close();
}
voidSwapRow(floatrowA[],floatrowB[])//交换两行的子程序
{
inti;
floattmp;
for(i=0;i
{
tmp=rowA[i];
rowA[i]=rowB[i];
rowB[i]=tmp;
}
}
voidGaussElimination(floatmatrix[ROW][COL])//高斯列主元消去法主程序
{
introw,col;
inttmpRow,maxRow,tmpCol;
floatratio;
floatroot[ROW];
//Calculatetheechelonmatrix
for(row=0;row{
col=row;
maxRow=row;
for(tmpRow=row;tmpRow{
if(fabs(matrix[tmpRow][col])>fabs(matrix[maxRow][col]))
maxRow=tmpRow;
}
SwapRow(matrix[row],matrix[maxRow]);
for(tmpRow=row+1;tmpRow{
ratio=matrix[tmpRow][col]/matrix[row][col];
if(ratio==0)
continue;
for(tmpCol=0;tmpCol
{
matrix[tmpRow][tmpCol]-=matrix[row][tmpCol]*ratio;
}
}
}
printf("TheEchelonMatrixis:
\n");
PrintMatrix(AugmentedMatrix);
PrintDelimiter();
//Substitutebackforthefinalroots.
root[ROW-1]=AugmentedMatrix[ROW-1][COL-1]/AugmentedMatrix[ROW-1][ROW-1];
floatsum;
intcnt;
for(tmpRow=ROW-2;tmpRow>=0;tmpRow--)
{
sum=0;
for(cnt=tmpRow+1;cntsum+=AugmentedMatrix[tmpRow][cnt]*root[cnt];
root[tmpRow]=(AugmentedMatrix[tmpRow][COL-1]-sum)/AugmentedMatrix[tmpRow][tmpRow];
}
printf("TheSolutionVectoris:
\n[");
for(cnt=0;cnt{
printf("%8.5f",root[cnt]);
}
printf("]\n");
}
voidmain()//主程序
{
ReadData("AugmentedMatrix.txt",AugmentedMatrix);
GaussElimination(AugmentedMatrix);
}
题目四
14题目
37.(上机题)3次样条插值函数
(1)编制求第一型3次样条插值函数的通用程序;
(2)已知汽车门曲线型值点的数据如下:
0
1
2
3
4
5
6
7
8
9
10
0
1
2
3
4
5
6
7
8
9
10
2.51
3.30
4.04
4.70
5.22
5.54
5.78
5.40
5.57
5.70
5.80
端点条件为
,
,用所编制程序求车门的3次样条插值函数
并打印出
。
15计算结果展示
通过运行程序,求得S(X)的表达式后,分别求各个点的值,得到如下结果:
图4-1运行结果
即有:
S(0.5)=2.90856425
S(1.5)=3.67842873
S(2.5)=4.38147084
S(3.5)=4.98818792
S(4.5)=5.38327747
S(5.5)=5.72370219
S(6.5)=5.59441375
S(7.5)=5.4298928
S(8.5)=5.65976505
S(9.5)=5.73229699
16完整代码
#include
#include
#include
usingnamespacestd;
cons
|
|
|
|
|
|
|
|
|
|