劉軒宇，羅 鯤，王 皓

(四川大學(xué)數(shù)學(xué)學(xué)院，成都 610064)

1 引 言

考慮下面的拋物型積分微分方程：

其中Ω?Rd(d=2,3) 是一個(gè)有界凸多邊形(多面體)，其邊界記為 ?Ω.在物理和工程領(lǐng)域，同時(shí)帶時(shí)間變量和空間變量的系統(tǒng)非常常見.在研究其中一些問(wèn)題時(shí)，必須考慮之前的狀態(tài)對(duì)現(xiàn)在狀態(tài)的影響.拋物型積分微分方程應(yīng)運(yùn)而生.這類方程在熱傳導(dǎo)、核反應(yīng)堆動(dòng)力學(xué)、熱彈性力學(xué)等領(lǐng)域有廣泛的應(yīng)用[1-24].針對(duì)這類方程，近幾十年來(lái)出現(xiàn)了很多數(shù)值方法，如混合有限元方法[7,13,18]，譜方法[8]，全離散自適應(yīng)有限元方法[17]，間斷Galerkin方法[14,15]，最小二乘Galerkin方法[9]，緊致的交替隱式差分格式[16]等.

弱Galerkin (Weak Galerkin,簡(jiǎn)稱WG)有限元法是近幾年發(fā)展起來(lái)的用于求解偏微分方程的一類有限元方法[25-29]，最早由Wang和Ye針對(duì)二階橢圓方程提出[21,23].WG方法通過(guò)引入對(duì)間斷函數(shù)適用的弱微分(弱梯度，弱散度)算子，使相應(yīng)變分方程容許間斷的逼近函數(shù).此方法既有間斷有限元法網(wǎng)格剖分靈活的特點(diǎn)，也有局部消去的特性，減少了變量個(gè)數(shù).近年來(lái)，WG方法更是被推廣應(yīng)用到多種偏微分方程的數(shù)值求解，如線彈性問(wèn)題[3]，Stokes方程[2,20,22,25,27]，Ossen方程[12]，Biot固結(jié)問(wèn)題[5]，對(duì)流擴(kuò)散反應(yīng)方程[1]，重調(diào)和方程[1]等.

針對(duì)二維拋物型積分微分方程，文獻(xiàn)[29]提出了一類WG有限元方法，給出了半離散格式與全離散格式的解的存在唯一性，導(dǎo)出了相應(yīng)的誤差估計(jì)結(jié)果.本文與該文獻(xiàn)所提出的格式的不同之處在于：

·本文的分析同時(shí)適用于二維與三維情形；

·本文除了適用于單純形網(wǎng)格，還適用于任意多邊形(多面體)網(wǎng)格，文獻(xiàn)[29]只針對(duì)二維的三角形網(wǎng)格；

·本文采用Crack-Nicolson格式，時(shí)間方向上精度為二階，文獻(xiàn)[29]對(duì)時(shí)間的離散采用了向后Euler差分格式，時(shí)間方向上精度為一階.

本文的剩余內(nèi)容安排如下： 第2節(jié)給出本文的記號(hào)和預(yù)備性結(jié)果；第3節(jié)給出半離散格式以及相應(yīng)的誤差估計(jì)；第4節(jié)給出了全離散格式以及相應(yīng)的誤差分析；第5節(jié)給出數(shù)值結(jié)果.

2 預(yù)備知識(shí)

2.1 記號(hào)

引入空間

L(0,T;Hs(D)):={v:(0,T]→Hs(D),v是可測(cè)函數(shù)，且<}

及相應(yīng)的范數(shù)

‖u‖L(0,T;Hs(D))

在本文中，除特殊說(shuō)明外，我們用C表示與網(wǎng)格尺度h和時(shí)間步長(zhǎng)Δt(稍后定義)無(wú)關(guān)的正常數(shù).

2.2 預(yù)備結(jié)果

Gronwall不等式： 設(shè)u(t),β(t),α(t)是定義在區(qū)間[a,b]上的實(shí)值連續(xù)函數(shù)，其中β(t)≥0對(duì)任意t∈[a,b],若

則有

?t∈[a,b]

(2)

離散的Gronwall不等式: 設(shè)(kn)n≥0,(pn)n≥0均為非負(fù)數(shù)列.若數(shù)列(φn)n≥0滿足

則

(3)

右矩形公式誤差估計(jì)： 設(shè)實(shí)數(shù)a,b滿足a

(4)

梯形公式誤差估計(jì)： 設(shè)實(shí)數(shù)a,b滿足a

(5)

則

(6)

3 半離散WG有限元方法

3.1 網(wǎng)格剖分

3.2 離散弱算子

根據(jù)文獻(xiàn)[21]，我們來(lái)定義離散的弱梯度算子.給定Ω上的一個(gè)劃分Th,對(duì)于T∈Th，用V(T)表示T上的弱函數(shù)空間，其定義如下：

V(T)={v={v0,vb}:v0∈L2(T),

vb∈H1/2(?T)}

(7)

選擇一個(gè)有限維向量空間G(T)?H(div,T),定義離散弱算子w,G,T:V(T)→G(T).

定義3.1對(duì)于任意的v∈V(T)，相應(yīng)的離散弱梯度算子w,G,Tv∈G(T) 滿足下列方程：

〈vb,τ·nT〉?T?τ∈G(T)

(8)

特別地，如果G(T)=[Pr(T)]d，我們直接記w,rw,G.

3.3 半離散格式

vhb|e∈Pk-1(e),?T∈Th,?e∈εh},

(9)

初值條件為uh(x,0)=Ehu0(定義見后文式(29))，這里

a(uh,vh)=(w,k-1uh,w,k-1vh),

as(uh,vh)=a(uh,vh)+s(uh,vh).

我們引入一種半范數(shù)：對(duì)于任意vh∈Vh,定義

‖|vh|‖2:=as(uh,uh)=‖wvh‖2+

(10)

引理3.2[21]對(duì)于qh={vh0,vhb}∈Vh，T∈Th，wvh=0 在T上成立當(dāng)且僅當(dāng)vh0=vhb=常數(shù).

引理3.3[3]設(shè)整數(shù)m滿足1≤m≤j+1.則對(duì)于任意T∈Th,e∈∈h,有

?v∈Hm(T),

?v∈Hm(T),0≤s≤m,

?v∈Hm(T),0≤s+1≤m.

設(shè)uh0(t)|T=Φ0α(t)，uhb(t)|T=Φbβ(t)，wuh(t)|T=Φγ(t)，其中給出一些矩陣的定義如下：

在式(8)中取τ=φj,j=1,2…r3，可得

(11)

(12)

在式(9)中取vh={0,vhb}={0,φjb},j=1,2…r2,可得

(13)

設(shè)

則

(14)

(15)

(16)

(17)

(18)

其中，

在方程(17)和(18)兩端對(duì)時(shí)間變量求導(dǎo)，則有

3.4 誤差估計(jì)

設(shè)

0≤t≤T,

(19)

引理3.5[28]設(shè)u∈L(0,T;Hk+1(Ω))為問(wèn)題(1)的真解. 記ρ{ρ0,ρb}=Ehu-Qhu.則當(dāng)h足夠小時(shí)，有

注1文獻(xiàn)[28]中只證明了二維RT元的情形.事實(shí)上，對(duì)(k,k-1,k-1)型元以及三維的情況可類似證明.

對(duì)(19)兩邊關(guān)于t求導(dǎo)，再由類似引理3.5的方式我們可以得到ρt,ρtt的估計(jì)如下.

引理3.6假設(shè)引理3.5的條件均成立，且ut∈L(0,T;Hk+1(Ω)),則當(dāng)h足夠小時(shí)，有

進(jìn)一步，若utt∈L(0,T;Hk+1(Ω)),則當(dāng)h足夠小時(shí)，有

引理3.7[21]假設(shè)v∈Hk+1(Ω),其中k≥1.則‖v-w(Qhv)‖≤Chk‖v‖k+1,其中

(20)

(21)

(22)

由Ehu的定義可知

(23)

結(jié)合(9),(22)和(23)式可得

((ξ0)t,vh0)+as(ξ,vh)=

(24)

取vh=ξt.由于

則

(25)

對(duì)時(shí)間變量從0到t積分，考慮到ξ(x,0)=0,則有

(26)

對(duì)于右端第二項(xiàng)，由分部積分公式可得

結(jié)合‖|·|‖的定義，有

利用Gronwall不等式，結(jié)合引理3.6,可得如下結(jié)果：

則式(20)可通過(guò)引理3.5,引理3.7和三角不等式得到.

在式(24)中取v=ξ，有

C‖(ρ0)t‖2+‖ξ0‖2+

易得

對(duì)時(shí)間變量從0到t積分，可得

利用Gronwall不等式，考慮到ξ(x,0)=0,結(jié)合引理3.5,可得如下結(jié)果

結(jié)合引理3.5,引理3.4和三角不等式，我們得出式(21).證畢.

4 全離散WG有限元方法

4.1 全離散格式

vh)=(f(·,tn+1/2),vh0)

(27)

定理4.1問(wèn)題(1)的的全離散格式(27)有唯一解.

4.2 誤差估計(jì)

在開始分析誤差之前，我們先給出一些逼近誤差的定義和估計(jì).當(dāng)n=0,1…N-1時(shí)，記

Eh0ut(·,tn+1/2).

當(dāng)n=1,2…,N-1時(shí)，記

(28)

(29)

(30)

(31)

根據(jù)Eh0u(·,tn+1)和Eh0u(·,tn)在tn+1/2處的Taylor展開式，可得

由復(fù)合梯形公式(7)和右矩形公式(4)的誤差估計(jì)，得

由梯形公式的誤差估計(jì)(5),有

(32)

由引理3.6,引理3.4,3.7和三角不等式，可得

(33)

類似可得

進(jìn)而可得出式(28)～(31).

注2證明過(guò)程中利用到了時(shí)間與空間的獨(dú)立性，即(wEhu)tt=wEhutt,(wEhu)t=wEhut,(Eh0u)tt=Eh0utt.

定理4.3設(shè)u是問(wèn)題(1)的解，且滿足u∈L(0,T;Hk+1(Ω))，utt1,2…N)是全離散格式(26)的解，則h足夠小時(shí)，有

(34)

(35)

證明 記

Qhu(x,tn+1/2)，n=0,1…,N-1.

(36)

在式(36)中取t=tn+1/2,結(jié)合式(27),(19),可得

(37)

取vh=ξn+1+ξn.我們有

(38)

其中

對(duì)于式(38)左端第一項(xiàng)，通過(guò)計(jì)算可得

(39)

下面估計(jì)式(38)的右端項(xiàng).利用Young不等式,有

(40)

(41)

利用Young不等式，可得

進(jìn)而有

(42)

結(jié)合式(38)～(42)以及|‖·‖|的定義，可得

上式兩邊同時(shí)乘上Δt,對(duì)n從0 取到m-1求和，由

再結(jié)合ξ0=0可得

其中

利用離散的Gronwall不等式可得

由引理3.5可得

再結(jié)合引理4.2有

結(jié)合引理3.5,引理3.7與三角不等式，我們得出式(34).

在式(37)中取vh=?tξn+1.對(duì)于 0≤n≤N-1,有

(43)

其中

對(duì)于方程(43)的左端項(xiàng)，有

as(ξn+1+ξn,ξn+1-ξn)=|‖wξn+1‖|2-

(44)

對(duì)于其右端項(xiàng)，我們有

(45)

同時(shí)，由

(46)

結(jié)合式(43)～(45)可得

(47)

上式對(duì)n從 0 取到m-1求和，再兩邊同時(shí)乘Δt，可得

(48)

現(xiàn)分別估計(jì)上式右端的三項(xiàng).我們有

對(duì)于第三項(xiàng)，采用與式(42)同樣方式可得

因此

(49)

結(jié)合式(48),(49)和 |‖·‖|的定義，有

其中

利用離散的 Gronwall 不等式，可得

由引理3.5和引理4.2,我們有

Hk+1(Ω))].

最后，根據(jù)引理3.5,引理3.4和三角不等式，即得式(35).證畢.

5 數(shù)值實(shí)驗(yàn)

本節(jié)將用數(shù)值算例來(lái)驗(yàn)證定理4.3中的理論結(jié)果.在所有算例中，取Ω=(0,1)×(0,1),T=1,選擇右端項(xiàng)f使初值條件使真解為u(x,t)=t2x1(1-x1)x2(1-x2).

由表1和表2可以看出，k=1時(shí)，u的逼近在時(shí)間精度上有2階，在空間精度上有2階，與定理4.3相符;u的逼近在時(shí)間精度上有2階，在空間精度上有1階，也與定理4.3相符.

表1 正方形網(wǎng)格，k=1Tab.1 Squaredgnid，k=1

表2 正方形網(wǎng)格，k=1Tab.2 Squaredgnid，k=1