王劉彭, 易年余
(湘潭大學 數(shù)學與計算科學學院,科學工程計算與數(shù)值仿真湖南省重點實驗室,湖南 湘潭 411105)
有限元方法是求解偏微分方程的一種行之有效的數(shù)值方法,廣泛應用于科學與工程計算的各領(lǐng)域.它的基本思想是分片函數(shù)(多項式)逼近與變分原理.隨著泛函分析在偏微分方程逼近理論中的成熟應用,有限元方法與偏微分方程逼近完美結(jié)合,有限元方法的數(shù)學理論得到了蓬勃發(fā)展.自20世紀60年代開始,逐漸建立了有限元方法的先驗誤差估計理論[1-3]和后驗誤差估計理論[4-7].
有限元誤差估計理論主要是基于誤差方程和Galerkin正交性.設Th是區(qū)域Ω的一個正則剖分,Vh是定義在Th上的有限元空間.在能量范數(shù)意義下,注意到有限元解是有限元空間Vh中偏微分方程解的最佳逼近,即
特別地,取vh為u在Vh上的插值uI,可將有限元誤差能量范數(shù)估計轉(zhuǎn)化為插值誤差估計
‖u-uh‖E≤‖u-uI‖E.
關(guān)于有限元誤差的L2范數(shù)估計,一般是通過引入對偶問題,利用Aubin-Nitsche技巧將誤差的L2范數(shù)轉(zhuǎn)化為誤差的能量范數(shù)控制,進而得到有限元誤差L2范數(shù)的一個上界
‖u-uh‖≤Ch‖u-uh‖E≤Ch‖u-uI‖E.
設Ω?Rd(d=1,2,3)為d維空間中的一有界區(qū)域,其邊界?Ω是Lipschitz連續(xù)的.記Th為區(qū)域Ω上的一正則剖分,Vh?H1(Ω)是定義在網(wǎng)格Th上的連續(xù)分片線性有限元空間.
對任意的函數(shù)u∈L2(Ω),其在有限元空間Vh上的L2投影記為Pu∈Vh,滿足:
(u-Pu,vh)=0,?vh∈Vh.
(1)
對于投影Pu,有
‖u‖2=‖Pu‖2+‖u-Pu‖2,?u∈L2(Ω).
(2)
易知
‖Pu‖≤‖u‖,?u∈L2(Ω),
(3)
當且僅當u∈Vh時,等式成立.
基于估計式(2),可以合理假設:存在常數(shù)0<θ<1,使得
‖Pu‖≤θ‖u‖,?u?Vh.
(4)
顯然,在L2范數(shù)意義下,Pu是u在有限元空間Vh中的最佳逼近,即
(5)
因此,投影誤差可由誤差‖u-vh‖來估計.下面定理給出投影誤差‖u-Pu‖與誤差‖u-vh‖的比較結(jié)果.
定理1設u∈L2(Ω),Pu∈Vh為u的L2投影且滿足假設(4),則對任一vh∈Vh,
(6)
從而有
(7)
證明由三角不等式和(4),對任意的vh∈Vh,有
‖u-vh‖≤‖u-Pu‖+‖Pu-vh‖=‖u-Pu‖+‖P(u-vh)‖≤
‖u-Pu‖+θ‖u-vh‖.
因此,(1-θ)‖u-vh‖≤‖u-Pu‖.即得(6).再結(jié)合(5)可得(7).
將定理1的結(jié)果應用到有限元的誤差估計中,分別得到有限元誤差L2范數(shù)和H1范數(shù)的估計.基于估計式(7),可以得到有限元誤差與插值誤差在L2范數(shù)下的等價性,進一步基于插值誤差的漸近展開式和Hessian重構(gòu)方法,構(gòu)造了一個后驗誤差估計子.
考慮如下模型問題:
(8)
(9)
a(uh,vh)=f(vh),?vh∈Vh.
(10)
記uI∈Vh為u的分片線性插值,在定理1中分別取vh=uh和vh=uI可得
由上述兩式,得到有限元誤差‖u-uh‖和插值誤差‖u-uI‖的等價性
‖u-uh‖≈‖u-uI‖.
(11)
因此,
(12)
可作為有限元誤差L2范數(shù)的一個可靠且有效的誤差指示子.對于問題(8),定義能量范數(shù)如下:
c|u-uh|1≤‖u-uh‖E≤‖u-Pu‖E≤C|u-Pu|1.
(13)
由三角不等式,Vh空間上的逆估計以及‖u-uh‖與‖u-Pu‖的等價性,有
再利用對偶論證可得‖u-uh‖≤Ch|u-uh|1.結(jié)合上述兩式,有|u-Pu|1≤C|u-uh|1,即得
‖u-Pu‖E≤C‖u-uh‖E.
(14)
由(13)和(14),即得有限元誤差‖u-uh‖E與投影誤差‖u-Pu‖E的等價性:
‖u-uh‖E≈‖u-Pu‖E.
(15)
由插值誤差估計理論,我們有
(16)
其中,(ξi,ηi),1≤i≤3介于點(x,y)和點 (xi,yi)所成線段上,φi(x,y)為頂點(xi,yi)對應的線性元基函數(shù).
基于誤差展開式(16),可以得到單元誤差估計的一個近似
其中,Hh(uh)表示由有限元解uh得到的Hessian重構(gòu),是D2u的分片常數(shù)近似.
(p(zi)-uh(zi))2.
則單元τ上重構(gòu)的Hessian定義為
Hτ(uh):=2pτ.
基于以上的Hessian重構(gòu),我們可以構(gòu)造相應的重構(gòu)型后驗誤差估計子ητ
(17)
它是關(guān)于有限元誤差的L2范數(shù)的一個可靠且有效的后驗誤差估計子.
我們測試了大量的數(shù)值算例來檢驗有限元L2范數(shù)誤差的重構(gòu)型誤差估計子(17)的有效性.考慮如下橢圓偏微分方程
(18)
其中,?Ω=ΓD∪ΓN,n為邊界上的單位外法向量.我們應用自適應有限元方法求解上述偏微分方程,其中誤差估計子采用(17),網(wǎng)格加密方法分別考慮了二分法[15]和CfCVDT[16]方法.所選取的算例主要來自[17],包括解含有大梯度問題,區(qū)域幾何形狀引起的奇性問題,界面奇性問題等.限于篇幅,我們僅列出了一個例子的結(jié)果.
例1考慮問題(18),其中A=I是一個2×2的單位陣,b=0,Ω=(0,1)2,ΓD=?Ω,真解
圖1畫出了解的圖形,這個解是光滑的,但是有一個陡峭的內(nèi)層,其中參數(shù)S反應內(nèi)層斜坡的陡峭程度,本例中取S=60.圖2畫出了自適應網(wǎng)格,誤差‖u-uh‖及誤差估計子η,誤差有效因子η/‖u-uh‖ 等圖形.可以看到,誤差估計子能夠引導網(wǎng)格在內(nèi)層處加密,且誤差估計子是漸近準確的,誤差‖u-uh‖=O(N-1),達到了最優(yōu)下降速度.
[1]CIARLET P G. The finite element method for elliptic problems[M].Amsterdam: North-Holland, 1997.
[2]王烈衡,許學軍.有限元方法的數(shù)學基礎(chǔ)[M].北京:科學出版社,1996.
[3]陳傳淼,黃云清.有限元高精度理論[M].長沙:湖南科學技術(shù)出版社,1995.
[4]AINSWORTH M,ODEN J T. A posteriori error estimation in finite element analysis[M].New York:Wiley Interscience,2000.
[5]BABUSKA I,STROUBOULIS T. The finite element method and its reliability[M].London:Oxford University Press,2001.
[6]REPIN S. A posteriori estimates for partial differential equations[M].Berlin:Walter de Gruyter, 2008.
[7]VERFüRTH R.A review of a posteriori error estimation and adaptive mesh-refinement techniques[M].Chichester:Wiley-Teubner,1996.
[8]CIARLET P G,WAGSCHAL C. Multipoint Taylor formulas and applications to the finite element method[J].Numer Math,1997,17:84-100.
[9]黃云清,楊偉,易年余.基于顯式多項式恢復的后驗誤差估計[J].湘潭大學自然科學學報,2011,33(3):1-12.
[10]HUANG Y,YI N.The superconvergent cluster recovery method[J].J Sci Comput,2010,44:301-322.
[11]易年余.基于梯度重構(gòu)的后驗誤差估計及自適應有限元方法[D].湘潭:湘潭大學,2011.
[12]AGOUZAL A,VASSILEVSKI Y.A discrete Hessian recovery forp1finite elements[J].J Numer Math,2002,10:1-12.
[13]GUO H,ZHANG Z,ZHAO R.Hessian recovery for finite element methods[J].Math Comput,2017,86:1671-1692.
[14]PICASSO M,ALAUZET F,BOROUCHAKI H,et al. A numerical study of some Hessian recovery techniques on isotropic and anisotropic meshes[J].SIAM J Sci Comput,2011,33:1058-1076.
[15]CHEN L.iFEM: An integrated finite element methods package in MATLAB[R].Technical report, University of California at Irvine, 2010.
[16]JU L,GUNZBURGER M,ZHAO W.Adaptive finite element methods for elliptic PDEs based on conforming centroidal Voronoi elaunay triangulations[J].SIAM J Sci Comput,2006,28:2023-2053.
[17]MITCHELL W F.A collection of 2D elliptic problems for testing adaptive algorithms[R].NISTIR 7668, 2010.