張建松, 張?jiān)轮?朱 江, 羊丹平

(1.中國(guó)石油大學(xué)(華東)應(yīng)用數(shù)學(xué)系,山東青島266580;2.巴西國(guó)家科學(xué)計(jì)算實(shí)驗(yàn)室,巴西里約熱內(nèi)盧25651-075;3.華東師范大學(xué)數(shù)學(xué)系,上海200062)

對(duì)流占優(yōu)擴(kuò)散方程的分裂特征混合有限元方法

張建松1, 張?jiān)轮?,朱 江2, 羊丹平3

(1.中國(guó)石油大學(xué)(華東)應(yīng)用數(shù)學(xué)系,山東青島266580;2.巴西國(guó)家科學(xué)計(jì)算實(shí)驗(yàn)室,巴西里約熱內(nèi)盧25651-075;3.華東師范大學(xué)數(shù)學(xué)系,上海200062)

利用修正的特征線方法,構(gòu)建一類求解對(duì)流占優(yōu)擴(kuò)散方程的分裂特征混合有限元算法.在新的算法中,混合系統(tǒng)的系數(shù)矩陣對(duì)稱正定,且原未知函數(shù)u與流函數(shù)σ=?ε?u可分離求解.推導(dǎo)了加權(quán)能量模意義下的最優(yōu)階誤差估計(jì),并給出數(shù)值算例驗(yàn)證理論上的分析結(jié)果.

混合有限元;特征線方法;分裂解;對(duì)稱正定系統(tǒng);對(duì)流占優(yōu)擴(kuò)散方程

§1 引 言

眾所周知,對(duì)流-擴(kuò)散方程在科學(xué)工程中有著廣泛的應(yīng)用.例如,它描述了地下水中溶質(zhì)的傳輸,油藏問(wèn)題中的石油驅(qū)動(dòng),空氣動(dòng)力學(xué)中空氣的運(yùn)移以及混溶驅(qū)動(dòng)流的流動(dòng)等等.

特征線方法是一種能夠較好地求解對(duì)流擴(kuò)散問(wèn)題的數(shù)值方法.許多科研人員對(duì)此方法做了大量的研究工作,見(jiàn)參考文獻(xiàn)[1-14].其中,Douglas和Russell在文獻(xiàn)[1]中首次提出了特征有限元和特征有限差分法求解對(duì)流擴(kuò)散方程.之后Suli[2]和Pironneau[3-4]將此方法推廣應(yīng)用到Navier-Stokes方程;Russell在文獻(xiàn)[5]中用此方法求解多孔介質(zhì)中混溶驅(qū)動(dòng)問(wèn)題;朱江研究了KdV方程和RLW方程的特征有限元方法.為了得到質(zhì)量守恒的特征有限元格式,Douglas,Huang和Pereira在文獻(xiàn)[10]中給出了一種對(duì)流校正的修正特征有限元算法,Wang在[11]中研究了此類方法一致收斂性.對(duì)于特征混合有限元算法,也有大量的文獻(xiàn)可供參考,見(jiàn)[15-18].然而,經(jīng)典的混合有限元方法通常導(dǎo)致鞍點(diǎn)問(wèn)題,這使得混合系統(tǒng)因失去對(duì)稱正定性而給計(jì)算帶來(lái)許多困難.

本文的主要目的就是建立一類新型分裂特征混合有限元方法求解對(duì)流占優(yōu)擴(kuò)散方程.算法將分裂混合有限元技巧(見(jiàn)文獻(xiàn)[13,19])與修正的特征線方法(MMOC)及對(duì)流校正的修正特征線方法(MMOCAA)相結(jié)合,得到的新混合系統(tǒng)系數(shù)矩陣對(duì)稱正定,且原未知函數(shù)u與新引入函數(shù)σ可分開(kāi)求解.研究了算法的收斂性并給出了ε一致最優(yōu)誤差估計(jì).

本文結(jié)構(gòu)安排如下:首先在§2形成新型分裂特征混合有限元算法求解對(duì)流占優(yōu)擴(kuò)散方程;之后在§3中分析算法的收斂性并給出相應(yīng)的誤差估計(jì);最后在§4中給出數(shù)值算例驗(yàn)證理論上的分析結(jié)果.

§2 算法形成

通篇文章,將用到Sobolev空間中的一些常規(guī)定義,符號(hào),見(jiàn)文獻(xiàn)[20].K無(wú)論是否帶有下標(biāo),在此文章中都代表某個(gè)常數(shù),不同地方出現(xiàn)時(shí)取值可能不同.

考慮如下的對(duì)流占優(yōu)擴(kuò)散問(wèn)題:

其中?是R2中的一有界區(qū)域,邊界記為??,n為邊界??的單位外法向;f(x,t)表示某個(gè)已知函數(shù);b(x)=(b1(x),b2(x))為速度函數(shù)并假設(shè)它是不可壓縮的,即?·b(x)=0,當(dāng)擴(kuò)散系數(shù)0＜ε0≤ε?1,方程(??a)則對(duì)流嚴(yán)格占優(yōu).為了分析方便，這里假定系數(shù)b(x)充分光滑,擴(kuò)散系數(shù)ε為常數(shù).

下面引入所謂的特征線方法.令

令σ=?ε?u.方程(??)重寫(xiě)為

設(shè)Nt為一正整數(shù),△t=T/Nt為時(shí)間步長(zhǎng).對(duì)任意的函數(shù)w用wn表示其在時(shí)刻t=tn=n△t,(n=0,...,Nt)處的值.設(shè)Thu和Thσ為區(qū)域?的兩個(gè)擬正則剖分類,他們可能相同也可能不同,且剖分單元直徑分別以hu和hσ為界.設(shè)Mhu?H1(?)和Vhσ? H(div;?)分別為定義在剖分Thu和Thσ上的有限元空間.

下面分別給出兩種分裂特征混合有限元算法.

2.1 SMMOC-MFE算法

2.2 SMMOCAA-MFE算法

文獻(xiàn)[10]中的MMOCAA算法主要目的是消除文獻(xiàn)[1]中MMOC算法引起的質(zhì)量損失.對(duì)于方程(??a),運(yùn)用修正的特征線∫方 法逼近解∫可得質(zhì)量守恒∫方程

從上面的算法可以看到,混合系統(tǒng)的系數(shù)矩陣對(duì)稱正定,并且原函數(shù)與引入的流函數(shù)可以分開(kāi)求解.顯然經(jīng)典混合有限元方法中要求的LBB條件,這里不再是必須的條件.從計(jì)算的角度,空間Vhσ可以選擇一般的連續(xù)有限元空間.在接下來(lái)的分析中，為了理論分析的方便以及推導(dǎo)算法的最優(yōu)階誤差估計(jì),仍然選擇Vhσ為經(jīng)典的混合有限元空間.

定理2.1SMMOC-MFE算法和SMMOCAA-MFE算法存在唯一解.

其中0＜ K0≤ min(1,△tε).從而可知,矩陣A和B對(duì)稱正定,故而系統(tǒng)(??)有唯一解.進(jìn)而SMMOC-MFE算法有唯一解.同理可證SMMOCAA-MFE有唯一解.

§3 收斂性分析及誤差估計(jì)

假設(shè)有限元空間Vhσ和Mhu有如下的逆估計(jì)性和逼近性質(zhì)(見(jiàn)文獻(xiàn)[21]):存在某些常數(shù)r,r1,k＞0,使得對(duì)1≤q≤∞和任意的ω∈H(div;?)∩[Wr+1,q(?)]2,

為了分析算法的收斂性,假定Vhσ為文獻(xiàn)[21]中提到的任一經(jīng)典混合有限元空間,并引入兩個(gè)投影算子.眾所周知，在經(jīng)典的混合有限空間中,存在一個(gè)投影算子Πh:H(div;?)→Vhσ使得對(duì)任意的1≤q≤∞,滿足

同時(shí)定義橢圓投影算子PM:H1(?)→Mhu滿足

定理3.1假設(shè)問(wèn)題(??)的解(u,σ)有如下的正則性:

則對(duì)于SMMOC-MFE算法,有先驗(yàn)誤差估計(jì)

圖1 b=[x2,?x1],Nt=N時(shí)的誤差估計(jì)‖σ?σh‖L2.

圖2 b=[x2?x1,x2],Nt=N時(shí)的誤差估計(jì)‖σ?σh‖L2.

§4 數(shù)值結(jié)果

在實(shí)際計(jì)算中,可以選取(??a)(或(??a))求解uh.(??a)(或者(??a))為經(jīng)典的特征有限元算法,因此只需要給出數(shù)值結(jié)果去檢驗(yàn)算法(??b)(或(??b))即可.考慮二維邊值問(wèn)題,并取?=[0,1]×[0,1],時(shí)間區(qū)間為(0,T]=(0,1].

圖3b=[x2,?x1],Nt=N時(shí)的誤差估計(jì)|‖σ?σh‖|?L(0,T;H(div)).

在這個(gè)算例中初始邊界條件由真解u=e?tsin2(πx1)sin2(πx2)給出.令σ= ?ε?u.首先將(0,T]分差Nt等份,△t=T/Nt.然后將區(qū)域?分成N×N一致矩形單元,令hu=hσ=1/N,再將每個(gè)矩形單元分成兩個(gè)三角形單元從而獲得一致三角網(wǎng)格剖分.基于此三角網(wǎng)格剖分選取分片線性多項(xiàng)式空間作為有限元空間.通過(guò)選取不同的網(wǎng)格參數(shù)以及對(duì)流項(xiàng)系數(shù)b和擴(kuò)散系數(shù)ε,得到相應(yīng)的收斂性估計(jì),見(jiàn)圖1,圖2和圖3.這些數(shù)值結(jié)果顯示本文的算法是穩(wěn)定的和收斂的.

致謝在此由衷地感謝各位評(píng)審老師,正是你們真摯的意見(jiàn)和建議,大大地提升了本文的質(zhì)量.

[1] Douglas Jr J,Russell T F.Numerical methods for convection-dominated di ff usion problems based on combining the method of characteristics with fi nite element or fi nite di ff erence procedures[J].SIAM J Numer Anal,1982,19(5):871-885.

[2] Suli E.Convergence and nonlinear stability of the Lagrange-Galerkin method for the Navier-Stokes equations[J].Numer Math,1988,53(4):459-483.

[3] Pironneau O.On the transport-di ff usion algorithm,its application to the Navier-Stokes equations[J].Numer Math,1982,38(3):309-332.

[4] Pironneau O.Finite Element Methods for Fluids[M].John-Wiley Sons,Chichester,1989.

[5] Russell T F.Time stepping along characteristics with incomplete iteration for a Galerkin approximate of miscible displacement in porous media[J].SIAM J Numer Anal,1985,22(5):970-1013.

[6] 朱江.KdV方程的特征數(shù)值方法[J].高等學(xué)校計(jì)算數(shù)學(xué)學(xué)報(bào),1988,10(1):11-27.

[7] 朱江.非線性RLW方程的特征數(shù)值方法[J].應(yīng)用數(shù)學(xué)學(xué)報(bào),1990,13(1):64-73.

[8] Pironneau O,Liou J,Tezduyar T.Characteristic-Galerkin and Galerkin/least-squares spacetime formulations for the advection-di ff usion equations with time dependent domains[J].Computer Methods Appl Mech Engrg,1993,27(3):205-227.

[9] Douglas Jr J,Furtado F,Pereira L F.On the numerical simulation of water fl ooding of heterogeneous petroleum reservoirs[J].Comput Geosciences,1997,1(2):155-190.

[10]Douglas Jr J,Huang Chieh-Sen,Pereira L F.The modi fi ed method of characteristics with adjusted advection[J].Numer Math,1999,83(3):353-369.

[11]Wang Kaixin,Wang Hong. A uniform estimate for the MMOC for two-dimensional advection-di ff usion equations[J].Numer Meth Part Di ff er Eq,2010,26(26):1054-1069.

[12]張建松,朱江,郭會(huì),付紅斐.對(duì)流擴(kuò)散反應(yīng)方程的特征分裂最小二乘方法[J].高等學(xué)校計(jì)算數(shù)學(xué)學(xué)報(bào),2012,34(4):289-299.

[13]Zhang Jiansong,Yang Danping,Shen Shuqian,et al.A new MMOCAA-MFE method for compressible miscible displacement in porous media[J].Applied Numerical Mathematics,2014 80(3):65-80.

[14]Zhang Jiansong,Guo Hui.A split least-squares characteristic mixed element method for nonlinear nonstationary convection-di ff usion problem[J].International Journal of Computer Mathematics,2012,89(7):932-943.

[15]Arbogast T,Wheeler M F.A characteristics-mixed fi nite element method for advectiondominated transport problems[J].SIAM J Numer Anal,1995,32(2):404-424.

[16]Ewing R E,Russell T F,Wheeler M F.Convergence analysis of an approximation of miscible displacement in porous media by mixed fi nite elements and a modi fi ed method of characteristics[J].Computer Methods Appl Mech Engrg,1984,47(1-2):23-92.

[17]Ewing R E,Yuan Yirang,Li Gang.Time-stepping along characteristics for a mixed fi niteelement approximation for compressible fl ow of contamination from nuclear waste in porous media[J].SIAM J Numer Anal,1989,26(6):1513-1524.

[18]Yang Daoqi. A characteristic mixed method with dynamic fi nite element space for convection-dominated di ff usion problems[J].J Comput Appl Math,1992,43(3):343-353.

[19]Zhang Jiansong,Yang Danping.A splitting positive de fi nite mixed element method for second order hyperbolic equations[J].Numer Meth Part Di ff er Eqs,2009,25(3):622-636.

[20]Adams R A.Sobolev spaces[M].New York:Academic,1975.

[21]Ciarlet P G.The fi nite element methods for elliptic problems[M].New York:North-Holland,1978.

Split characteristic mixed fi nite element methods for advection-dominated di ff usion equation

ZHANG Jian-song1,ZHANG Yue-zhi1,ZHU Jiang2,YANG Dan-ping3
(1.Department of Applied Mathematics,China University of Petroleum,Qingdao 266580,China;2.Laborat′orio Nacional de Computa?co Cient′? fi ca,25651-075 Petr′opolis,RJ,Brazil;3.Department of Mathematics,East China Normal University,Shanghai 200062,China)

A split modi fi ed method of characteristics mixed fi nite element(SMMOC-MFE)and a split modi fi ed method of characteristics with adjusted advection mixed fi nite element(SMMOCAAMFE)are proposed for solving advection-dominated di ff usion equations,in which the mixed element systems are symmetric positive de fi nite,and the original variable u and the di ff usive fl ux σ = ?ε?u can be solved separately.The optimal-order error estimates in weighted energy norm are derived and some numerical implementations are given to con fi rm the convergence results.

advection-dominated di ff usion equations;mixed element;the method of characteristics;split solution;symmetric positive de fi nite system

65M25;65M60;65M12;65M15

O175.14

A

:1000-4424(2016)03-0338-13

2016-04-10

國(guó)家自然科學(xué)基金(11126084;11401588);山東省自然科學(xué)基金(ZR2014AQ005);中央高?；A(chǔ)研究專項(xiàng)基金(R1510063A)