4. 偏微分方程的数值解法 下载本文

内容发布更新时间 : 2024/6/2 13:36:40星期一 下面是文章的全部内容请认真阅读。

§4 偏微分方程的数值解法

一、 差分法

差分法是常用的一种数值解法.它是在微分方程中用差商代替偏导数,得到相应的差分方程,通过解差分方程得到微分方程解的近似值. 1. 网格与差商

在平面 (x,y)上的一以S为边界的有界区域D上考虑定解问题.为了用差分法求解,分别作平行于x轴和y轴的直线族.

?x?xi?ih (i,j=0,?1,?2,…,?n) ?y?y?jhi?作成一个正方形网格,这里h为事先指定的正数,称为步

长;网格的交点称为节点,简记为(i,j).取一些与边界S接近的网格节点,用它们连成折线Sh,Sh所围成的区域记作Dh.称Dh内的节点为内节点,位于Sh上的节点称为边界节点(图14.7).下面都在网格Dh + Sh上考虑问题:寻求各个节点上解的近似值.在边界节点上取与它最接近的边界点上的边值作为解的近似值,而在内节点上,用以下的差商代替偏导数:

?u1??u?x?h,y??u?x,y???xh?u1??u?x,y?h??u?x,y???yh

图14.7

?2u1?u?x?h,y??2u?x,y??u?x?h,y????x2h2?2u1?u?x,y?h??2u?x,y??u?x,y?h????y2h2

?2u1?2?u?x,y?h??u?x?h,y??u(x,y?h)?u?x,y???x?yh11 注意, 1? 式中的差商?u?x?h,y??u?x,y??称为向后差商,而?u?x,y??u?x?h,y??称为向

hh1前差商,?u?x?h,y??u?x?h,y??称为中心差商.也可用向前差商或中心差商代替一阶偏导数.

2h 2? x轴与y轴也可分别采用不同的步长h,l,即用直线族

?x?xi?ih (i,j=0, ±1, ±2,? ) ?y?y?jhj?作一个矩形网格.

2. 椭圆型方程的差分方法

[五点格式] 考虑拉普拉斯方程的第一边值问题

??2u?2u?x,y?D???0???x2?y2 ??u???x,y?S??式中?(x,y)为定义在D的边界S上的已知函数.

采用正方形网格,记u(xi,yj)=uij ,在节点(i,j)上分别用差商

ui?1,j?2uij?ui?1,jui,j?1?2uij?ui,j?1 ,h2h2?2u?2u代替2,2,对应的差分方程为

?x?yui?1,j?2uij?ui?1,jui,j?1?2uij?ui,j?1??0 (1) 22hh或

1uij?ui?1,j?ui?1,j?ui,j?1?ui,j?1

4即任一节点(i,j)上uij的值等于周围相邻节点上解的值的算术平均,这种形式的差分方程称为五点格式,在边界节点上取

uij??xi*,y*j??i,j??Sh? (2)

??式中(xi*,yj*)是与节点(i,j)最接近的S上的点.于是得到了以所有内节点上的uij值为未知量的若干个

线性代数方程,由于每一个节点都可列出一个方程,所以未知量的个数与方程的个数都等于节点的总数,于是,可用通常的方法(如高斯消去法)解此线性代数方程组,但当步长不很大时,用高斯消去法将会遇到很大困难,可用下面介绍的其他方法求解.

若h?0时,差分方程的解收敛于微分方程的解,则称差分方程为收敛的.

在计算过程中,由于进行四则运算引起舍入误差,每一步计算的舍入误差都会影响以后的计算结果,如果这种影响所产生的计算偏差可以控制,而不至于随着计算次数的增加而无限增大,则称差分方程是稳定的.

[迭代法解差分方程] 在五点格式的差分方程中,任意取一组初值{uij},只要求它们在边界节点(i,j)上取以已知值?(xi*,yj*),然后用逐次逼近法(也称迭代法)解五点格式:

1n??n?1??n??n??n??n?0,1,2,?? uij?ui??1,j?ui?1,j?ui,j?1?ui,j?14逐次求出{uij(n)}.当(i+1,j),(i-1,j),(i,j-1),(i,j+1)中有一点是边界节点时,每次迭代时,都要在这一点上取最接近的边界点的值.当n?∞时,uij(n)收敛于差分方程的解,因此n充分大时,{uij(n)}可作差分方程的近似解,迭代次数越多,近似解越接近差分方程的解.

[用调节余数法求节点上解的近似值] 以差商代替Δu时,用节点(i+1,j),(i-1,j),(i,j+1),(i,j-1)上u的近似值来表示u在节点(i,j)的值将产生的误差,称此误差为余数Rij,即

u?xi?h,yj??u?xi?h,yj??u?xi,yj?h??u?xi,yj?h??4u?xi,yj??Rij

?? 设在(i,j)上给uij以改变量?uij,从上式可见Rij将减少4?uij,而其余含有u(xi,yj)的差分方程中的余数将增加?uij,多次调整?uij的值就可将余数调整到许可的有效数字的范围内,这样可获得各节点上u(x,y)的近似值.这种方法比较简单,特别在对称区域中计算更简捷.

例 求Δu=0在内节点A,B,C,D上解的近似值.设在边界节点1,2,3,4上分别取值为1,2,3,4(图14.8)

解 记u(A)=uA,点A,B,C,D的余数分别为

图14.8

-4uA+ uB+ uc +5=RA uA-4 uB + uD+7=RB

u -4 uc+ uD+3=RC A

uB+ uc-4uD+5=RD

以边界节点的边值的算术平均值作为初次近似值,即

uA(0)=uB(0)=uC(0)=uD(0)=2.5

则相应的余数为:

RA=0, RB=2, RC= -2, RD=0

最大余数为±2.先用δuC=-0.5把RC缩减为零,uC相应地变为2,这时RA, RD也同时缩减(-0.5),新余数是RA=-0.5,RB=2,RC?0, RD=-0.5.类似地再变更δuB=0.5,从而 uB变为3,则得新余数为RA?RB?RC?RD?0.这样便可消去各节点的余数,于是u在各节点的近似值为:

uA=2.5, uB=3, uC=2, uD=2.5

现将各次近似值及余数列表如下:

次调 整 值 数 0 1 δuC = -0.5 2 δuB = 0.5 结果近似值

uA 2.5 2.5 2.5 2.5 RA 0 -0.5 0 uB 2.5 2.5 3 3 第n次近似值及余数 RB uC RC 2 2.5 -2 2 2 0 0 2 0 2 uD 2.5 2.5 2.5 2.5 RD 0 -0.5 0 [解重调和方程的差分方法] 在矩形D(x0≤x≤x0+a,y0≤y≤y0+a)中考虑重调和方程

?4u?4u?4u2?u?4?222?4?0

?x?x?y?ya取步长h?,引直线族

n?x?x0?ih (i, j = 0, 1, 2,?, n) ??y?y0?jh作成一个正方形网格.用差商代替偏导数

1u?x,y???8?u?x?h,y??u?x?h,y??u?x,y?h??u?x,y?h??20?2?u?x?h,y?h??u?x?h,y?h??u?x?h,y?h??u?x?h,y?h?? ??u?x?2h,y??u?x?2h,y??u?x,y?2h??u?x,y?2h???上式表明了以(x,y)为中心时,u(x,y)的函数值与周围各点函数值的关系,但对于邻近边界节点的点(x,y),如图14.9中的A,就不能直接使用上式,此时将划分网格的直线族延伸,在延伸线上定出与边界距离为h的点,称这些点为外邻边界节点,如图14.9以A为中心时,点E,C为边界节点,点J,K为E,C的外邻边界节点,用下法补充定义外邻边界节点J处函数的近似值uJ,便可应用上面的公式.

1? 边界条件为

?uuS??1?P?,??2?P??P?S? 图14.9

?xS