魏雪丹,李夢軍,戴厚平

(吉首大學(xué)數(shù)學(xué)與統(tǒng)計學(xué)院,湖南 吉首 416000)

對流是由流體粒子的整體運動和流體分子的隨機運動形成的熱量、質(zhì)量和動量傳遞[1],經(jīng)典的對流擴散問題在不考慮擴散過程時是純對流問題.目前,對流方程廣泛應(yīng)用于泥沙運移、污染物輸送等實際問題中[2-5].但是在復(fù)雜不均勻介質(zhì)中,由于不同介質(zhì)之間的相互作用,傳統(tǒng)的整數(shù)階對流方程已不能準(zhǔn)確刻畫復(fù)雜流體的演化過程.隨著分?jǐn)?shù)階微分理論的發(fā)展,人們逐漸發(fā)現(xiàn)分?jǐn)?shù)階微積分算子的非局部性能夠很好地刻畫具有時間記憶和全局依賴特點的演化過程,于是分?jǐn)?shù)階對流方程應(yīng)運而生.然而,也正是分?jǐn)?shù)階微分算子的全局相關(guān)性,使得大多數(shù)分?jǐn)?shù)階微分方程的精確解不能顯式表達,從而數(shù)值解逼近精確解在分?jǐn)?shù)階微分方程求解方面越來越重要[6-11].對于空間Caputo型分?jǐn)?shù)階微分方程,陳雪娟等[12]運用二次多項式樣條函數(shù)數(shù)值求解了一類空間分?jǐn)?shù)階Fisher方程.為了進一步拓展空間Caputo型分?jǐn)?shù)階微分方程的數(shù)值解法,筆者擬考慮一類分?jǐn)?shù)階對流方程,利用積分中值定理和線性插值方法將分?jǐn)?shù)階對流方程轉(zhuǎn)化為標(biāo)準(zhǔn)對流方程的形式,運用格子Boltzmann方法(Lattice Boltzmann Method,LBM)進行數(shù)值求解.

1 分?jǐn)?shù)階對流方程的格子Boltzmann方法

1.1 分?jǐn)?shù)階對流方程

考慮如下形式的分?jǐn)?shù)階對流方程：

(1)

事實上,由于Caputo型分?jǐn)?shù)階導(dǎo)數(shù)與Riemann-Liouville分?jǐn)?shù)階導(dǎo)數(shù)二者可以相互轉(zhuǎn)化,因此方程(1)也可以表示為

1.2 積分離散化

1.3 格子Boltzmann方法

分布函數(shù)采用如下形式的演化方程：

(2)

(3)

根據(jù)Chapman-Enskog多尺度展開分析,引入二階時間尺度和一階空間尺度,并將g(X,t)展開到二階,可得

(4)

(5)

g(X,t)=ε2g(2)(X,t).

(6)

在小Knudsen數(shù)ε的假設(shè)下展開fi(X,t),可得

(7)

將方程(3)～(7)代入方程(2),可得

(8)

比較方程(8)兩邊ε的各階系數(shù),ε的系數(shù)

(9)

ε2的系數(shù)

(10)

(11)

將方程(11)代入方程(10),可得

(12)

(13)

(14)

(15)

將方程(9)和方程(12)在i方向上求和,并恢復(fù)到原尺度,可得

對于一維分?jǐn)?shù)階對流方程,采用D1Q3 LB模型求解,由方程(13)～(15)可求得平衡態(tài)分布函數(shù)的表達式為

對于二維分?jǐn)?shù)階對流方程,采用D2Q5 LB模型求解,由方程(13)～(15)可求得平衡態(tài)分布函數(shù)的表達式為

2 數(shù)值算例

為了量化格子Boltzmann模型的精確性,引入誤差檢驗公式,全局相對誤差(Global Relative Error，用EGRE表示)和最大誤差(Maximum Error，用EME表示)分別定義為

問題1在有限區(qū)間[0,1]上考慮方程(1)的一維初邊值問題，其初邊值條件為

源項為

驗證該一維初邊值問題的精確解u(x,t)=x2(1-x)2et.

取參數(shù)k=0.001,τ=1.95,Δt=0.001,N=64,此時全局相對誤差見表1.由表1可知：相同時間T,隨著分?jǐn)?shù)階α的增加,全局相對誤差逐漸增加；相同分?jǐn)?shù)階α,隨著時間T的增加,由于誤差的累計,全局相對誤差逐漸增加,但誤差范圍基本保持在10-3以內(nèi).

表1 問題1在N=64時的全局相對誤差

當(dāng)T=1,Δt=0.001時，不同分?jǐn)?shù)階α對應(yīng)的空間最大誤差及其收斂階見表2.由表2可知,LBM在空間方向收斂.

表2 問題1在T=1,Δt=0.001時的最大誤差及其收斂階

圖1示出了α=0.4時不同時刻的數(shù)值解與精確解.由圖1可見,數(shù)值解與精確解吻合較好,且可以長時間保持穩(wěn)定.

圖1 α=0.4時的數(shù)值解與精確解

問題2在區(qū)間[0,1]×[0,1]上考慮方程(1)的二維初邊值問題,其初邊值條件為

源項為

驗證該二維初邊值問題精確解u(x,y,t)=x2y2sint.

取參數(shù)k1=k2=0.001,τ=1.25,Δt=10-4,T=0.3,N=100,當(dāng)α=0.5時數(shù)值解與精確解如圖2所示,數(shù)值解與精確解的絕對誤差如圖3所示.由圖2,3可見,數(shù)值解與精確解的全局相對誤差為7.888 8×10-4,最大誤差為6.188 1×10-4,二者在任意時刻和位置的誤差達到10-4數(shù)量級.

圖2 α=0.5時的數(shù)值解與精確解

圖3 數(shù)值解與精確解的絕對誤差

當(dāng)T=0.1時,不同α和網(wǎng)格數(shù)下的全局相對誤差見表3.由表3可知：相同α?xí)r,隨著網(wǎng)格數(shù)的增加,全局相對誤差逐漸減小,說明LBM在空間上收斂；相同網(wǎng)格數(shù)下,隨著α的增加,全局相對誤差逐漸增加,但控制在10-2以內(nèi).

表3 問題2在T=0.1時的全局相對誤差

當(dāng)網(wǎng)格數(shù)N=32時,不同時刻和α下的全局相對誤差見表4.由表4可知，數(shù)值解能較好逼近精確解,進一步說明了LBM的有效性.

表4 問題2在不同時刻和α下的全局相對誤差

3 結(jié)語

運用構(gòu)建的格子Boltzmann模型數(shù)值求解了一類空間Caputo型分?jǐn)?shù)階對流方程.運用積分中值定理和線性插值方法將分?jǐn)?shù)階對流方程轉(zhuǎn)化為標(biāo)準(zhǔn)對流方程,使得分?jǐn)?shù)階算子產(chǎn)生的記憶部分保留在源項中,極大地簡化了數(shù)值計算過程.數(shù)值解與精確解的實例比較結(jié)果表明,格子Boltzmann模型可以有效地求解Caputo型分?jǐn)?shù)階對流方程.本研究可以推廣到空間三維問題及時間分?jǐn)?shù)階對流問題,進而拓展LBM在分?jǐn)?shù)階微分方程數(shù)值解求解方面的應(yīng)用.