数值分析上机题

1、数值分析上机题习题117.上机题舍入误差与有效数N设Sn1,其精确值为(1)编制按从大到小的顺序Sn(2)编制按从小到大的顺序Sn12211N14L32112(N1)1-,计算Sn的通用程序.11,计算Sn的通用程序.(3)按两种顺序分别计算002,S,S106,并指出有效位数.编制程序时用单精度通过本上机题你明白了什么?按从大到小的顺序计算Sn的通用程序为:#include<iostream.h>floatsum(floatN)floatj,s,sum=0;for(j=2;j<=N;j+)s=1/(j*j-1);sum+=s;)returnsum;按从小到大的顺序计算Sn的

2、通用程序为:#include<iostream.h>floatsum(floatN)floatj,s,sum=0;for(j=N;j>=2;j-)s=1/(j*j-1);sum+=s;)returnsum;)从大到小的顺序的值从小到大的顺序的值精确值后效位数从大到小从小到大S1020.7400490.740050.74004965S1040.7498520.74990.749944S1060.7498520.7499990.74999936)通过本上机题,看出按两种不同的顺序计算的结果是不相同的,按从大到小的顺序计算的值与精确值有较大的误差,计算得到的结果的有效位数少.致计算

3、结果的精度有所降低,的算法,其结果的相对误差较小.而按从小到大的顺序计算的值与精确值吻合.从大到小的顺序计算机在进行数值计算时会出现“大数吃小数的现象,导我们在计算机中进行同号数的加法时,采用绝对值较小者先加20.(上机题)Newton迭代法(1)给定初值Xo及容许误差,编制Newton法解方程f(x)0根的通用程序.给定方程f(x)x3/3x0,易知其有三个根X1J3,x20,x3J3.1.由Newton方法的局部收敛性可知存在0,当比(,)时,Newton迭代序列收敛于根x2.试确定尽可能大的2.试取假设干初始值,观察当x0(1),(1,),(,),(,1),(1,)时floatx0,x1

4、,a;intk=0;cout<<"请输入初值x0:"cin>>x0;doa=-f(x0)/df(x0);x1=x0+a;k+;x0=x1;while(fabs(a)>eps);cout<<k<<'t'<<x0;输出迭代的次数和根值Newton序列是否收敛以及收敛于哪一个根.(3)通过本上机题,你明白了什么？解：(1)编制的通用程序：#include<iostream.h>#include<math.h>#defineeps0.000001/给定容许误差floatf(fl

5、oatx)定义函数f(x)floatf;f=x*x*x/3-x;/f(x)的表达式；return(f);floatdf(floatx)定义函数df(x),计算f(x)的导函数floatdf;df=x*x-1;/f(x)导函数的表达式；return(df);voidmain(void)(2)计算迭代序列收敛于根x2的尽可能大的#include<iostream.h>#include<math.h>voiddelay(intn)/定义延时函数for(n=10000;n>0;n-);#defineeps0.000001floatf(floatx)定义函数f(x)floa

6、tf;f=x*x*x/3-x;/f(x)的表达式；的函数为：return(f);floatdf(floatx)定义函数df(x),计算f(x)的导函数floatdf;df=x*x-1;/f(x)导函数的表达式；return(df);=0.774,即在区间-0.774,0.774内迭代序列*对于不同得初始值U敛于不同的根,Xo在-8,-1内收敛于Xi,在-0.774,0.774内收敛于x2,在1,+°°内收敛于x3,但在内0.774,1和一1,0.774均可能收敛于*一一.一、一一.X1和x3.X,x2,x3分别为万程的精确解.分析：对于不同的初值,迭代序列会收敛于不同的根,

7、所以在某个区间内求根对于初值的选取有很大的关系.产生上述结果的原因是区间不满足大范围收敛的条件.intjudgement(floatz)(intcount=5;floatx0,x1,type,type1;x0=z;while(count->0)(x1=x0-f(x0)/df(x0);type=fabs(x1);type1=fabs(x1-x0);调试值用cout<<"count="<<count<<'t'<<"type="<<type<<'t'&

8、lt;<"type1="<<type1<<'n'if(fabs(x1-x0)<eps)return1;x0=x1;delay(30000);调试值用当步长为0.001时,程序计算出的8的为8收敛于0.return0;voidmain(void)(floatdelta=0;intflag=1;while(flag=1)(cout<<方程的根为："<<'n'delta+=eps;flag=judgement(delta);cout<<"输出方程根收敛的区间

9、彳K:n"cout<<delta-eps;输出收敛的区间值39.(上机题)列主元Gauss消去法对于某电路的分析,归结为求解线性方程组RI=Vo其中,R=31-13000-10000-1335-90-1100000-931-100000000-1079-30000-9000-3057-70-500000-747-300000000-3041000000-50027-2000-9000-229V=(-15,27,-23,0,-20,12,-7,7,10)#include<iostream.h>cout<<'n'#include<

10、math.h>输出数组voidmain(void)cout<<<<'n'/进行第一行和第一inti,j,n,k,q;/floata1011,s10,s110;intt=1;cout<<"请输入n的值:"for(i=1;i<=n;i+)cin>>n;si=ai1;cout<<"输入数组a:"<<endl;floatmax=fabs(s1);for(i=1;i<=n;i+)for(i=2;i<=n;i+)for(j=1;j<=(n+1);j+

11、)if(fabs(si)>max)cin>>aij;给矩阵a赋值for(i=1;i<=n;i+)max=fabs(si);t=i;for(j=1;j<=(n+1);j+)cout<<aij<<'t'for(j=1;j<=(n+1);j+)(1)编制解n阶线性方程组Ax=b的列主元三角分解法的通用程序；(2)用所编制的程序解线性方程组RI=V,并打印出解向量,保存五位有效数;(3)本编程之中,你提升了哪些编程水平？程序为：aU元素的求取'<<'n'/进行第k步分解'/a加尸b；

12、)/进行第一列主元互换for(i=2;i<=n;i+)ai1=ai1/max;/第一列除以a11for(i=1;i<=n;i+)(for(j=1;j<=(n+1);j+)cout<<aij<<'t'cout<<'n')/输出进行第一步变换的数组acout<<"for(k=2;k<=n;k+)(for(i=k;i<=n;i+)(floatsum=0;for(q=1;q<k;q+)sum+=aiq*aqk;s1i=aik卜sum;)intl=k;floatm=fabs(s1

13、k);for(i=k;i<=n;i+)比拟第k步分解的第k列值的大小(if(fabs(s1i)>m)(m=fabs(s1i);l=i;返回行值)for(j=1;j<=n+1;j+)交换两行元素(floats2=akj;akj=alj;结果：方程的解为：x1=-0.28923,x2=0.34544,x3=-0.71281,x4=-0.22061,x5=-0.43040,x6=alj=s2;)for(j=k;j<=n+1;j+)算出第k行行元素的值(floatsum1=0;for(q=1;q<k;q+)sum1+=akq*aqj;akj=akj-sum1;)for(i

14、=k+1;i<=n;i+)/算出第k列列元素的值(floatsum2=0;for(q=1;q<k;q+)sum2+=aiq*aqk;aik=(aik-sum2)/(akk);)第k步分解结束for(i=1;i<=n;i+)(for(j=1;j<=(n+1);j+)cout<<aij<<'t'cout<<'n')输出改变后的数组/输出解/floatx10;for(i=n-1;i>=1;i-)(xn=ann+1/ann;floatsum3=0;for(j=i+1;j<=n;j+)sum3+=ai

15、j*xj;xi=(ain+1-sum3)/aii;)回代过程for(i=1;i<=n;i+)(cout<<'x'<<i<<'='<<xi<<endl;/输出解向量)0.15431,x7=-0.057823,x8=0.20215,x9=0.29023.分析：我感觉是提升了查错误点的水平和编循环语句的水平,即利用很规整的迭代公式进行编程.另外列主元三角分解法的阶梯步骤有了更深的了解37.(上机题)3次样条插值函数(1)编制求第一型3次样条插值函数的通用程序(2)汽车曲线型值点的数据如下：xi01234

16、5678910Yi2.513.304.044.705.225.545.785.405.575.705.80端点条件为y0=0.8,y10=0.2.用所编制程序求车门的3次样条插值函数S(x),并打印出S(i+0.5)(i=0,1,9)o解：通用程序：#include<iostream.h>voidmain(void)floatx11;/存放数组xjfloaty11;/存放数组yjfloath11;/存放数组hjfloatu11;/存放数组ujfloatv11;/存放数组vjfloatd11;/存放数组djfloatM11;/存放数组Mjfloatb11;/存放数组bjfloatt1

17、1,l11,yy11,s4,aa1,aa2,aa3,aa4;floats110;inti,j,n;floatxx;/x为区间值将初值初始化cout<<"请输入n的值:n"cin>>n;cout<<"输入数组x:n"for(i=0;i<=n;i+)cin>>xi;cout<<"输入数组y:n"for(i=0;i<=n;i+)cin>>yi;输入端点值floatdf2;cout<<"输入两个端点值：n"for(i=0;i&l

18、t;2;i+)cin>>dfi;利用书本上的算法求出所需要的值求出hj的值for(j=0;j<=n-1;j+)hj=xj+1-xj;cout<<'h'<<''<<j<<''<<'='<<hj<<'t'cout<<endl;求出uj和vj的初值v0=1;un=1;for(j=1;j<=n-1;j+)uj=hj-1/(hj-1+hj);vj=hj/(hj-1+hj);求出dj的值for(j=1;j&l

19、t;n;j+)dj=6*(yj+1-yj)/hj-(yj-yj-1)/hj-1)/(hj+hj-1);d0=6*(y1-y0)/h0-df0)/h0;dn=6*(df1-(yn-yn-1)/hn-1)/hn-1;for(j=1;j<=n;j+)cout<<'u'<<''<<j<<T<<'='<<uj<<'t'cout<<endl;for(j=0;j<n;j+)cout<<'v'<<&#

20、39;'<<j<<T<<'='<<vj<<'t'cout<<endl;for(j=0;j<=n;j+)(cout<<'d'<<''<<j<<''<<'='<<dj<<'t')cout<<endl;利用书本上的追赶法求解方程组for(i=0;i<=n;i+)bi=2;cout<<endl;

21、t0=b0;yy0=d0;消元过程for(i=1;i<=n;i+)li=ui/ti-1;ti=bi-li*vi-1;yyi=di-li*yyi-1;回代过程Mn=yyn/tn;for(i=n-1;i>=0;i-)Mi=(yyi卜vi*Mi+1)/ti;(2)编制的程序求车门的3次样条插值函数将Mj的值输出for(i=0;i<=n;i+)cout<<'M'<<''<<i<<''<<'='<<Mi<<endl;输出插值多项式的系数fo

22、r(j=0;j<n;j+)s0=yj;s1=(yj+1-yj)/hj-(Mj/3+Mj+1/6)*hj;s2=Mj/2;s3=(Mj+1-Mj)/(6*hj);cout<<"当x的值在区间"<<'x'<<''<<j<<''<<"至u"<<'x'<<',<<(j+1)<<''<<"时,输出插值多项式的系数：n"fo

23、r(intk=0;k<4;k+)cout<<'s'<<''<<k<<''<<'='<<sk<<'t'cout<<endl;S(x):x属于区间0,1时;S(x)=2.51+0.8(x)-0.0014861(x)(x)-0.00851395(x)(x)(x)x属于区间1,2时;S(x)=3.3+0.771486(x-1)-0.027028(x-1)(x-1)-0.00445799(x-1)(x-1)(x-1)x属于区

24、间2,3时；S(x)=4.04+0.704056(x-2)-0.0404019(x-2)(x-2)-0.0036543(x-2)(x-2)(x-2)x属于区间3,4时；S(x)=4.7+0.612289(x-3)-0.0513648(x-3)(x-3)-0.0409245(x-3)(x-3)(x-3)x属于区间4,5时；S(x)=5.22+0.386786(x-4)-0.174138(x-4)(x-4)+0.107352(x-4)(x-4)(x-4)x属于区间5,6时；S(x)=5.54+0.360567(x-5)+0.147919(x-5)(x-5)-0.268485(x-5)(x-5)(x

25、-5)x属于区间6,7时；S(x)=5.78-0.149051(x-6)-0.657537(x-6)(x-6)+0.426588(x-6)(x-6)(x-6)x属于区间7,8时；S(x)=5.4-0.184361(x-7)+0.622227(x-7)(x-7)-0.267865(x-7)(x-7)(x-7)x属于区间8,9时；S(x)=5.57+0.256496(x-8)-0.181369(x-8)(x-8)+0.0548728(x-8)(x-8)(x-8)x属于区间9,10时；S(x)=5.7+0.058376(x-9)-0.0167508(x-9)(x-9)+0.0583752(x-9)(

26、x-9)(x-9)S(0.5)=2.90856S(3.5)=4.98819S(6.5)=5.59441S(9.5)=5.7323S(1.5)=3.67843S(4.5)=5.38328S(7.5)=5.42989S(2.5)=4.38147S(5.5)=5.7237S(8.5)=5.6597623.(上机题)常微分方程初值问题数值解(1)编制RK4方法的通用程序；(2)编制AB4方法的通用程序(由RK4提供初值)；(3)编制AB4-AM4预测校正方法的通用程序(由RK4提供初值)；(4)编制带改良的AB4-AM4预测校正方法的通用程序(由RK4提供初值)(5)对于初值问题Hx2y2(0x1.5

27、)ty(0)3取步长h0.1,应用(1)(4)中的四种方法进行计算,并将计算结果和精确解3y(x)3/(1x)作比拟；(6)通过本上机题,你能得到哪些结论？解：程序为：#include<iostream.h>#include<fstream.h>#include<stdlib.h>#include<math.h>ofstreamoutfile("data.txt");/此处定义函数f(x,y)的表达式用户可以自己设定所需要求得函数表达式doublef1(doublex,doubley)doublef1;f1=(-1)*x*x*

28、y*y;returnf1;/此处定义求函数精确解的函数表达式doublef2(doublex)doublef2;f2=3/(1+x*x*x);returnf2;/此处为精确求函数解的通用程序voidaccurate(doublea,doubleb,doubleh)doublex100,accurate100;x0=a;inti=0;outfile<<"输出函数准确值的程序结果:n"doxi=x0+i*h;accuratei=f2(xi);outfile<<"accurate"<<i<<"=&quo

29、t;<<accuratei<<'n'i+;while(i<(b-a)/h+1);此处为经典Runge-Kutta公式的通用程序voidRK4(doublea,doubleb,doubleh,doublec)inti=0;doublek1,k2,k3,k4;doublex100,y100;y0=c;x0=a;outfile<<"输出经典Runge-Kutta公式的程序结果:n"do(xi=x0+i*h;k1=f1(xi,yi);k2=f1(xi+h/2),(yi+h*k1/2);k3=f1(xi+h/2),(yi+h*

30、k2/2);k4=f1(xi+h),(yi+h*k3);yi+1=yi+h*(k1+2*k2+2*k3+k4)/6;outfile<<"y"<<""<<i<<"="<<yi<<'n'i+;while(i<(b-a)/h+1);/此处为4阶Adams显式方法的通用程序voidAB4(doublea,doubleb,doubleh,doublec)(doublex100,y100,y1100;doublek1,k2,k3,k4;y0=c;x0=a

31、;outfile<<"输出4阶Adams显式方法的程序结果:n"for(inti=0;i<=2;i+)(xi=x0+i*h;k1=f1(xi,yi);k2=f1(xi+h/2),(yi+h*k1/2);k3=f1(xi+h/2),(yi+h*k2/2);k4=f1(xi+h),(yi+h*k3);yi+1=yi+h*(k1+2*k2+2*k3+k4)/6;intj=3;y10=y0;y11=y1;y12=y2;y13=y3;do(xj=x0+j*h;y1j+1=y1j+(55*f1(xj,y1j)-59*f1(xj-1,y1j-1)+37*f1(xj-2,

32、y1j-2)-9*f1(xj-3,y1j-3)*h/24;n'j+;while(j<(b-a)/h+1);/主函数voidmain(void)(doublea,b,h,c;cout<<"输入上下区间、步长和初始值:n"cin>>a>>b>>h>>c;accurate(a,b,h);RK4(a,b,h,c);AB4(a,b,h,c);结果为：y(xi)和精确值和由经典Runge-Kutta公式得出的结果列在下面的表格中,以及精确值数值解的误差：iXiyiy(xi)|y(xi)-yi|0033010.12

33、.9972.9971.87138e-00720.22.976192.976193.91665e-00730.32.921132.921137.58342e-00740.42.819552.819551.61101e-00650.52.666662.666673.17735e-00660.62.46712.467115.00551e-00670.72.23382.23385.77233e-00680.81.984121.984134.12954e-00690.91.735111.735111.15554e-007101.01.500011.55.80668e-006111.11.287011.2

34、871.13075e-005121.21.099721.099711.54242e-005131.30.9383970.938381.77272e-005141.40.80130.8012821.83754e-005151.50.6857320.6857141.78e-005由AB4方法得出的结果为:Y10=3y15=2.66467y110=1.50219y115=0.685335y11=2.997y16=2.4652y111=1.28876y12=2.97619y17=2.23308y112=1.10072y13=2.92113y18=1.98495y113=0.93871y14=2.818

35、39y19=1.73704y114=0.801135通过本上机题我明白了各种求微分方程的数值方法,经典Runge-Kutta公式,AB4方法以及Runge-Kutta公式的精度,四AB4-AM4预测校正方法求解公式的精度是不同的.其中经典阶Adams显式方法(AB4)具有4阶精度.10.抛物线方程Crank-Nicolson格式(1)编制用Crank-Nicolson格式求抛物线方程22u/t-au/x=f(x,t)(0<x<l,0tT)<u(x,0)=(x)(0x1)二u(0,t)=x,u(1,t)=t(0<t1)数值解的通用程序.(2)就a=1,f(x,t)=0,x

36、=exp(x),t=exp(t),t=exp(1+t),M=40,N=40,输入点(0.2,1.0),(0.4,1.0),(0.6,1.0),(.8,1.0)4点处u(x,t)的近似值.(3)所给方程的精确解为u(x,t)=exp(x+t),将步长反复二分,从(0.2,1.0),(0.4,1.0),(0.6,1.0),(0.8,1.0)4点处精确解与数值解的误差观察当步长缩小一半时,误差以什么规律缩小.算法:此题选择空间步长1.h=0.025,时间步长I=0.025置初值；2. 通过点(i-1,k),(i-1,k+1),(i,k),(i+1,k),(i+1,k+1)的值；3. 输出结果求解点(I,k+1)原程序：#include<iostream>#include<math.h>floath=0.025,k=0.025;intm=40;intn=40;floaty4040,r=a*k/(h*h);voidInput()inti,j;cout<<"LoadingInputData."<<endl;for(i=0;i<m;i+)for(j=0;j<n;j+)if(i=j)a皿=1+r;for(j=0;j&