有限差分法-20161116_第1页
有限差分法-20161116_第2页
有限差分法-20161116_第3页
有限差分法-20161116_第4页
有限差分法-20161116_第5页
已阅读5页,还剩19页未读 继续免费阅读

下载本文档

版权说明:本文档由用户提供并上传,收益归属内容提供方,若内容存在侵权,请进行举报或认领

文档简介

1、第4章有限差分法目录第4章有限差分法14.1 有限差分法的原理24.1.1 定解区域离散化34.1.2 微分方程离散化34.1.3 初边值条件离散化44.1.4 解差分方程组44.2 椭圆型方程的有限差分法54.2.1 微分方程的离散化61、五点差分网格62、九点差分网格74.2.2 边值条件的离散化81、第一边界条件离散化82、第三边界条件离散化94.3 抛物型方程的有限差分法104.3.1 一维抛物型方程的差分格式101、向前差分格式112、向后差分格式113、六点对称格式124、Richardson格式124.3.2 高维方程的差分格式131、二维问题差分方法132、分裂格式算法144.

2、4 双曲型方程的有限差分法164.4.1 一阶线性双曲型方程的差分格式171、一阶线性双曲型方程初值问题172、一阶线性双曲型方程初边值问题184.4.2 二阶线性双曲型方程的差分格式191、波动方程及其特征192、二阶方程的差分格式203、定解条件的处理21参考文献23第4章有限差分法刚刚介绍的变分法和加权余量法,都是将实际问题的控制方程转化为积分形式进行求解,而有限差分方法不同,它是用网格覆盖求解区域和时间域,用差分近似替代控制方程中的微分,把求解微分方程的问题改换成求解代数方程问题,进行近似求解。有限差分法是一种经典数值方法,当网格细化或节点较多时,近似解的精度可以得到改进,对于几何形状

3、规则且边界条件简单的问题,差分法有很高的求解精度和收敛性,但用于求解几何形状和边界条件复杂的问题时,它的精度将会降低,甚至产生困难。有限差分法在流体力学和爆炸力学中得到广泛应用,也曾用以求解弹性力学边值问题。目前,弹性静力学一般都不用差分法进行求解了,但是对于弹性动力学等初值问题,仍然要借助于有限差分法予以解决。有限差分法的关键是构造差分格式,保证解的收敛性和解的稳定性。这里简要介绍椭圆型方程、抛物型方程和双曲型方程的差分格式构造方法。4.1 有限差分法的原理差分法是微分方程的一种近似数值解法。它不是去寻求函数式的解答,而是去求出函数在一些网格结点上的数值。具体地讲,差分法就是先将求解域离散化

4、,在这些离散点上用差商近似地代替导数,从而把基本方程和边界条件近似地改用差分方程来表示,把求解微分方程的问题改换为求解代数方程的问题。考虑一维热传导方程的第一类初边值问题Fu一=au(x,0)u(0,t)-2二u-+f(X):x2=(X)=u(1,t)(X,t)'JXW(0,1)(4-1)t0,T其中,复=(x,t)|xW0,1,tW0,T,a为正常数。假定f(X)和9(X)在相应区域光滑,并且在X=0,l满足相容条件,使上述问题有唯充分光滑的解。用差分法求解上述问题时,通常情况下,首先采用矩形网格对求解区域进行离散化,接着,利用差商或Taylor展式对微分方程和初边值条件等连续性方程

5、,在网格节点处做离散化处理,最后解差分方程组,获得所有网格节点处的解的近似值,用近似研究微分方程初边值问题的解。4.1.1 定解区域离散化由于有限差分法求解偏微分方程问题必须把连续问题进行离散化,因此首先要对求解区域进行网格剖分。对于不同的求解问题,求解区域也各不相同。对于问题(4-1)而言,如4-1所示,设N是一给定的正整数,取空间步长Ax=l/N,再取一定的时间步长At=T/M,做两族平行线,用平行直线族把区域分成若干个小矩形。位于区域D内部的网格节点(x:,tJjK.z(简记为(j,k)称为内节点,所有内节点用Dh表示。位于边界上的节点称为界点,所有界点记为rh。图4-1差分法矩形网格剖

6、分4.1.2 微分方程离散化由微分学知道,函数的导数是函数的增量与自变量增量之比的极限,即u(xLX)-u(x)u(x)-u(x-.:x)u二hm二lim;.x-Qx-x-Dxu(x.:x)-u(x)limxDxlim一.x0::xux.:x)-u(x)xxu(x)-u(xxx=lim.x-Dux:,x)-2u(x)ux-x)Txp当IAx很小时,u可以近似地用差商u(x:x)u(x)u(x)u(x-:x)代替,u”可以近似地用二阶差商u(x-:,x)一2u(x)u(x一,x)(3x7代替,从而一个微分方程可近似地用一个差分方程来代替。因此,对于问题(4-1)而言,其一维热传导方程2:uc2;

7、:u一二a-2:t?x可以用方程u(x,t+At)u(x,t)_a2u(x+Ax,t)2u(x,t)+uxAx,t)比(:x)2近似代替,其截断误差为O(At)+0(Ax)2)。差分格式(4-2)称为解热传导方程(4-1)的古典显格式,它是一个两层4点格式,所用到的节点图式如图4-2所示。(j,k+1)(j-1,k)(j,k)(j+1,k)图4-2古典显式差分格式节点图在有限差分法中,为简单起见,通常将节点(占上k)处的函数值i(x,t)简记为u:,依次,节点(xj+1,tk)处的函数值u(xj+Ax,t)简记为uj,节点Oj,tk+1)处的函数值.k+1u(xj,t+&)简记为u:1

8、。这样的话,式(4-2)就可以写成k-1kuj-Ujkkkuj1-2ujuj-1h2k+fj(4-3)4.1.3 初边值条件离散化对于问题(4-1)中的初边值条件,在网格节点处取值,有u:二(xj),u;=u:=。,u;=u:=0,k=1,2,M.4.1.4 解差分方程组对于问题(4-1),将微分方程和初边值条件的离散化方程联立,得差分方程组kk41krkckk-IkkUj十=Uj+rUj2Uj+Uj甲十Tfjkckk£k=rUj1+(1-2r)Uj+rUj,1+if:,11-jj*j(4-4)°kk八Uj=Xj),U°=Un=0,r=a2T/h2,j=1,2,,

9、N_1,k=1,2,,M_1这样,利用k=0层上的数值,即可逐点算出k=1层上的全部离散点处解的近似值,再利用k层上的节点值算出k+1层上的节点值。这样逐层逐点地计算数值解是非常方便的。同样地,若给定的是初值问题,仍采用前面的古典显格式和初值离散方法,也容易写出其差分方程组,但其计算求解范围与初边值问题情形有所不同。此时,利用k=0时间层上的节点近似值只能计算出k=1时间层上关于j=1,2,N-1处的节点值近似值,依次下去,纯初值问题的差分解被局限于图4-1中0A1上的网格节点处。4.2椭圆型方程的有限差分法这里以Possion方程边值问题为例,讨论椭圆型方程的差分格式建立,边界条件处理,以及

10、差分格式的收敛性问题的讨论。考虑Possion方程-2-2U=丹+丹=f(x,y)(x,y)WD(4-5):xy其中D是x-y平面内一有界区域,其边界用r表示,假定它是分段光滑的曲线所组成,如图4-3所示,对求解域D剖分矩形网格,在平面上取一固定点(x0,y0),设x轴方向的步长为h,y轴方向的步长为T,过该点做垂直于x轴和y轴的平行线x=x。+ih,i=0,1,N1,y=y。j,j=0,1,N2,称网格线的交点的心)为网格点。位于D内部的网格点称为内部节点,简称内点,它们的全体记为网格线与区域边界的交点称为边界节点,简称边界点,它们的全体记为h。为了便于后面讨论,将内点分为两类。若某内点的网

11、格线上四个相邻的网格点都是内点,则称它为正则内点,否则称为非正则内点。4.2.1微分方程的离散化1、五点差分网格现在我们来构造Possion方程的五点差分格式。为简单起见,把点(xi,yj)记为(i,j)参见图4-4所示,我们这里采用Taylor展开方法给出差分格式。设(i,j)为内节点,则有Uj1-2uj11+Uj-12J:u;2.xi24jh1fui6uj124:x3606.x+0(h6)(4-6)3uj+1-2duj-132d-y224jh2FUi46jh2u;36077+0(h;)(4-7)44-5)中的J;x利用式(4-6)和(4-7)中沿x和y方向的二阶中心差商直接代替方程(二2u

12、和二。,就得到y1ijj1i+1ii-1i(u;+-2u;+u;j)+(U;2u;+u;)=f;(4-8)h1h2由于在建立差分方程(4-5)时用到了节点(iJ)及其相邻的四个节点,所以称式(4-8)为五点差分格式。差分方程相当于在节点(iJ)分别沿x、y方向的二阶中心差商代替二阶导数的结果。由式(4-6)和(4-7)可知五点差分格式的截断误差为Rj=Auj-uj=0(h2)其中,h=max(h1,h2)。截断误差反映了差分算子对微分算子的相容逼近,也反映了差分解的局部误差状况。(i,j+1)O(i-1,j)iIO41(i,j)(i+1,j)O(i,j-1)图4-4五点差分格式节点图如果h1=

13、h2,即用正方形网格,则差分方程化为i1jjii+1i-1.ihU=7T(4Ui+Uii+Ui+1+Ui+Ui)=fi(4-9)h2一形式特别简单,若f=0,则1U=+Ui刊+Ui,i中+Uij)(4-10)差分解在节点(iJ)的值等于其周围四点值的平均。2、九点差分网格为了提高差分格式的精度,考虑九点差分格式,为推导简单起见,令编号如图4-5所示。这里我们按Taylor展开推导差分格式,得Ju,-Ju";+ui-1l+U,:1+U,1-4u,i2i+1i+1i-1i-1i2h1244jUU,622x二x二yh4+3606Quuu+15+15a6cx4cy24、UU七,-6”二U24

14、二x2二y4(4-11)O(h6)+U+1+U:1-4uj(4-12)360c6Uc6u0(h6)算子1w和0分别表示式中所涉及的节点除外分别是以为中心的正方形和菱形的顶点。从式(4-11)和(4-12)分别得到逼近Possion方程(4-5)的差分方程Wu.=fj(4-13)冤=fj(4-14)它们的截断误差都是qh2)。与五点差分格式(4-14)相比,差分格式(4-13)并无什么优点,但我们可以综合式(4-13)和(4-14),记二uj2 j.1j=组+-WUi(4-15)3 i3二uj,jh22i=Ui+12U上360+O(h6)(4-16)(4-17)略去上式右端的Q(h6)项,得到微

15、分方程(4-5)在节点(i,j)处的又一种差分方程二ujh43.:3f+2360其截断误差为O(h6)。在建立差分方程时共用到九个节点,故称其为九点差分格式。需要说明的是,如果求解域不规则时,要想在处理边界条件时也能得到截断误差为O(h6)的方程存在困难,因此,九点差分格式特别适用于比较规则的区域。4.2.2边值条件的离散化上节讨论了椭圆型方程差分格式的构造方法,在对边界条件进行离散化时,需要注意两个问题:第一,边界点处的离散方程的截断误差要和内节点处离散方程的截断误差匹配;第二,边值条件中含有外法向导数,其方向指向边界外侧。由于曲边界和矩形网格的原因,边界点处的差分方程不可能像正则内点差分格

16、式那样规则,必须小心处理,要特别注意外法线矢量的方向,否则会造成错误。1、第一边界条件离散化第一类边界条件又称为Dirichlet边界条件。下面这里介绍三种处理办法。(1)简单转移若边界点E最靠近非正则内点(i,j),取其截断误差为Qh)。Uij=cp(E)(4-18)(2)不等距格式。如图4-6所示,设P是非正则内点,与它相邻的四个节点记为S、T,其中任何一个点都可以是边界点。Q、R、ru(T)k-ok+Q)(4-19)k、图4-6up(3)线性插值。如图4-7所示,可以用T,Q两点作线性插值,设TP=6<k,则有图4-7Upk+U;Q)(4-20)k、联立所有边界点差分方程和内节点处

17、的差分方程,就得到差分方程组,Poisson边值问题在节点处的近似解。即可解得原2、第三边界条件离散化对于第二类边界条件£U.:n=q(x,y),(x,y)和第三类边界条件(x,y)wT(4-21).:u一+;(x,y)u=,:(x,y),.n其中n是区域C的边界的外法线方向,。之0,但至少在的一部分上仃0。在这样的假设下,边值问题有唯一解。当。三0,式(4-21)便成为第二类边界条件,因此第二类边界条件是第三类边界条件的特例,这里只讨论第三类边界条件。(1)边界点是节点对于边界点是节点的情形,如图4-8所示,边界条件(4-21)可离散为(4-22)1,、1,、(uP-uQ)cos(

18、uP-uR)cosI二PUP=9Ph1PQh2PRPPP当外法线处于轴向位置时,上式可做相应改动。图4-8(2)边界点不是节点对于边界点不是网格节点的情形,如图4-9所示,边值条件可离散为1,(4-23)MS(uM-us)+仃MuM='M其中,可利用内插方法获得。图4-9差分解的误差估计、提高差差分方法中的难点和技巧性强的地方在于边值条件的处理、分解精度及差分方程的高效率求解等问题。4.3抛物型方程的有限差分法4.3.1 一维抛物型方程的差分格式4-1)问题4-1就是一种最简单的一维抛物型方程,按照定解条件的不同给法,可将(的定解问题分为两类:(1)初值问题(Cauchy问题):求具有

19、所需次数偏微商的函数u(x,t),满足方程(4-1)和初始条件:(x,0)=5(x),xw(g,)(4-24)(2)初边值问题(也称混合问题):求具有所需次数偏微商的函数u(x,t),满足方程(4-1)和初始条件:(x,0)=x),xw(0,l)(4-25)及边值条件(0,t)=u(l,t)=0,tW0,T(4-26)在4-1节我们给出了一维热传导方程的混合初边值问题的古典差分格式,这里我们进一步给出该问题的其他差分格式。1、向前差分格式加(x:,tk)在节点(j,k)上用优1,tk)在t方向的向前差商代替j-,用(xtk)在jktjkx方向的二阶中心差商代替2方u(xrtk):x2,即得向前

20、差分格式"k书kkck,kuj_ujuj书2uj+uj-1.f=a+fjhhfk=f(xj,tk)uj0=*=9y;=uN=0其中,j=1,2,,M1.以r=aJh2表示网比,将上式改写成便于计算的形式,可得u;*=ru;书+(1-2r)ujk+ru;+计/(4-27)无需求解线性代数方程组,通过第k层值可以将k+1的值直接表示出来,这种差分格式称为显格式。2、向后差分格式_(4/+)1在节点(j,k+1)上用(xj,tk+)在t方向的向后差商代替一,用jFt(Xi,tk+1)在x方向的二阶中心差商代替jk+1/;:2u(Xj,tk+1),即得向后差分格式kk_1kk+1-k+1k+

21、1Uj中_UjUj邛_2Uj+uj1-a-5七h2fjk=f(Xj,tk)Uj0=*=%xj)Uo=UN0fj其中,j=1,2,,M1.以r=ai/h2表示网比,将上式改写成便于计算的形式,可得.111-叫:十(1+2r)%叫:%+zfj(4-28)令k=0,1,2,则可利用u0和边值条件确定u;,利用uj1和边值条件确定U:3。现在第k+1层的值不能通过第k层值直接表示出来,而是需要求解线性代数方程组(4-28)来确定,这种差分格式称为隐格式。3、六点对称格式将向前差分格式和向后差分格式做算术平均,即得六点对称格式k书k-k+1、k+1.k+1uj-uja棋书2uj+uj-1t"2

22、h"力k=f(Xj,tk)Uj0=*=%Xj)kk八0=UN=0L+f:可将上式改写为k1rk+1(1+r)Uj-万Ujrr=Uj由+(1-r)Uj+Uj+Tfj(4-29)利用Uj°和边值条件便可逐层求到U;。六点对称格式是隐格式,由第k层值计算第(k+1)层值时需要求解线性代数方程组,它联系的网点分布如图4-6所示。4、Richardson格式上述3种格式都是双层格式,这里介绍一种简单的三层格式。如对u(Xj,tk)用(j,k)点的中心差商来代替ft(Xj,tk),得Richardson格式:uk1ukuk2ukukuj一ujuj1-2ujuj-1k=a2+fj(4-3

23、0)2h2j其联系网点如图4-7所示。它与双层格式不同的是,在求第(k+1)层上的值时,要用到前两层的值。由于在x和t方向上都用中心差商代替,所以这种格式的截断误差显然为qJ+h2)。除了上述四种差分格式之外,还可以构造出多种差分格式,但并非每一种差分格式都是可用的。衡量一个差分格式是否有效,主要看这种差分格式的收敛性和收敛速度,稳定性等几方面的因素,当然也希望计算简单,等等。4.3.2高维方程的差分格式1、二维问题差分方法前面介绍的一维常系数抛物型方程的差分求解方法,相关的概念、差分构造方法、稳定性分析等都可简单地推广到多维情况,这里以二维抛物型方程为例来说明这种推广的简单性。相应的区域划分

24、参数记为hph2t,还记r1=a2T/h,2,r2=a2jt/h。对于二维模型问题u=a什2-2-2uuuu+3例)(x,y,t)u(x,y,0)=/x,y)u(0,y,t)=u(1,y,t)=u(x,0,t)xW(0,1)(4-31)=u(x,1,t)=0t0,T可以仿照一维的相应情况,采用直接差分化方法得到对应于上节的各种差分格式。为了书写简便起见,我们用6:和丁分别表示关于xyx和y的二阶中心差分,即、.2ukxiju:1,j也-uu-i这里k可以不是整数。对于三维问题,同样用6;表示关于z的二阶中心差分。例如,用一阶向前差商代替,用二阶中心差商ft、:Ujk/h;=u:j_2u;u:1

25、j/h;5yUjk/h;=(Uik,j_2Ujk+Uik,j.)/h;.2.2分别代替二U2和则得到向前差分格式;x汽1k.1k.122k.22k(U>Uikj=agxWj/h1+6yUikj/h2)(4-32)将Uik,j扩展定义为k时间层上的函数U:(x1,x2)并拓展到整个x1-x2平面上,就可将向前显格式改写为Ujk*(x1,x2)=(1-2r2r2)心人)kk+r(1-h1,x2)Uh(x一"+r(2Uk(x1,x2h2)Uk(x1,x2-h2)虽然,二维问题的差分方法可由一维情形简单推广得到,但是它也有自己的特殊性。二维显格式计算简单,但稳定性条件比一维相应情形要苛

26、刻得多。如一维显格式的稳定性条件是网比rE以,而n维显格式的稳定性条件将是rE以。所以,要维持多维22n显格式的稳定性需要付出较大的计算代价。二维隐格式的最大优点是绝对稳定。在实际计算中,步长选取不受限制,这是非常方便实用的,例如,二维Crank-Nicholson格式是一个截断误差为Qt2+h2)的绝对稳定格式。但也应该看到,解二维隐格式一般要求解五对角方程组,其计算量比求一维隐格式所对应的三对角方程组的计算量要大得多,维数越高,这一现象就越突出。因此,对高维问题,构造绝对稳定且计算量又小的差分格式一直是抛物型方程数值解法中重要的研究课题。2、分裂格式算法从上述二维差分格式的介绍可以看出,由

27、于高维问题隐格式的计算不能采用追赶法,即使采用各种最佳迭代格式,其计算量也大大高于显格式。因此,隐格式在多维情况下并不比显格式优越,需要引进新的计算量较小而又无条件稳定的格式。这里介绍的交替方向就是这样一类绝对稳定的,而每层计算能分解成若干隐格式的计算格式。交替方向隐格式(),简称ADI格式,有时也称为分裂格式,因为通过对一个高维差分格式引进不同的中间变量,可分裂得到不同的ADI格式。由一维交替显隐格式可知,几个简单格式可以复合成一个复杂算法,由此可以得到启发,若能将高维差分格式分解为几个可简单求解的低维,特别是一维的差分格式,就一定能有效地降低高维差分格式的求解计算量。观察式()的Crank

28、-Nicholson格式1uk1-uk.1一二2uk1+2uk12u.k.+2u.k.I,jI,j-22|;x1,j、yI,j'XI,j'yI,j.J将其变形为uk1-2uk1-2uk1-Uk.-2Uk.-2Uk.ui,j2、xui,j2-yui'J-uij2-xui,j2-yuij两边分别加上1r26;6轴丁和1r26;6:u1kj,且分解因式,就有4y,J4y,Ju/rr2"(r2"k1+6;1+6;uk(4-33)T2xJI2y值得注意的是,在某些差分格式中适当增添数值项,或在微分方程中适当增添某些连续项后再做差分离散,是抛物型和双曲型差分方法

29、中常常使用的方法之一。但一定要注意,添加某些项之后,一定要保证新的差分格式对原问题的相容性,不能降低原格式的截断误差和稳定性,要便于数值计算。有关这方面的内容请参见相关专著。式(4-20)能分解为两个简单的一维格式,下面介绍几种可由式生成的分裂算法。1) P-R(Peaceman-Rachfoud)格式作为模型,考虑二维热传导方程的初边值问题,取X,y方向的空间步长均为h=1/N,时间步长为七。最早的ADI格式是Peaceman-Rachfoud(1955)提出的,他们在第k层到第k+1层中引入一个过渡层,即K+1/2层,从第k层到第K+1/2层时,对在过渡层,对在第k层上用二阶中心差商逼近,

30、然后从第k+1/2层到第k+1层,对仍在过渡层,对在第k+1层上用二阶差商逼近,于是得P-R格式。r-21-6X2Ik1/2ui,jr21+6y.2jkui,juik::1/2uQt2+h2)的绝i,jP-R格式消去中间变量u1k12就得式(4-),因此P-R格式是截断误差为i,j对稳定格式。P-R格式是分别在x和y方向交替使用一维隐、显格式的结果,每个时间层只需解两个三对角方程组,计算量小,只相当于解二维问题Crank-Nicholson格式的计算量的1/7o不过P-R格式在计算形式上不能推广到三维情形。2) Douglas格式将式(4-)变形为1-S2.2Ik1u=i,j+r(V+6y2)

31、k,j上式可分裂为k.1/2ui,jk一ui,j=r整理得到Douglas格式:2k1-ui乃2/k甲/25y(ui,jk.1/2,jk+1/2-ui,j2a2.2k=h2、x+yui,jkk.1/2一ui,j=ui,j2a2k.1=声'yui,j2k',yui,jk_ui,j2k,yui,jDouglas格式与P-R格式一样,分别在x和y方向交替使用一维隐、显格式的结果,它的计算量、截断误差和稳定性都同于例如,对于三维模型问题P-R格式,但它可以推广到三维情形。-u2=a:t%2u02-2m二u+歹f(x,y,z,t)其初始条件与边界条件可与问题()式为中类似给出,此略。该问

32、题相应的Douglas差分格1一匚、2xujmFm,iii,m尸r仅:+E+H)uikj,m+*霁,"mk1i,j,m-uk+2/31,j,m2一曳、2uk.2/3.uk2yui,j,mui,j,mh2=2',:uikj1m-uikjmh2zi,j,mi,j,m还可以通过适当引入中间量的方法分裂式()得到一些新的格式,但它们的精度与稳定性及计算特征都类似,在此就不一一罗列了。值得指出的是,若在本节中介绍的抛物型方程中固定时间t=tk,那么它就成为椭圆型方程。这表明,椭圆型方程的离散化方法,可被有效地引进到抛物型方程的离散化处理中,求解椭圆型方程的高效率高精度算法也可以被有效地

33、引进到求解抛物型方程中,本书就不详细讨论这些问题了,读者可根据需要自行拓广应用。4.4双曲型方程的有限差分法然后再介绍二阶线性双曲型方我们先介绍简单的一阶线性双曲型方程差分格式的构造,程的差分格式。4.4.1 一阶线性双曲型方程的差分格式1、一阶线性双曲型方程初值问题首先考虑一阶常系数线性方程的初值问题fu二u_一a=0xR,t.0tAex(4-34)u(x,0)=q:(x)xeR其中a为给定常数,这是最简单的双曲型方程,一般称其为对流方程。由于双曲型方程差分格式的性质与定解问题解析解的性质之间有着密切的关系,因此我们先来研究解析解的一些特性。一般情况下,对于各种双曲型方程,我们做如下的代换:

34、=xat=x-at可以得到方程的通解形式u(x,t)=f(xat)f便-at)其中fl、f2都是二次连续可微函数,其具体形式需要通过定解条件确定。上式有明确的物理意义,以式(4-34)表示的一维双曲型方程为例,如图4-10所示,f仅+at)代表一个以速度a沿x轴负方向传播的行波,称为左行波。同样,f%-at)代表一个以速度a沿x轴正方向传播的行波,称为右行波。区间,x-at,x+atI称为点(x,t)的依赖区间。在x-t平面上斜率为士1的两族直线(x士at)=常数,称为一维波动的特征线。a1)迎风格式迎风格式在实际计算中引起了普遍的重视,从而产生了很多好的方法和技巧。所谓迎风格式的基本思想简单

35、,就是双曲型方程中关于空间偏导数用在特征性方向一侧的单边差商来代替,(4-34)式的迎风格式是a>0(4-35)Ujn17Ujn1-UjnTa<0(4-36)其中t、h分别为时间步长和空间步长。如果差分格式(所用的网格点)与微分方程的特征线走向一致,那么网格比满足一定条件想是稳定的,否则,差分格式不稳定。迎风格式的节点分布如图4-11所示。图4-11迎风差分格式2)Courant-Friedrichs-Lewy条件先分析差分格式的解的依赖区域,然后从差分格式解的依赖区域和对流方程初值问题的依赖区域的关系推导出差分格式收敛的一个必要条件,这个条件称为Courant-Friedrich

36、s-Lewy条件,简称CFL条件。2、一阶线性双曲型方程初边值问题这里考虑下面的一阶常系数线性方程的初边值问题:t«u(x,0)u(0,t)::uca=0;x=(x)=yt)0:t_T,0_x:二(4-37)其与x轴正向的夹角为锐角,因此只a<0,特征线与x轴正向的夹角为钝角,uu;uru:;a=0,k二二0,1,2,;j1,2,(4-38)其中a假定为给定正常数,方程的特征线的斜率为正,能在x的变化区域的左边界上给出边界条件。假若只能在X的变化区域的右边界上给出边界条件,否则将导致定解问题的不适定。1)最简单隐格式用向后差商代替对时间的偏导数,得.a令r=丁,将上式改写为hk

37、Uj(4-39)再由初边值条件q°=*,u;=/,j=1,2,;k=1,2,可计算初区域内所有网格点上的函数近似值,而且由于左边界上的需中已知,因此可计算出u1k+,琮平,无需求解方程组。它的截断误差为Oe+h)。2) Wendroff格式将迎风格式()和()做加权平均,得.k1.kj一jh/u;1.u3(1一)ujk-ujk_1=0(4-40)般地,格式(4-27)的截断误差仍为Qt十h)。为了得到更高精度的格式,通过Taylor112展开,可知当为0=-时截断误差提高为二阶,即Of?22rh2)。此时,差分格式称为Wendroff格式UjkUjk1- rkk+11T7M心(4-4

38、1)同样由于左边界条件已知,此格式实际上也是显式的。可以证明,Wendroff是绝对稳定的。4.4.2二阶线性双曲型方程的差分格式1、波动方程及其特征最简单的线性二阶双曲型偏微分方程是一维波动方程1u:t2=a2四(4-42)二x其中,常数a>0是常数。根据二阶偏微分方程理论,与()式相应的特征方程为dx2-a2dt2=0,利用特征方向dx/dt=均,即可求得两族特征线在研究波动方程的各种定解问题时,它沿x,t的偏导数,那么特征起着重要作用,如果我们用u沿特征的偏导数表示2=a于是方程()化为从而得到()的通解为如果u在x轴的初值为u(x,0)2u2.xu(x,t)22,-u2fu-=&

39、quot;-N-2二1:1:2-2.u-二0,1,22fu-2:2=fi(1)f2(2)=fi(x-at)f/xat)(4-43)=Q(x),ut(x,0)=Q(x),g<xc8则可确定u(x,t)1,、=2广0(xat)1x-at°(x-at)2ax童")d(4-44)这就是大家熟悉的Dalembert(达朗贝尔)公式。区间x0公式(4-31)告诉我们,解u在点(x0,t0)(t0>0)的值仅依赖于初始函数邛0(x),-at0,x0+at0上的值,与区间外的值无关,因此称区间"x0-at0,x0十at0为点(x0,t0)的依存域。其实,区间x0一at

40、0,x0+at0上的初值不只确定u(x0,t0),而且确定了u再以(x0at0,0),(x0+at0,0),(x0,t0)为顶点的三角形域的值,因此称此三角形域为区间x0-at0,x0+at0I的决定域。从公式(4-31)还可以看出,过x轴上点(x0,0)的特征线形成的角形域,为点(x0,0)的影响域,如图4-12所示。图4-122、二阶方程的差分格式1)显格式现在构造()的差分格式,取空间步长h和时间步长t,用两族平行线离散求解域,仿照抛物型方程的直接差分化方法,可以构造二阶方程()的显格式为1k.1kk1-Uj-2UjUj-Th2kuj1这是一个截断误差为q/+h2)的三层显格式。2)隐格式为了得到绝对稳定的差分格式,用第k-1层、k层、k+1层的中心差商的权平均去逼近uxx,得到下列差分格式1k平ckkJa,k+1ck+1k+1(uj_2q*%)=彳f(uj_1-2uj*q中)hh-111+(1-29)(u-2ujk+uM/6(u;:-2u;1+u;J其中0M8M1是参数。当时就是显格

温馨提示

  • 1. 本站所有资源如无特殊说明,都需要本地电脑安装OFFICE2007和PDF阅读器。图纸软件为CAD,CAXA,PROE,UG,SolidWorks等.压缩文件请下载最新的WinRAR软件解压。
  • 2. 本站的文档不包含任何第三方提供的附件图纸等,如果需要附件,请联系上传者。文件的所有权益归上传用户所有。
  • 3. 本站RAR压缩包中若带图纸,网页内容里面会有图纸预览,若没有图纸预览就没有图纸。
  • 4. 未经权益所有人同意不得将文件中的内容挪作商业或盈利用途。
  • 5. 人人文库网仅提供信息存储空间,仅对用户上传内容的表现方式做保护处理,对用户上传分享的文档内容本身不做任何修改或编辑,并不能对任何下载内容负责。
  • 6. 下载文件中如有侵权或不适当内容,请与我们联系,我们立即纠正。
  • 7. 本站不保证下载资源的准确性、安全性和完整性, 同时也不承担用户因使用这些下载资源对自己和他人造成任何形式的伤害或损失。

评论

0/150

提交评论