無(wú)斜率選擇的分子束外延方程兩階能量穩(wěn)定的線性格式

潘晗霜，羅智文，王淑芬

(復(fù)旦大學(xué) 數(shù)學(xué)科學(xué)學(xué)院，上海 200433)

1 基礎(chǔ)知識(shí)

本文考慮無(wú)斜率選擇的分子束外延模型[1-3].設(shè)Ω是2上的有界區(qū)域，u和w滿足如下方程：

?tu-ε2Δw+

(1)

w+Δu=0.

(2)

并考慮邊值條件：

?nu=?nw=0 在區(qū)域?Ω×[0,T]上，

(3)

和初值條件：

u(x,0)=u0(x),

(4)

容易驗(yàn)證u，w具有守恒性質(zhì)，這里(·)表示L2內(nèi)積：

(u,1)=(w,1)=0.

(5)

對(duì)于上述方程的能量泛函定義為

(6)

將能量中的第一項(xiàng)稱為Ehrlich-Schwoebel能量項(xiàng)，該項(xiàng)具有非常重要的物理意義.

文獻(xiàn)[4]中通過(guò)擾動(dòng)分析的方法對(duì)方程進(jìn)行分析，得出了對(duì)應(yīng)方程初邊值問(wèn)題的適定性.除此之外，還得出了很多有趣的粗化性質(zhì)，例如Li等[5]證明了當(dāng)ε→0時(shí)，對(duì)于長(zhǎng)時(shí)間t，能量下界被O(-lnt)所限制.由于梯度流粗化過(guò)程的長(zhǎng)時(shí)間性，對(duì)外延方程提出的數(shù)值方法需要具有長(zhǎng)時(shí)間穩(wěn)定以及數(shù)值精度高的性質(zhì).

關(guān)于能量穩(wěn)定的數(shù)值格式，受Eyre[6]所做工作的影響，在研究中涌現(xiàn)出了一類重要的數(shù)值處理方法——凹凸分離的方法，文獻(xiàn)[7]中首次對(duì)外延生長(zhǎng)模型提出了時(shí)間一階無(wú)條件能量穩(wěn)定的凹凸分離格式，但格式的非線性導(dǎo)致數(shù)值求解比較具有挑戰(zhàn).針對(duì)這個(gè)問(wèn)題，文獻(xiàn)[8]提出了一種基于凹凸分離的線性格式.通過(guò)對(duì)能量巧妙地分離，使得非線性項(xiàng)對(duì)應(yīng)的能量位于凸部，并且通過(guò)顯式處理，使得得到的格式為線性格式，極大地簡(jiǎn)化了數(shù)值計(jì)算難度.

為了實(shí)現(xiàn)兩階時(shí)間精度，文獻(xiàn)[9]中基于同樣的凹凸分離方法，通過(guò)使用修正的Crank-Nicolson格式得到兩階數(shù)值格式，并且格式保持無(wú)條件能量穩(wěn)定.類似地，還可以使用后向差分(Backward Differentiation Formula, BDF)格式構(gòu)造兩階精度凹凸分離格式.但是，這兩個(gè)格式都是非線性的，使得數(shù)值實(shí)現(xiàn)變得比較困難.為了處理非線性的問(wèn)題，文獻(xiàn)[10]中提出了一種線性迭代的能量穩(wěn)定格式，通過(guò)迭代算法來(lái)有效地處理高度非線性的格式.事實(shí)上，對(duì)MBE的能量穩(wěn)定格式還有很多研究[11-14].

文獻(xiàn)[15]中提出了一種無(wú)條件能量穩(wěn)定的線性格式，其格式設(shè)計(jì)為：

Aτ(vh)=0,

(7)

(8)

使用后向差分格式處理對(duì)時(shí)間的導(dǎo)數(shù)項(xiàng)，而對(duì)非線性項(xiàng)使用顯式外推格式，同時(shí)，為了保證格式能量穩(wěn)定，還人為添加了一個(gè)兩階時(shí)間精度O(τ2)的正則項(xiàng)AτΔ2(un+1-un).

本文在文獻(xiàn)[15]提出的格式基礎(chǔ)上，考慮了新的正則項(xiàng)設(shè)計(jì)以及非線性項(xiàng)外推，從而得到了一種新的兩階數(shù)值格式，這個(gè)數(shù)值格式同樣可以線性求解且能量穩(wěn)定，本文將給出相應(yīng)的穩(wěn)定性和收斂性分析，最后使用數(shù)值實(shí)驗(yàn)驗(yàn)證了理論分析中得出的各項(xiàng)性質(zhì).

2 新的兩階數(shù)值格式

首先，使用有限元方法對(duì)方程進(jìn)行空間上的離散： 用Th={K}表示Ω上的擬一致三角剖分，這里h表示網(wǎng)格劃分大小的離散參數(shù).記P=Pq(K)為在單元網(wǎng)格K上定義的不高于q階的多項(xiàng)式函數(shù)集合，并定義分片多項(xiàng)式函數(shù)空間Xh≡{v∈X∩C0(Ω)|v|K∈P(K)，?K∈Th}?X.這里X表示在Ω上積分為零的Sobolev函數(shù)空間X={v∈Hq(Ω)|(v,1)=0}.

((Rhφ-φ),χ)=0, (Rhφ-φ,1)=0 ?χ∈Xh.

(9)

(10)

(11)

對(duì)于n≥1時(shí)的中間步格式設(shè)置，我們將文獻(xiàn)[15]的格式中的正則項(xiàng)由Aτ(vh)修改為(vh)，對(duì)于非線性項(xiàng)外推，我們把其由形式修改為從而得到新的中間步格式： 給定和尋找使得

(12)

(13)

注1對(duì)正則項(xiàng)修改是非常自然的，使得兩個(gè)時(shí)間方向離散的處理統(tǒng)一起來(lái)，對(duì)于非線性項(xiàng)，現(xiàn)在的處理方式也是非常自然的，它是對(duì)非線性項(xiàng)函數(shù)做外推，來(lái)取代對(duì)原來(lái)解的外推.

3 穩(wěn)定性分析

首先，定義離散的能量泛函：

(14)

在穩(wěn)定性證明中，我們會(huì)用到文獻(xiàn)[15]中的如下引理.

引理1給定函數(shù)φ1,φ2∈X和v∈X，我們定義函數(shù)Φφ1,φ2,v(s): [0,1]→，

(15)

則Φφ1,φ2,v滿足

(16)

對(duì)數(shù)值格式(12)，(13)有如下能量穩(wěn)定結(jié)果.

定理1給定A≥1，數(shù)值格式(12)，(13)具有能量下降的性質(zhì)，即

(17)

證 首先從式(13)可以推出

(18)

(19)

接下來(lái)我們依次估計(jì)式(19)右邊的項(xiàng).對(duì)于第1項(xiàng)，由Cauchy-Schwarz不等式可得

(20)

對(duì)于式(19)右邊第2項(xiàng)，由式(13)可得

(21)

對(duì)式(19)右邊第3項(xiàng)存在估計(jì)：

(22)

直接引用文獻(xiàn)[15]中的結(jié)果，對(duì)Ehrlich-Schwoebel能量項(xiàng)存在估計(jì)：

(23)

由此成立下面的不等式：

(24)

下面我們依次估計(jì)式(24)右邊積分中的3項(xiàng)S1，S2和S3.首先對(duì)于S1，利用引理1的結(jié)論以及Cauchy-Schwarz不等式：

-|

(25)

類似于對(duì)S1的估計(jì)，我們可以得到：

(26)

而對(duì)于S3，注意到對(duì)任意正實(shí)數(shù)a，b，成立a(a+b)<(1+a2)(1+b2)，因此，

(27)

綜合上面3項(xiàng)的估計(jì)，可得

(28)

對(duì)于式(28)中的梯度項(xiàng)，可以利用Cauchy-Schwarz不等式做放縮：

(29)

(30)

綜合式(19),(20),(21),(22),(28)可得到下式：

(31)

(32)

從而對(duì)修正的能量格式成立能量下降的性質(zhì).

4 收斂性分析

首先，通過(guò)引用文獻(xiàn)[16]的結(jié)論，對(duì)于任意φ∈Hq+1(Ω)∩X的Ritz投影Rh，成立如下估計(jì)：

(33)

這里我們將證明格式是兩階收斂的.用(u,w)表示方程(1)，(2)的真實(shí)解，為了使后面分析成立，需要對(duì)(u,w)作出如下正則性假設(shè)：

u∈L∞(0,T;Hq+1)∩H1(0,T;Hq+1)∩H2(0,T;H1)∩W2,∞(0,T;L2)∩H3(0,T;L2),
w∈L∞(0,t;Hq+1)∩H1(0,T;H1).

(34)

下面將證明格式的誤差收斂階.

這里C2為任意大于0的常數(shù)，則有如下的誤差估計(jì)：

證 首先，定義如下的誤差函數(shù)：

(35)

(36)

對(duì)于中間步,將方程和數(shù)值格式(12),(13)相減可得誤差方程:

(37)

其中：

根據(jù)映射Rh的定義可知： (χ)=0,?χ∈Xh，由此式(37)可改寫(xiě)為： 當(dāng)n≥1時(shí)，對(duì)任意vh∈Xh,φh∈Xh，成立

(38)

(39)

(40)

這里用到了一個(gè)變換：

(41)

運(yùn)用這個(gè)記號(hào)，式(40)中左邊第1項(xiàng)可以寫(xiě)為：

(42)

(43)

接下來(lái)估計(jì)式(40)的右邊項(xiàng).對(duì)于前兩項(xiàng)，運(yùn)用Ritz投影的性質(zhì)(33)以及Cauchy-Schwarz不等式可得

(44)

以及

(45)

對(duì)于第3項(xiàng)和第4項(xiàng)的分析，同樣可參考文獻(xiàn)[15]的結(jié)論：

(46)

(47)

為了估計(jì)式(40)右邊第5項(xiàng)，注意到：

(48)

以及

(49)

根據(jù)Taylor展開(kāi)和Cauchy-Schwarz不等式，結(jié)合式(48)和(49)，對(duì)于式(40)右邊第5 項(xiàng)成立：

(50)

(Nn+1，

(Nn+1,

(51)

將式(42)～(51)帶回誤差方程(40)，得到:

(52)

對(duì)式(52)做從n=1到n=m求和，再在兩邊乘上2τ，可得

(53)

容易驗(yàn)證

根據(jù)式(33)和(36)的結(jié)果可得

(54)

(55)

結(jié)合式(55)以及式(14)的結(jié)果即可得到定理2的結(jié)論.

5 數(shù)值結(jié)果

這一節(jié)中，我們將通過(guò)幾個(gè)數(shù)值實(shí)驗(yàn)來(lái)驗(yàn)證前面提到的一些性質(zhì)，首先，我們驗(yàn)證格式的收斂階，在下面的例子中我們都是使用文章第1節(jié)提出的數(shù)值格式，這里設(shè)定區(qū)域Ω=[0,1]2和ε2=0.05，方程的真值解為

ue(x,y,t)=cos(πx)cos(πy)e-t.

(56)

為了滿足方程以及邊界條件，我們?cè)诜匠痰挠疫吿砑由贤饬?xiàng)：

?tue+

這里，

表1 空間收斂性

表2 時(shí)間收斂性

下面，我們將進(jìn)行數(shù)值模擬來(lái)驗(yàn)證結(jié)果是否具有能量下降等性質(zhì).對(duì)于方程本身的參數(shù)，我們?nèi)《é?=0.005，而將方程的初值選取為每個(gè)點(diǎn)均為-0.05至0.05之間的隨機(jī)數(shù)，設(shè)定空間區(qū)域Ω為[0，12.8]×[0，12.8]的正方形區(qū)域.對(duì)于模型的參數(shù)設(shè)定，首先就空間的離散，考慮到計(jì)算的效率我們?cè)O(shè)定N=128，為了觀察方程解的長(zhǎng)時(shí)間性質(zhì)，取定T=20000，而對(duì)時(shí)間步長(zhǎng)，為了快速計(jì)算，我們?cè)诓煌臅r(shí)間區(qū)間選取了不同的時(shí)間步長(zhǎng)： 當(dāng)t≤200時(shí)，時(shí)間步長(zhǎng)τ=0.004；當(dāng)2002000時(shí)，取時(shí)間步長(zhǎng)為0.15.

我們首先觀察該數(shù)值格式解出的解在不同時(shí)刻的形態(tài)，圖1中所示為解在t=1,500,1000,5000,10000,20000時(shí)刻的形態(tài).

接下來(lái)，我們來(lái)驗(yàn)證解是否具有能量下降的性質(zhì).圖2展示了能量隨時(shí)間的變化情況，從該圖中我們可以明顯觀察到能量是隨時(shí)間下降的，也就從數(shù)值上驗(yàn)證了方法的能量穩(wěn)定性.

定義表面粗度R(t)以及平均梯度W(t)分別為:

(57)

從文獻(xiàn)[1，5] 中我們已經(jīng)知道，對(duì)于沒(méi)有斜率選擇的分子束外延方程，E(t),R(t)和W(t)分別滿足如下性質(zhì)：

(58)

以下，我們將對(duì)數(shù)值結(jié)果分別驗(yàn)證其是否滿足上面提到的3條性質(zhì).為此，圖3(a)顯示了能量E(t)對(duì)lnt的變化情況，圖中紅色虛線表示線性擬合的結(jié)果，這里我們擬合出的結(jié)果為E(t)=-38.1890 lnt-60.1026，而對(duì)于表面粗度R(t)和時(shí)間的關(guān)系，圖3(b)顯示了ln(R(t))相對(duì)于lnt的變化情況，同樣這里給出了線性擬合的結(jié)果，即R(t)=0.4454t0.482 4，最后對(duì)于平均梯度W(t)，可以查看圖3(c)的結(jié)果，這里擬合出的結(jié)果為W(t)=2.3245t0.243 8.從這些結(jié)果我們可以發(fā)現(xiàn)，該格式計(jì)算出的數(shù)值結(jié)果與文獻(xiàn)[1，5]中的結(jié)論相符.