反常次擴散問題的有限元逼近

吳春紅,劉青霞

(1.廈門理工學院應用數學學院,福建 廈門 361024;2.廈門大學數學科學學院,福建 廈門 361005)

1 預備知識

自然界與工程應用中存在很多擴散現象,如地下水溶質運移、污染物在土壤中的遷移、核廢料的擴散、藥物在身體中的擴散、石油滲流、湍流等,這些擴散現象一般不滿足經典的Fick梯度擴散律,通常稱為“反常”擴散．“反常”擴散現象很難用整數階微分方程表征,近年來研究者們廣泛采用分數階微分方程來描述反常擴散過程,并求助數值方法進行模擬．文獻[1]最早采用譜方法處理時間分數階擴散方程,隨后進行推廣,將譜方法應用到空間時間分數階擴散方程[2]．文獻[3]采用有限元方法考慮空間分數階對流-擴散方程．更多的文獻均采用有限差分方法處理分數階次擴散問題[4-6]．文獻[7]采用Adomian分裂方法給出時間分數階Klein-Gordon方程的解析近似解．本文嘗試有限元方法處理分數階導數施加于右側擴散項的次擴散問題．

本文討論下面帶有初邊值條件

u(x,0)=ω(x),x∈Ω=[0,L],

(1)

u(0,t)=u(L,t)=0,t∈[0,T]

(2)

的反常次擴散問題(ASFE):

(3)

其中Γ(·)是Gamma函數．

本文采用有限元方法求解上述ASFE問題,分別給出了其時間半離散、時間空間全離散形式,討論了兩種形式的穩定性、收斂性,并且給出數值例子與有限差分方法[5]進行比較,得到了較好的結果．

2 ASFE半離散格式

定義函數空間:

Hm(Ω)={v∈L2(Ω),

(4)

(5)

2.1 ASFE半離散格式

(6)

其中v(x,t)=u(x,t+τ)-u(x,t)．

為了方便,引入記號:

bs=(s+1)α-sα,s=0,1,2,…,n,

(7)

并對式(6)右端后兩項分別進行如下逼近:

(8)

|Rn+1|≤C*(τ2+bnτα+1) (C*=

max{C1,C2})．

(9)

下面引入關于系數性質的引理[5]:

引理1 式(7)定義的系數bs滿足:

(i)b0=1,bn>0,n=0,1,2,…;

(ii)bn>bn+1,n=0,1,2,…;

由引理1,截斷誤差式(9)滿足:

(10)

可知C″與τ無關．

假設un為u(x,tn)的逼近解,ASFE (1)～(3) 的隱式差分格式為:

(11)

(12)

2.2 ASFE半離散格式的穩定性

定理1 對任意的τ,離散格式(12)是無條件穩定的．

2.3 ASFE半離散格式的誤差分析

取v=ηn+1,有:

(13)

綜上得到半離散問題(12)的誤差分析:

接下來考慮ASFE(1)～(3) 的變分問題．引入記號:

(14)

以及F(v)=(f,v),弱形式(12) 可寫為:

A(un+1,v)=F(v)．

(15)

由Lax-Milgram定理可得下面定理:

3 ASFE的全離散問題

(16)

(17)

由式(5)和逼近性質可得

(18)

(19)

(20)

(21)

(22)

由(21)和(22)得:

下面利用Nitsche技巧[8]得到L2(Ω)中的最優估計．假設w為下列變分問題的解:

(23)

(24)

由式(12)和(16)可得

再由式(18)、定理4有

(25)

此處常數C與式(18)和定理4中的常數有關,與μ、τ和h均無關．

(26)

由定理3以及式(26)可得下面定理:

4 數值例子

例1 利用有限元法近似下列問題:

x∈[0,π],t>0,

u(0,t)=u(π,t)=0,t>0,u(x,0)=0,

x∈[0,π],

(27)

表1給出了二階有限差分方法和采用P1(K)的有限元方法的誤差以及收斂率,從中可以看出,雖然有限差分方法可以得到二階格式,但是顯然沒有有限元方法的誤差精度高．從圖3 可以看出,只要選取到合適的基函數,有限元方法可以達到更高階的精度．圖1是誤差與時間步長的關系圖,圖2、3是誤差與空間步長的關系圖．可以看出,針對不同的階數α,格式對時間步長具有一階精度,對空間步長分別具有二階和三階精度．

表有效元法與有限差分方法的誤差、收斂率Tab.1 The error and convergence of finite element method and finite difference method when

圖1 固定采用P1(K)的數值誤差Fig.1 The numerical errors by use of P1(K) when

圖2 固定采用P1(K)的數值誤差Fig.2 The numerical errors by use of P1(K) when

圖3 固定采用P2(K)的數值誤差Fig.3 The numerical errors by use of P2 (K) when

例2 考慮帶有初邊值條件:

u(x,0)=

u(0,t)=u(2,t)=0,t>0

的問題:

(28)

圖4描述了α=0.7時式(28)解的次擴散情況,可以看出,隨著時間的增加,源在擴散．圖5描述了式(28)不同階數α的次擴散情形,可以看出隨著分數階導數階數的增加,擴散速度在加快．

圖4 α=0.7時式(28)的數值解Fig.4 The numerical solution of (28) when α=0.7

圖5 α∈[0,1],t=1.0時式(28)的數值解Fig.5 The numerical solution of (28) when α∈[0,1]

5 結 論

本文針對分數階導數施加于右側擴散項的次擴散問題,采用有限元法分別給出了其時間半離散、時間空間全離散形式,并且討論了兩種形式的穩定性、收斂性,數值例子中與有限差分方法進行比較,得到了較好的結果．通過本文的研究,不僅探索了次擴散問題有限元方法的處理方式和技巧,同時為求解其他類型的次擴散問題提供了思路．

[1]LinYumin,XuChuanju．Finitedifference/spectralapproximationsforthetime-fractionaldiffusioneuqation[J]．JComputPhys,2007,225(2):1533-1552．

[2]LiXianjuan,XuChuanju．Aspace-timespectralmethodforthetimefractionaldiffusionequation[J]．SIAMJNumerAnal,2009,47(3):2108-2131．

[3]ZhengYunying,LiChangpin,ZhaoZhengang．Anoteonthefiniteelementmethodforthespace-fractionaladvectiondiffusionequation[J]．ComputersandMathematicswithApplications,2010,59(5):1718-1726．

[4]LiuF,YangC,BurrageK．Numericalmethodandanalyticaltechniqueofthemodifiedanomaloussubdiffusionequationwithanonlinearsourceterm[J]．JComputApplMath,2009,231(1):160-176．

[5]ZhuangP,LiuF,AnhV,etal．Newsolutionandanalyticaltechniquesoftheimplicitnumericalmethodfortheanomaloussubdiffusionequation[J]．SIAMJNumerAnal,2008,46(2):1079-1095．

[6]YuQ,LiuF,AnhV,etal．Solvinglinearandnon-linearspace-timefractionalreaction-diffusionequationsbytheAdomiandecompositionmethod[J]．InternationalJournalforNumericalMethodsinEnginerring,2008,74(1):138-158．

[7] 郭鵬,王藝紅,陶春興,等．時間分數階Klein-Gordon型方程的解析近似解[J]．科技導報,2009,27(2):47-52．

[8]QuarteroniA,ValliA．Numericalapproximationofpartialdifferentialequations[M]．NewYork：Springer,1997．