数值分析自主上机题

2016级数值分析上机实践报告

机械工程学院 10 吴臻

标准题：迭代格式的比较

设方程f(x)=x 3

- 3x –1=0 有三个实根 x *1=1.8793 , x *2

=-0.34727 ,x *3

=-1.53209现采用下面三种不同计算格式，求 f(x)=0的根 x *1

或x *2

1、 x = 2

13x x + 2、 x = 31

3-x

3、 x =

3

13+x

数学原理：

简单迭代法是根据f(x)=0这个方程，对其进行等价变换化为x= ψ(x)并由该式产生逼近解x *的迭代数列{x k }，这就是简单迭代法的思想。对于同一f(x)=0可以等价变换出不同的x= ψ(x)而且它们各自的收敛性不同。 程序设计：

#include #include #include using namespace 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 main() {

double x,x0;

int N;//最大迭代次数 int k;

cout<<"输入初解:"; cin>>x0;

cout<<"输入最大迭代次数:"; cin>>N;

for(k=1;k<=N;k++) {

x=f(x0); if(fabs(x-x0)<0.0000001) { cout<<"迭代次数:"<

cout<<"已达到最大迭代次数:"<

实验结果：

程序运行结果讨论和分析：

对于第一种迭代格式，收敛区间[-8.2 -0.4]，在该收敛区间内迭代收敛于-1.53209，只能求得方程的一个根； 对于第二种迭代格式，收敛区间[-1.5 1.8]，在该收敛区间内迭代收敛于-0.34730，同样只能求得方程的一个根；

对于第三种迭代格式，收敛区间[-0.3 +∞)，在该收敛区间内迭代收敛于1.87937，只能求得方程的一个根； 由以上结果很容易发现，初值的选取对迭代敛散性有很大影响。以第一种迭代格式为例，当初值大于等于-0.3时，迭代格式发散；当初值小于等于-8.3时，迭代格式也发散；只有初值在-0.3和-8.3之间时，迭代格式才收敛于—1.53209。其他迭代格式也有这样的性质，即收敛于某个数值区间，超出这个区间迭代格式就是发散的，这就是所谓迭代格式的收敛性。 自主题：机械运动的数值仿真

背景: 描述物理学里把物体位置的变化叫机械运动。 如我们所知，力的作用效果有： 改变物体的运动状态 改变物体的形状 改变物的运动状态大多会引起物体的位置变化，引起机械运动。 改变物体的形状而不改变它的运动状态就叫是非机械运动中的一种。在工程实际与生产生活中我们常常需要对特定对象的机械运动进行研究分析得出其机械运动的规律，然后将该规律应用于对人们有益的方面。

实例：如下图1（左）所示，假设有一烟花火箭，其初始条件为零。将其放在地方然后点火，该烟花火箭的初始质量为0120m g =，其中粉末燃料占70g 。经过实验得知，燃料的持续时间为 2.0c t s =。燃料所产生的恒定推力为 5.2T N =。这也说明燃料的消耗率恒定。空气产生的阻力和烟花火箭的速度的平方成正比：2422, 4.010(/)

R kv k Ns m -= =⋅。这里，要求选择一种数值方法对其运动过程进行仿真并且其截断误差为4

()O h 或者更高。要求计算出该烟花火箭的最高高度，同时计算出从燃料消耗到该烟花火箭运动到最高点的时间延迟。

数学原理及数学模型：该实际问题要求其截断误差要求大于或等于4()

O h ，这就

使得较为简单的欧拉法，中点法不适合本例。龙格-库塔法以其优异的数值特性成为解决本问题的首选。

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

加速上升过程：

2

11

111100010001203512035(0)0,(0)0dv T kv g dt t t dh v dt v h ⎧=--⎪

--⎪⎪⎨= ⎪⎪

= = ⎪⎩

（1）

减速上升过程：

2

222

21211000120352(0)(2)(0)(2)dh v dt dv kv g dt

h h v v ⎧= ⎪⎪

⎪=--⎨

-⨯⎪= ⎪⎪

=

⎩ （2）

式中，h 为上升的高度，v 为上升过程的速度，g 为重力加速度。1(2)h 表示加速上升过程的最终高度，1(2)v 表示加速上升过程的最终速度。

使用龙格-库塔法求解如上的微分方程组。该系统的数值仿真结果如图2所示。图2（左）是该烟花火箭的上升过程高度的数值仿真，图2（右）是其上升过程速度的数值仿真。同时亦可以得到烟花火箭上升的最大高度和问题中所需的时间延迟：

从如上的分析和仿真可知，使用数值方法进行机械运动的数值仿真，可以简化本身复杂变化的物理运动过程。同时注意到，即使改变机械运动的初始条件或者系统的某些特征，其对应的仿真只需要进行简单的调整。而且其仿真精度也可以得以预见。

图2烟花火箭的上升过程高度（左）和速度（右）的数值仿真 编程过程： 加速上升过程

f[x_,y_]:=5200/(120-35t)-0.4x^2/(120-35t)-10; g[x_,y_]:=x; {x,y}={0,0}; h=0.1;