潘悅悅, 楊曉忠
(1.華北電力大學 控制與計算機工程學院, 北京 102206;2.華北電力大學 數(shù)理學院 信息與計算研究所, 北京 102206)
KdV-Burgers方程是一種帶阻尼的KdV方程模型,也是一個典型的非線性波動方程.其一方面,可以描述同時含有色散、非線性水平對流、黏性流體的淺水長波的傳播過程[1-2];另一方面,可以用來描述同時含有色散和耗散的等離子物理中的一些重要物理現(xiàn)象[3-5].KdV-Burgers方程求解的快速算法一直是計算數(shù)學研究的熱點問題,其高效數(shù)值解法具有重要的科學意義和應用價值.
隨著越來越多大規(guī)??茖W計算問題的出現(xiàn),人們開始對串行計算程序的合理化和并行化進行深入研究.近年來的大多數(shù)并行差分方法是條件穩(wěn)定的,或者是無條件穩(wěn)定但精度只有一階.為了得到穩(wěn)定性條件更寬松、計算精度更高的并行方法,實現(xiàn)串行差分方法的合理并行化,一種新的高效并行方法是必不可少的.針對拋物型方程,Evans和Abdullah[6]提出了分組顯式(group explicit,GE)的思想,設計出交替分組顯式(AGE)差分格式.通過計算網(wǎng)格區(qū)域中的兩個相鄰點,將得到的隱式方程轉(zhuǎn)化為顯式方程.交替使用這種方法,截斷誤差項可以部分抵消,從而提高了算法的精度.
隱格式具有良好的穩(wěn)定性,但不適合并行計算.受到構(gòu)造AGE方法的啟示,張寶琳等[7]提出了用Saul’yev非對稱格式構(gòu)造分段隱式的思想,并且恰當?shù)厥褂媒惶婕夹g(shù)建立了多種交替分段顯-隱式(ASE-I)和交替分段Crank-Nicolson(ASC-N)并行方法,得到了穩(wěn)定性和并行性兼顧的研究成果.周毓麟院士[8]將最一般的拋物型方程的顯隱混合格式稱之為具有并行本性的差分格式,對其差分解的存在性、唯一性、收斂性和穩(wěn)定性等理論問題進行了研究,建立了拋物方程并行本性差分方法的基本理論.現(xiàn)在,并行本性差分方法已經(jīng)推廣到很多發(fā)展方程的數(shù)值求解中[9-10].
王文洽[11]使用一組非對稱格式構(gòu)造了具有本性并行的ASC-N差分格式,求解了KdV方程,給出了線性絕對穩(wěn)定性的理論分析.Yuan等[12]針對非線性拋物型方程組提出了一類無條件穩(wěn)定,空間二階精度的并行差分格式.Borhanifar和Abazari[13]應用無條件穩(wěn)定的差分格式求解了電報方程.數(shù)值算例表明了新方法的穩(wěn)定性和準確性.針對四階拋物型方程,Guo和Liu[14]給出了經(jīng)典的顯式、隱式和C-N格式,主要分析了格式的無條件穩(wěn)定性.Namjoo等[15]提出了廣義Burgers-Fisher方程的3種格式:兩種精確有限差分格式和一種非標準有限差分格式.通過誤差分析,驗證了該格式的精確性和有效性.Xue和Feng[16]建立了求解拋物型方程的有效算法.在第n層使用兩種區(qū)域分解法,兩種方法的解取平均值,得到第n+1層的解.該方法保持了并行性和穩(wěn)定性.并行本性方法既適用于并行計算,又是無條件穩(wěn)定或線性穩(wěn)定的.然而,這些方法誤差稍大,精度較低[17-18].
受到上述差分方法的啟發(fā),本文結(jié)合經(jīng)典的顯式、隱式和C-N格式,提出了一種求解KdV-Burgers方程的新并行差分方法.MASC-N方法具有更精確的結(jié)果和更寬松的穩(wěn)定性條件.理論分析和數(shù)值算例表明,MASC-N格式線性絕對穩(wěn)定,收斂階為二階.本文格式的計算效率遠高于C-N格式的計算效率,新方法對于求解KdV-Burgers方程是高效合理的.
KdV-Burgers方程的一般形式為[1-2]
ut+αuux+βuxx+εuxxx=0,L1 (1) 其中α>0,β<0,ε>0.初始條件和邊界條件為 u(x,0)=f(x),L1≤x≤L2, (2) u(L1,t)=g1(t),u(L2,t)=g2(t), 0 (3) 其中L1,L2都是適當大的數(shù). 以下是方程(1)的顯格式、隱格式和C-N格式: (4) (5) (6) 其中 現(xiàn)在討論KdV-Burgers方程的交替分段并行差分格式,本文的MASC-N格式(圖1)將結(jié)合以上幾種差分格式使用. 圖1 MASC-N格式構(gòu)造示意圖 (7) 式中 Un,F1,F2是M-3維向量,I是M-3階單位矩陣,A1和A2均是M-3階對角矩陣,且滿足A1+A2=I,A1=diag(θ1,θ2,…,θM-4,θM-3), 引理1(Kellogg引理)[19]設θ>0,矩陣G+GT是非負定的,則(θI+G)-1存在,并且有 ‖(θI+G)-1‖2≤θ-1. 引理2 MASC-N格式的A1G和A2G是非負定矩陣. 其中c<0.顯然A1G+(A1G)T是弱對角占優(yōu)矩陣,對角線上的元素是非負實數(shù).因此A1G+(A1G)T是非負定矩陣,A2G+(A2G)T也是非負定矩陣,則A1G和A2G是非負定矩陣. 由初始條件得,U0已知,假設U2n已知,求U2n+1.由式(7)得到求解U2n+1的公式為 (I+A1G)U2n+1=(I-A2G)U2n+F1. (8) 顯然,等式右邊已知,(I+A1G)-1存在,式(8)有唯一解.同樣地,使用MASC-N格式計算2n+2層的點: (I+A2G)U2n+2=(I-A1G)U2n+1+F2. (9) 同理可得式(9)有唯一解,定理得證. 定理1 KdV-Burgers方程MASC-N格式(7)的解是存在且唯一的. 引理3[19]設θ>0,矩陣G+GT是非負定的,則‖(θI-G)(θI+G)-1‖2≤1. 定理2 KdV-Burgers方程的MASC-N格式(7)是線性絕對穩(wěn)定的. 古典顯格式為 古典隱格式為 由圖1可知,在每段的“內(nèi)點”使用C-N格式,其計算精度為二階.在每段的“內(nèi)邊界點”使用古典顯式和古典隱式,計算精度為O(τ+h2τ+h2).在T1(τ,h)和T2(τ,h)的表達式中,分別包含了絕對值相同但符號相反的項.古典顯式和古典隱式在每一時間層交替使用,部分項將相互抵消.因此,MASC-N格式的計算精度近似為二階. 定理3 KdV-Burgers方程的MASC-N格式(7)的計算精度為O(h2+τ2),內(nèi)邊界點的計算精度為O(τ+h2τ+h2). ‖en+2‖≤‖(I+G2)-1(I-G1)(I+G1)-1(I-G2)en‖+ ‖(I+G2)-1(I-G1)(I+G1)-1Rn+1‖+‖(I+G2)-1Rn+2‖≤ C′n+2(τ2+h2). 定理4 KdV-Burgers方程的MASC-N格式(7)是收斂的. 為了驗證理論分析,給出數(shù)值算例證明MASC-N格式的計算精度、穩(wěn)定性、收斂性和并行性.同時給出MASC-N格式與ASE-I格式和ASC-N格式的比較試驗、MASC-N格式擴展到mKdV-Burgers方程的數(shù)值試驗. 例1 令ε=1,γ=1,μ=1,考慮以下KdV-Burgers方程[20-21]: ut+εuux-γuxx+μuxxx=0, -50≤x≤50, 0≤t≤1. (10) 初始條件為 邊界條件為 解析解為 取時間層N=500,空間層M=500,圖2和3分別給出了式(10)的孤立波解從時間t=0到t=1的精確解和MASC-N格式解的波形圖.從圖中可以看到,數(shù)值結(jié)果穩(wěn)定可靠,保持了與精確解幾乎完全一致的波形. 圖2 精確解的波形 圖3 MASC-N格式解的波形 表1 兩種格式解相對于解析解的絕對誤差 以下給出C-N格式和MASC-N格式收斂階的數(shù)值試驗.在點(xi,tn)處,定義誤差E∞(h,τ)、空間收斂階S和時間收斂階T[22-23]: 令M=25,50,100,200,400,N=500和N=100,200,400,800,1 600,M=100,表2和3分別給出了C-N格式和MASC-N格式的空間和時間精度.由表中數(shù)據(jù)可知,兩種格式的空間和時間精度為二階.因此本文MASC-N 格式和已有的C-N格式有相同的階數(shù),數(shù)值試驗結(jié)果與理論分析相符. 表2 兩種格式的空間收斂階 表3 兩種格式的時間收斂階 定義加速比Sp=T1/Tp和效率Ep=Sp/p[24],T1和Tp分別是C-N和MASC-N格式的運行時間,p是處理器的數(shù)量.使用八核處理器,令N=1 000,M=600,1 100,1 600,2 100,2 600,計算結(jié)果見表4.表4表明當M增加時,兩種格式的運行時間都增加.但與C-N格式相比,MASC-N格式的運行時間增加緩慢.Sp和Ep都逐漸增加.隨著M的增大,MASC-N格式的效率大幅度提高.當M>600后,MASC-N 格式的運行時間比C-N格式節(jié)約56%. 表4 兩種格式運行時間的比較 圖4給出了C-N和MASC-N格式的計算時間.可以看到C-N格式的計算時間迅速增加,變化率大.本文構(gòu)造的MASC-N格式的計算時間更少,且變化率小.C-N格式是串行格式,通過LU分解計算5對角方程.但使用MASC-N 格式計算時,全局問題被分解成多個不相關(guān)的問題同步進行計算.程序循環(huán)有效地提高了計算效率,縮短了計算時間.從時間的角度來看,MASC-N格式具有明顯的優(yōu)勢. 圖4 C-N和MASC-N格式的計算時間 圖5 不同CPU核數(shù)的MASC-N格式計算時間 圖5給出了MASC-N格式在四核和八核CPU下的計算時間.可以看到隨著空間網(wǎng)格數(shù)的增加,八核CPU的計算時間比四核CPU更少,MASC-N格式并行性明顯. 為了驗證MASC-N格式的穩(wěn)定性和精度,給出隨時間和空間變化的相對誤差.定義時間層相對誤差和(the sum of relative error for every time level,SRET,ΔSRET)與誤差能量(difference total energy,DTE,ΔDTE): 取N=200,M=200,圖6表明當N增加時,兩種格式解的ΔSRET迅速減小,且趨向于0.因此,MASC-N方法求解KdV-Burgers方程是計算穩(wěn)定的. 圖6 C-N和MASC-N格式解的SRET 圖7 C-N和MASC-N格式解的DTE 取N=500,M=500,從圖7可以看出兩種格式ΔDTE的最大值是5×10-3.表明兩種格式數(shù)值解與解析解非常接近,誤差很小,MASC-N方法有良好的計算精度. 例2 令L=100,T=1,考慮以下KdV-Burgers方程[25]: ut+uux-uxx+uxxx=0. (11) 初邊界值由解析解給出,解析解為 比較3種并行格式的計算時間.固定時間層N=10 000,計算結(jié)果見圖8和表5.ASC-N格式有最長的運行時間.ASE-I和MASC-N格式運行時間相似且MASC-N格式計算時間略短. 表5 3種并行格式的計算時間 圖8 3種并行格式計算時間比較 比較解析解和ASE-I[9]、ASC-N[11]、MASC-N格式解.令N=210,M=210,t=0.5,表6給出了計算結(jié)果.表中數(shù)據(jù)表明3種格式的數(shù)值解相似,且MASC-N格式解最接近解析解. 表6 數(shù)值解和解析解 最后,取計算精度的權(quán)重為0.6,計算時間的權(quán)重為0.4.根據(jù)計算精度和時間比較3種并行格式,見表7.可以看到,MASC-N格式的加權(quán)排序總和數(shù)最小,可知其具有最好的性能和實用性. 表7 3種并行格式的比較 例3 令β=-1,ε=-1,考慮以下mKdV-Burgers方程[26-27]: ut+u2ux+βuxx+εuxxx=0, -50≤x≤50, 0≤t≤1. (12) 初邊界值由解析解給出,解析解為 令N=210,M=210,t=0.5,表8為mKdV-Burgers方程解析解和兩種格式數(shù)值解的比較.可以看出C-N格式解和MASC-N格式解幾乎相同.圖9和10分別是C-N和MASC-N格式的絕對誤差,對比發(fā)現(xiàn)兩種差分格式絕對誤差最大值的級數(shù)為10-2. 表8 數(shù)值解和解析解 圖9 C-N格式解的絕對誤差 圖10 MASC-N格式解的絕對誤差 令時間層N=1 000,圖11是C-N和MASC-N格式的運行時間.可以看到C-N格式的運行時間遠大于MASC-N格式的運行時間.與C-N格式相比,MASC-N格式求解mKdV-Burgers方程更為高效. 圖11 C-N和MASC-N格式的運行時間 KdV-Burgers方程的MASC-N并行差分格式是線性絕對穩(wěn)定的,格式解存在唯一,與C-N格式的收斂階同為二階.在保證計算精度的前提下,MASC-N格式比C-N格式有更高的效率.與ASE-I和ASC-N格式相比,MASC-N并行差分格式對求解KdV-Burgers方程有更好的性能.同時,將MASC-N并行差分格式擴展到求解mKdV-Burgers方程,本文的數(shù)值試驗結(jié)果表明,MASC-N格式的計算效率明顯高于C-N格式的計算效率,MASC-N 格式求解各種類型的并行計算系統(tǒng)是有效的. 致謝本文作者衷心感謝德國阿爾弗雷德韋格納研究所Sergey Danilov教授和王強博士在寫作過程中提出的寶貴建議和意見;同時衷心感謝華北電力大學國內(nèi)外聯(lián)合培養(yǎng)博士生項目(2020)對本文的資助.1.2 MASC-N并行差分格式的構(gòu)造
2 MASC-N并行差分方法的數(shù)值分析
2.1 MASC-N格式解的存在唯一性
2.2 MASC-N格式的線性絕對穩(wěn)定性
2.3 MASC-N格式的收斂性
3 數(shù) 值 試 驗
4 結(jié) 論