第十章 偏微分方程数值解法
一、 典型的偏微分方程介绍
1．椭圆型方程
科学技术中经常遇到一些重要的、典型的偏微分方程。在研究有热源稳定状态下的热传导，有固定外力作用下薄膜的平衡问题时，都会遇到Poisson方程
（10.1）
其中D表示平面区域。特别在没有热源或没有外力时，就得到Laplace方程
（10.2）
此外，当研究不可压缩理想流体无旋流动的速度势以及静电场的电位等，也会遇到（10.1）或（10.2）类型的方程。
2．抛物型方程
在研究热传导过程、气体扩散现象、电磁场的传播等问题中以及在统计物理、概率论和重子力学中，经常遇到抛物型方程。这类方程中最简单、最典型的是热传导方程。
（10.3）
其中a是常数。它表示长度为L的细杆内，物体温度分布的规律。
3．双曲型方程
在研究波的传播、物体的振动时，常遇到双曲型方程。这类方程中最简单、最典型的是波动方程
（10.4）
它表示长度为L的弦振动的规律。
二、定解问题
偏微分方程（10.1）~（10.4）是描述物理过程的普遍规律的。要使它们刻划某一特定的物理过程，必须给出附加条件。把决定方程唯一解所必须给定的初始条件和边界条件叫做定解条件。定解条件由实际问题提出。对方程（10.3）来说，初始条件的提法应为
，其中f (x)为已知函数，它表示物体在初始状态下温度分布是已知的。边界条件的提法应为物体在端点的温度分布为已知，即
（10.5）
其中j（t）和y（t）为已知函数。
对（10.4）来说，边界条件的提法和（10.5）形式一样，它表示弦在两端振动规律为已知。初始条件的提法为
其中f (x)和g (x)为已知函数。它表示在初始时刻弦振动的规律和振动的速度。对于方程（10.1）和（10.2）来说，因为它们反应稳定状态的情况，与时间无关，所以不需要提初始条件。边界条件的提法为：
其中j (x, y)为已知边界，s是区域D的边界。
由偏微分方程和定解条件所构成的问题叫定解问题。许多实际问题提出的定解问题无法求出解析解，少数问题即使求得解析解，计算过程、解的表达式也可能很复杂。因此必须寻求方程简单、可以在计算机上计算的数值方法。本章将对三类典型的偏微分方程定解问题给出数值计算方法。
本章主要针对几个典型的微分方程介绍常用的差分方法和有限元方法。这些方法基本思想是：把一个连续问题离散化，通过各种手法化成有限形式的线性方程组，然后求其解。
§1 差分法简介
差分法是求偏微分方程数值解的重要方法之一，它的主要做法是把偏微分方程中所有偏导数分别用差商代替，从而得到一代数方程组——差分方程，然后对差分方程求解，并以所得的解作为偏微分方程数值解。
为此，必须对区域进行剖分，用网格点来 u=0
代替连续区域，因此差分法亦称“网格法”。 1
我们用一个简单例子来说明差分法的基本思想和具体要求。 U=0 取一边长为1的正方形均匀薄板， 
上下侧面绝热，四周保持恒温（如图10 .1），求板内各点的稳定温定分布。
这个总是如在数学物理方程中所知，它可
以化为拉普拉斯方程第一边值问题： 0 u=0 1
(10.1)
一般的来说对这类问题我们无法求出解的解析表达式，有的即便能求出也是很复杂的，在实际问题中往往也并不需要求出u在区域W内每点值，实际上能求出在区域内某些点的近似值也就满足需要了。
在图10.1中作平行于坐标轴间隔为
的两族直线，我们求u在网格点（落在W内两族直线的交点）上的值，并且以后采用下列记号：
我们利用u在这些点满足主程
（10.2）
求出u在网格点上的值，（10.2）中( )_表示_u在（i, k）点上的值。
从方程（10.2）中是无法直接求出u值，而我们求出u在网格点上的近似值也就可以了，为此，和常微分方程的差分方法一样，将（10.2）中偏导数用差商代替，则有
（10.3）
（10.4）
于是就得到u (i, k)的近似u，所满足的线性代数方程组：
（10.5）
其中
，
用迭代法来解方程组（10.5）。首先将方程组（10.5）化成迭代形式
（10.6）
然后用下边方法取初始值
先用线性插值，注意边界条件
给定区域内部的四个网格点的值（表10.1），然后再用（10.6）算出其五个网格点的值，则得到初始值，如表10.2。
计算
时可将（10.6）与成简单迭代形式：
然后，用（表10.2）中初始值进行计算
| 表10.1 | 表10.2 | |||||||||
|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 0 | 0 | 0 | 0 | k=4 | 0 | 0 | 0 | 0 | 0 |
| 0 | 0 | 0.354 | 0 | 0.707 | k=3 | 0 | 0.151 | 0.354 | 0.453 | 0.707 |
| 0 | 0.25 | 0 | 0.75 | 1 | k=2 | 0 | 0.25 | 0.427 | 0.75 | 1 |
| 0 | 0 | 0.354 | 0 | 0.707 | k=1 | 0 | 0.151 | 0.354 | 0.453 | 0.707 |
| 0 | 0 | 0 | 0 | 0 | k=0 | 0 | 0 | 0 | 0 | 0 |
| i=0 | i=1 | i=2 | i=3 | i=4 | u(0) |
| 表10.3 | 表10.4 | |||||||||
|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | |
| 0 | 0.151 | 0.258 | 0.453 | 0.707 | 0 | 0.134 | 0.243 | 0381 | 0.707 | |
| 0 | 0.182 | 0.427 | 0.583 | 1 | 0 | 0.182 | 0.386 | 0.573 | 1 | |
| 0 | 0.151 | 0.258 | 0.453 | 0.707 | 0 | 0.151 | 0.258 | 0.453 | 0.707 | |
| 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | |
| u(1) (简单迭代) | u(1) (采德尔迭代) |
当我们计算时只要将周围四个点加起来除以4，将所得的值填表10.3 (_i_,_ k_)位置上，这样就得到表10.3，再用这个方法由表10.3计算出，如此下去算到满足所需要的精度为止。<br /> 同样我们也可以用采德尔迭代法来解上述方程组，作法可由左到右，由下到上，从图10.2可知_k_小的先作；对固定_k_，_i_小的先作，于是便有下述迭代公式：<br /> <br />计算时初始值仍为表10.2，先由表10.2中的<br />值计算出并计入表10.4中位置（1.1）上。然后用 _ i-1,k _ _i,k i_+1_,k_<br />表10.2中（_i_+1，_k_），（_i_，_k_+1）位置值和表10.4<br />中（_i_-1，_k_），（_i_，_k_-1）位置上的值相加除以4，将所 图 10.2<br />得的值填入表10.4中（_i_, _k_）的位置上，得出表10.4。如此继续下去就可能计算出，，……直到所需要的精度为止。<br />由上面的例子可以看出用差分法解椭圆型方程需要考虑三个问题：<br />1．选用网格，将微分方程离散化为差分方程。<br /> 2．当网格步长_h_ ®0时差分方程的准确解是否收敛于微分方程的解？<br /> 3．如何解相应的代数方程组？<br /> 关于第3个问题，在第三章中已经讨论，这里就不再重复，下面就第1，第2个问题进行讨论。<br /> <br />**§2 椭圆型方程的差分解法**<br /> 椭圆型方程最简单的典型问题就是拉普拉斯方程<br /> <br />和泊松方程<br /> <br /> 下面以泊松方程第一边值问题为例，来建立差分方程。<br /> 考虑泊松方程第一边值问题：<br /> <br /> <br />**（一） 矩形网格**<br /> 设W为_xy_平面上一有界区域，¶W为其边界，是分段光滑曲线。<br /> <br /> <br /> <br /> <br /> <br /> <br /> <br /> <br /> <br /> <br />图10.3<br /> 取定沿_x_轴和_y_轴方向的步长分别为_h_1和_h_2。<br />。作为坐标轴平行的两族直线：<br /> <br /> <br /> 两族直线的交点（_ih_1, _kh_）称为网点或节点，记为（_x_，_y_）或（_i_，_k_）。以表示所有属于W内部节点集合，并称此类节点为**内点**。以_¶ _W_h_表示网线_x = x_或_y = y_与_¶ _W的交点集合，并称此类的点为**边界点**。令，则就是代替连续区域的网点集合。若两个节点之间距离等于一个步长称此两个节点为邻点。若内点（_x_，_y_）的相邻点都属于W_h_，就称为**正则内点**；否则就称做**非正则内点**。图10.3中打“〇”号的点为正则内点，打“×”号的点为非正则内点，打“• ”号的点为界点。<br />** **<br />** （二）五点差分格式**<br /> 现在假设（_i_，_k_）为正则内点。沿着_x_，_y_轴方向分别用二阶中心差商代替_u_，_u_，则得<br />  （10.8）<br />称（10.8）为差分方程。式中_u_表示节点（_i_，_k_）上的网函数。若以_u_，_f_表示网函数，<br />，则差分方程（10.8）可简写成：<br />  （10.9）<br />利用_Taylor_展式<br /> <br /> <br /><br /><br /> 这四个式子两两相加便有：<br />  （10.10）<br />  （10.11）<br />于是可得差分方程（10.9）的截断误差<br /> <br />其中_u_是方程（10.7）的光滑解。<br /> 由于差分方程（10.8）中只出现_u_在（_i_，_k_）及其相邻四个点上的值，故称之为五点差分格式，其图式如图10.2。<br /> 特别当取正方形网格：<br /> <br />则差分方程（10.8）简化为<br />  （10.12）<br />若_f_ º0，则有<br />  （10.13）<br />对每一个正则点，都可以得到一个这样方程，而对非正则点，一般不用上述方法列方程。<br />** （三）边值条件的处理**<br /> 这里我们只讨论第一边值条件<br />  （10.14）<br />以表示非正则内点集合，¶W_h_表示边界点集合。所谓处理边值条件就是利用（10.14）列出中点的补充方程。通常用下述三种方法解决这个问题：<br />** （1）直接转移法**<br /> 对（_x_, _y_）Î，我们用边界上距离这点最近的点的值作为（_x_, _y_）的值，即<br />  （10.15）<br />** （2）线性插值法**<br /> 图10.4中1点属于，对此点我们取2、4两点沿_x_轴方向作线性插值，1点与4点的距离为_a_，则_u_在这些点上的值有近似关系：<br /> <br /> <br /> <br /> <br /> <br /> <br /> <br />图10.4<br /> <br />于是得到：<br />  （10.16）<br /> 这样在每个边界点上都可列出一个方程，而且不用引进新的未知数，有多少个边界点就列出多少方程，与（10.13）联系起来就得到方程个数与未知数相同的线性代数方程组。<br />** （3）列不等距差分方程**<br /> 把非正则点看成和正则点一样列不等距差分方程逼近泊松方程。由图10.4所示点1Î，如果点1按正则点列差分方程就要用到区域外点6的_u_值，列方程时不用这个邻点上的_u_值而改用这个方向上网格直线与边界的交点4上的_u_值。仍然用中心差商代替偏导数。于是在节点1的不等距方程为：<br />  （10.17）<br />其中，_f_1为_f_在1点的值。这样就可以得到一个方程个数和未知数个数相同的方程组。<br />  （10.18）<br />  (10.19)<br />** （四）差分方程解的存在唯一性问题**<br />** 定理1 （极值原理）**<br /> 假设（1）_u_是定义在网格点上一组值；<br /> （2）_u_≡常数；<br /> （3）<br />则_u_不可能在内部节点上达到正的最大值。<br /> 类似地，如果第（3）个条件改为，则_u_不可能在内部节节点上达到负的最小值。<br /> **证明 ** 用反证法。<br /> 假定_u_在内点取正的最大值。因为_u_≡常数，故必存在如此的内点（_i_0, _k_）ÎW_h_：在（_i_0, _k_）上_u_取正的最大值，且至少有一个邻点上的值_u _< _u __k _。因此<br /> <br />这与假设相相矛盾，所以_u_不可能在W_h_内取正的最大值，即<br /> <br />定理的第二部分可类似证明。<br />** 定理2 差分方程边值问题（10.18），（10.19）的解存在且唯一。**<br /> **证明** 这只需证明相应的齐次问题<br />  （10.20）<br />只有零解就可以了。<br /> 利用极值原理，易证明（10.18）只有零解。<br /> 事实上，若_u_ º常数，则_c = _0，若_u_≡常数，则由极值原理第一部分可知，_u_只能在边界上取到正的最大值，但是在_¶_W_h_上_u=_0，因此，。<br /> 再应用极值原理第二部分可知，_u_只能在边界上取到负的最小值，但是在_¶_W_h_上_u=_0，所以必有<br /> <br />综合上面两个结果，我们有：<br />  证完。<br />** （五）差分方程的收敛性与误差估计**<br />** 定理3 （比较定理）**<br /> 设（1）_V_，_U_是定义在上的两个网函数；<br /> （2） （10.21）<br /> （3） （10.22）<br />则在网格区域（即上），以下不等式成立：<br />  （10.23）<br /> **证明:** 因在W_h_上，等价于<br /> <br />于是便有<br />  （10.24）<br />因为在_¶_W_h_上有<br />  （10.25）<br />即<br />  <br />由极值定理可知，在网格区域上处处有<br /> <br />即不等式<br />  证完。<br /> 现在考虑对于任意一个内部节点<br /> <br />当_h_ ®0时是否有_u_ ®_u_?<br /> 由于微分方程的解_u_ (_x_, _y_)在W_h_内满足方程<br /> （10.26）<br />而<br /> <br />设<br /> <br />则<br />  （10.27）<br />现在令<br /> <br />则有 <br /> 为讨论简便，考虑边界无误差的情形，即<br />  （10.28）<br />根据对_R_ (_u_)的估计式（10.27），我们有<br />  （10.29）<br />为了能利用比较定理证明收敛性，我们构造函数使之满足差分方程：<br /> <br />  （10.30）<br />于是有<br /> <br />又（10.30）式第二式<br /> <br />由比较定理可知：<br /> <br />因此只要对_Q_作出估计，就可以得到误差为的估计。<br />因为<br /> <br />所以希望构造的函数_Q = Q_ (_x, y_)是一个二次曲面，_Q_ (_x_, _y_) ³ 0，并且覆盖_W = W_ (_x_,_ y_)。<br /> 作辅助函数<br />  （10.31）<br />它与_xy_平面交于一个圆<br /> <br />以_r_为半径，（_x_0, _y_）为圆心。我们知道对于任意具有四阶连续导数的函数_F_ (_x_, _y_)，都有<br /> <br />因为（10.31）所给出的函数是二次的，所以它的四阶偏导数为零，即_R_ (_u_) = 0，于是<br />  （10.32）<br />因此，函数一定满足不等式<br /> <br />（因由（10.27），（10.28）式而上面不等式成立）。<br />而<br /> <br />故有<br />  (10.33)<br /> 上面的讨论可归纳为如下的定理：<br /> **定理4 **若(10.7)，(10.7)的解则五点差分格式收敛，且有估计式<br /> 以上定理证明了对任意一个固定节点，当_h_ ®0时且有误差估计<br /> <br />但一般是不知道的，这个估计只是告诉我们：差分方程当_h_ ®0时，（）是收敛的，这样估计为事前估计。<br /> 为了要在实际计算中估计误差，我们要采用事后估计的办法，现叙述如下：设<br /> <br /> 式中表示步长为_h_在节点（_i_, _k_）上_u_的近似值与准确值之差，而_c_是与_h_，_n_无关的常数。再用_h_/2作为步长算得，如此有<br /> <br />所以<br /> <br />即<br />  （10.34）<br />这就是在实际计算过程中估计误差的公式，其中和分别表示用_h_和_h_/2作为步长所算得的差分方程的解，那么<br /> <br />就是误差的近似值。
若有收获，就点个赞吧
0 人点赞
