《第9章 微分方程问题的求解》由会员分享，可在线阅读，更多相关《第9章 微分方程问题的求解（22页珍藏版）》请在人人文库网上搜索。
1、第9章 微分方程问题的求解,9.1 单步方法 9.2 线性多步法 9.3 一阶微分方程组和高阶微分方程组 9.4 边值问题的求解 9.5 偏微分方程的数值解法 9.6 实例解析,本章目标：求 的数值解,9.1 单步方法,为求解前述初值问题的数值解，我们采用离散化方法。在求解区间上取一组节点： 称 为步长，为简单起见，仅考虑等距步长。
离散化方法的基本特点是依照某一递推公式，按节点从左至右的顺序依次求出y(xk)的近似值yk，若在计算第k+1的近似值时，只用到它前一步xk,xk+1,yk的信息，那么称此方法为单步方法。单步法主要有Euler方法和Runge-Kutta方法等。,一、 欧拉方法,
2、欧拉法的建立及其几何意义：,1差商方法,2积分方法,对 在区间x, x+h上积分得：,特别地，当x=xn时，有,再用yn代替y(xn)，便有欧拉法公式（9.1）.,3欧拉法的几何意义,亦称为欧拉折线法,（9.1）,据矩形公式, 隐式欧拉法和二步欧拉法：, 隐式欧拉法,由于未知数 yi+1 同时出现在等式的两边，不能直接得到，故称为隐式 欧拉公式，而前者称为显式
欧拉公式。,用此公式计算 yi+1 时要用到前两步的信息 yi-1 , yi ，故称为二步欧拉法（公式）。, 二步欧拉法（中点欧拉公式）,一般先用显式 计算一个初值， 再迭代求解。,需要2个初值 y0和 y1来启动递推 过程，这样的算法
3、称为双步法 ， 而前面的二种算法都是单步法 。, 梯形公式, 显、隐式两种算法的平均,注：的确有局部截断误差 ， 即梯形公式具有2 阶精度，比欧拉方法有了进步。但注意到该公式是隐式公式，计算时不得不用到迭代法，其迭代收敛性与欧拉公式相似。,为提高精度，我们再从积分近似的角度来看欧拉公式的改进。,对 在区间xi, xi+h上积分得：,用梯形公式求积计算积分，可得,隐式,
改进的欧拉法,注：此法亦称为预测-校正法 。可以证明该算法具有 2 阶精度，同时可以看到它是个单步递推格式，比隐式公式的迭代求解过程简单。后面将看到，它的稳定性高于显式欧拉法。,改进的欧拉法的“嵌套”表示形式,改进的欧拉法也可以
) 点出发，以某一斜率沿直线达到 ( xi+1 , yi+1 ) 点。欧拉法及其各种变形所能达到的最高精度为2阶。,建立高精度的单步递推格式。,斜率 一定取K1 K2 的平均值吗？,步长一定是一个h 吗？,二、龙格
yi+1 与 y( xi+1 ) 在 xi 点的泰勒展开作比较,要求 ，则必须有：,这里有 个未知数， 个方程。,3,2,存在无穷多个解。所有满足上式的格式统称为2阶龙格 - 库塔格式。,注意到， 就是改进的欧拉法。,若
6、取 则得到另一个龙格库塔法 的计算公式：,此即中点方式。此公式称作变形的欧拉公式。,为了进一步提高精度，设除xi+p外再考虑一点：,xi+q=xi+qh (pq1),用xi , xi+p xi+q点的斜率值K1、K2、K3 加权平均作为平均斜率的近似值，这时计算格式具有形式：,yi+1= yi+h(1K1+2K2+3K3
7、2) ),仿二阶情形的推导，只要把K1,K2,K3代入yn+1的表达式中在( xn,yn)处作泰勒展开, 再与y(xn+1)在xn处的泰勒展开式比较,欲使公式的局部截断误差=O(h4),可得含参数的方程组,类似知有无穷多解；得到的公式称为三阶龙格库塔公式。,一个比较简单而又重要的三阶龙格库塔公式，称为库塔公式，其为：,继续上述这个过程,可以进一步讨论四阶龙格库塔公式。,
最常用的是4阶经典龙格-库塔法：,9.2 线性多步法,用若干节点处的 y 及 y 值的线性组合来近似y(xi+1)。,其通式可写为：,当 10 时，为隐式公式; 1=0 则为显式公式。, 基于数值积分的构造法, 亚当姆斯显式公
9、d final value yi+1,外推技术 /* extrapolation */,9.3 一阶微分方程组和高阶微分方程组, 一阶微分方程组,IVP的一般形式为：,前述所有公式皆适用于向量形式。, 高阶微分方程,化作一阶微分方程组求解。,引入新变量,初值条件为：,9.4 边值问题的求解,2 阶常微分方程边值问题, 打靶法,先猜测一个初始斜率 y (a) =
s，通过解初值问题,找出s*使得(s*) = ，即把问题转化为求方程 (s) = 0 的根。,每计算一个(s) 都必须解一个ODE.,9.5 偏微分方程的数值解法,一、几类特殊的偏微分方程 椭圆型偏微分方程 抛物型偏微分方程 双曲型偏微
10、分方程 特征值型偏微分方程（椭圆型偏微分方程的一个特例）,二、偏微分方程的界面求解 （1）绘制求解区域 选择Options|Axes Limits菜单，将打开一个设置x,y轴坐标范围的对话框，这里设置x,y轴的范围均为-0.1,1.1。 单击矩形按钮 ，画一个边长为1的正方形区域 ，该正方形自动命名为R1，并显示在区域上方的公式栏（Set
formula）中，为了定位光标，可选择Options|Grid菜单，这时将在区域内画网格。如果要精确地确定正方形的位置，可以在正方形内部双击鼠标左键。在出现的对话框中输入正方形的左边界（Left），下边界（Bottom），宽度（Width）和高度（Heig
11、ht）的数值。另外也可以在Name栏中修改正方形的名称。之后选择Options|Application|heat transfer菜单或直接在工具栏右侧的列表框中选择以确定偏微分方程的类型。 （2）设定边界条件 单击按钮 ，进入边界模式，这时区域由灰色变成白色，而边界变成红色，选择菜单Boundary|Show Edge
Labels给四条边界标上序号1、2、3、4。边界条件的默认设置是在所有的边界上，根据题意，边界1和2的边界条件需要修改，双击边界1，,在出现对话框的边界条件类型（Condition type）中选择Dirichlet（默认），这时描述边界的方程（Boundary condi
PDE）中选择椭圆型（Elliptic）（默认），这时方程的形式为-div(k*grad(T)=Q+h*(Text-T),只要取k=1,Q=0,h=0,Text=0即可。 （4）建立网格 单击按钮 ，将区域划分为三角形网格，为了能得到更高的精度，再单击按钮 ，将得到更小的三角形网格。 （5）输出图形 单击按钮 ，得到用颜色表示的解的分布图，为了得到更好的视觉效果，可单击按钮
，选择Hight(3-D plot)和Plot x-y grid，在这里Plot x-y grid是用矩形网格画图，速度要快一些，而Show mesh在是用三角形网格作图，速度稍微慢一些。最后单击对话框中Plot按钮即得到效果图。,
