霍冠澤, 林亭秀, 王林君

(1. 吉林大學(xué) 數(shù)學(xué)學(xué)院, 長春 130012； 2. 江蘇大學(xué) 理學(xué)院, 江蘇 鎮(zhèn)江 212013)

1 預(yù)備知識

1.1 解的猝滅現(xiàn)象

半線性拋物方程解的猝滅問題在工業(yè)生產(chǎn)、 彈性體的斷裂過程研究[1]及化學(xué)催化反應(yīng)和酶反應(yīng)的動力學(xué)描述[2]等領(lǐng)域應(yīng)用廣泛, 它反映燃燒過程以及火箭燃料在爆炸過程中的穩(wěn)態(tài)與不穩(wěn)態(tài)之間的轉(zhuǎn)換現(xiàn)象, 對于猝滅條件的研究有利于構(gòu)建最佳的燃燒環(huán)境, 了解更精確的燃燒過程[3].

本文考慮如下拋物問題猝滅解的數(shù)值計算方法：

(1)

定義1令u=u(t,x)是方程(1)定義在t>0,x∈Ω上初值問題的解,Ω?m是有界區(qū)域, 如果‖ut‖∞在有限時間T處趨于無窮, 則稱解u具有猝滅現(xiàn)象, 此時T定義為猝滅時間.

由定義1可知, 方程的解是否表現(xiàn)出猝滅現(xiàn)象的決定條件是u對t導(dǎo)數(shù)的最大模是否在猝滅時間處趨近于+∞, 而并不是u是否在猝滅時間處趨近于0, 例如問題(1)的解在時間接近猝滅時間時, 解u→1.

1.2 B方法

解發(fā)展方程數(shù)值解的B方法[12]用來解決帶有爆破現(xiàn)象的半線性發(fā)展方程. 傳統(tǒng)的時間迭代方法, 如方程ut=Δu+F(u)向后Euler格式：

un+1-un=h(Δun+1+F(un+1)),

其中:un+1表示u(tn+1);h表示時間步長. 其數(shù)值解的收斂速度依賴于u和u的導(dǎo)數(shù)大小, 而這兩項在爆破時間附近時至少有一個非常大, 因此數(shù)值解的收斂速度很慢, 有時對于給定的時間步長h甚至不收斂, 而B方法很好地克服了上述問題.

假設(shè)常微分方程ut=f(u)的精確解已知, 考慮擾動方程：

ut=f(u)+εl(u),

(2)

其中εl(u)為空間微分項, 即擾動項. B方法是在解非線性常微分方程變異系數(shù)法的基礎(chǔ)上進(jìn)行修正得來的, 適用于解方程(2)這種半線性或擬線性的拋物方程.

首先, 假設(shè)U=U(t,C)為方程Ut=f(U)的通解, 其中C為積分常數(shù)； 其次, 設(shè)問題(2)的一個解u滿足u=U(t,C(t)), 這里C隨時間t變化, 此時u對t求導(dǎo), 得

ut=Ut+UCCt=f(u)+εl(u).

(3)

從方程u=U(t,C(t))中反求C得C=U-1(t,u), 則式(3)可轉(zhuǎn)化為

(4)

該方法稱為變異常數(shù)的向前Euler法(VCFE), 類似地可推導(dǎo)出變異常數(shù)的向后Euler法(VCBE), 這類方法即稱為B方法. B方法在求解爆破解和爆破時間時效果較好[12], 本文將B方法用于研究解的猝滅現(xiàn)象.

2 主要結(jié)果

2.1 猝滅解的B方法

下面將B方法的思想應(yīng)用于帶有猝滅現(xiàn)象的拋物方程, 其與帶有爆破解的方程根本差別在于非線性項不同, 爆破解方程的非線性項一般為指數(shù)大于1的增長函數(shù), 而猝滅解方程的指數(shù)一般為小于1的函數(shù), 這也決定了具有猝滅解的方程的解隨著時間變化會趨向于某一固定值, 但是其解u對于時間的導(dǎo)數(shù)卻具有爆破現(xiàn)象.

考慮如下方程：

(5)

將B方法應(yīng)用于方程(5), 由式(3)得

(6)

已知方程解具有猝滅現(xiàn)象的條件[13]為

其中0

2)f=-u-p(p>1).

由g-1=G得

g(u)=C+δt?C=g(u)-δt,

以上述結(jié)果代入式(6)得

(7)

式(7)可視為由變異常數(shù)的向前Euler方法(4)得到的, 簡記為VEFE. 式(7)即為B方法在具有猝滅現(xiàn)象拋物方程中的推導(dǎo)應(yīng)用. 類似地可運(yùn)用B方法分別寫出向后Euler法(VCBE)、 梯形法(VCTR)和中點法(VCMR)的計算格式：

為使G(y)有意義, 其定義域和值域分別為(-∞,1/2]和(-∞,1]. 又由于g和G互為反函數(shù), 易得f與g的定義域為u<1, 其VCFE格式為

2.2 解的存在性

證明： 定義Au∶=-Δu, 則式(7)可寫成Aun+1=F(x,un+1), 其中

為了應(yīng)用定理1, 需要證明F(x,0)>0, 即

由例1得f(0)=1,g(0)=0, 且g(un)在滿足定義域的條件下大于零, 故F(x,0)>0得證.

下面證明常數(shù)Bn是一個上解, 即下列不等式成立：

由猝滅解的定義知, 在方程(1)中函數(shù)f(u),g(u)在u=1處存在奇異點, 由例1知0≤u<1, 因此f(un)≥0, 且g(un)是單調(diào)遞增函數(shù), 往證g(Bn)≥g(un)+δh. 由于G函數(shù)此時是單調(diào)遞增的, 因此Bn≥G(g(un)+δh)>0. 因此, 取常數(shù)Bn=G(g(‖un‖∞)+δh)即為一個滿足定理1條件的上解, 由定理1知方程組(1)的解存在.

3 數(shù)值實驗

下面用B方法(VCFE)數(shù)值研究方程組(1)的猝滅解及猝滅時間的近似模擬, 在方程組(1)中取l=5, 空間剖分步長Δx=0.05.

圖1為B方法下計算的數(shù)值解u隨時間t的變化值, 圖2為B方法下計算數(shù)值解u對t的導(dǎo)數(shù)隨時間t的變化值. 由圖1和圖2可見： 隨著時間的增長,u及u對t的導(dǎo)數(shù)值也相應(yīng)增長； 當(dāng)時間趨近于猝滅時間時, 方程組的解u逐漸趨近于1,u對t的導(dǎo)數(shù)趨近于∞. 表明將B方法應(yīng)用于解具有猝滅解的拋物方程時可以很好地刻畫解的猝滅過程.

圖1 B方法下方程解u隨時間t的變化趨勢Fig.1 Variation tendency of solution u with time t of B-method

圖2 B方法下解u對t的導(dǎo)數(shù)隨時間t的變化趨勢Fig.2 Variation tendency of derivative of solution u to time t with time of B-method

考慮問題

在滿足Dirichlet邊值條件下, 初值u0=5×10-5時的解. 在分別取區(qū)間長度l=2,π,5,10的條件下, 對區(qū)間在空間上均勻剖分, 用Δx表示每個剖分單元長度, 分別取Δx=0.1,0.05,0.02,0.01. 表1列出了不同剖分下向前Euler法和B方法的猝滅時間, 這里根據(jù)猝滅時間的定義記u未達(dá)到1的最后一個剖分節(jié)點對應(yīng)的時間為猝滅時間. 表1中TE表示運(yùn)用經(jīng)典向前Euler法所得到的猝滅時間,TB表示運(yùn)用B方法得到的猝滅時間. 由表1可見, 運(yùn)用B方法求得的猝滅時間與經(jīng)典的向前Euler法求得的猝滅時間幾乎相等, 隨著剖分單元長度Δx逐漸減小, 猝滅時間也逐漸減小, 與文獻(xiàn)[3,10]的數(shù)值結(jié)果相近, 表明B方法在模擬猝滅時間上具有與經(jīng)典Euler方法近似的誤差, 因而B方法適用于解具有猝滅解的拋物方程.

表1 問題(1)在不同單元長度空間離散下的猝滅時間

[1] Chan C Y, Ke L. Damage Models of Elastic Materials [J]. Dyn Contin Discrete Impul Syst, Ser A: Math Anal, 2001, 8(1): 1-14.

[2] Phillips D. Existence of Solutions of Quenching Problem [J]. Appl Anal, 1987, 24(4): 253-264.

[3] Sheng Q, Khaliq A Q M. A Compound Adaptive Approach to Degenerate Nonlinear Quenching Problems [J]. Numer Methods Partial Differential Equations, 1999, 15(1): 29-47.

[5] Acker A, Walter W. On the Global Existence of Solutions of Parabolic Differential Equations with a Singular Nonlinear Term [J]. Nonlinear Anal, 1978, 2(4): 499-504.

[6] Levine H A. Quenching and Beyond: A Survey of Recent Results [J]. GAKUTO Internat Ser Math Sci Appl, 1993, 2: 501-512.

[7] GUO Jongsheng, HU Bei. Quenching Problem for a Singular Semilinear Heat Equation [J]. Nonlinear Anal, 1997, 30(2): 905-910.

[8] GUO Jongsheng. On the Quenching Behavior of the Solution of a Semilinear Parabolic Equation [J]. J Math Anal Appl, 1990, 151(1): 58-79.

[9] Fila M, Levine H A. Quenching on the Boundary [J]. Nonlinear Anal, 1993, 21(10): 795-802.

[10] Liang K W, Lin P, Tan R C E. Numerical Solution of Quenching Problems Using Mesh-Dependent Variable Temporal Steps [J]. Appl Numer Math, 2007, 57(5/6/7): 791-800.

[11] Sheng Q, Cheng H. An Adaptive Grid Method for Degenerate Semilinear Quenching Problems [J]. Comput Math Appl, 2000, 39(9/10): 57-71.

[12] Beck M, Gander M J, Kwok F. B-Methods for the Numerical Solution of Evolution Problems with Blow-Up Solutions Part Ⅰ: Variation of the Constant [J]. SIAM J Comput, 2015, 37(6): A2998-A3029.

[13] Nabongo D, Boni T K. Quenching for Semidiscretizations of a Heat Equation with a Singular Boundary Condition [J]. Asymptot Anal, 2008, 59(1/2): 27-38.

[14] Amann H. On the Existence of Positive Solutions of Nonlinear Elliptic Boundary Value Problems [J]. Indiana Univ Math J, 1971, 21: 125-146.