翻译——抛物问题混合有限元逼近的后验误差估计

抛物问题混合有限元逼近的后验误差估计

摘要 对于抛物问题的混合形式，我们在空间上利用Raviart-Thomas-Nedelec元，在

时间上采用向后欧拉法，得到了一个基于后验误差估计的残差量。这个误差范数是通过通量的能量范数在整个时间上的积分所定义的。为了得到一个最优上界，我们利用了逐单元后处理方法。对于椭圆问题来说这是一个很常用的技巧。最后的误差上界包含了空间离散化误差和时间离散化误差等几项。

1.引言

具有通量高精度的守恒方法通常利用的就是RTN元来减少计算量。在多孔介质流中，即一个典型的抛物压力方程再加上一个波动方程，用到的就是这个方法。这个通量来自于压力方程，然后在波动方程中变成了一个对流场。因此通量的误差需要在全局范数上来控制。

准备工作 对于椭圆问题的混合方法，其后验误差估计已经得到了[4,6,7]。在文献[6]中，得到了一个关于通量的Hdiv,范数。这个Hdiv,范数可以通过分部积分直接计算得到。当要估计通量的L2范数时，通量σ所在的空间要比位移量u所用的空间要大，例如RTN元空间，这是总所周知的困难所在。这个原因就是如果通量空间比位移量的空间要大，那么通过通量和位移量的梯度相关的方程a1σu0而产生的自然残差量会变大。

在文献[15]中，Lovadina和Stenberg得到了一个关于通量的L2范数的后验误差估计，这是基于RTN元的方法，其中利用了一个典型的对于u的近似解的后处理。这个证明基于一个用到后处理近似了的等效法的后验误差估计。在最近的文献[13]中，得到了关于一簇元的后验误差估计。这里的估计和[15]中Lovadina和Stenberg得到的估计很接近，

但是这个证明更一般化，且表明了一个事实，即在残量计算的时候可以利用任何分片多项式来对于位移量逼近。

关于抛物问题混合有限元法的文献并不是太多。有关参考文献包括了书[21]和随后的工作[9]，其中得到了关于热传导方程的先验误差估计。在最近的文献[8]中，得到了一个通量的散度在弱形式下的后验误差估计。这个后验误差估计的技巧揭示了一个固定问题如在文献[13,15]中我们没有把知识拓展到抛物问题上来。另外一方面，对于热传导的标准形式的后验误差估计，有广泛的文献如[10,22]。随后椭圆重构的方法的被运用到来得到抛物问题的后验误差估计上，见文献[2,11,12,16]。在这一类工作中，对于抛物问题来说，与相应的椭圆问题的有关的一个后验误差估计可以运用到误差导数的估计上来。这里我们也运用了这个技巧。

新的贡献 对于抛物问题我们在空间上利用RTN元，在时间上采用向后欧拉法，得到了一个关于通量在能量范数下的后验误差估计。在抛物方程误差界的推导中，我们将椭圆重构的思想运用到了混合形式上，且利用了一个后验误差界，对于相对应的椭圆问题，这是众所周知的最优阶。

本文架构 我们通过一个模型问题开始，在第二部分中介绍有限元法。在第三部分中我们给出了一个后验误差估计。

2.连续问题和有限元法

2.1 连续问题

我们考虑具有齐次Neumann边界条件的热传导方程的混合形式：求u和σ满足

. fin,ua1u0in, in,uu0 n0 on,t0t 0 (2.1) t0t0d1,2,3，这里Rd，是具有边界的多项式空间；fx,tL2，满足fx,tdx0；

Wx。满足u0dx0；满足0a0ax，aW1,，u0H1，a0R，

1,是由函数和其导数几乎处处有界的函数构成的空间。此外，我们按如下规则令为一个衡量问题正则性的参数。令1/21，则相应的椭圆问题为，求v且vdx0, 使得下

式成立

av nv0xx (2.2)

其中L2，且满足dx0。这个问题有解vH1，这里Hs是具有光滑度s的Sobolev空间，见文献[1]，且有

vH1CL2 (2.3)

这意味着，特别的，如果是凸集，则在给定的关于a和f的假定下我们有1。

在这些假定下，对于方程（2.1）存在一个唯一的弱形式满足udx0，使得

σL2[0,T],H0div,，这里

21. ，且H0d,ivvHd,iv:nv0onuL[0,T],H，uL[0,T],H1，这些空间的具体定义见文献[5]。 2.2 连续问题的弱形式

为了得到方程（2.1）的弱形式我们将（2.1）的第一个方程两端都乘以检验函数

wL2，第二个方程两端都乘以检验函数vH0div,，并积分，利用分部积分得到

可以抛物方程的变分方程：即求

u,L2H0div,

满足如下方程：

u2t,w,wf,w wL ,t01vHdiv,,t0 (2.4) a,v u,v0 xux,0u0其中,指的是L2空间的内积。

2.3 离散空间和逼近性质

我们将时间区间0,T离散0t0t1tNT，对应的时间步长ntntn1，

n1,N。对于每一时间层，对空间进行离散。我们令hK是的剖分，K是

剖分单元，hKdiamK是单元的长度，dK是单元K的内含球的最大直径。假设剖分单元是规则的，即存在常数C使得K，hK/dKC。为了计算简单，我们假设所有时间层的区域剖分是相同的。在注释3.3中讨论了怎么将网格拓展到变网格上去。

定义离散的Raviart-Thomas-Nedelec有限元空间和差值函数。对于任意的k0，

Kh，有限元空间RTkK定义为

~dkKkKxkK (2.5)

~其中kK是K上的k阶多项式，PkK是K上的k阶齐次多项式。kK上的差值函数定义如下，见文献[3]第10页的引理3.2：存在唯一的KvkKvHT，

d使得下式成立

如果k1，还成立

FiKvnipkdsvnipkds pkkFi (2.6)

FiTKvpk-1dxvpk-1dx pk-1k-1Td (2.7)

T其中Fi是单元K的边界曲面，1/21，ni是曲面的外法向量。由插值函数的定义可以得到插值的最优误差估计式12：vHK，有

dv-KvL2KChKvHK (2.8)

其中C0是常数。

我们还定义了L2K到PkK的L2K投影。对于pHmT

p-PKpL(K)2mChKpHm(K) (2.9)

其中C0是常数，0mk1。

我们也定义了一个相对应的全局空间和一个全局插值。我们令

kvHdi,v:v|KkK,Kh (2.10) 和投影:Hdiv,KhHKk，通过v|KKv,Kh所定义。对于标

d量我们定义

kpL2:p|KkK,K (2.11) 和投影P:L2k，通过Pu|KPKu,Kh所定义。最后我们给出一个很有用的等式

vPv (2.12) 对于所有的vHdiv,，见文献[3]。对于这些空间更完善的解释，包括inf-sup条件，我们指的是文献[5]中Brezzi-Fortin提到的。但是，我们的后验误差分析并没有显示的用到inf-sup条件。 2.4 有限元方法

我们构造了一个在时间上用向后欧拉法和在空间上用混合元法的格式。为了简便起见，我们令vn是定义在时间tn上的一个函数且tvn是时间上向后欧拉差分算子

tvnvn-vn-1n， n1,,N (2.13)

nn用这些记号我们有如下的数值格式：求uhk和hRTk，对于任意的n1,,N，满

足

u,w,wfnth1nhn,w wPk (2.14)

a和初始条件

ntσh,vhn,vfn,v vRTk (2.15) 0u,wu,w　　　　　　w (2.16)

0hka10,vu0,v=0 vRTk (2.17)

hnn，我们利用线性差值将其拓展到,uh通过离散格式可以求得在tn时刻的近似解σh

tn1,tn中的任意时间上

uhttn1nnuhtntnn1，huhttn1nnhtntnn1 (2.18) httn1,tn，n0,1,N。

3.后验误差估计

基于椭圆重构我们得到了对于误差a0T1/2hL的时空能量范数。这里引入一个

22中间函数w，这个总误差h被分为两部分hwwh，这个w选取为椭圆重构在时间上的离散解。椭圆重构w按照很接近h的w的Galerkin投影来构造。误差的第二部分wh可以利用椭圆问题的后验误差估计来估计。近些年已经得出了这些估计，见文献[13,15,20]。为了估计误差的第一部分w，我们注意到它满足一个具有右端项的抛物方程，则可以再一次的利用椭圆问题的后验误差估计。对于抛物问题，运用这个稳定估计我们得到了一个关于w的后验误差。

注意到对于抛物方程来说这个逼近和标准先验误差分析很相似。在那里是把精确解到有限元空间的Ritz投影作为中间函数而不是椭圆重构的逼近解。总结一下这个思想就是误差的一部分可以通过Ritz投影的先验误差估计来得到，第二部分可以注意到它满足一个具有右端项的抛物方程，则可以再一次的利用Ritz投影的先验误差，参考文献[14]了解更多的细节。 3.1 uh的后处理

为了得到后验误差估计的上界，我们需要对压力进行后处理。我们利用文献[15]中

n,nLovadina和Stenberg提出的后处理方法，记uh为uh的后处理结果，n0,...,N，我们

按照如下定义方式

n,n,定义3.1 （后处理方法） 求uh使得uh|Kk1K满足

n,n PKuhuh|K (3.1)

和

n,n,va1h,v v(I-K uhP)k1(K) (3.2) KKn,

为了计算uh我们只需在每个单元上求解一个小的问题即可，因此总的计算量是非

*常少的。我们利用线性插值定义在0,T上的uh

* uhttn1nn,*uhtntnn1,* (3.3) uh对ttn1,tn，n0,1,N。我们可以注意到对于所有的vk，都有

u,vu,v成立。

*hh3.2 椭圆重构

nn的椭圆重构ωn,n定义为下面问题的解：求ωn,n 满足,uh数值解σhndx0， 下式成立

nωnh1nnaω0nωn0in in  (3.4) on nn由通量的边界条件和Green公式知hdxnhds0。利用这个事实，我们知

nn,uh道上述问题在Hdiv,L2上是适定的。根据椭圆重构的定义知h是ωn,n的

有限元近似解。将tn时刻的解拓展到tn1,tn中的任意时间上

httn1nnhtntnn1, ωhhttn1nωnhtntn1ωnh (3.5)

ttn1,tn，n0,1,N。

3.3 椭圆重构的一个后验误差估计

我们通过对于椭圆重构离散解和离散解之间的差得到一个后验误差估计开始分析。我们通过有界正数C来表示网格大小和时间步长之间的性。

nn,uh引理3.1 令h是方程(2.14)-(2.17)的解，ωn,n是由(3.4)式定义的{hn,uhn}的椭圆

重构，则下面的后验误差估计式成立：

a－1/2ω－nnh2L22L2CKKh2Tnhn,,uh (3.6)

n,tnuh22nn,ChTTth,tuh (3.7) TThn1,,N，其中

2 Tv,wa－1v－w2L2T－1wL2T (3.8) hK2n,*1/21是根据方程(2.3)而来的相应的椭圆重构的正则性参数。这里uh是按照定义

3.1所定义的。

nn证明： 由于h是wn,n的Galerkin投影，我们可以利用文献[13]中对于椭圆问题,uh的混合有限元逼近的技巧来证明所期望的估计。

n首先证明(3.6)式，由ωnh0和nωnhn0我们有

　　　　　　　,ω　　　　　　　u,ωnnnhn,*hnnh　n,ωnhnn,ωnhnn,nωnhn (3.9)

利用等式a1ωnn和(3.9)，我们可以得到

ω＝　aω,ω＝,ωa,ω＝u,ωa,ω＝u,ωaσ,ωua－1/2nnhL2－1nnhnnhnnnh－1nhnnhn,*hnnh－1nhnnhn,*hnKKhnhK－1nhnnhKKKhKKKh2

n,*hn,nωnhKnn,*n＝a1h-uh,ωnhuhn,*,nωnhnKnωnhKKhnn,a1huhLK2nωnhLK2KKhh1/2Tn,uhLK2L2K由于a是有界的，所以

v1/2CavL2L2 (3.10)

从而

a－1/2ω－nnh2L2Ca－u－1nhKKhn,2hL2Kh－1Tun,h2L2K (3.11)

(3.6)式得证。

下面我们来证明(3.7)式。我们引入对偶问题：求,满足

x1 a0x (3.12)

n0x其中H1，H，L2，满足dx0。根据对偶问题的定义，下

d面的等式成立

　tmn,　atmn,a1atmn,tmωn,　　　　　　　tmωn,tmσn,tmσn,a1 (3.13) 　　　　　　　a1tmσn,n,*将(3.12)的第一个方程与tnuh做L2内积，并且利用等式(3.13)得：

u,＝u,ntn,*hntn,*hnn,*＝a1th,tuh,nn,*＝a1th,Ttuh，Tnn,*　　a1th,Ttuh,T (3.14)

利用等式(2.14)知

a1nnth,Kfn,Kh,K (3.15)

n,*再根据uh的定义和等式(2.14)，我们可以得到

nn,*a1th,Ktuh,K0 (3.16)

由Green公式得

tKKhnn,*uh,KKKKKnn,*n,*＝a1th,Ktuh,Ktuh,nKnn,*n,*＝a1thtuh,Ktuh,nTKKhKKhnn,*Cta1htuhKKh (3.17)

L2KTKL2KCKKhhn,*hT1/2tuhL2KhKHL2KhK(K)L2K利用差值误差估计式(2.8)和稳定性估计式(2.3)得到

KLK2KChKL2K (3.18)

通过误差估计式(2.12)我们可以得到下面的结论

 hKKLK2ChKL2K (3.19)

这里我们利用了hKChK，对于所有的0hKh0都成立，这里h0是最大的网格剖分参

n,*数。最后我们令tnuh，得到

tunn,*h2L2nn,*Cta1htuhKKhLK2hKL2K　　　　　　　C由Holder不等式，有

n,tn-uh2L2KKhh1/2Kun,*thL2KhK (3.20)

L2K21nn,ChKtahuhKKh21n,huKthL2K2 (3.21)

L2K从而引理4.1得证。 □

现在我们准备给出主要定理。

定理3.1 令,u是方程(2.1)的解，h,uh是有限元格式(2.14)-(2.17)的解，uh是uhn按

照定义3.1所定义的后处理，则

Tt0a1/2h2L20,dtu0uhN2L22ChT2Th0,uh0,TTh2L2Tn,n ChhT2tuhhfnn1TThN2 ChThn,uhn,

n1TThN2 ChhT2Tthn,tuhn,n1TThN2L22L2 (3.22)

3n Chthn1Ntn Cn1tn1ffndt其中1/21是根据方程(2.3)而来的相应的椭圆重构的正则性参数，

2T(v,w)a－1v－w2L2(T)－1hK[w]L2(T)。

2注释3.1 定理3.1中的给出的误差上界包含了7项。第一项和第二项衡量的是初始条件和它的逼近间的误差的影响。接下来的三项衡量的是空间离散化的影响，第六项衡量的是时间离散化的影响，最后一项衡量的是右端项f和它在离散点间的差在时间上的影响。这个上界的形式和文献[12]中的上界类似，在[12]中，对于热传导方程用到了相同的时间离散和一个标准的空间离散。

证明：我们首先将误差分成两部分

** euuhuuhe1e2 (3.23)

 hωωh12 (3.24) 由三角不等式知

2 a1/L2/2a11L2a21/22III (3.25)

L接下来我们将分别估计这两项。

I. e1和ε1满足下面的方程ttn1,tn，n1,,N

,,v,ve,vω,vf,vn,*nn,*nuhhfn,vu,vh,vffn,v,vωh,v (3.26)

令ve1，我们对(3.26)左端两项关于时间从0到Tt积分再利用Green公式得到

1,e1dt0e0TtTt1d21122e1dtdxe1TtL2e10L2 (3.27) 2dt22

,edt,edt,adt1011011011TtTtTtTt0a1/212L2dt (3.28)

由(3.26)-(3.28)式我们可以得到

Tt2211/2e1Tt2a12dtLL02Ntn21n,*nhe10L2uhfn,e1dtt2n1n1 (3.29)

Ntntn1n1Ntntn1ff,e1dtnn1Ntntn1,e1dtun,*hn1,edtnh1I.II.III.IIII.IVI.V

下面我们开始估计I.II.V各项.

0,*

再利用三角不等式得 I.I. 对该项加一项减一项uh

12u00L222110,*0,*2u0uhu 　　 (3.30) h0L2L22210,*22200,　　u0uhChK,uh2KhL2TThI.I其中最后一步用到了引理3.1中m0时的结论。

warz不等式和差值估计式(3-2-8)得 I.II. 由于e1L2，知PTe1PT。利用Cauchy-SchII　n1NNtntn1tnuun,*hnhfn,e1dtnhfn,e1PTe1dt

　n1Ntn1tnn,*h (3.31)

Cn1tn1KKhhKn,*ntuhhfnL2Ke1L2KL2Kdt这里我们用到了插值估计式e1Pke1L2KChKe1。由于e1a11，我们得到

Ntnn,*nII　ChTtuhhfntn1n1TTh

L2TTt0a1/2ε12L2Te1L2Tdt (3.32)

1/21/2I.III. 因为e1dx0，利用Poincare-Friedrich不等式和a的有界性，得到

 e1L()2Ce1L2Ca112Ca1/212 (3.33)

LL利用Cauchy-Schwarz不等式得 NIIItnn,e1dtn1tn1ffN

tnffne1n1tn1L2L2dt 21/21/2NCtnffndtNt2n1tn1L2nen1tn11L2dtI.IV. 利用Cauchy-Schwarz不等式我们可以得到 NI.IVtntun,*n1h,e1dt

n1N tn n,*h n1tn1tuL2e1L2dt利用引理3.1的第二个不等式和不等式(3.33)，则有

1/2/2 I.IVCN22n1ε/22hhTTt,hun,th2n1KKdthTt0a1L1 I.V. 利用椭圆重构的定义(3.4)得到

I.VNtnnh,e1n1tn1dtNtnnhh,e1dt

n1tn1Nt nnhdtn1tn1h,e11/2CNtTt1/2nn2hhn1tnL2dt10aε21L2dt1/2我们注意到ttn1,tn

(3.34)

(3.35)

(3.36)

(3.37)

n hhtntnnhn1htntt n h (3.38)

所以

n1Ntntn1hnh2L23ndtCnthn1N2L2 (3.39)

得到第I.V项的估计式

N3n I.VCnthn12dtL21/2Tt0a1/2ε12L2dt1/2 (3.40)

利用I.II.V的估计式(3.30)、(3.32)、(3.34)、(3.36)和(3.40)，还有选取适当的的Young不等式aba2/2b2/2我们得到

0,　Iu0uhN2LChK2K2h,u0h0,KKh C CNNn1KKh2hKhuf2hKthn,tuhn,2L22L2n,thnhn2L2K2hK (3.41)

n1KKh3n Chthn1Ntn Cn1tn1ffndtII. 由(3.3)和(3.5)的第二个等式知

2ωhttn1nωnhntntnωn1n1h (3.42)

利用引理4.1的第一个不等式得，ttn1,tn，n1,,N

a1/22

L2na1/2ωnhL2n1a1/2(ωn1h)1/2L21/22nn,CT,uhhKKh2n1n1,CT,uhhKKh (3.43)

因此第II项的上界是

IICNn1KKh2nn,,u (3.44) nThh

则通过不等式(3.25)、(3.41)和(3.44)我们得到了定理的证明。 □

注释3.2 在方程(2.17)中如果将a替换成与时间有关的atn。误差分析将会改变，第二项在最终的估计上的数值震荡将会增加，那是因为a和atn在每个时间区间上的差是成比例的。

注释3.3 在本文中我们考虑的是一个静态网格，在时间上没有变化。很自然的拓展便是允许在不同的时间层有不同的网格剖分。特别的，为了提高逼近解得质量而构造自适应算法时这是非常有用的。混合元方法可以很容易的拓展到变网格上来，即将上一个时间层的网格上的解作为当前时间层的新网格的投影。在误差分析中会引入一个新的项来衡量这个插值的误差。这部分的分析在文献[12]中的第3.7节可以看到。