1、Chapter5 常微分方程边值问题的数值解法5.1 数学模型二阶常微分方程的边值问题如下：控制方程：* MERGEFORMAT (,)()yfxyaxb(1.1)其边界条件（Bounding Condition, B.C）可分为如下三类：（1）第一边值条件: （1.1.1）(),()yay（2）第二边值条件： （1.1.2） b（3）第三边值条件： （1.1.3） 0101()(),()+,yyb其中:第一边值条件给定了在 边界处的函数值；第二边值条件给定了在边界处的导数值；第三边值条件为混合边界条件， 给出了边界处的函数值和 导数值之间的关系。当方程* MERGEFORMAT (1.1)
2、为线性方程时，其一般形式为：* ()()ypxqyfxabMERGEFORMAT (1.2)设其系数 、 及自由项 均为 区间内 的连续已知函数，故方()pxqf(,)b程* MERGEFORMAT (1.2)为变系数的线性方程， 为消除其中的一 阶项，处理方法如下：将* MERGEFORMAT ( 1.2)式两端乘以 ，即得：()pxde* ()()() pxd pxdpxdyeyqyfe MERGEFORMAT (1.3)由求导公式 ，则上式可改写为：()uv(1.3a)()()()pxdpxdpxdyeqeyfe令： ，则方程（1.3a）可以新变量 表为如下形式：-()pxdtet* M
3、ERGEFORMAT (1.4)2()dyQtRt此形式仍为二阶变系数线性常微分方程，但其中已不 显含一 阶导数 项，而 边值条y件的形式不变。这表明，对于一般形式的二阶常微分方程，总可以通过上述变量变换的处理将其转化成不含一阶导数项的形式。 为此，在下面考 虑线 性方程的情况时，可以假设方程中不再显含 项。y对于方程* MERGEFORMAT (1.4)在多数情况下，无法求得其解析解，故只能采用数值方法求解。关于边值问题的数值解法有多种，本章只介 绍其中比较常用的两种方法：差分法和试射法。5.2 解线性边值问题的差分方法5.2.1 差分方程的建立现对线性方程的情况来讨论之， 设方程的形式为：
4、* (),()0,)yqxrqxaxbMERGEFORMAT (2.1)1. 将区间离散化：取等间距 将区间 划分为 等分， 诸结 点为：(bahn,bn0(0,12.),i nxhixab如此，我们将原在 上求解方程 * MERGEFORMAT (2.1)解 的问题转,a ()yx化为求各结点 上近似解值 问题。为此，我 们首先(012.)ixn(,.)iy应将方程* MERGEFORMAT (2.1)中的变量进行离散化处理，其具体作法是：2. 构造差分方程对内部结点 ，利用二 阶中心差商来近似原二 阶微商，由二 阶中(,.)ix心差商的关系式：* MERG2 2 (4)112() ,.;,
5、iiiiix iiiiydhy ynxEFORMAT (2.2)其中 ，将它代入方程* MERGEFORMAT (2.1)中，得在 结点 处()iiyx ix所满足的关系式：ii 2(4)112iiiiiiyhqyrh其中， (),()iiiiqxr在上式中，若略去截断误差项 （它是间距 的二阶小量）， 则 可得原微分方程* MERGEFORMAT (2.1)的近似差分方程为：* 112(1,2.)iiiiiyqyrnhMERGEFORMAT (2.3)这是含有 共 个未知数的线性方程组，而方程的个数为(0,iyn个，欲使此方程组有唯一解，还需由两个边值条件补 充两个方程。1n对于第一边值条件
6、，可直接由 给定的边值条件给出两个补 充方程，即：* MERGEFORMAT (2.4)0ny将方程* MERGEFORMAT (2.3)与* MERGEFORMAT (2.4)联立即构成了第一边值问题完整的差分方程组。对于第二和第三边值条件，由于两者均 给出了边界处的一 阶导数信息，相 应于差分方程的形式，对于边界的导数 值或表达式，我们亦须用差商来近似表示之。因 为我们无法利用区间 之外结点的信息，所以在引进两个边界 导数,ab的差分近似表达式时，就不能再利用中心差商公式。0()()nyy和若要求的截断误差为一阶的，即 ，则我们可利用简单 的前、后差商表达式来近()Oh似原导数边值条件，即
7、有：* MERGEFORMAT 100nyh;（ 前 差 ）（ 后 差 ）(2.5)若要求的截断误差为二阶的，即 ，则原导数边值条件的差商近似表示需利用2()ONewton 等距插值公式（补图）：前插公式：* MERGEF 2210000431()()yyayhhORMAT (2.6a)后插公式：(2.6b) 2121()()nnnnnyyyby 其中 为向前差分符号。1iii将方程* MERGEFORMAT (2.6a) 、（2.6b）与方程组* MERGEFORMAT (2.3)联立，即构成了第二边值问题完整的差分方程组。对于（1.1.3）式给出的第三边值条件，利用 * MERGEFORM
8、AT (2.6a)和(2.6b)式，其边值条件可以差分方程近似表为：* MERGEFORMAT (2.210143yyh7a)(2.7b)12014nnnyh将方程* MERGEFORMAT (2.7a) 、(2.7b)与方程组* MERGEFORMAT (2.3)联立，即构成了第三边值问题完整的差分方程组。至此，我们通过离散化的处理，将原微分方程的边值问题近似 转化成为一个差分方程组（线性方程组）的边值问题了。5.2.2 差分方程组的解法追赶法由上述建立的差分方程组，其系数矩 阵为三对角形的，故通常利用求解线性方程组的追赶方法求解。现以第一边值问题的差分方程组求解为例： 0112(1,2.)
9、iiiiinyqyrnh分别将 代入内部结点 和 的差分方程中，消去未知数 ，0,y1xn 0,ny则上差分方程组可表为：* MERGEFORMAT (2.8.1)2211()qhyr（2.8.2）223()qhyr.（2.8.n-2）223 1()nnnnyh（2.8.n-1）21qhyr此方程组的系数矩阵是三对角形，且 （即主对角元素均 为主元素），故可按自0i然顺序利用追赶法求解，其计 算步骤如下：先从方程（2.8.1）中解出 1y* 211 1212hryUyVqMERGEFORMAT (2.9.1)将上式代入方程（2.8.2）中解出（2.9.2）212323211()()VhryyU
10、yVqhUq一般可设有下列递推关系式：（2.9.i）1(,.)iiiyn（2.9.i+1）121iiiyV为了求得 的递推公式，可将（ 2.8.i）式代入（2.8.i+1）式中，解出 为：iiUV和 1iy211221()()iii ii iVhryyqhUqU从而有 的递推关系式为：UV与* MERGEFORM121(+)(0,12.)iiiiqhnrAT (2.10)即由上式，从 开始，可逐次求出 之值，首先需要 之值。0i,iUV0,UV由（2.9.i）式可看出，当 时，由 B.C，则有：i010y从而应有： ,V这样便可利用递推公式* MERGEFORMAT (2.10)按下标从小到大
11、的顺序逐个求得系数 ，此 过程称为追的过程。,(1,2.)iUVn而求 的过程正好相反。因 为已知 ，于是倒 过来，我们先iy ny从（2.9.i）中最后一个方程中直接求出 ，即：1ny1nUV将此式代入（2.9.i）中的倒数第二个方程中解得 ，依次往上逐次求解，即可求得2y。求 的过程中下标由大到小，故此 过程称为赶的过程。3421,.,nyy i对于系数矩阵呈三对角形式的线性方程组，追赶法是一种十分有效的方法，它在微分方程边值问题的数值解法中有着广泛的应用。例 1 （计算方法武大版 ）235P (01)(),yx解:取步长 ,则结点 ,按(2.8)式,其差分10h0(,2.0)iih方程为
12、:2 212 22 292(10)000.()1().10()y 其解为: 1234567890.78364.5910.4.8160.5yyy例 2 (常微分方程数值解法 南大 )137P(1)0),yx解：取 ，结点0.2h.2(,.5)ixihi其差分方程为：（a）0115()0(,34)()iiiyi边 界 点 内 结 点边 界 点由（2.9.i）式，将上方程组化为如下形式：（b）1(,2)iiiyUV其中： （c）1,02i ii由于 ，由（b）式得 ，由 递推式（c）逐次解得：0y012341234, 05UUV再由 和（b）式逐次求得：5y4321,5yy5.3 差分算法的稳定性与
13、解的收敛性稳定性和收敛性是反映某一算法好坏的两个重要的概念，两者是不同的两个概念。本节将对差分算法的稳定性和解的收敛性进行简单的讨论。5.3.1 算法的稳定性所谓算法的稳定性简单地讲就是指在计算过程中的某一扰动误差不会在后续的计算中被逐渐放大传播。稳定性是反映某一计算步骤中出现 的误差对计算结果的影响。追赶法在解上述差分方程的过程中，具有其 计算误差在传 播中不增长的优点，即：在第 i 个节点上对 的 计算误差 ，在左右相 邻的两节点上该误差不,iiUVy,iiiUVy会增长。下面对此进行定量的分析证明：设：在计算系数 中产生一误差 ，由此影响到计算 中产生了误差 。i i 1i1iU现证：
14、1iiU事实上，由于 和 均满足递推关系式* MERGEFORMAT (2.10)iiU1iiU中的第一式，即有：* 121()()i ii iiUqhMERGEFORMAT (2.11)将* MERGEFORMAT (2.11)与* MERGEFORMAT (2.10)中第一式相减，则有：* 1112 21 12 2121()()()=()ii i iiiiiiii UUii iiUiiUiiiUqhqh MERGEFORMAT (2.12)我们用归纳法可证 10i设： 01iU2,i iiqqh由* MERGEFORMAT (2.10)式， 可知：121()iiiUqhU10i再由* ME
15、RGEFORMAT (2.12)式可知： 1ii再看求 的过程，从* MERGEFORMAT (2.10)式中第二式和（ 2.9.i）式出发，iVy和同上类似推导，可得如下关系式： 1i iiiVVyyU于是，由 和 可知： ， 。01iU1i1iiV1iiy此结果表明：追赶法在计算过程中某一计算误差不会在计算过程中被放大传播，即该算法是稳定的。5.3.2 差分方程解的收敛性所谓收敛性是反映所用计算格式（公式）自身的截断误差对计算结果的影响。设： 是边值问题的精确解在结点 上的数值， 是差分方程在 结点 上的数()iyxixiyix值解，则两者误差： ()(0,12.)iiiEyn利用极值原理
16、和推论以及数学上的推导处理（可参见计算方法武大版 ）.可最23P终获得误差的上限估计式为:* 2241()()96iiiyxhMbaMERGEFORMAT (2.13)其中， ， 为 的四阶导数。(4)4maxbMy4(从此估计式可见，误差的上限是与 三项 因子成正比的，当 边值问224,)hba题给定之后，因子 均为定值，在计算中是不 变的。从而易见，当步长24,()，误差 。这表明：当步长 逐 渐缩短时，差分方程的解 将0h0iiiEyxiy逐渐收敛于原微分方程边值问题的精确解。5.4 试射法（Shooting Method）基本作法：将边值问题转化为初值问题，并通 过给定和修正初始斜率来
17、求解。现仅以第一边值问题为例介绍试射方法。设二阶方程的第一边值问题：* MERGEFO :(,)yfxy控 制 方 程RMAT (4.1)* MERGEFORMA.:(),()BCabT (4.2)试射法的出发点就是设法确定初值点 处的斜率 之值 ，使方程* a()yamMERGEFORMAT (4.1)满 足初值问题 的解亦能满足边值问题* (),yMERGEFORMAT (4.2)，即 。也即意味着要从微分方程* bMERGEFORMAT (4.1)的解曲 线经过点 而具有不同初值斜率的积分曲线中，去(,)a寻找一条经过点 的曲线。(,)b问题的关键是如何选取一个能够满足上述要求的初值斜率
18、值。试射法顾名思义就是：根据经验，或对方程进行定性分析，或按照实际问题中存在的运 动规律，选取一个初值斜率 ，以这个斜率进行试算，即解如下初值问题：1m* MERGEFORMAT 1(,)yfxyam(4.3)获得相应的一个解曲线 。1()yx若 ，或 （允许误差）， 则 即为 原边值问题* 1()ybb1()yxMERGEFORMAT (4.1)和* MERGEFORMAT (4.2)的解。否 则，我们可依据与 的差距来适当地将初值斜率 值进行修正，修正公式：1() 1m* MERGEFORMAT (4.4)21再以 作为初始斜率，解如下初 值问题：2m* MERGEFORMAT 2(,)y
19、fxyam(4.5)解上方程可得另一解 ，同前，若 ，或 满足 ，则2()yx2()yb2()yb即为所求之解，否则再 对初始斜率 值进行适当修正，最 简单的修正办法是同时2()yx利用 和 的信息，通过线 性插值方法来求出新的初始斜率 值，即：1m2 3m* MERGEFORMAT 2131()m(4.6)以新的初始斜率 重复上述求解过程，直到 获得满足精度要求的解 为止。3可见，试射法就类似于射击校正 问题，它是通 过不断调整修正初射的角度（斜率），使弹着点最终能落在给定的目标点 或该点邻近处（补图 ）。(,)b从理论上讲，试射法对于高阶 的微分方程或方程组也是适用的，然而在实际使用中，可
20、能遇到不只一个未知初值需做选择，此时尝试的工作量将会大大的增加。因此，该法通常仅适用于只有一个未知初值的情况。例 用试射法解边值问题： 310,1()()yx给定精度 310解：取步长 ，取在初始点 处的初始斜率 ，解如下初 值问题：.h(0,)1m 311,0()()yx用龙格-库塔法解得相应解 的终点值1x，它与给定值 相差尚远，11()0.64,0.64yb故应对初始斜率进行修正，由修正公式* MERGEFORMAT (4.4)得：211.50.64m以 作为初射斜率，再重新解如下初 值问题：2m 321,1(0)().5yxm得 3222()1.6,.60yb故再应对初射斜率 值进行修
21、正，由修正公式 * MERGEFORMAT (4.6)，则得：2m2131 1.5()(0.64)064.以 为初射斜率，解初始问题 ：3m 3310,1()()4.yxm得 ，显然它离目标点值 又近了一些。33()10.685yb重复上述过程，当以 为初射斜率时，解得 ，43m88()10.9yb此时有 。8.91则近似解 即为原边值问题的满足精度要求的数值解。有关计算结果见下表：()yx12.534.4m5 7m84.3x()()()yx() ()yx0 0 0 0 0 0 0 0 00.1 0.1 0.155 -0.441 -0.4820.2 0.2 0.309 -0.869 -0.94
22、8 0.9 0.667 0.726 0.256 0.4971.0 0.646 0.610 0.685 0.9995.5 非线性差分方程及迭代解法5.5.1 非线性差分方程设二阶非线性第一类边值问题为：* (,)()yfxyaxbMERGEFORMAT (5.1)* MERGEFORMAT (5.2)(),()其中 为 的非线性函数。 (,)yfxy(,xy（1（ 将积分区间 离散化：)ab取步长为 ，结点为hn(0,12.)ixhin（2（ 在诸内点 处以一阶和二阶中心差商近似方程* (1,2.)inMERGEFORMAT (5.1)中的一阶微商 和二阶微商 ，即有：dyx2dyx* 1112
23、(,)(1,.)2ii iiiyf nhhMERGEFORMAT (5.3)（3（ 由边条* MERGEFORMAT (5.2) 列补充方程：* MERGEFORMAT (5.4)0ny将* MERGEFORMAT (5.3)和* MERGEFORMAT (5.4)联立组成了一个关于的完整的方程组，但由于* MERGEFORMAT (5.3)中的右端项01,.ny为变量 的非线性函数，故完整方程组为非1()iiiyfxh(1,2.)iyn线性的方程，所以不能直接求解，通常须利用迭代法求其近似解。在此我们介绍两种迭代方法：Jacobi 法和 Newton-Raphson 法。5.5.2 Jaco
24、bi 法最简单的迭代方法即为 Jacobi 法，其具体 实施方法如下：先将初始值（人为给定）近似解 代入方程组* (0),12.)iynMERGEFORMAT (5.3)的右端项，即有 中使之成为已知(0)(0)()1,2i iifxyyh量，此时完整方程组* MERGEFORMAT (5.3)、* MERGEFORMAT (5.4)便成为关于未知量 的线性方程组，可利用通常的求解线性方程组的各种方法直(0,12.)iyn接求解。设求得的新解为 ，然后再以 为初始(1)0,2.)iyn(1)0,2.)iyn解，重复上述迭代过程求出更好的近似解。其具体迭代求解格式为：* (1)(1)()()()
25、12()01 , ,.1)2(0,.)kkkkiii iiiknfxyhhy MERGEFORMAT (5.5)可以证明：在一定条件下，此种迭代 过程是收敛的，且当步长 ，其数值（差分）0h解收敛于原边值问题* MERGEFORMAT (5.1)、* MERGEFORMAT (5.2)的精确解。5.5.3 Newton-Raphson 法设所考虑的边值问题的控制方程中不含一阶导数项，即 为* MERGEFORMAT (5.6)(,)yfxab其中 为 的非线性函数，且满足：(,)fx(,)1. 在 中， 为 的连续 函数；,xy(,)fxy(,)2. 对于 ，任意两个函数 值 ，存在着一个常数
26、 ，使得：ab*L成立；*(,)(,)fxyfLy3. 偏导数 连续，且 。,yfxf(,)0yfx在此情况下，我们可利用 Newton-Raphson 法迭代方法求解边值问题* MERGEFORMAT (5.6)的数值解。首先，同前处理，将积分区间 离散化为（n+1）个结点,ab，在每一内点 构造一个差分方程：(0,12.)ixahin(1,2.)ixn* 201(0(1,2.)ii iiiyhffinMERGEFORMAT (5.7)Chapter6 椭圆 型偏微分方程的数值解法6.1 概述6.1.1 方程椭圆型偏微分方程是数学、物理及传热学、 电磁学、弹性力学、流体力学等工程科学技术问题
27、中常遇到的一类重要的方程。在平面域 G 上，椭圆型二阶偏微分方程的一般形式为：* ()()(,) (,)UpqUfxyGxyMERGEFORMAT (1.1)式中： 为方程* MERGEFORMAT (1.1)待求的解函数，而,Uy均为平面域 G 上自变量 的连续已()0(,)0,()pxqxyfxy(,)xy知函数，且与解函数 U 无关。若 .，则方程* MERGEFORMAT (1.1)可表为：Const（1.1.a）22(,)1(,) (,)qxyUfxyGxp若 ，则（1.1.a）成为 ：(,)0qxy（1.1.b）;22(,) (,)fxyx其中 ;1(,)(,)fxyfyp方程（1
28、.1.b）称为 Poisson（泊松）方程。它描述的是平面域 G 内部有源的静电磁场或稳定温度场的场分布规律。若方程（1.1.b）中的右端项函数 ，则（1.1.b）成为：;(,)0fxy（1.1.c）220 (,)UxyGxy此方程称为 Laplace（拉普拉斯）方程，它描述的是平面域 G 内无源的静电磁场或稳定温度场的分布规律，或等截面直杆的扭 转时杆截面的翘 曲分布。引进 Laplace 算子符号： 2xy则泊松方程（1.1.b）和拉普拉斯方程（1.1.c）可分别记为： 2(,) (,)0 UfxyG6.1.2 边界条件（B.C）( 补图)在实际应用中，通常需要求解函数 在区间 G 内满
29、足椭圆型偏微分方程，并(,)xy在 G 的边界上满足给定的边 界条件的解，即所 谓的边值问题。如同常微分方程的边值问题， 椭圆型二阶偏微分方程的边界条件（B.C ）通常有以下三 类：（1（ 第一类 B.C：* MERGEFORMAT (1.2)U|=,)xy（其中 为定义在 上的已知函数，即在平面域 G 的全部边界上给出了函数 U 的,)xy（分布规律。在数理方程中，将寻求满足控制方程 * MERGEFORMAT (1.1)和 B.C * MERGEFORMAT (1.2)的定解问题称为 Dirichlet（狄利赫莱）问题。（2（ 第二类 B.C：* MERGEFORMAT (1.3)U=,)
30、xyn（其中 n 表示边界曲线 的“外法线”方向。此 B.C 给出了函数 U 在所有边界点上其外法线方向上的梯度（导数）分布规律。在数理方程中，将寻求满足控制方程 * MERGEFORMAT (1.1)和 B.C * MERGEFORMAT (1.3)的定解问题称为 Neumann（牛曼）问题。（3（ 第三类 B.C：* MERGEFORMAT (1.4)U(,)=,)xyxyn（其中 、 均为定义在 上的已知函数，且 。(,),)（ (,)0xy此 B.C 给出在所有边界 上函数 U 之值与其外法线方向导数 之间的关系式。Un在数理方程中，将寻求满足控制方程 * MERGEFORMAT (1
31、.1)和 B.C * MERGEFORMAT (1.4)的定解问题称为 Robin（诺宾）问题。偏微分方程边值问题的求解是一个较为复杂的问题，除少数经典问题外，要求得 问题的精确解一般是十分困难的。因此，人们不得不着眼于各种近似的数值解法，其中工程上最常用的数值解法有两种：有限差分法和有限单元法。本章将着重介绍前者。6.2 差分方程构建现以平面域 G 上的二阶泊松方程 为例， 讨论其差分方程的构建及三种差分格式。控制方程：* 22(,) (,)UfxyGxMERGEFORMAT (2.1)其差分方程的构建主要由三部分步骤组成：1）将求解区域 G 离散化；2）将控制微分方程以数值差商的形式表出；
32、3）由 B.C 构造补充方程（将在下一 节讨论）。（1（ 区域 G 离散化如图（补图），我们用平行于 和 轴的两组平行线族，将 问题 的平面域 G 离散化：xyx0y0xGhl两族相交直线的任一交点（网格）（ ）的坐标为：,ixy0(,12.)ijhyjl其中 分别为网格沿 轴的步长。0,hl,x显然，若 时，为矩形网格；当 时，为正方形网格。hl所有属于平面域 G 的网格点称 为内结点，如图中所有划“o”之点。（2（ 微分方程的差分逼近取 正方形网格。hl设： 是区域 G 的内结 点，简记为（ ），写出控制方程在点（ ）处的表达(,)ixy,ij ,ij式，即有：* 222()()()iji
33、jijijUfxyMERGEFORMAT (2.2)现我们用二阶中心差商近似代替方程* MERGEFORMAT (2.2)中的二阶微商，即：* 21,1,2,2()ijijijijijijijijxhUy;MERGEFORMAT (2.3)将* MERGEFORMAT (2.3)式代入方程* MERGEFORMAT (2.2)中，则得到 的ijU近似解 所满足的差分组为：ijU* 1,1,2( 4)ijijijijijijUUfhMERGEFORMAT (2.4)由于方程* MERGEFORMAT (2.4)中只有结点 及其四个菱形端点上的 值(,)ij U出现（补图），故通常称为“五点菱形格
34、式”。仿照拉氏算子的定 义，将此种格式定义为一种差分算子，用符号 表示，即：* 1,1,2( 4)ijijijijijijUUhMERGEFORMAT (2.5)从而差分方程* MERGEFORMAT (2.4) 亦可用五点菱形差分算子简写为：（2.4a）ijijf它与原微分方程* MERGEFORMAT (2.1) 的逼近误差为：* MERGE22(4) (4)1 2(),)1ijij ijxijyijhEUUhUxh FORMAT (2.6)其中： ； ， 分别表示 U 对 和 的四阶偏导数。12,(4)x()y若以 表示 在区域 G 内所有四阶偏导数绝对值 的上确界， 则可得逼近误4M(
35、,)差 的上限估计式为：ijE* MERGEFORMAT (2.7)246ijhEM此时表明：“五点菱形差分格式”的截断误差为二阶的， 记为 。2()Oh* MERGEFORMAT (2.5)式表明：差分逼近 是网格点 及其四个邻点上ijU,ij的 U 值的线性组合。由此产生这样的想法： 能否可用网格点 及其另外四2()ij()个矩形结点的 U 值去作线性 组合呢？这一问题的回答是肯定的，事 实上，我们只要简单地将网格旋转 ，即可得 * MERGEFORMAT (2.5)式的另一种差分逼近格式的方程45为：* 1,1,1,1,2( 4)ijijijijijijUUfhMERGEFORMAT (
36、2.8)同前，我们可定义“五点矩形格式 ”的差分算子为：* 1,1,1,1,2( 4)ijijijijijijh;MERGEFORMAT (2.9)则差分方程* MERGEFORMAT (2.8) 可简记为：（2.8a）ijijUf;其逼近误差为：* 2*2(4)(4)(4)123()6ijijijxxyyhEUQQ ;MERGEFORMAT (2.10)式中： ；点 均在闭正方形域： 之4(4)2)xyU123,Q,i ixhyh中。令 为 所有四阶偏导数绝对值在 G 内的上确界。4M(,)xy则逼近误差 的上限估计式可表为：*ijE* MERGEFORMAT (2.11)*243ijEhM
37、此式表明：“五点矩形差分格式”的截断误差亦为二阶的， 记为 。2()Oh欲提高差分逼近的精确度，则须 引入更多的结点信息。下面直接给出逼近拉氏算子的“九点格式”的差分算子，它是“五点菱形差分算子”和“五点矩形差分算子”的线性组合算子，其具体表达式为：* 1,1,1,+1,2,364( )20ijijijijijijijijUUhU ;MERGEFORMAT (2.12)该差分算子的逼近误差的上限估计为：* 68*2503!()ijijijMhEMERGEFORMAT (2.13)式中： 是 在区域 G 上的所有 8 阶偏导数绝对值的上确界。8M,Uxy这表明“九点格式”的差分算子其截断误差为六
38、阶的，即 ，其逼近精度将大大6()Oh地提高。6.3 边界条件的处理由前构造的差分方程是针对区域 G 内各结点的，在靠近边界 的各点的方程中，将需要边界 之外有关结点的函数 U 之信息，如前图中的“ ”点。而这些点上的 U 值亦是未知的。这样，由内结点列出的差分方程的数目将少于未知数的个数，这些均须根据给定的边界条件建立补充的差分方程，使差分方程 组中方程的数目与未知数个数一致。6.3.1 第一类 B.C 的处理* MERGEFORMAT (1.2)U|=,)xy（1（ 若区域 G 的边界 恰为规则的矩形， 则我们进行网格划分时，就将边界结点设置在边界线 上，如下图中的 A、B点。此时， 对于
39、任何五点格式或九点格式，我们都可以根据给定的第一 B.C 之值， 简单地给出边界 上各个 节点的函数值：A B CG，这即为补充方程。(),(),.U（2（ 对于一般情形，区域 G 的边界 是由一般曲线组成，它们难得包含网格的一个结点，如下图所示（补图）。为了保证计算的精确度，就需根据 边界曲线 与网格节点之 间所处的位置以及问题的精确度要求，采用不同的近似处理方法。通常有以下三种 处理方法。1） 零次插值（简单迁移）即对每一靠近边界曲线 的内结点，如 图中的 点，选定一个与之最近的边界PG上点 ，使：Q* MERGEFORMAT (3.1)()UQ其中 Q 为 上离 P 的某一最近点，为简单
40、起见亦可取：x 方向上最近点()E或 y 方向上最近点()N上式即为补充方程。由于将 用 之值代替，可以看成是在 y 方向、在 P 点处用 处之值()UP) N的零次多项式插值的结果。同理 可视为是在 x 方向用 处之值的零次多()N()UEE项式插值的结果，故此近似处 理方法称为零次插值方法。对于 为光滑函数的情况，由微分学的中值定理可知：当步 长 时， 这种(,)xy 0h插值产生的局部误差为一阶的，即 。()Oh2） 一次插值这是一种比较准确的处理方法。它不是直接取定靠近边界 之 P 点的 之值，()U而是沿着水平或铅垂网格线，利用 线性插值方法确定 P 点之值 。()参见下图（补图），
41、设 U 沿 x 方向呈线性分布，则由比例关系：()EUwPh从而线性插值公式为：（x 方向）* ()()EhPMERGEFORMAT (3.2.1)同理： （y 方向） （3.2.2）()()NUSh上式即为补充方程。其中， 和 为定值； 、 为未知数，即将未知数 用一定值和()UE)N()wS()UP一未知数的线性插值形式表出，并没有增加任何新的未知数。若函数 为充分光滑的，则这种线性插值产生的局部误差为二阶的， 。(,)xy 2()Oh3） 二次内插即在节点 P 上列不等距的差分方程。在点 P 处用以下二阶差商公式近似代替方程* MERGEFORMAT (2.1)中二阶微商，即有： 22(
42、)()() 2()()() 2EPNPUUPwhhxPShhUy;（ ）（ ）代入方程* MERGEFORMAT (2.1) 中，由此导出在 P 点的差分方程为：* 2()()2()()2)()PENENNEWUSUfhhhhMERGEFORMAT (3.3)上式即为补充方程。这种二次内插处理相对较为复杂，但精确度并不高，由它产生的局部误差为 。()Oh此外，它还有一个很大的缺点：破坏了差分方程 组的对称正定性。6.3.2 第二、三 B.C 的处理这两种边值条件可统一写成如下形式：* MERGEFORMAT (,)Uxyn(3.4)当 时，* MERGEFORMAT (3.4)式即为第二 B.
43、C；当 时，即 为第三 B.C。0 0由于其边界条件中含有 U 的法向 导数，故它的近似处理要比第一 B.C 复杂些，下面分两种情况讨论之。（1）矩形网格的节点 P 刚好落在边界曲线 上此时，就在 P 点直接由 B.C 列出一个差分方程，实现这个问题的关键是如何用差分近似代替法向导数 。Un边界上 P 点的外法向 与坐标轴平行（补图）此时，显然有：* (P)-UQ()PUnxh前 插 公 式MERGEFORMAT (3.5)代入 B.C * MERGEFORMAT (3.4)中， 则得到 U 在 P 点应近似满足的差分方程，即补充方程为：* ()()UPQPhMERGEFORMAT (3.6.
44、a)边界上 P 点的外法向 与坐标轴不平行（补图）n此时，显然有：* cos(,)cos(,) ()cos)P PUUnxnynQWhhMERGEFORMAT (3.7)代入 B.C * MERGEFORMAT (3.4)中， 则得 U 在 P 点应近似满足的差分方程，即补充方程为：（3.6b）()()coscos()UPQPhh（2）网格的边界点 P 不在边界曲线 上（补图）此时仍按（1）中的情况在 P 点列差分方程，但需将 上与 P 点相近的某点 的*外法向作为 P 点的外法向 ，并近似令（ 3.6b）中的 。n *()至此，我们对正方形网格的情况其各种 边界条件进行了近似 处理，从而 获
45、得了补充方程，将原泊松方程的边值问题转 化为一个完整的线性代数方程 组求解。数学上可以证明：该差分方程组的解存在且唯一，并当步 长 ，其解收 敛于原0h微分方程的精确解。证明可参 见武大、山大编计算方法 。2475P例 （高应才编数学物理方程及其数值解法 ）189求解拉氏方程的第一 B.C 问题：无源稳定的温度场的分布。220(4,03)PUxyxy其中 的取值如图所示（补图）。解：用方形网格，取 ，共有 6 个内结点，14 个边界结点，所有内结点的函数 U 值1h是未知的，所有边界结点的函数 U 之值均已给定。计算结果如下：ij0.867 .41 0.2961358Chapter7 加权残值
46、方法简介本章参考文献：1.徐次达，加权残值法和有限元混合法， 应用数学和力学第 32 期讲义2.徐文焕，陈虬，加 权余量法在结构分析中的应用， 铁道版，19853.邱吉宝，加权残值法的理论与应用，宇航版， 19914.徐次达，固体力学加权残值法，同济版， 19877.1 方法的基本概念大量的应用科学和工程学问题最终都归结为求解微分方程（组）的边值问题或初值问题， 这个微分方程（组）可以是常微分方程（组）或偏微分方程（组），方程（ 组）可以是线性的，亦可以是非线性的。众所周知，对于微分方程（组）的 边值问题或初值问题，除极少数经典问题之外，解析方法是无能为力的，只能依 赖于各种数值方法，如前几章学 习的常微分方程初
