清华大学谢金星数学实验作业

1、实验二 常微分方程的数值解土木系 结23 李会平【实验目的】1、 掌握用matlab解数值微分方程2、 了解龙格-库塔方法的基本原理3、 用这些手段解决一些实际的问题【实验内容】4-5 核废料问题l 首先列出问题的运动学方程，由牛顿第二定律，md2s/dt2=G-F-f其中m=G/g, f=kv,由于我们熟悉的单位是公制的，所以在定义函数的时候进行了单位转换。以下是函数的定义：这其实是一个关于s的二阶常微分方程，需要定义两个变量x(1),x(2)将其化为一阶微分方程组，此处x(1),x(2)实际上分别代表速度和深度，相应的微分方程如代码中所示：dx=(G-F-k*x(1)/m;x(1)（如课堂

2、提醒，中间应该是分号，这点容易出错。）l 接写来进行m文件的命令编写，如下所示：（1）执行pause之前的代码，得到的结果如下图所示：该图代表了假设水无限深的情况下，物体速度的变化情况，容易看出速度是有一个上限的，这也符合直观的感受，因为速度不可能无限大，否则阻力无限大，物体将无法继续运动。（2）红线为题给的速度阈值，可见初步判断速度是有可能超过该阈值的，需要进一步判断沉底的时候是否超过该速度，于是需要做出s和v的关系曲线，即x(2),x(1)的关系曲线，如下所示：从图像中看出，s趋于时候，v也趋于极限值，题给的smax=300*0.3048m,在编程的时候，如果要绘制该水平直线，需要将其转化

3、为数组才能进行绘图，否则会出现错误。（3）但是从这幅图像中，由于s非常大，红线触底，无法直接看出s,v的相交情况，需要进一步限定坐标轴范围进行细化，如下所示：从图中明显读出，当s达到smax时候，其相应的vvmax,所以从题给的情况看，工程师们的说法是更有道理的。在这个例题中，加深了对微分方程数值解的理解，同时在用matlab绘图时也有了更多的心得。绘图时候尤其要注意为数组和数组的对应！另外在通过查询相关知识知道，如果要绘制直接过（a,b）和（c,d）两点的直线，可以直接用plot(a,c,b,d)命令快速实现，方便了一些作图过程。4-6小船渡河问题l 首先给出问题的理论解法：建立小船运行的数

4、学模型，用极坐标，以B为原点进行计算。取定BA为极轴。小船位置设为C，可以得到方程组： ddt=v1sin-v2ddt=v1cos/ 消去td=v1sin-v2v1cosd积分得到：ln=ln1cos12kln1+sin1-sin+Ck=v1v2可以解得：=Lcos(1+sin1-sin)-12kl 接下来用程序进行数值解和理论解答的分析：（1）用matlab运算时，采用普通坐标更加方便，所以以A为原点建立直角坐标系，容易得到 dxdt=v1-v2xx2+d-y2dydt=v2(d-y)x2+d-y2 据此定义函数m文件如下所示：（2）相应的主程序为：其中pause前的代码表示进行v1=1，v

5、2=2时候的数值运算求解绘图。pause后的代码表示分别绘制出v1=0,0.5,1,1.5,2时候的运动轨迹并实现题目（2）中要求的数值解与理论解的对比情况，运行的结果如下所示：v1=1m/s,v2=2m/s时候的轨迹图：（2）然后是v1=0,0.5m/s,1m/s,1.5m/s,2m/s的理论曲线，并实现v1=1m/s时候的理论曲线和实际曲线的对比：该图中有两点需要注意的地方：1）标注均在曲线的右方，而v1=0时候的曲线实际上和y轴重合，所以绘出的曲线没有显现出来；2）注意到v1=1m/s时候的曲线，实际上此处还有一条蓝色曲线，图上可以看到蓝色曲线的痕迹，这是v1=1m/s时候的数值解曲线，可以看出二者的重合度是非常高的，几乎完全重合，这说明之前的数值解的结果是可信的。另外我们注意到1）当v1=2m/s的时候，船已经不可能正好到达河的正对岸了，从物理规律也可以直观的理解，因为当v1v2时候，v2的任何一个分量都不可能大于v1，这样一定会随着水流向下游走一定的距离，而不可能完全到达对岸；2）当水流速度逐渐增大的时候，路径越来越凸向右侧，即要到达正对岸会更加费力气，这也是符合生活规律的。【实验总结】本次实验总的来说思路较为清晰，通过实验进一步加深了对matlab的了解，巩固了课堂所学的常微分方程的一些知识。在matlab的实现过程中，有很多需要注意的