抛物型方程的计算方法

1、分类号：O241.82侠、力冲展丸学SHAANXINORMALUNIVERSITY本科生毕业论文（设计）题目：一类抛物型方程的计算方法作者单位数学与信息科学学院作者姓名专业班级2011级数学与应用数学创新2班指导教师论文完成时间二Q五年四月一类抛物型方程的数值计算方法（数学与信息科学学院数学与应用数学专业2011级创新2班）指导教师摘要：抛物型方程数值求解常用方法有差分方法、有限元方法等。差分方法是一种对方程直接进行离散化后得到的差分计算格式，有限元方法是基于抛物型方程的变分形式给出的数值计算格式.本文首先给出抛物型方程的差分计算方法，并分析了相应差分格式的收敛性、稳定性等基本理论问题.然后，

2、给出抛物型方程的有限元计算方法及理论分析.关键词：差分方法，有限元方法，收敛性，稳定性NumericalcomputationmethodsforaparabolicequationYanqian(Class2,Grade2011,CollegeofMathematicsandInformationScience)Advisor:NiehuaAbstract:Thecommonmethodstosolveparabolicequationsincludedifferentialmethod,finiteelementmethodetc.Themainideaofdifferentialmeth

3、odistoconstructdifferentialschemesbydiscretizingdifferentialequationsdirectly.Finiteelementschemeisbasedonthevariationalmethodofparabolicequations.Inthisarticle,wegivesomedifferentialschemesforaparabolicequationandanalyzetheirconvergenceandstability.Moreover,thefiniteelementmethodandthecorresponding

4、theoreticalanalysisforparabolicequationareestablished.Keywords:differentialmethod,finiteelementmethod,convergence,stability1绪论1.1 引言自然界里中热的传播，溶质在液体中弥散，多孔介质中渗流等随时间发展的现象和过程，都可以用抛物型方程来描述.因此，抛物型方程是刻画自然界的一类很重要的方程.然而，很多的方程我们并不能求出它的精解确，或者表达式过于复杂，所以需要采用数值方法去计算它们的近似解.抛物型方程最基本的计算方法当属有限差分法1,通过离散化便可得到计算格式，该方法构造

5、简单，易于操作.但是在处理一些复杂的边值问题时计算会很复杂，因此我们需要探讨一些新的处理手段.有限元计算方法起源于椭圆型方程的计算，它将求解椭圆型方程的解转换为求解其变分形式的解1,从而极大地丰富了偏微分方程的计算手段.正式由于其在椭圆型方程计算中的巨大优势，以及抛物型方程与椭圆型方程的密切联系，所以该方法很自然的被推广到了抛物型方程初边值问题的计算上4.本文系统的总结了一类抛物型方程的计算方法，包括有限差分法和有限元方法.并且通过数值算例给出了两类方法的一个比较.为此,本文需要先给出一些基本的分析知识作为研究该问题的基础6,7,下来就给出了抛物型方程的变分形式，这个是构造有限元计算格式的基础

6、，在此基础上，给出了有限元计算格式并讨论了其收敛性和稳定性.1.2 准备知识(1.1.1)抛物型偏微分方程是一类典型的发展方程，其一般形式如下:其中u(x,t)是空间自变量线性椭圆型微分算子，即-L(u)=f(x).tX=(X1.Xn)和时间t的未知函数，L是关于空间变量的(n/2n内L三士工aj+£bi,i,jCXiGXjicXi其系数aij,bijd右端项f为自变量X=(*.Xn)的实函数，且在方程(1.1.1)的定义域jRn中满足椭圆性条件nn2'aij(x)ij-：(x)vi(1.1.2)i,j1i1«(x)>0,V(.1)wRZt,Vx-Q当L是非线

7、性椭圆型微分算子或者f是u的非线性函数时，则称相应的抛物型方程为非线性的.(1.1.3)(1.1.4)(1.1.5)(1.1.6)下面给出抛物型方程的定解条件：初值条件，不妨设初始时刻t=0,则u(x,0)=u°(x),-xc第一类边值条件：u(x,t)=uD(x,t),-x,/：.:,-t0第二类边值条件：:u,(x,t)=g(x,t),一x:.，一t0.:v第三类边值条件：.二u(一：u)(x,t)=g(x,t),_x-:,_t0ft其中uD,g,a是u(x,t)的已知函数40,且至少在一部J&aJfe叫，v为EC的单位外法向量.2,有限差分法本章将给出抛物型方程最基本的

8、计算方法一有限差分法。我们以一维热传导方程为例，给出其差分格式并讨论其收敛性，稳定性等基本问题.本章内容主要引用文献1.用差分法计算抛物型方程的初边值问题时，可以先考虑在区域。上引入空间网格，例如在直角坐标系中采用平行于坐标轴的等距离直线族形成的矩形网格，其次，将定义在QXR+上的函数u(x,t)替换成定义在空间网格节点集上的离散函数U(t);然后，用适当的差分格式将微分算子L替换成差分算子Lh,这一过程称为半离散化.对由半离散化得到的常微分方程初值问题，再进一步对时间离散化，选用适当的求解常微分方程初值问题的数值方法，就得到求解抛物型方程的初边值问题的全离散化格式.接下来，将按照这一处理思路

9、对热传导方程的差分计算格式进彳T探讨.2.1 差分格式考虑一维热传导方程：二二2四=a±4+f(x),0<t<T(2.1.1)t:x其中a是正常数，f(x)连续。下面给出两类定解条件：第一，初值问题：求可微函数u(x,t),满足(1.1.1)和初始条件：第二，初边值问题：求可微函数u(x,t),满足方程(1.1.1)和初始条件：u(x,0)=#(x),0<x<l(2.1.3)以及边值条件u(Qt)=u(l,t)=0,0<t<T(2.1.4)现在考虑边值问题(1.2.1),(1.2.3),(1.2.4)的差分格式.取步长空间h=%和时间步长丁二%,其

10、中J,N都是自然数.用两族平行直线x=xj=jh(j=0,1.J)和t=tn=”6=0,1一.3)将矩形域G=0<x<l;0<t<T分割成矩形网格，网格节点为(xj,tn).以Gh表示网格内点集合，即位于开矩形G的网点集合；Gh表示所有位于闭矩形G的网点集合；Th=GhGh是网格界点集合.其次，用u；表示定义在网点(xj,tn)上的函数0WjWJQWnWN,用适当的差商代替方程(1.2.1)中相应的偏微商，便得到以下几种最简单的差分格式.2.1.1 向前差分格式考虑n-1nnnnUj-ujUj2Uj+Uj，(2.1.5)=afj，i其中fj=f(xj)，j=1,2.JT

11、u0=j=(xj)；=u：=0n=0,1N-1以二a%2表示网比，将(2.2.5)整理成易于计算的形式，使得第n层值，即上标为n在等式右边，第n+1层值在等式左边，则可得到u；41=run七十(12r)un+run+小(2.1.6)这样的话，又(2.1.6)取n=0,利用初值条件/=%和边值条件口01=山=0可计算出u1.再将n=1的值带入计算，从而就可逐次迭代计算出所以的u；,并且视其为精确解u(xj,tn)的近似，由于第(n+1)层的值通过第n层值明显表示，无需求解线性代数方程组，如此差分格式称为显示格式.下来给出这种计算格式的误差分析：；：2uLu-a2t；x(1)nLhUjn：；1nn

12、nnuj-ujuj1-2ujujja72h11二2c显然截断误差(2.1.7)2220(.2h2)=0(.h2)2.1.2 向后差分格式将上式改写为n1nuj一uj一ah2uj=j=(Xj),u；=uJ=0,fj，(2.1.8J=1,2.J-1,N=0,1.N-1一ru：+(1+2r)u：4ru：=u：+fj.(2.1.9)显然，第（n+1）层的值不能用第n层值明显表示，而是由线性代数方程组（2.1.9）确定，这样的差分格式称为隐格式.n1nUj-uj-aTnFn1n"1uj2ujh2则截断误差为R：(u)=L(；)u(Xj,tn)-Lu；=-A"(2.1.10)0(2h2

13、)=0(h2).此外，还有六点差分格式以及Richardson格式，具体可以参见文献1,都是简单的抛物型方程差分格式.2.2差分格式的稳定性与收敛性差分格式的稳定性概念见文献1，此处本文只给出相关的稳定性定理及实例分析.2.2.1 判别稳定性的直接估计法(矩阵法)命题11(必要条件)以P(C)表示矩阵C(e)的谱半径，则差分格式稳定性的必要条件是存在与T无关的常数M使P(C)<1+MT(P(C)<1+O(z)(2.2.1)命题2(充分条件)若C(e)是正规矩阵，及C和它的共腕转置C成绩可交换：CCC,则(2.2.1)也是差分格式稳定的充分条件.推论1若S是对称矩阵，C(。是矩阵S的

14、实系数有理函数：C(t)=R(S),则差分格式稳定的充要条件是myKl)M1+MT,其中*是S的特征值。(只需注意R(S)是实数和矩阵S的四则运算)下面引出两个例子，来具体分析有限差分计算的稳定性判定：例1对向前差分格式(以下设(2.2.4)中的l=1),则C=(1-2r)I+rS,*=12r+2rcosjnh=14rsin2jh,为使<1+Mr或1-MrE*=1-4rsin2MhE1十Mr当且仅当jj24rsin2圮M2+MEj=1,2,J-1，从而4r<2,r<1.所以向前差分格式当r<2221 时稳定，当r_1时不稳定.22.2.2收敛性与敛速估计如最简差分格式，

15、考虑热传导方程的初边值问题：UcLu=-a.2二u二f(x),0二x：l,0<t<T,(2.2.2u(x,0)=(x),u(0,t)=u(l,t)=0.相应的差分格式为Lun=fj,j=1,2,.,J-1,n=0,1,.,N-1,J(2.2.3)u0"：(x)U=un=0.其向量形式如U"=CUn+叭¥,其中C为增长矩阵.那么差分逼近的截断误差R：(u)=Lhu(Xj,tn)-Lu；(2.2.4)u(x,t)是0MxMl,0SWT上的任一充分光滑函数，称差分算子Lh是边值问题(2.2.2)的相容逼近，如果相容条件成立，其中Rn是分量为Rn（u）的向量，

16、II是RJ中的范数.先对差分解作出某种估计：将（2.2.3）的解分解为nnnUjtVjWj其中v；满足零初值和非齐右端方程:Vn1=C(.)VnAF(V°=0)而w；满足非零初值和齐右端方程:Wn1=C(.)Wn(W0=U°)其中Vn,Wn依次为以V；,w；为分量的向量。若差分格式按初值稳定，则亦按右端稳定，于是有常数Ki,K2,使Vn<KiF,Wn<K2W0=K2U0从而(2.2.6)Un|<K(|U0|+|F|),K=max(Ki，K2)现在估计差分解的误差：设u（x,t）是热传导方程（2.2.2）的解，u；是差分方程（2.2.3）的解，误差e；=u（

17、Xj,tn）-u；=u；-u；那么R；(u)=Lhu；-Lu=Lhu；-fj=Lhu；-Lhu；=Lhe；,即误差e；满足差分方程：Lhe；=R；(u),e：=0,j=1,2,.,J-1,其向量形式为En1=C()EnARn,这里En,Rn依次为以e；,Rjn为分量的向量.由估计式（2.2.5）得若相容条件(2.2.4)成立，则-un|=0,其中un表小以U(Xj,tn)为分量的向量.证明了如下：定理11:若差分方程满足相容条件，且按初值稳定，则差分解收敛到热传导方程的解，且有误差估计式(2.2.7).推论1川：当网比r414时，向前差分格式的解有收敛阶O(e+h2)。对任何网比r>0,

18、向后差分格式的解有收敛阶O(i+h2),六点对称格式的解有收敛阶0(2h2).3,有限元计算有限元计算方法产生于椭圆型方程的计算网，其优越的计算性能使得很多学者开始探索将其用于发展方程的计算之中，文献3给出了这方面的具体研究.本章给出热传导方程的有限元计算格式，首先给出有限元计算的基本理论，之后建立热传导方程的变分形式，从而在此基础上给出有限元计算格式.3.1 基本理论本节先给出有限元计算的基本数学理论，包括索伯列夫空间初步和初边值条件下解的存在性与正则性1,6,7, Sobolev空间构造如下Sobolev空间:Wm'p(J)=vLloc。)：/|m,DvLp(c)赋范

19、数忖皿(。)=|DPvP)；p,1<p<«(3.1.1):加p,ZesssupDmv,p=o0(3.1.2)p宝Qq定理3.17上述赋范数的Sobolev空间是Banach空间.特别的，当p=2,Wm,2(C)记为Hm(C),引入下述内积(u,v)mQ=£(D%,D°v)(3.1.3)：'.m定理3.2口Hm(C)是Hilbert空间.另外，用记号C(0,T;X)表示映射族v(t):(0,T)TX),其中任一v关于tw(0,T)按空间X的度量是连续的.类似的记号还有L2(0,T;X),U°(0,T;X).3.1.2 解的存在性，正则性

20、先给出如下定理9,假定u0wHstfwL2(0,T；Hs),s之0为整数，则初边值问题(2.2.1)-(2.2.4)存在唯一的解u(x,t),满足uL2(0,T；Hs2),UtL2(0,T;Hs)t20|2(£|d%(,7)|dT<C0(|u0|Hs+I"Hs和估计式2f(d)|sdT).(3.1.4)Hs3.2 有限元计算格式本节讨论初边值问题(1.1.1)-(1.1.4)的有限元近似.首先，给出其变分形式.记H0=vwH1:v=0当xw!7.设fwC(0,T;H),用函数v(x)wH；与(1.1.(1) 1.1.4)的两端做内积，有(型,v)十(Au,v)=(f,

21、v)利用格林公式，；:tdd.dd-(Au,v)=(二:ajcuv)dx-.(.-ag-)vdx,1 iijd二为二为二yj二二xj因v=0当x-T,故上式右端第二项为零，引进双线性泛函dd_a(u,v)=(-a。cuv)dx,1yjw%K那么可得，:u1(,v)+a(u,v)=(f,v),VvH0ftu(；0)=u°()(3.1.5)则称(3.1.5)为问题(1.1.1)-(1.1.4)的变分形式.有限元方法的第一步，就是将求解区域C刨分成有限个互不重叠的子区域，称其为单元.用hmax表示刨分中单元的最大直径，记相应的刨分为I=Uh,0<h<hmax)，其代表了一个刨分

22、族.以P(K)和仃(K)表示h中单元的外接圆和内切圆的直径.如果存在不依赖于h的常数C,使得P(K)/<r(K)<C,则称刨分族I='h,0<h<hmax是正则的.第二步是构造H0的一族有限维子空间Sh,0<h<h。,要求它具有如下逼近性质：对于某一整数r之2,有infv-Vh|hv-Vh1<Chsvs,1<s<r,Vh-Sh对任意VWHscH0.(3.1.6)一般，&是通过在刨分上作分片多项式差值的方法去构造的.由索伯列夫空间插值理论7,当刨分族h=九0<卜<欧为正则时，由T；上所有属于C(C)的,次数Er-1

23、的分片多项式组成的子空间Sh满足条件.对于给定有限元空间&UH0后，初边值问题(1.1.1)-(1.1.4)的有限元近似定义为：求映射Uh(t):=0,TT&,它满足(-,W)a(Uh,vh)=(f,/),一叫S7、ct(3.1.7)Uh(0)=u0其中u；wSh是函数u0(x)的某个近似.这里，可以看出(3.1.7)是变分问题(3.1.5)的一个近似.设%(x)N%是空间6的一个基底，则近似问题(3.1.5)又可以表示为：求解函数表达式NhUh(x,t)-%：j(t)j(x)j1即确定其中系数j(t)NL使得NhNh工;j(t)(j,i)':j(t)a(j,i)=(f

24、,i)j1j1i=1,2.Nh,(3.1.8)%(0)=小j=1,2.Nh,.diNh其中叫=,勺为/(x)=£j*(x)的系数.dtJjT所以，(3.1.8)是以%(t)NZ为未知函数的一个一阶常微分方程组.由于此处时间t仍然是一个连续变量，所以说(3.1.8)是问题(1.1.1)-(1.1.4)的一个半离散格式.引入一些记号：M=(mj)Nh5网=(一,j),B=(bj)NhNhM=a(i,)a(t)=(t),.二Nh(t)TF=(."),=(f,i),r=(r1.rNh)T,其中M是一个Gram矩阵，它是非退化的，从而可以将(3.1.8)改写为矩阵形式，即«

25、(t)+MB"(t)=MF,0<t<T,(3.1.9)：(0)=r.由常微分方程基本理论可知，初值问题(3.1.9)对于任意F和r存在唯一的解口(t),从而近似问题存在唯一解Uh(X,t).3.3收敛性分析和误差估计本节介绍抛物问题有限元的理论分析，将通过能量估计证明有限元法的收敛性和近似解的误差估计.首先给出方程dd_-'、(aij(x)c(x)w=f,x凸jT为因(3.1.10)w=0,x三.利用上一节的有限元空间ShUH0,可以得到(3.1.10)的一个有限元近似：求Wh亡Sh,使得a(Wh,Vh)=(f,Vh),VVheSh(3.1.11)其中a(；)是上

26、节中所定义的的双线性泛函.因为a(,)正定，通过Lax-Milgram定理6可以得到，有限元方程(3.1.10)存在唯一的解Wh-Sh.下面，本文将不加证明的给出关于抛物型方程有限元计算方法的收敛性分析和误差估计，具体证明过程可以参考文献3,4.定理1假定空间3具有逼近性质(3.1.6)，边值问题(3.1.10)的解wWHscH；,则由(3.1.11)所定义的近似解Wh是收敛的，并且满足|wh-w|+h|wh-w|1<Chs|Ws,1«s«r.(3.1.12)定理2(Gronwall不等式)设y(t)于0,T)连续且满足(3.1.13)ty(t)-y0,()y()d,

27、0其中九之0且九(7产L1(0,T),则ty(t)<y0exp(z(x)dT)(3.1.14)-0对于半离散非齐次方程(3.1.7),此处将给出其解所满足的一个先验估计式.定理3半离散方程(3.1.7)的解5(t),0MtwT满足ttNh(t)|+.|UhI*-Uh(0)|+*f巫.(3.1.15)估计式(3.1.15)给出了半离散问题(3.1.7)的适定性.特别的，当f=0时，由上面定理可知：|Uh(t)|4/(0)|,口>0,即半离散齐次方程的解在L2。)范数意义下是稳定的.下面给出非齐次半离散问题(3.1.7)解的误差估计：定理4假定空间Sh，0<h<h。具有逼近

28、性质(3.1.6)，并且近似初值U0满足u0-u0|<Chr|u0|r,(3.1.16)则半离散问题(3.1.7)的解uh(t)满足uh(t)-u(t)|<Chr|u0|r+J0|ut(T)|rdi(3.1.17)至此，本节就给出了抛物型问题有限元法的收敛性和误差估计，也就意味着本论文的目标，抛物型方程(热传导)的有限元计算方法完整的给出.4,总结将有限元方法推广到抛物型方程，可以丰富这类方程的计算手段，同时也可以有效的解决传统的差分格式存在的缺陷.近年来，自适应有限元计算方法成为这一领域的研究前沿9-12,它在传统的有限元方法的基础上增加了“智能”的因素，使得该方法可以更好的刻画

29、客观世界.本文详细的给出了抛物型方程的简单差分计算格式和有限元计算格式.通过构造计算格式，得出数值计算的矩阵方程，以及给出了收敛性，稳定性，误差估计等一些基本的理论，比较全面的介绍了抛物型方程的一些基本计算方法.更为详尽的资料可以参考文献4.当然，本文只给出了抛物型方程计算的半离散格式，并没有对时间变量t进行离散化，如何给出抛物型方程的全离散格式，以及探讨抛物型方程的其他计算方法，包括偏微分方程的应用，例如现在有学者通过凸几何的相关内容如Brunn-Minkowski不等式来研究偏微分方程的凸性等12,14,这些都是自己今后学习的方向.参考文献1 李荣华.微分方程数值解法M.北京：高等教育出版社,20022 李志平.偏微分方程数值解讲义M.北京：北京大学出版社,20093 黄明