《计算方法与实习》实验报告

课程:计算方法与实习

学期：

2010-2011学年第三学期

学院：电气工程学院 学号：双号 姓名：XXX 2011年5月26日

习题一：

10000用两种不同的顺序计算

?nn?1?2?1.644834，分析其误差的变化。

思路分析

用一个循环语句，对n?2从1到10000进行叠加，两种不同顺序指从1叠加到10000和

从10000叠加到1，每隔一定的叠加次数就比较一次误差。

用C++语言编程

（1）从1叠加到10000源代码如下： #include #include #include using namespace std; int main(){ double N=10000,i=0; int a;

double n=0,S=1.644834; for(i=1;i<=N;i++){

n+=1/(i*i); a=i; if(aP0==0)cout<

} return 0;

运行结果如下：

迭代500次时，和为：S=1.642936 误差为：i=0.001897934

迭代1000次时，和为：S=1.643935 误差为：i=0.0008994333

迭代1500次时，和为：S=1.644268 误差为：i=0.0005663776

迭代2000次时，和为：S=1.644434 误差为：i=0.0003998082

迭代2500次时，和为：S=1.644534 误差为：i=0.0002998532

迭代3000次时，和为：S=1.644601 误差为：i=0.0002332109

计算方法与实习

迭代3500次时，和为：S=1.644648 误差为：i=0.0001856066

迭代4000次时，和为：S=1.644684 误差为：i=0.0001499019

迭代4500次时，和为：S=1.644712 误差为：i=0.0001221307

迭代5000次时，和为：S=1.644734 误差为：i=9.991315e-005

迭代5500次时，和为：S=1.644752 误差为：i=8.173481e-005

迭代6000次时，和为：S=1.644767 误差为：i=6.658593e-005

迭代6500次时，和为：S=1.64478 误差为：i=5.376747e-005

迭代7000次时，和为：S=1.644791 误差为：i=4.278009e-005

迭代7500次时，和为：S=1.644801 误差为：i=3.32576e-005

迭代8000次时，和为：S=1.644809 误差为：i=2.492534e-005

迭代8500次时，和为：S=1.644816 误差为：i=1.757329e-005

迭代9000次时，和为：S=1.644823 误差为：i=1.103809e-005

迭代9500次时，和为：S=1.644829 误差为：i=5.19077e-006

迭代10000次时，和为：S=1.644834 误差为：i=7.184807e-008

Press any key to continue

（2）从10000叠加到1的源代码如下：

#include #include #include using namespace std; int main(){ double N=10000,i=0; int a;

double n=0,S=1.644834; for(i=N;i>=1;i--){ n+=1/(i*i);

a=i;

if(a==N)cout<<\以下为n大于200时每隔500所积累的和：\

if(aP0==0&&a>=200)cout<

3

：

计算方法与实习

S=\误差为：i=\

if(a==200)cout<

return 0; }

运行结果如下：

以下为n大于200时每隔500所积累的和：

n=10000时，和为：S=1e-008 误差为：i=1.644834

n=9500时，和为：S=5.273698e-006 误差为：i=1.644829

n=9000时，和为：S=1.112228e-005 误差为：i=1.644823

n=8500时，和为：S=1.765898e-005 误差为：i=1.644816

n=8000时，和为：S=2.501281e-005 误差为：i=1.644809

n=7500时，和为：S=3.334722e-005 误差为：i=1.644801

n=7000时，和为：S=4.287235e-005 误差为：i=1.644791

n=6500时，和为：S=5.386299e-005 误差为：i=1.64478

n=6000时，和为：S=6.668556e-005 误差为：i=1.644767

n=5500时，和为：S=8.183971e-005 误差为：i=1.644752

n=5000时，和为：S=0.000100025 误差为：i=1.644734

n=4500时，和为：S=0.0001222519 误差为：i=1.644712

n=4000时，和为：S=0.0001500363 误差为：i=1.644684

n=3500时，和为：S=0.0001857601 误差为：i=1.644648

n=3000时，和为：S=0.0002333939 误差为：i=1.644601

n=2500时，和为：S=0.000300085 误差为：i=1.644534

4

计算方法与实习

n=2000时，和为：S=0.00040013 误差为：i=1.644434

n=1500时，和为：S=0.0005668939 误差为：i=1.644267

n=1000时，和为：S=0.0009005052 误差为：i=1.643933

n=500时，和为：S=0.001902006 误差为：i=1.642932

以下为n小于200时每隔20所积累的和

n=180时，和为：S=0.005471021 误差为：i=1.639363

n=160时，和为：S=0.006169577 误差为：i=1.638664

n=140时，和为：S=0.007068433 误差为：i=1.637766

n=120时，和为：S=0.008268157 误差为：i=1.636566

n=100时，和为：S=0.009950172 误差为：i=1.634884

n=80时，和为：S=0.01247846 误差为：i=1.632356

n=60时，和为：S=0.01670633 误差为：i=1.628128

n=40时，和为：S=0.02521511 误差为：i=1.619619

n=20时，和为：S=0.05117083 误差为：i=1.593663

n=10时，和为：S=0.1050663 误差为：i=1.539768

n=9时，和为：S=0.117412 误差为：i=1.527422

n=8时，和为：S=0.133037 误差为：i=1.511797

n=7时，和为：S=0.1534452 误差为：i=1.491389

n=6时，和为：S=0.181223 误差为：i=1.463611

n=5时，和为：S=0.221223 误差为：i=1.423611

n=4时，和为：S=0.283723 误差为：i=1.361111

n=3时，和为：S=0.3948341 误差为：i=1.25

5

计算方法与实习

n=2时，和为：S=0.6448341 误差为：i=0.9999999

n=1时，和为：S=1.644834 误差为：i=7.184806e-008

Press any key to continue

结果分析

从这两种不同顺序的结果来看，随着叠加次数的增加，误差是逐渐变小的，当全部叠加完之后，其结果达到了预期值1.644834，而且两种方法出来的最终误差是相同的

习题二：

用牛顿法求下列方程的根：（1）x2?ex?0；

思路分析

1.给定初始值

x0，e为根的容许误差，N为迭代次数的容许值 '2.若

f(x0)?0或迭代次数大于N，则算法失败，结束，否则转向3 xf(x0)1?x0?3.计算f'(x0)

4.若x1?x0?e则输出x1，否则令

x0?x1，转向2

C++语言编程源代码如下：

#include #include #include using namespace std; int main(){ double x0=1,x1,e=0.00001,u; int i,N=100;

cout<<\迭代过程数据如下：\

cout<<\误差\ for(i=1;;i++){ x1=x0-(x0*x0-exp(x0))/(2*x0-exp(x0)); u=fabs(x1-x0);

cout<N)break;

else x0=x1;

6

计算方法与实习

}

}

if(i>N)cout<<\迭代次数超出限制，算法失败\

if(u

实验结果如下：

迭代过程数据如下：

k x0 x1 误差 1 1 -1.39221 2.39221 2 -1.39221 -0.835088 0.557124 3 -0.835088 -0.709834 0.125253 4 -0.709834 -0.703483 0.00635069 5 -0.703483 -0.703467 1.59816e-005 6 -0.703467 -0.703467 1.0107e-010 算法成功，迭代次数为：n=6 Press any key to continue

第三章第1题（1）：

编写用追赶法解三对角线性方程组的程序，并解下列方程组：

??4?1?????????1?411?4?1?1??41???27?????15??????15??，b10?1???

??????15?1?????4?????15??Ax?b，其中A10?10

思路分析

1.用一个二维数组存储an,bn,cn,xi,dn

ai?ai?aibi?bibi?1?bi?1?aibi?cidi?1?di?1?aibi?di2.用一个循环语句迭代使矩阵变换为上三角矩阵

3.求出xn

3.用一个循环语句完成回代法,

和，，

xi?di?cixi?1

C++语言编程源代码如下：

#include

7

计算方法与实习

#include #include using namespace std; int main(){ double s[5][10]; int i,j,n=5,m=10; for(i=0;i<=n-1;i++){ for(j=0;j<=m-1;j++){ if(i==0)s[i][j]=1;//an赋值 else if(i==1)s[i][j]=-4;//bn赋值 else if(i==2)s[i][j]=1;//cn赋值 else if(i==3)s[i][j]=0;//xn赋值

else if(i==4&&j==0)s[i][j]=-27;//dn赋值 else if(i==4&&j!=0)s[i][j]=-15;

}

}

for(i=0;i<=m-2;i++){ s[4][i+1]=s[4][i+1]-(s[0][i]/s[1][i])*s[4][i]; s[1][i+1]=s[1][i+1]-(s[0][i]/s[1][i])*s[2][i]; s[0][i]=0;

}

s[0][9]=0;

s[3][9]=s[4][9]/s[1][9]; for(i=m-2;i>=0;i--){ s[3][i]=(s[4][i]-s[2][i]*s[3][i+1])/s[1][i]; }

cout<<\方程的解如下：\

for(i=0;i<=m-1;i++){ cout<<\ } return 0;

}

实验结果如下：

方程的解如下： x0=8.70576 x1=7.82303 x2=7.58637 x3=7.52245 x4=7.50344 x5=7.49131 x6=7.46179 x7=7.35584 x8=6.96156 x9=5.49039

8

计算方法与实习

Press any key to continue

习题4、分别用雅可比迭代法与高斯塞德尔迭代法解下列方程组:

?31??13??0??0，其中R??0??0?0??0??0?1335?90000000?9?31?100000000?1079?3000000?110?305770?50?10000?747?300000000?3041000000?50027?20??0?0???9?0??0?0???2??29?RI?V，

??15???27????23???0??V???20?

??12????7????7?????10?一、雅可比迭代法

思路分析

1.用一个二维数组存储方程组的各项系数，两个数组x0,x1分别存储x和x第二层和第三层进行雅克比迭代

3.输出方程的解

kk?1

2.用三个嵌套的循环语句，第一层把上一次循环的x1赋给x0,并判断是否符合方程解,

C++语言编程源代码如下：

//雅克比迭代式 #include #include #include using namespace std;

int main(){ double

s[9][10]={{ 31,-13,0,0,0,-10,0,0,0,-15},{-13,35,-9,0,-11,0,0,0,0,27},{0,-9,31,-10,0,0,0,0,0,-23},{0

9

计算方法与实习

,0,-10,79,-30,0,0,0,-9,0},{

0,0,0,-30,57,-7,0,-5,0,-20},{0,0,0,0,7,47,-30,0,0,12},{0,0,0,0,0,-30,41,0,0,-7},{0,0,0,0,-5,0,0,27,-2,7},{0,0,0,0,0,0,0,-2,29,-10}} ,x1[9],x0[9],e=0.00001,u,s2=0; int i,j,h,N=0,n=0; for(i=0;i<=8;i++){

x1[i]=x0[i]=1; }

cout<<\迭代过程数据如下：\ for(i=0;;i++){ for(j=0;j<=8;j++){ x0[j]=x1[j];

}

for(j=0;j<=8;j++){ for(h=0;h<=8;h++){ if(h!=j)s2+=(s[j][h]*x0[h]); }

x1[j]=(-s2+s[j][9])/s[j][j];

s2=0; }

for(j=0;j<=8;j++){ u=fabs(x1[j]-x0[j]); if(u

}

if(N==9)break;

cout<<\ cout<<\ cout<<\ cout<

n++;

}

cout<<\算法成功,迭代次数为：n=\ return 0;

}

实验结果如下：

迭代过程数据如下：

x1=0.258065 x2=1.71429 x3=-0.129032 x4=0.620253 x5=0.385965 x6=0.744681

10

计算方法与实习

x7=0.560976 x8=0.518519 x9=-0.275862

x1=0.475243 x2=0.955405 x3=-0.044158 x4=0.0988084 x5=0.112508 x6=0.555905 x7=0.374157 x8=0.3103 x9=-0.309068

x1=0.0961066 x2=0.971952 x3=-0.432686 x4=0.00192468 x5=-0.203384 x6=0.477386 x7=0.236028 x8=0.2572 x9=-0.323428

x1=0.0777173 x4=-0.168851 x7=0.178575

x1=-0.0781318 x4=-0.197411 x7=0.148488

x1=-0.101146 x4=-0.255381 x7=0.12877

x1=-0.15659 x4=-0.265276 x7=0.125632

x1=-0.167924 x4=-0.284944 x7=0.11854

x1=-0.185996 x4=-0.288647 x7=0.120638

x1=-0.190742 x4=-0.295165 x7=0.118057

x1=-0.19611 x4=-0.296589 x7=0.120185

x1=-0.197987 x4=-0.298657

x2=0.631942 x5=-0.268676 x8=0.197638 x2=0.597791 x5=-0.368833 x8=0.185276 x2=0.468877 x5=-0.388258 x8=0.166424 x2=0.449305 x5=-0.420949 x8=0.162763 x2=0.404005 x5=-0.427669 x8=0.156613 x2=0.395401 x5=-0.438207 x8=0.15535 x2=0.380363 x5=-0.4407 x8=0.153367 x2=0.376868 x5=-0.443948 x8=0.152899 x2=0.37219 x5=-0.444895 x3=-0.459135 x6=0.436266 x9=-0.32709 x3=-0.612936 x6=0.409319 x9=-0.331197 x3=-0.632064 x6=0.405031 x9=-0.33205 x3=-0.688191 x6=0.395338 x9=-0.33335 x3=-0.697065 x6=0.398205 x9=-0.333603 x3=-0.716561 x6=0.394678 x9=-0.334027 x3=-0.720254 x6=0.397587 x9=-0.334114 x3=-0.726722 x6=0.396311 x9=-0.334251 x3=-0.728196 x6=0.398153

11

计算方法与实习

x7=0.119252 x8=0.152287 x9=-0.334283

x1=-0.199355 x2=0.370816 x3=-0.730221 x4=-0.299207 x5=-0.445811 x6=0.397698 x7=0.1206 x8=0.15211 x9=-0.334325

x1=-0.200078 x2=0.369499 x3=-0.730798 x4=-0.299816 x5=-0.446171 x6=0.398695 x7=0.120267 x8=0.151937 x9=-0.334337

x1=-0.200308 x2=0.368969 x3=-0.731376 x4=-0.300027 x5=-0.446385 x6=0.398536 x7=0.120996 x8=0.151869 x9=-0.334349

x1=-0.200582 x2=0.368668 x3=-0.731598 x4=-0.300183 x5=-0.446521 x6=0.399034 x7=0.12088 x8=0.151829 x9=-0.334354

x1=-0.200548 x2=0.368466 x3=-0.731736 x4=-0.300263 x5=-0.446546 x6=0.39898 x7=0.121244 x8=0.151803 x9=-0.334357

x1=-0.20065 x2=0.368436 x3=-0.731821 x4=-0.30029 x5=-0.446597 x6=0.399216 x7=0.121205 x8=0.151798 x9=-0.334358

x1=-0.200586 x2=0.36836 x3=-0.731838 x4=-0.300321 x5=-0.446583 x6=0.399198 x7=0.121377 x8=0.151789 x9=-0.334359

算法成功,迭代次数为：n=19 Press any key to continue

二、高斯—赛德尔迭代法

思路分析

1.用一个二维数组存储方程组的各项系数，两个数组x0,x1分别存储x和x求,并判断是否符合精度要求,第二层和第三层进行高斯塞德尔迭代

3.输出方程的解

kk?1

2.用三个嵌套的循环语句，第一层把上一次循环的x0赋给x1用于判断是否符合精度要

C++语言编程源代码如下：

//高斯—赛德尔迭代法

12

计算方法与实习

#include #include #include using namespace std; int main(){

double

s[9][10]={{ 31,-13,0,0,0,-10,0,0,0,-15},{-13,35,-9,0,-11,0,0,0,0,27},{0,-9,31,-10,0,0,0,0,0,-23},{0,0,-10,79,-30,0,0,0,-9,0},{

0,0,0,-30,57,-7,0,-5,0,-20},{0,0,0,0,7,47,-30,0,0,12},{0,0,0,0,0,-30,41,0,0,-7},{0,0,0,0,-5,0,0,27,-2,7},{0,0,0,0,0,0,0,-2,29,-10}} ,x1[9],x0[9],e=0.00001,u,s2=0;

int i,j,h,N=0,n=0; for(i=0;i<=8;i++){ }

x1[i]=x0[i]=1;

cout<<\迭代过程数据如下：\for(i=0;;i++){

for(j=0;j<=8;j++){ x1[j]=x0[j];

}

for(j=0;j<=8;j++){ for(h=0;h<=8;h++){

if(h!=j)s2+=(s[j][h]*x0[h]);

}

x0[j]=(-s2+s[j][9])/s[j][j]; s2=0;

}

for(j=0;j<=8;j++){ }

if(N==9)break;

cout<<\cout<<\cout<<\cout<

u=fabs(x0[j]-x1[j]); if(u

13

计算方法与实习

}

}

cout<<\算法成功,迭代次数为：n=\return 0;

实验结果如下：

迭代过程数据如下：

x1=0.258065 x4=0.49346 x7=0.470126

x1=0.401988 x4=-0.0297297 x7=0.260482

x1=0.107902 x4=-0.195084 x7=0.177436

x1=-0.0529445 x4=-0.260029 x7=0.144144

x1=-0.138131 x4=-0.284608 x7=0.130703

x1=-0.174984 x4=-0.294117 x7=0.125252

x1=-0.190247 x4=-0.29785 x7=0.123034

x1=-0.196416 x4=-0.29933 x7=0.122128

x1=-0.19889 x4=-0.29992 x7=0.121758

x2=1.43871 x5=0.119365 x8=0.355438 x2=0.957825 x5=-0.227786 x8=0.19335 x2=0.661571 x5=-0.364219 x8=0.167256 x2=0.493434 x5=-0.414628 x8=0.157788 x2=0.419683 x5=-0.433983 x8=0.154155 x2=0.389018 x5=-0.441562 x8=0.152733 x2=0.376639 x5=-0.444566 x8=0.15217 x2=0.371691 x5=-0.445766 x8=0.151945 x2=0.369716 x5=-0.446249 x8=0.151854 x3=-0.00166493 x6=0.875839 x9=-0.320315 x3=-0.304677 x6=0.589325 x9=-0.331493 x3=-0.559457 x6=0.47583 x9=-0.333293 x3=-0.661611 x6=0.430329 x9=-0.333946 x3=-0.703972 x6=0.411961 x9=-0.334196 x3=-0.720804 x6=0.404511 x9=-0.334294 x3=-0.727465 x6=0.401479 x9=-0.334333 x3=-0.730106 x6=0.400242 x9=-0.334349 x3=-0.731157 x6=0.399736 x9=-0.334355

14

计算方法与实习

x1=-0.199882 x2=0.368925 x3=-0.731577 x4=-0.300157 x5=-0.446444 x6=0.399529 x7=0.121606 x8=0.151817 x9=-0.334357

x1=-0.20028 x2=0.368608 x3=-0.731745 x4=-0.300253 x5=-0.446523 x6=0.399444 x7=0.121544 x8=0.151803 x9=-0.334358

x1=-0.200441 x2=0.36848 x3=-0.731813 x4=-0.300292 x5=-0.446555 x6=0.399409 x7=0.121519 x8=0.151797 x9=-0.334359

x1=-0.200505 x2=0.368429 x3=-0.731841 x4=-0.300307 x5=-0.446568 x6=0.399394 x7=0.121508 x8=0.151794 x9=-0.334359

x1=-0.200532 x2=0.368408 x3=-0.731852 x4=-0.300314 x5=-0.446573 x6=0.399388 x7=0.121504 x8=0.151793 x9=-0.334359

x1=-0.200542 x2=0.368399 x3=-0.731856 x4=-0.300316 x5=-0.446576 x6=0.399386 x7=0.121502 x8=0.151793 x9=-0.334359

算法成功,迭代次数为：n=15 Press any key to continue

实验分析

从两种迭代法的实验结果来看，高斯塞德尔迭代法比雅可比迭代法的收敛速度更快，可以用更少的迭代次数来得到方程的解。

习题5、按下列数据

xi 0.30 0.42 0.50 0.58 0.66 0.72 yi 1.04403 1.08462 1.11803 1.15603 1.19817 1.23223 作五次插值，并求x1?0.46，x2?0.55，x3?0.60时的函数近似值。

思路分析

可用拉格朗日差值多项式作五次插值，用一个双层嵌套的循环语句来计算拉格朗日多项

15

计算方法与实习

nli??j?0,j?ix?xjxi?xjnLn(x)?式，内层计算

，外层计算

?li?0i(x)y（x）i

C++语言编程源代码如下：

#include #include #include

using namespace std; double chazhi(double x){

double

xi[6]={0.30,0.42,0.50,0.58,0.66,0.72},yi[6]={1.04403,1.08462,1.11803,1.15603,1.19817,1.23223},li[6]={1,1,1,1,1,1},Ln=0; int i,j;

for(i=0;i<=5;i++){ for(j=0;j<=5;j++){ }

if(j!=i)li[i]*=(x-xi[j])/(xi[i]-xi[j]); }

li[i]*=yi[i]; Ln+=li[i];

return Ln; } int main(){ cout<<\五次插值的结果如下：\ cout<<\时,y=\ }

cout<<\时,y=\cout<<\时,y=\return 0;

实验结果如下：

五次插值的结果如下：

x=0.46时,y=1.100724 x=0.55时,y=1.141271 x=0.60时,y=1.166194 Press any key to continue

习题6、试分别用抛物线

y?a?bx?cx2bx和指数曲线

y?ae拟合下列数据

16

计算方法与实习

xi yi 1 1.5 2 122.65 2.5 159.05 3 189.15 3.5 214.15 4 238.65 4.5 252.50 33.4 79.50 xi 5 5.5 6 6.5 7 7.5 8 yi 267.55 280.50 296.65 301.40 310.40 318.15 325.15

比较两个拟合函数的优劣。

一、抛物线拟合

思路分析

套用最小二乘法的解题思路，先求取正规方程组，再用高斯消元法解得正规方程组的解即为抛物线拟合的系数，a，b，c。

C++语言编程源代码如下：

#include #include #include using namespace std; int main(){

double

x[15]={1,1.5,2,2.5,3,3.5,4,4.5,5,5.5,6,6.5,7,7.5,8},y[15]={33.4,79.50,122.65,159.05,189.15,214.15,238.65,252.50,

267.55,280.50,296.65,301.40,310.40,318.15,325.15},s[3][4]={0},k,p,a,b,c; int i,j,n=14; for(i=0;i<=14;i++){ s[0][0]+=1;

s[1][0]=s[0][1]+=x[i];

s[1][1]=s[2][0]=s[0][2]+=x[i]*x[i]; s[1][2]=s[2][1]+=x[i]*x[i]*x[i]; s[2][2]+=x[i]*x[i]*x[i]*x[i]; s[0][3]+=y[i];

s[1][3]+=x[i]*y[i]; s[2][3]+=x[i]*x[i]*y[i]; }

for(i=0;i<=1;i++){

for(j=0;j<=3;j++){

if(j==0){ k=s[i+1][i]/s[i][i]; }

17

计算方法与实习

}

s[i+1][j]-=k*s[i][j]; if(i==0){ if(j==0){ p=s[i+2][i]/s[i][i]; } } }

s[i+2][j]-=p*s[i][j];

c=s[2][3]/s[2][2];

b=(s[1][3]-c*s[1][2])/s[1][1];

a=(s[0][3]-c*s[0][2]-b*s[0][1])/s[0][0];

cout<<\抛物线拟合的结果为：y=\ cout<<\与原数据相比，拟合结果如下：\

cout<

cout<}

return 0;

实验结果如下：

抛物线拟合的结果为：y=-45.3333+94.2302x-6.13161x^2 与原数据相比，拟合结果如下：

x y y' y-y’ 1 33.4 42.76529 9.365294 1.5 79.5 82.21588 2.715882 2 122.65 118.6007 4.049334 2.5 159.05 151.9196 7.130356 3 189.15 182.1728 6.977182 3.5 214.15 209.3602 4.789813 4 238.65 233.4818 5.168248 4.5 252.5 254.5375 2.037511 5 267.55 272.5275 4.977466 5.5 280.5 287.4516 6.951616 6 296.65 299.31 2.659961 6.5 301.4 308.1025 6.702502 7 310.4 313.8292 3.429237 7.5 318.15 316.4902 1.659832 8 325.15 316.0853 9.064706 Press any key to continue

18

计算方法与实习

二、指数曲线拟合

思路分析

上述数据用x[i],y[i]数组存储，对y[i]取自然对数，用最小二乘法进行一次函数拟合，转化为正规方程组后用高斯消元法求得b与a’,e即为a。

a'C++语言编程源代码如下：

#include #include #include using namespace std; int main(){

double

x[15]={1,1.5,2,2.5,3,3.5,4,4.5,5,5.5,6,6.5,7,7.5,8},y[15]={33.4,79.50,122.65,159.05,189.15,214.15,238.65,252.50, 267.55,280.50,296.65,301.40,310.40,318.15,325.15},s[2][3]={0},k,a,b; int i,j,n=14;

cout<for(i=0;i<=n;i++){

y[i]=log(y[i]); }

for(i=0;i<=n;i++){ s[0][0]+=1; s[0][1]=s[1][0]+=x[i];

s[1][1]+=x[i]*x[i]; s[0][2]+=y[i];

s[1][2]+=x[i]*y[i]; }

for(j=0;j<=2;j++){

if(j==0){k=s[1][0]/s[0][0];} s[1][j]-=k*s[0][j];

}

b=s[1][2]/s[1][1];

a=(s[0][2]-b*s[0][1])/s[0][0]; a=exp(a);

cout<<\抛物线拟合的结果为：y=\ cout<<\与原数据相比，拟合结果如下：\

cout<

19

计算方法与实习

}

return 0;

实验结果如下：

抛物线拟合的结果为：y=67.4026e^(0.23896x) 与原数据相比，拟合结果如下：

x y y' y-y' 1 33.4 85.5964 52.1964 1.5 79.5 96.4595 16.9595 2 122.65 108.701 13.9488 2.5 159.05 122.497 36.5534 3 189.15 138.043 51.1073 3.5 214.15 155.562 58.5883 4 238.65 175.304 63.3458 4.5 252.5 197.552 54.9479 5 267.55 222.624 44.9264 5.5 280.5 250.877 29.6232 6 296.65 282.716 13.9342 6.5 301.4 318.595 17.1954 7 310.4 359.028 48.6285 7.5 318.15 404.593 86.443 8 325.15 455.94 130.79 Press any key to continue

实验分析

从以上两种方法的拟合来看，抛物线拟合可以使误差变得更小，而指数拟合只能实现一部分曲线与数据的误差比较小，其他点的误差则较大，说明这组数据更适合用抛物线拟合来进行拟合计算而不适合用指数拟合来计算。

?习题7、用复化梯形和复化辛卜生公式计算积分

?I1(f)??201?cos2xdx和

I2(f)??4tanxx0dx，观察n为多少时，所得近似值具有6位有效数字？

?一、

I1(f)??201?cos2xdx

思路分析

20

计算方法与实习

用两层嵌套循环语句，内层套用复化梯形公式和复化辛卜生公式，外层判断是否满足精度要求，若不满足，增加n，若满足，跳出循环。

C++语言编程源代码如下：

#include #include #include using namespace std;

int main(){ double f1=0,f2=0,p=0,q=0,e=0.000001,a,b,h,x1,x2,x3;

int i,j=1,n1=1,n2=1;

cout<<\计算结果如下：\

cout<<\以下为用复化梯形公式的计算结果\

cout<

h=(b-a)/n1; for(j=1;j<=n1;j++){

x1=a+(j-1)*h; x2=a+j*h;

f1+=h/2*(sqrt(1+cos(x1)*cos(x1))+sqrt(1+cos(x2)*cos(x2))); if(fabs(f1-p)<0.000001)break;

}

cout<

p=f1; n1++; f1=0;

}

cout<

n=\ cout<<\以下为用复化辛卜生公式的计算结果\

cout<

b=3.14159/2; h=(b-a)/n2; for(j=1;j<=n2;j++){ x1=a+(j-1)*h;

21

计算方法与实习

x2=a+j*h;

x3=a+(j-0.5)*h;

f2+=h/6*(sqrt(1+cos(x1)*cos(x1))+4*sqrt(1+cos(x3)*cos(x3))+sqrt(1+cos(x2)*cos(x2))); if(fabs(f2-q)<0.000001)break;

}

cout<

}

cout<

n=\ }

return 0;

实验结果如下：

计算结果如下：

以下为用复化梯形公式的计算结果 f1 前后项差值 1.896117 1.896117 1.90997 0.01385312 1.910096 0.0001251545 1.910098 1.937996e-006

1.910098 4.239134e-008

用复化梯形公式计算得f1=1.910098,迭代次数为：n=5

以下为用复化辛卜生公式的计算结果 f1 前后项差值 1.914588 1.914588 1.91014 0.004448252

1.910098 4.165611e-005 1.910098 6.414982e-007

用复化辛卜生公式计算得f1=1.910098,迭代次数为：n=4 Press any key to continue

?二、I2(f)??4tanxx0dx

思路分析

用两层嵌套循环语句，内层套用复化梯形公式和复化辛卜生公式，外层判断是否满足精

22

计算方法与实习

tan(x)度要求，若不满足，增加n，若满足，跳出循环。需要注意的是，由于有意义，因此可以把初始值设为0.00000001,这对结果是没有影响的。

x在x=0处没

C++语言编程源代码如下：

#include #include #include using namespace std;

int main(){ double f1=0,f2=0,p=0,q=0,e=0.000001,a,b,h,x1,x2,x3;

int i,j=1,n1=1,n2=1;

cout<<\计算结果如下：\

cout<<\以下为用复化梯形公式的计算结果\

cout<

h=(b-a)/n1; for(j=1;j<=n1;j++){

x1=a+(j-1)*h; x2=a+j*h;

f1+=h/2*(tan(x1)/x1+tan(x2)/x2); if(fabs(f1-p)<0.000001)break;

}

cout<

p=f1; n1++; f1=0;

}

cout<

n=\ cout<<\以下为用复化辛卜生公式的计算结果\

cout<

b=3.14159/4; h=(b-a)/n2; for(j=1;j<=n2;j++){ x1=a+(j-1)*h;

23

计算方法与实习

}

x2=a+j*h;

x3=a+(j-0.5)*h;

f2+=h/6*(tan(x1)/x1+4*tan(x3)/x3+tan(x2)/x2); if(fabs(f2-q)<0.000001)break;

cout<

q=f2; n2++;

f2=0; } cout<

return 0;

实验结果如下：

计算结果如下：

以下为用复化梯形公式的计算结果 f1 前后项差值 0.892698 0.892698 0.860562 0.0321359 0.85419 0.00637239 0.851919 0.00227055 0.850861 0.00105857 0.850284 0.000577125 0.849935 0.000348715 0.849708 0.000226625 0.849553 0.000155509 0.849441 0.000111302 0.849359 8.23878e-005 0.849296 6.26836e-005 0.849248 4.87951e-005 0.849209 3.87252e-005 0.849178 3.12466e-005 0.849152 2.55764e-005 0.849131 2.11994e-005 0.849113 1.77669e-005 0.849098 1.50372e-005 0.849085 1.28394e-005 0.849074 1.10498e-005 0.849065 9.57802e-006

24

计算方法与实习

0.849056 8.35643e-006

0.849049 0.849042 0.849037 0.849032

7.33406e-006 6.47191e-006 5.73977e-006 5.11403e-006

0.849027 4.57605e-006 0.849023 4.11096e-006 0.849019 3.70682e-006 0.849016 3.35397e-006 0.849013 3.04452e-006 0.84901 2.772e-006 0.849007 2.53106e-006 0.849005 2.31725e-006 0.849003 2.12686e-006 0.849001 1.95677e-006 0.848999 1.80435e-006 0.848998 1.66735e-006 0.848996 1.54388e-006 0.848995 1.43231e-006 0.848993 1.33123e-006 0.848992 1.23945e-006

0.848991 1.15591e-006 0.84899 1.07972e-006 0.848989 1.01007e-006 0.848988 9.46293e-007

用复化梯形公式计算得f1=0.848988,迭代次数为：n=47

以下为用复化辛卜生公式的计算结果 f1 前后项差值 0.84985 0.84985 0.849038 0.000811965 0.848981 5.68013e-005 0.848971 1.02701e-005 0.848968 2.89011e-006

0.848967 1.0536e-006 0.848967 4.5551e-007

用复化辛卜生公式计算得f1=0.848967,迭代次数为：n=7 Press any key to continue

实验分析

从实验中可以发现，利用复化辛卜生公式进行迭代比用复化梯形公式迭代有更大的收敛速度，达到同样的精度要求，用复化辛卜生公式比复化梯形公式需要迭代的次数更少，当复化梯形公式需要迭代的次数很多时，复化辛卜生公式的优势尤其明显。实验结果中，计算

25

计算方法与实习

?I2(f)??4tanxx0dx时，用复化辛卜生公式与用复化梯形公式计算的结果略有出入，这可能

是由于复化梯形公式经过多次迭代之后误差累计造成的。

26

本文来源：https://www.bwwdw.com/article/0or7.html