数值分析自主上机题

1、2016级数值分析上机实践报告机械工程学院2016021910吴臻标准题：迭代格式的比较3设方程 f(x)=x - 3x -=0 有三个实根 x1 =1.8793 , x2 =-0.34727 ,x3=-1.53209现采用下面三种不同计算格式，求f(x)=0的根x1或x23x 121、x = xX312、x =33、x = 33x 1数学原理：简单迭代法是根据f(x)=0这个方程，对其进行等价变换化为x=书(x)并由该 式产生逼近解x*的迭代数列xk，这就是简单迭代法的思想。对于同一 f(x)=0可 以等价变换出不同的x=书(x)而且它们各自的收敛性不同。程序设计：#in clude<

2、iostream>#in clude<cmath>#in clude<cstdlib> using n amespace std;double f(double i)/外调函数f(x),每次更新新的函数/以第一种迭代方式为例子double k,m,sum;k=3*i+1;m=pow(i,2.0);sum=k/m; return sum;int mai n()double x,x0;int N;最大迭代次数int k;cout«"输入初解:"cin> >x0;coutvv"输入最大迭代次数:"cin>

3、;>N;for(k=1;k<=N;k+)x=f(xO);if(fabs(x-x0)<0.0000001)cout«"迭代次数:"vvkvvendl; cout«"输出得到的解:"<<x<<endl; system("pause");return 0;else x0=x;coutvv"已达到最大迭代次数:"<<N<<endl; cout«"输出得到的解:"<<x<<endl; s

4、ystem("pause");return 0;实验结果：揄入最大迭代枕数 陆代枚数输出得到的解= -1-53209 陳袴任意毎靈续愉人初解：-乩* 输入嘉大迭代枕数迓 已达到最犬迭代次数応回 输出得到的解:-输入初解= -9-3 输入最大迭代枚数 迭代初数泊输出得至|的解=-0-3472?6 请按任童键继缮输人初解=15输入最大迭代次数;5回 迭代枕数:输出得到的M = -0-347296 请按任盍键継缭-输入初=输入最大迭代初数;別 迭代枕数:解轴出得到的解= 1-87929 请無任意键输入初i?=2-0驚入最夭迭祀枕数汚0 逹代枕数七输岀得1的解= 1-87939 谙

5、转任竟键继续.程序运行结果讨论和分析：对于第一种迭代格式，收敛区间-8.2 -0.4，在该收敛区间内迭代收敛于-1.53209,只能求得方程的一个根；对于第二种迭代格式，收敛区间-1.51.8，在该收敛区间内迭代收敛于-0.34730,同样只能求得方程的一个根；对于第三种迭代格式，收敛区间-0.3 + g），在该收敛区间内迭代收敛于1.87937,只能求得方程的一个根；由以上结果很容易发现，初值的选取对迭代敛散性有很大影响。以第一种迭代格式为例，当初值大于等于-0.3时，迭代格式发散；当初值小于等于-8.3时，迭代 格式也发散；只有初值在-0.3和-8.3之间时，迭代格式才收敛于一1.5320

6、9。其他 迭代格式也有这样的性质，即收敛于某个数值区间，超出这个区间迭代格式就是 发散的，这就是所谓迭代格式的收敛性。自主题：机械运动的数值仿真背景：描述物理学里把物体位置的变化叫机械运动。如我们所知，力的作用效果有：改变物体的运动状态 改变物体的形状 改变物的运动状态大多会引起物 体的位置变化，引起机械运动。改变物体的形状而不改变它的运动状态就叫是非机械运动中的一种。在工程实际与生产生活中我们常常需要对特定对象的机械 运动进行研究分析得出其机械运动的规律，然后将该规律应用于对人们有益的方面。实例：如下图1 （左）所示，假设有一烟花火箭，其初始条件为零。将其放在地 方然后点火，该烟花火箭的初始

7、质量为 m° 120g，其中粉末燃料占70g。经过实 验得知，燃料的持续时间为tc 2.0s。燃料所产生的恒定推力为T 5.2N。这也说 明燃料的消耗率恒定。空气产生的阻力和烟花火箭的速度的平方成正比：2422R kv,k 4.010 （Ns /m）。这里，要求选择一种数值方法对其运动过程进行仿真4并且其截断误差为°（h）或者更高。要求计算出该烟花火箭的最高高度，同时计 算出从燃料消耗到该烟花火箭运动到最高点的时间延迟。4数学原理及数学模型：该实际问题要求其截断误差要求大于或等于O（h ），这就使得较为简单的欧拉法，中点法不适合本例。龙格-库塔法以其优异的数值特性成为解决本

8、问题的首选。图1烟花火箭的机械运动的数值仿真（左）和结构力学问题的数值求解（右）很显然，该问题属于变质量的运动学问题， 在该运动过程中，其前两秒是在驱动 力和阻力的共同作用下加速上升的， 而后的时间内，该烟花火箭是在空气的阻力 下减速上升的，同时注意到空气的阻力和速度的平方成正比。 为了对该运动过程 进行数值仿真，那么必须建立相应的微分方程组。 分析该运动过程可知，应该将 该运动过程分为两部分：加速上升过程和减速上升过程。从而得到相应的微分方 程组。加速上升过程：dv11000T1000kv2dt120 35tg120 35tdgvdt1Vi(O)0, h1(0)0减速上升过程:dh2V2图2

9、烟花火箭的上升过程高度（左）和速度（右）的数值仿真dtdv21000kv2gdt 120 35 2h2（0） hi（2）V2（0） Vi（2）（2）式中，h为上升的高度，v为上升过程的速度，g为重力加速度。h（2）表示加速 上升过程的最终高度，v表示加速上升过程的最终速度。使用龙格-库塔法求解如上的微分方程组。该系统的数值仿真结果如图2所示。图2 （左）是该烟花火箭的上升过程高度的数值仿真，图2 （右）是其上升过程速度的数值仿真。同时亦可以得到烟花火箭上升的最大高度和问题中所需的时间 延迟：hmax198.462mt delay6.185（ S）从如上的分析和仿真可知，使用数值方法进行机械运动

10、的数值仿真，可以简化本 身复杂变化的物理运动过程。同时注意到，即使改变机械运动的初始条件或者系 统的某些特征，其对应的仿真只需要进行简单的调整。 而且其仿真精度也可以得 以预见。编程过程：加速上升过程fx_,y_：=5200/(120-35t)-0.4xA2/(120-35t)-10;gx_,y_：二x；x,y=0,0;h=0.1;t=0.1;xx=Table0,i,1,70;yy=Table0,i,1,70;tt=Table0,i,1,70;Doa=fx,y;xa=x+h (a+fx+h,y+h*a)/2; b=gx,y;ya=y+h (b+gx+h,y+h*b)/2;Prin tk,&qu

11、ot;",t,"",xa,"",ya;t,x,y,xxk,yyk,ttk=t+h,xa,ya,xa,ya,t+h,k,1,20£7B9101112131-4IS lfiIT1320O1 , 20.30,4 aS0.60.7a.eo甲9i.1.1le 21.31- 41-S1. fi1.71. S1_ *2.:a 4 3 =.2 . 00E7.060920.35352mml0®74414,esse2.1515713.64143.621837*7712631314510 474335,"4EE13.E工口740.260

12、317190344-843121.221349 47S-3ZSIOCS4,139230.32ss t eoss3 6 亠 05223.4石号$41.9684S 313T2.591255.130377.01113945fil .2 90*77CK100C65< 397 77EZ34Ttt0.6,1.1,1.6,2.1,2.6,3.1,3.6,4.1,4.6,5.1,5.6,6.1,6.6,7.1,7.6xx17.3174,37.9564,61.2759,85.0712,104.674,113.404,109.842,75.5189,8.73997,-59.2801,-96.8113,-11

13、1.485,-117.825,-121.341,-123.752yy 0.125,8.90869,28.0119,58.7749,101.435,153.898,210.725,265.771,303.655,308.15,278.635,230.354,174.737,115.949,55.4041ListPlotTabletti,xxi,i,1,20,PlotStylePoi ntSize0.02ListPlotTable tti,xxi,i,1,20,PlotStylePointSize0.02180L-60一-40一20*1 1 1 1 1 1 1i1aH1I1111iHi0.51.0

14、1.52.0ListPlotTabletti,xxi,i,1,20,PlotJoi nedTrue减速上升过程fx_,y_：=-0.4xA2心20-35*2)-10;gx_,y_：二x；x,y=85.3977,78.2347;h=0.1;t=2.1;Doa=fx,y;xa=x+h (a+fx+h,y+h*a)/2; b=gx,y;ya=y+h (b+gx+h,y+h*b)/2;Prin tk,"",t,"",xa,"",ya;t,x,y,xx20+k,yy20+k,tt20+k=t+h,xa,ya,xa,ya,t+h,k,1,502.

15、2 72. 13494.百ML2.3 67-3S95L01+3062.4 62.751106.-&5E丐58.5956X14+93L54.S44412C.7952-751.433612.2552,843-3132131*4332 汽45.4421.眦，迪3,讹3.14C-31S4145.1023,236,014"145,133.3S5.B5Efi152.9453 + 433.E24215,5353,531.9062159,9233.6MoosesK3-H93*72S,362S146-1333.32-7161-62.3743 + 3£5+li36ii,esi4,23.

16、f55B174,174.122-15717-E.5334.22C,79L4le.7c2313.4433】麹.S4744IS.13994.5W,e?S2154.£1S123456189I口II1213U1$IE1"1 =空20212Z232425107-377189.327Z23293031325334箔3537IS3940鼻42434445箕47454350生首4$i.l14.44534.B13.250712.133511.0197.9217191,S79192,9365.2e .a 421.5193.9335.37-77089194-S7Z5*4.729351954655,333s.e4空适19孔迹S.73.«4SS?157-379&.£L3-.7435.9l.«32L6198.018乞(L 629£»4ISE.iefi5.2fi.36.45.7-l,37O5fi-2.37136-3*37627-4,35513-S*40016.42307-7-45556-B.4994419C.254196,222138-09197«S5£137.525197.09