1、练习 16 常微分方程（一）知识背景求解常微分方程在工程技术中运用比较广泛，在水力学、电工电子技术、环境模型、生态模型以及经济管理等各个领域都有重要的运用。由于常微分方程在工程技术中运用的广泛性、重要性，我们分两个练习来重点讲解和练习。求解常微分方程，有时候会遇到无法求出解析解的情况，即不能通过等式左右同时积分法不定积分求出。但是我们有时候在给定的条件下，可以求出常微分的数值解来逼近数值解，以完成对模型的研究。求解常微分方程数值解的方法通常有欧拉法和龙格贝塔法，对于微分方程，欧拉方法解法的基本思想是在小区间 上用差商代替方程的导1,nX数，并通过把小区间不断地划分求极限，从而最终得到数值解。龙
2、格贝塔方法的基本思想也是用差商来代替导数，只不过龙格贝塔在小区间是运用了微分中值定理，在小区间内多取点而取加权的平均值来构造精度更高的计算式。而在 MATLAB 中，主要采用龙格贝塔法来计算常微分方程的数值解。主要内容【本练习考查知识点】本练习主要通过几个解一阶常微分方程或其模型来辅导读者掌握解决常微分方程数值的基本方法。本练习主要用到的函数是 ode23 和 ode45，其中ode45(阶及阶)是首选的方法，而 ode23（阶及阶）是比 ode45 低阶的函数。其主要表达式是：t,x=ode23(f,t0,tf,x0,TOL)t,x=ode45(f,t0,tf,x0,TOL)f 是将所求解的
3、常微分方程写成的函数并存于 m 文件中， to、tf 分别为自变量的处始值和最终值，而 x0 则是函数的初始值，是误差范围（缺省值是0.000001） 。这些知识点请读者熟记。练习过程（1） 先以一个简单的例子来说明解常微分方程的解法：例：求解常微分方程： 。3x其中初始值为 1)0(x这是一个较简单的微分方程，可以通过在等式两边求积分从而求出其表达式，为了后面研究方程的数值解，这里我们采用求数值解的方法，在求解中我们只给出在区间0,1 的数值解及图形，函数 yfun.m 文件：function exer=yfun(t,x)exer=-x.3程序如下：t,x=ode(yfun,0,1,1)pl
4、ot(t,x,k+,t,x,r)结果如图 16-1、图 16-2、图 16-3 所示。本例仅通过一个简单的练习让读者熟悉和理解求解微分方程数值解时所用的基本方法和基本函数，下面我们通过在生态学中的一个重要的捕食模型来继续微分方程数值解的研究。（2） 求解二阶以上的常微分方程数值解，通常是将表达式迭代的过程，即使将方程在函数中叠代表达，从而在命令行中编程实现。我们以一个常见的物理模型为例，来说明二阶以上常微分方程的一般解法在单摆模型中，以垂直方向为平衡位置，即 x=0 时，以右边为正方向建立摆角x，小球所受重力沿运动轨迹方向的分力为-mgsinx ，利用牛顿运动第二定律有：xmgxlsin设初始
5、偏离角度为 a，此时无初速，则方程的初始条件为：0)(,)0(xa。由于方程没有解析解，故我们求其数值解。初始角度为 0.1745（弧度）时，求数值解如下：先将所列的方程存入一个函数文件中，然后再在命令区下计算调用函数文件。函数文件如下：high.mfunction xexer=high(t,x)g=9.8length=25xexer=x(2);g/length*sin(x(1);在命令区中输入的程序如下：ts=0,10;x0=0.1745,0t,x=ode(high,ts,x0)运行结果为： t x0 0.17450.0007 0.17450.0015 0.17450.0022 0.1745
6、0.0030 0.17450.0066 0.17450.0103 0.17450.0140 0.17450.0177 0.17450.0362 0.17450.0546 0.17440.0731 0.17430.0915 0.17420.1838 0.17440.2761 0.17190.3683 0.16990.4606 0.16730.7106 0.15760.9606 0.14401.2106 0.12691.4606 0.10681.7106 0.08401.9606 0.05912.2106 0.05912.4606 0.00582.7106 -0.02142.9606 -0.048
7、23.2106 -0.07373.4606 -0.09743.7106 -0.11873.9606 -0.13724.2106 -0.15234.4606 -0.16374.7106 -0.17114.9606 -0.17435.2106 -0.17345.4606 -0.16815.7106 -0.15885.9606 -0.14576.2106 -0.12896.4606 -0.10916.7106 -0.08666.9606 -0.06197.2106 -0.03587.4606 -0.00877.7106 0.01857.9606 0.04538.2106 0.07108.4606 0
8、.09498.7106 0.11658.9606 0.13539.2106 0.15089.4606 0.16269.5955 0.16749.7303 0.17099.8652 0.173310.0000 0.1744【练习小结】本练习主要介绍了一阶以及二阶常微分方程的解法，并介绍了 MATLAB 中常微分方程常用的函数及其用法实例。通过实例介绍了求解二阶常微分方程的迭代算法，并介绍了其在 MATLAB 中的实现。【思考题】1请运用本练习介绍的算法以及函数来解常微分方程： 1)0(,2 x2求解二阶常微分方程：0)(,2 x3ode23 和 ode45 的精度有什么差别？4如何求解二阶以及二阶以上的常微分方程？【思考题求解】1函数文件 first.mfunction exer=first(t,x)exer=-x.2;程序如下：t,x=ode23(first,0,1,1);plot(t,x,k+,t,x,r)23ode45 为阶及阶是 MATLAB 首选的方法，而 ode23 为阶及阶是比ode45 低阶的函数。4求解二阶以上的常微分方程数值解，通常是将表达式迭代的过程，即使将方程在函数中叠代表达，从而在命令行中编程实现。
