馬亮亮， 劉冬兵

（攀枝花學(xué)院 數(shù)學(xué)與計算機學(xué)院，四川 攀枝花 617000）

變系數(shù)時間－空間分?jǐn)?shù)階對流－擴散方程的數(shù)值算法比較

馬亮亮， 劉冬兵

（攀枝花學(xué)院 數(shù)學(xué)與計算機學(xué)院，四川 攀枝花 617000）

文章分別采用顯式差分格式、隱式差分格式和Crank－Nicholson格式，數(shù)值求解變系數(shù)時間－空間分?jǐn)?shù)階對流－擴散方程，并從局部截斷誤差、穩(wěn)定性和計算量3個方面對數(shù)值算法進行了比較分析，通過數(shù)值算例驗證了分析結(jié)果。

對流－擴散方程；有限差分格式；穩(wěn)定性；收斂性；變系數(shù)

0 引 言

分?jǐn)?shù)階微分方程具有深刻的物理背景和豐富的理論內(nèi)涵，近年來特別引人注目。目前，分?jǐn)?shù)階導(dǎo)數(shù)和分?jǐn)?shù)階積分在物理、生物、化學(xué)等多個學(xué)科領(lǐng)域有著廣泛的應(yīng)用，如具有混沌力行為的動力系統(tǒng)、擬混沌動力系統(tǒng)、復(fù)雜物質(zhì)或者多孔介質(zhì)的動力學(xué)、具有記憶的隨機游走等。

目前相對應(yīng)用較多且較成熟的方法依然是有限差分法和級數(shù)法，主要是Adomian分解和變分迭代方法。理論分析工具主要有傅里葉方法、能量估計、矩陣方法（特征值）和數(shù)學(xué)歸納法等，還有一些其他的數(shù)值方法，但大都不能作為普適性的數(shù)值方法或缺乏相對較完善的理論分析。文獻［1－5］研究了分?jǐn)?shù)階微分方程數(shù)值算法，借助于一定條件下 Riemann－Liouville分?jǐn)?shù)階導(dǎo)數(shù)與Grünwald－Letnikov分?jǐn)?shù)階導(dǎo)數(shù)的等價性，用移位的 Grünwald－Letnikov技巧逼近 Riemann－Liouville分?jǐn)?shù)階導(dǎo)數(shù)，得到時間、空間、時間和空間分?jǐn)?shù)階擴散方程的有限差分離散近似格式，進而把相應(yīng)的離散格式解釋成時間、空間、時間和空間上的離散隨機游走模型；文獻［6］給出了變系數(shù)空間分?jǐn)?shù)階對流－擴散方程的有限差分格式，并給出了誤差分析；文獻［7］考慮了變系數(shù)空間分?jǐn)?shù)階擴散方程，借助于移位Grünwald－Letnikov技巧和Crank－Nicolson法，得到關(guān)于空間步長一階、時間步長二階收斂的無條件穩(wěn)定數(shù)值離散格式，并對空間變量采用外推技術(shù)，使空間步長的收斂階達到二階；文獻［8］考慮了時間分?jǐn)?shù)階對流－擴散方程，利用Mellin變換和拉普拉斯變換得到了此方程的基本解；文獻［9］將文獻［8］中的結(jié)論推廣到n維全空間和半空間上，得到了相應(yīng)的基本解，并考慮了空間－時間分?jǐn)?shù)階對流－擴散方程的解析解。

分?jǐn)?shù)階微分方程的數(shù)值算法研究起步不久，理論分析和對算法改進方面的研究有限。分?jǐn)?shù)階算子本身的非局部性這一特殊結(jié)構(gòu)，使得分?jǐn)?shù)階微分方程比整數(shù)階需要花費更多的計算時間和更高的存儲要求。

本文考慮與時間和空間都相關(guān)的變系數(shù)時間－空間分?jǐn)?shù)階對流－擴散方程，即

其中，0≤t≤T；0≤x≤L。

初始條件和邊界條件分別為：

其中，0＜α；β≤1；1＜γ≤2；v，d≥0。

（x，t）是 Riemann－Liouville空間分?jǐn)?shù)階導(dǎo)數(shù)［11］，即

1 差分格式

做網(wǎng)格剖分，令τ、h分別是時間和空間步長，即xi＝ih，i＝0，1，…，N，h＝L／N；tk＝kτ，k＝0，1，…，M，τ＝T／M。

Riemann－Liouville空間分?jǐn)?shù)階導(dǎo)數(shù)采用G算法［11］，由于0＜β≤1，1＜γ≤2，所以（x，t）和（x，t）的 G 算法中，最佳移位數(shù)分別是p＝0和p＝1，即采用G1和GS（1）算法離散為：

1.1 顯式差分格式

令為點u（xi，tk）的近似值，即有≈u（xi，tk），則（1）～（3）式的顯式差分離散可表示為：

其中，＝v（ih，kτ）；＝d（ih，kτ）；＝kτ）＝ταΓ（2－α），i＝1，2，…，N，k＝0，1，2，…，M。

初邊值條件為：

其中，i＝0，1，2，…，N；k＝0，1，2，…，M。

為了簡便起見，記＝rm，m＝1，2。當(dāng)r1－r2γ＜2－21－α＝1－b1時，顯式差分（9）式關(guān)于初值條件穩(wěn)定，且（1）～（3）式的數(shù)值解與精確解的誤差滿足［11］：

即（9）式也是條件收斂的，且收斂階為O（τ＋h）。

1.2 隱式差分格式

使用與顯式差分格式相同的時間-空間分?jǐn)?shù)階導(dǎo)數(shù)定義、時間與空間節(jié)點的表達方式以及步長的劃分方法。

采用（1）～（3）式，可得：

或改寫為：

（10）式是無條件收斂的，而且收斂階為O（τ＋h）。

1.3 Crank－Nicholson格式

Crank－Nicholson格式實質(zhì)上是一種隱式差分格式，它主要是在時間方向上運用第j層與第j＋1層值的加權(quán)表示空間分?jǐn)?shù)階導(dǎo)數(shù)項［12］。

使用與顯式差分格式相同的時間-空間分?jǐn)?shù)階導(dǎo)數(shù)定義、時間與空間節(jié)點的表達方式及步長劃分方法，運用時間－空間分?jǐn)?shù)階導(dǎo)數(shù)的Crank－Nicholson格式對（1）～（3）式進行離散，可得：

或改寫為：

顯然（12）式是無條件穩(wěn)定的，且方程的數(shù)值解與精確解的誤差滿足［7］：

即Crank－Nicholson差分格式（12）式無條件收斂，且收斂階為O（τ2＋h）。對（12）式中的空間變量采用外推技術(shù)，可進一步得到關(guān)于空間步長的二階收斂階離散方法。

1.4 3種差分格式的分析與比較

根據(jù)已建立的顯式差分格式、隱式差分格式和Crank－Nicholson格式，分析比較3種差分格式的局部截斷誤差、穩(wěn)定性和計算量，結(jié)果見表1所列。

從表1可以看出，對于變系數(shù)時間－空間分?jǐn)?shù)階對流－擴散方程，顯式差分格式是計算量最小的一種差分格式。

但由于顯式差分格式是條件穩(wěn)定的，所以計算精度較差。另外，利用顯式差分格式進行數(shù)值求解時，還需要注意其穩(wěn)定性。隱式差分格式無條件穩(wěn)定，計算量較小，且計算精度較好，是一種適中的差分方法。Crank－Nicholson格式在局部截斷誤差、穩(wěn)定性和收斂性方面都是三者之中最好的，但其計算量也最大。

表1 3種差分格式的對比

2 數(shù)值算例

分別采用顯式差分格式、隱式差分格式和Crank－Nicholson格式，對如下變系數(shù)時間－空間分?jǐn)?shù)階對流－擴散方程進行數(shù)值求解：

3種差分方法得到的數(shù)值解與精確解之間的相對誤差比較，見表2所列。

表2 3種差分格式的相對誤差比較（x＝0.50）

t顯式 隱式Crank－Nicholson 0.1 0.001 557 755 8 0.000 135 973 6 0.000 052 805 3 0.2 0.001 305 128 2 0.000 379 487 2 0.002 307 692 3 0.3 0.001 223 241 6 0.000 103 975 5 0.000 146 789 0 0.4 0.001 103 448 3 0.000 239 310 3 0.000 137 931 0 0.5 0.000 976 640 0 0.000 308 906 7 0.000 017 066 7 0.6 0.000 862 745 1 0.000 345 098 0 0.000 249 019 6 0.7 0.000 769 574 9 0.000 216 554 8 0.000 109 172 3 0.8 0.000 671 544 7 0.000 565 853 7 0.000 455 284 6 0.9 0.000 427 256 0 0.000 172 081 0 0.000 153 856 4

從表2可以看出，顯式差分格式的初始相對誤差較大，而隱式差分格式和Crank－Nicholson格式則不存在該問題。另外，顯式差分格式隨著時間節(jié)點的增加，初始相對誤差逐漸減小，但始終無法消除初始相對誤差較大的問題。顯式差分格式的相對誤差最大，Crank－Nicholson格式的相對誤差最小，與理論分析結(jié)果一致。

3 結(jié)束語

本文分析比較了數(shù)值求解變系數(shù)時間－空間分?jǐn)?shù)階對流－擴散方程的3種差分格式，通過比較發(fā)現(xiàn)，顯式差分格式的計算量最小，但相對誤差最大，且穩(wěn)定性欠佳；隱式差分格式無條件穩(wěn)定，但相對誤差也較大；Crank－Nicholson格式的相對誤差最小，且無條件穩(wěn)定，但計算量最大。

除本文列舉的差分方法之外，還有一些其他的近似求解方法，如有限元方法［13］、Adomian分解法［14］、配置法［15］及譜方法［16］等，但它們只能求解一些特定的對流－擴散方程，不能作為普適性的數(shù)值算法。因此，變系數(shù)時間－空間分?jǐn)?shù)階對流－擴散方程的數(shù)值計算問題還有待進一步

研究［12，17－23］。

［1］Gorenflo R，Mainardi F.Random walk models for spacefractional diffusion processes［J］.Fractional Calculus ＆Applied Analysis，1998，1（2）：167－191.

［2］Mainardi F，Gorenflo R.Approximation of Lévy－Feller diffusion by random walk ［J］.J Anal Appl，1999，18：231－246.

［3］Gorenflo R，F(xiàn)abritiis G D，Mainardi F.Discrete random walk models for symmetric Lévy－Feller diffusion processes［J］.Physcica A，1999，269：79－89.

［4］Mainardi F，Gorenflo R，Luehko Yu.Wright functions as scale－invariant solutions of the diffusion－wave equation［J］.J Comput Appl Math，2000，118：175－191.

［5］Gorenflo R，Mainardi F.Discrete random walk models for space－time fractional diffusion［J］.Chemical Physics，2002，284：521－541.

［6］Meerschaert M M，Tadjeran C.Finite difference approximations for fractional advection－dispersion flow equations［J］.J Comput Appl Math，2004，172：65－77.

［7］Tadjeran C，Meerschaert M M，Scheffler H P.A secondorder accurate numerical approximation for the fractional diffusion equation ［J］.J Comput Phys， 2006，213：205－213.

［8］Liu F，Anh V，Turner I，et al.Time fractional advectiondispersion equation［J］.J Appl Math Comp，2003，13（1／2）：233－246.

［9］Huang F，Liu F.The time fractional diffusion equation and the advection－dispersion equation［J］.Anziam J，2005，46：317－330.

［10］Podlubny I.Fractional differential equations［M］.San Diego：Academic Press，1999：50－78.

［11］郭柏靈，蒲學(xué)科，黃鳳輝.分?jǐn)?shù)階偏微分方程及其數(shù)值解［M］.北京：科學(xué)出版社，2011：230－247.

［12］孫洪廣，陳 文，蔡 行.空間分?jǐn)?shù)階導(dǎo)數(shù)“反?！睌U散方程數(shù)值算法的比較［J］.計算物理，2009，26（5）：719－724.

［13］Roop J C.Computational aspects of FEM approximation of fractional advection dispersion equations on bounded domain inR2［J］.J Comput Appl Math，2006，193：243－268.

［14］Momani S.Analytic and approximate solutions of the space－and time－fractional telegraph equations［J］.Appl Math Comput，2005，170：1126－1134.

［15］Rawashdeh E A.Numerical solution of fractional integrodifferential equations by collocation method ［J］.Appl Math Comput，2006，176：1－6.

［16］Liu Y，Xu C.Finite difference／spectral approximations for the time－fractional diffusion equation［J］.J Comput Phys，2007，225（2）：1533－1552.

［17］馬亮亮.時間分?jǐn)?shù)階擴散方程的數(shù)值解法［J］.數(shù)學(xué)的實踐與認(rèn)識，2013，43（10）：248－253.

［18］馬亮亮.一種Caputo分?jǐn)?shù)階反應(yīng)－擴散方程初邊值問題的隱式差分格式［J］.貴州師范大學(xué)學(xué)報：自然科學(xué)版，2013，31（1）：58－61.

［19］馬亮亮.一種時間分?jǐn)?shù)階對流擴散方程的隱式差分近似［J］.西北民族大學(xué)學(xué)報：自然科學(xué)版，2013，34（1）：7－12.

［20］馬亮亮.變系數(shù)階空間分?jǐn)?shù)階對流－擴散方程的有限差分解 法 ［J］.沈 陽 大 學(xué) 學(xué) 報：自 然 科 學(xué) 版，2013，25（4）：341－344.

［21］馬亮亮.變時間分?jǐn)?shù)階非定常對流擴散方程的數(shù)值分析［J］.遼東學(xué)院學(xué)報：自然科學(xué)版，2013，20（3）：220－223.

［22］馬亮亮，田富鵬.空間分?jǐn)?shù)階Edwards－Wilkinson方程的顯式差分近似［J］.沈陽大學(xué)學(xué)報：自然科學(xué)版，2013，25（3）：250－252.

［23］陳一鳴，劉麗麗，孫 璐，等.Legendre小波求解非線性分?jǐn)?shù)階積分微分方程數(shù)值解［J］.合肥工業(yè)大學(xué)學(xué)報：自然科學(xué)版，2013，36（8）：1019－1024.

Comparative study of numerical algorithms for time－space fractional convection－diffusion equation with variable coefficients

MA Liang－liang， LIU Dong－bing
（College of Mathematics and Computer，Panzhihua University，Panzhihua 617000，China）

In this paper，the explicit difference scheme，implicit difference scheme and Crank－Nicholson difference scheme are used respectively to solve the time－space fractional convection－diffusion equation with variable coefficients，and their performance is analyzed in terms of local truncation error，stability and computing expense.Finally，a numerical example is given to validate the results.

convection－diffusion equation；finite difference scheme； stability； convergence；variable coefficient

O241.82

A

1003－5060（2014）06－0757－04

10.3969／j.issn.1003－5060.2014.06.023

2013－06－28；

2013－09－28

國家自然科學(xué)基金資助項目（10671132；60673192）；四川省科技廳資助項目（2013JY0125）；攀枝花學(xué)院院級培育資助項目（2012PY08）；攀枝花學(xué)院校級科研資助項目（2013YB05）和攀枝花學(xué)院院級科研創(chuàng)新資助項目（Y2013－04）

馬亮亮（1986－），男，甘肅天水人，攀枝花學(xué)院講師.

（責(zé)任編輯 呂 杰）