劉忠敏,安 靜,陳 悅
(貴州師范大學 數(shù)學科學學院,貴州 貴陽550025)
二階橢圓方程問題是一類重要而又基礎的問題,在工程和科學領域中很多復雜的非線性問題最終都歸結為求解二階橢圓方程,如非線性Schr?dinger方程[1-3]、Allen-Cahn方程[4]和Navier-stokes方程[5-7]。關于二階橢圓方程的理論分析和數(shù)值計算已有很多成果[8-13]。這些數(shù)值方法主要基于低階的有限元法,對于一些正則性較低的擬線性二階橢圓方程,有限元法確實是一種較好的數(shù)值方法,尤其是與自適應網格相結合的有限元法。然而,對于正則性較高的二階橢圓方程,提出一些基于高階多項式逼近的有限元法也是有意義的,特別是對于一些帶有曲邊或曲面邊界的特殊區(qū)域(如圓域、球域等)上的二階變系數(shù)和非線性的橢圓方程。我們在用有限元求解這些特殊區(qū)域上的二階橢圓問題時,通常以直代曲去逼近這些邊界,從而會引入一些額外的誤差。譜方法是一種具有譜精度的高階數(shù)值方法[14-16]。目前還未有將譜方法應用于求解圓域上二階非線性橢圓方程的研究。
因此,本文提出了圓域上二階變系數(shù)橢圓方程的一種有效的譜方法。通過利用極坐標變換,將原問題轉化為極坐標下的一種等價形式。根據(jù)極條件、邊界條件以及方向的周期性,引入了適當?shù)腟obolev空間,建立了極坐標系下二階變系數(shù)橢圓方程的一種弱形式及其離散格式。然后,利用Lax-Milgram引理證明了弱解的存在唯一性。再由非一致帶權Sobolev空間中投影算子的逼近性質和傅里葉基函數(shù)的逼近性質,證明了逼近解的誤差估計。此外,還將提出的算法延伸到奇異非線性二階橢圓方程的計算,并給出了數(shù)值算例,數(shù)值結果表明本文算法是收斂的和高精度的。
考慮如下二階變系數(shù)橢圓方程:
(1)
其中:Ω={(x,y):x2+y2≤R2};α(x,y)為非負有界函數(shù);?Ω為Ω的邊界。
Δu(x,y)=Lσ(t,θ)=
(2)
由于極坐標變換引入了極點奇性,為了(2)式的適定性,σ(t,θ)需要滿足以下本質極條件:
?θσ(-1,θ)=0。
令f(t,θ)=g(x,y),β(t,θ)=α(x,y),則(1)式在極坐標下的等價形式為
(3)
令ω=t+1表示一個權函數(shù)。為了推導(3)式的弱形式,需要引入一些如下的帶權Sobolev空間:
相應的內積和范數(shù)分別為
和
|?θσ|2dtdθ<∞,
σ(1,θ)=?θσ(-1,θ)=0,
其內積和范數(shù)分別為
A(σ,δ)=F(δ),
(4)
其中
定義有限元空間
χMN=span{σmN(t)eimθ:σmN(t)∈QmN,|m|≤M},
其中QmN={q∈pN:mq(-1)=q(1)=0},pN為N次多項式空間。則(4)式的離散格式為:取σMN∈χMN,使得
A(σMN,δMN)=F(δMN),?δMN∈χMN。
(5)
我們先證明弱解的存在唯一性。為了敘述方便,用a?b表示存在一個正的常數(shù)c,使得a≤cb。
|A(σ,δ)|?‖σ‖1,*‖δ‖1,*,
證明由邊界條件和Cauchy-Schwarz不等式可得
則有
則由β(t,θ)的有界性和Cauchy-Schwarz不等式有
‖σ‖1,*‖δ‖1,*。
另一方面,由β(t,θ)的非負性有
A(σ,σ)=
證畢。
‖δ‖1,*,
引理2設σ和σMN分別是(4)和(5)式的解,則有
證明由(4)和(5)式有
(6)
A(σMN,δMN)=F(δMN),?δMN∈χMN。
(7)
在(6)式中取δ=δMN可得
A(σ,δMN)=F(δMN)。
(8)
將(8)和(7)式作差,由A(σ,δ)的雙線性性質可得
A(σ-σMN,δMN)=0。
(9)
由引理1有
A(σ-σMN,σ-σMN-δMN+δMN)=
A(σ-σMN,σ-δMN)?
‖σ-σMN‖1,*‖σ-δMN‖1,*。
由δMN的任意性有
證畢。
ρ(x+2π)=ρ(x)},
其內積和范數(shù)分別為
‖ρ(θ)-ρM(θ)‖k?Mk-s|ρ(θ)|s。
其內積和范數(shù)分別為
由于σ(t,θ)在θ方向上以2π為周期,則由傅里葉基函數(shù)展開和截斷可得
‖σ-σMN‖1,*?M1-s+N1-q。
證明由于
(1+t)|?t(σM-δMN)|2dtdθ+
則有
從而有
由引理3可得
由引理4可得
由引理2可得
M1-s+N1-q。
證畢。
我們對(5)式的算法進行詳細描述。令
Ψmi(t)=Li(t)-Li+2(t),(i=0,1,…,N-2),
其中Li(t)為i次Legendre多項式。則逼近空間χMN可表示為
χMN=span{Ψmi(t)eimθ:-M≤m≤M,
i=0,1,…,N-1-sign(|m|)}。
令
eimθe-inθdtdθ,
我們將尋找
(10)
將表達式(10)代入(5)式,讓δMN取遍χMN中的所有基函數(shù),則離散格式(5)可化為以下線性系統(tǒng):
(A+B+D)H=f。
其中:
H=[h-M,0,h-M,1,…,h-M,N-2,…,h00,h01,…,
h0,N-1,…,hM,0,hM,1,…,hM,N-2]T,
A=(akjnm),B=(bkjnm),
D=(dkjnm),f=(fkn)。
由Fourier基函數(shù)和Legendre多項式的正交性質可知,矩陣A、B都是稀疏的分塊帶狀矩陣,矩陣D的稀疏性依賴于變系數(shù)α(x,y)的性質。
在前文提出的算法基礎上,給出非線性二階橢圓方程問題的算法。我們考慮如下非線性二階橢圓方程問題:
(11)
類似(3)式的推導,(11)式在極坐標下的等價形式為
(12)
(13)
其中
則(13)式相應的離散格式為:取σMN∈χMN,使得
A(σMN,δMN)=F(δMN),?δMN∈χMN。
(14)
為了求解非線性方程(14),建立如下迭代格式:
因此,我們把非線性格式(14)轉化成變系數(shù)的迭代格式,從而可用第3節(jié)提出的算法有效地求解。
為了表明算法的收斂性和高精度,在MATLAB 2016b平臺上進行一系列數(shù)值實驗。設uMN(x,y)為精確解u(x,y)的逼近解,則由極坐標變換和變量替換有
定義逼近解uMN(x,y)與精確解u(x,y)間的誤差
e(u(x,y),uMN(x,y))=
‖u(x,y)-uMN(x,y)‖L∞(Ω)=
‖σ(t,θ)-σMN(t,θ)‖L∞(D)。
例1我們考慮二階變系數(shù)橢圓方程(1),取α(x,y)=ex+y+1,u(x,y)=(x2+y2-1)cosxy,Ω={(x,y):x2+y2≤1},顯然u(x,y)滿足方程(1)的邊界條件,再將精確解u(x,y)代入方程(1)可得到方程(1)右端的g(x,y)。我們用前文提出的算法求解問題(1)。對于不同的N和M,表1列出了逼近解與精確解之間的誤差。為進一步表明本文提出算法的收斂性和高精度,圖1和圖2給出了精確解與N=25、M=16時的逼近解對比圖與不同的N和M時逼近解與精確解間的誤差圖像。
表1 對于不同的M和N,逼近解與精確解間的誤差e(u(x,y),uMN(x,y))
圖1 精確解和N=25、M=16時的逼近解對比圖
圖2 當N=30、M=20和N=50、M=25時的逼近解與精確解之間的誤差圖像
從表1可以觀察到,當N≥15且M≥16時,逼近解達到了10-14左右精度。另外,從圖1和圖2也可觀察到本文算法是收斂的和高精度的。
例2考慮奇異非線性二階橢圓方程(11)。類似地,取Ω={(x,y):x2+y2≤1},α(x,y)=ex+y+1,u(x,y)=(x2+y2-1)ex+y,將u(x,y)代入方程(11)可得到右端的g(x,y)。對于不同的N和M,表2列出了逼近解與精確解間的誤差。圖3和圖4分別給出了精確解與N=35、M=15時逼近解的對比圖和不同的N和M時逼近解與精確解間的誤差圖像。
表2 對于不同的M和N,逼近解與精確解間的誤差e(u(x,y),uMN(x,y))
圖3 精確解和N=35、M=15時的逼近解對比圖
圖4 當N=15、M=12和N=30、M=14時逼近解與精確解之間的誤差圖像
從表2可以觀察到,當N≥15且M≥12時,逼近解達到了10-13左右精度。另外,從圖3和圖4也觀察到本文算法是收斂的和高精度的。
本文提出了圓域上二階變系數(shù)橢圓方程的一種高精度數(shù)值方法,并證明了弱解的存在唯一性和逼近解的誤差估計。將該算法應用于奇異非線性二階橢圓方程的計算,數(shù)值結果表明本文算法是收斂的和高精度的。但本文僅討論了二階橢圓方程,沒有對一些更復雜的高階非線性問題進行分析和討論,這將是未來的研究目標。