王新光， 毛枚良，， 鄧小剛， 涂國(guó)華

（1.中國(guó)空氣動(dòng)力研究與發(fā)展中心，四川 綿陽(yáng)621000；2.空氣動(dòng)力學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，四川 綿陽(yáng)621000）

0 引 言

高階精度差分算法同目前成熟的二階精度算法相比，普遍存在收斂速度慢等問(wèn)題，成為高精度算法推廣應(yīng)用的瓶頸之一。為了克服這一困難，參照多重網(wǎng)格算法的思想，提出了一種所謂的p－multigrid方法。眾所周知，多重網(wǎng)格方法（multigrid）的實(shí)質(zhì)是為了快速獲得精細(xì)網(wǎng)格上的計(jì)算結(jié)果，在迭代過(guò)程中，引入與精細(xì)網(wǎng)格存在特定關(guān)系的粗網(wǎng)格上的計(jì)算過(guò)程技術(shù)。類似地，p－multigrid方法的實(shí)質(zhì)是為了穩(wěn)定快速地得到高精度算法的計(jì)算結(jié)果，在迭代過(guò)程中，引入低階精度算法的計(jì)算過(guò)程技術(shù)。

p－multigrid方法是由R?nquist等首先提出的。隨后在間斷有限元中得到了廣泛的應(yīng)用，Maday等［1］分析了p－multigrid方法用于拉普拉斯方程的穩(wěn)定性，Bassi和 Rebay［2］率先將p－multigrid方法用于求解一維的歐拉方程。而Joseph等［3］建立了一套快速低存儲(chǔ)量的p－multigrid方法，并將其用于求解可壓縮歐拉方程。Krzysztof等［4］將p－multigrid方法和單元線光滑器（element line smoother）結(jié)合，建立了一個(gè)新的算法，并將其用于求解可壓縮的NS方程，結(jié)果表明這種方法可以明顯減小計(jì)算時(shí)間，加快收斂速度。以上這些都是p－multigrid方法在間斷有限元中的應(yīng)用。最近幾年，Kris Van Den Abeele等［5］將p－multigrid成功用于一維譜體積法（Spectral Volume Method）中，Aleksundar等［6］將其用于RPM（Recursive Projection Method）中，都可以減少計(jì)算時(shí)間加快收斂。以上這些都證明了p－multigrid方法在加速收斂等方面具有很大的潛力，但從檢索的文獻(xiàn)中還沒(méi)有發(fā)現(xiàn)p－multigrid算法用于有限差分中。

由于差分方法和間斷有限元方法的構(gòu)造思想不同：有限元方法是通過(guò)不同精度的基函數(shù)來(lái)達(dá)到不同的精度，各階導(dǎo)數(shù)是直接計(jì)算的，而有限差分方法中，并沒(méi)有這樣的基函數(shù)，而各階導(dǎo)數(shù)是由變量分布通過(guò)差商得到的。因此，差分方法的p－multigrid方法和間斷有限元方法的p－multigrid方法是不同的。在多重網(wǎng)格法中，使用一種虧損修正法（defect－correction），實(shí)現(xiàn)了將二階精度算法結(jié)果修正到四階精度。

WCNS－E－5格式為5階精度的非線性耗散緊致顯式格式，是鄧小剛研究員等構(gòu)造的高精度格式之一［7－8］，具有良好的耗散、色散特性，在 Euler、N－S和RANS方程等一系列的數(shù)值考核試驗(yàn)中，WCNS格式都表現(xiàn)了良好的性能［9－10］。本文的目標(biāo)是，基于WCNS－E－5格式和上述虧損修正法的基本思想，通過(guò)引入1階精度迎風(fēng)格式和3階精度非線性加權(quán)差分格式，設(shè)計(jì)了一類p－multigrid方法，通過(guò)典型算例，考察了不同方式的p－multigrid方法對(duì)計(jì)算過(guò)程收斂速度的影響。

1 基本理論

1.1 空間離散

考慮標(biāo)量雙曲型線性方程：

此處，c＞0為常數(shù)。對(duì)于網(wǎng)格間距h的均勻網(wǎng)格，式（1）的半離散方程：

式（2）中一階導(dǎo)數(shù)的離散采用3種精度的格式進(jìn)行離散，一階迎風(fēng)格式

最優(yōu)3階精度的加權(quán)格式

最優(yōu)5階精度的 WCNS－E－5格式其中，權(quán)值的定義為：

其中ε＝10－10是為避免分母為零而加的小量。ISk為模板的光滑性度量。

權(quán)值設(shè)計(jì)目標(biāo)是：在光滑的流動(dòng)區(qū)，權(quán)值應(yīng)逼近最優(yōu)值，即可以達(dá)到格式的最優(yōu)精度；在間斷附近，含有間斷的模板被賦予的權(quán)值接近于零，這樣就避免了跨間斷插值，消除了數(shù)值波動(dòng)。

1.2 P－multigrid

基于上述三種精度格式和p－multigrid方法的思想，類似于多重網(wǎng)格方法（即p－multigrid方法中最高精度格式相對(duì)于多重網(wǎng)格方法中的最密的網(wǎng)格），可以構(gòu)造算法為5階精度的多種迭代循環(huán)方法。圖1給出了典型的3種循環(huán)方法V、W、Saw，圖2和圖3分別給出了組合循環(huán)方法Pre＿V和FMG。圖中圓點(diǎn)表示使用了對(duì)應(yīng)精度的格式。

圖1 在三重精度上的V、W、S循環(huán)Fig.1 V，S ，W－circle based on three different schemes

圖2 在三重精度上的Pre＿V循環(huán)Fig.2 Pre＿V－circle based on three different schemes

圖3 在三重精度上的FMG循環(huán)Fig.3 FMG based on three different schemes

下面基于三階、五階格式的小V循環(huán)的計(jì)算過(guò)程敘述如下：

1）利用五階格式迭代ν1步，計(jì)算得到一個(gè)初始解u0；

2）利用三階和五階的加權(quán)格式分別計(jì)算得到RHS3（u0）和RHS5（u0），計(jì)算驅(qū)動(dòng)源項(xiàng)：F＝RHS5（u0）－RHS3（u0）；

3）將u＊0→u＊；

4）利用三階格式計(jì)算RHS3（u＊），再利用驅(qū)動(dòng)源項(xiàng)F計(jì)算：RHS＊＝F＋RHS3（u＊），并利用這個(gè)RHS＊得到一個(gè)新的u＊1；

5）重復(fù)4），迭代ν2步后得到的結(jié)果u＊ν2；

6）將u＊ν2→ui，回到1），完成了一個(gè)循環(huán)。

上述過(guò)程中，隨著殘差逐漸收斂到零，使得上述過(guò)程中ui＝u＊i＝u＊，則使得其中的第四步：RHS＊＝F＋RHS3（u＊）＝RHS5（u＊），最終結(jié)果得到的是一個(gè)五階精度算法的結(jié)果。因此，對(duì)于本文給出的p－multigrid方法，計(jì)算收斂后的結(jié)果是 WCNS－E－5格式的解。

2 算例考察

2.1 一維粘性Burgers方程

這個(gè)算例的主要目的是考察圖1中的V循環(huán)，在各個(gè)精度上選擇不同的步數(shù)，使用不同初值的計(jì)算收斂情況。對(duì)于非線性Burgers'方程：

其精確解為：

其中通過(guò)下式通過(guò)迭代求解：

計(jì)算條件：Re＝25，u｜x＝0，u｜x＝1＝－1。

表1給出了圖1中V循環(huán)在各層精度上使用不同步數(shù)的三種循環(huán)，圖4給出了四種不同的初值。圖5為由不同的初值得到的收斂結(jié)果，可見(jiàn)它們只與算法精度相關(guān)，表2進(jìn)一步定量地展示了計(jì)算結(jié)果同算法的關(guān)系，所采用的V循環(huán)計(jì)算結(jié)果同 WCNS－E－5格式的一致，但計(jì)算所消耗的CPU時(shí)間接近于3階精度格式，比5階精度格式減少約30%，圖6全面展示了算法和初值對(duì)收斂速度的影響，p－multigrid方法的收斂速度均比 WCNS－E－5格式的要快，同時(shí)，初值對(duì)收斂速度存在明顯的影響。

表1 三種不同的V循環(huán)Table 1 Three different V－circles

圖4 四種不同的初值Fig.4 Four different initial values

圖5 精確解和數(shù)值解的比較Fig.5 Comparison of exact solution and numerical solution

表2 Initial＿1的結(jié)果比較Table 2 Comparison of different initial＿1's solution

2.2 二維無(wú)粘圓柱

采用圖7所示的網(wǎng)格（O型網(wǎng)格：128×33，物面距離0.001，遠(yuǎn)場(chǎng)取到20倍直徑），計(jì)算條件為M＝0.3的無(wú)粘圓柱繞流。

在不可壓假設(shè)下，低速無(wú)粘圓柱繞流的壁面壓力系數(shù)的理論解為：

其中θ角表示圖7中以其中的粗線為基準(zhǔn)，逆時(shí)針旋轉(zhuǎn)的角度，前駐點(diǎn)處θ＝90°。圖8給出了流場(chǎng)的壓力等值線。圖9給出了數(shù)值解和理論解的比較，考慮到可壓縮性效應(yīng)，兩者一致性良好。

圖6 不同算法不同初值消耗的CPU時(shí)間的柱狀圖Fig.6 Bar chart of CPU time of different schemes and different initials

通過(guò)無(wú)粘圓柱來(lái)研究圖1和圖2所示的五種不同循環(huán)的效率以及脈動(dòng)范圍的大小。由于在無(wú)粘圓柱的計(jì)算中后駐點(diǎn)壓力收斂速度最慢，可以通過(guò)后駐點(diǎn)來(lái)觀察圓柱的收斂情況，圖10展示了各種算法給出的后駐點(diǎn)壓力收斂情況，使用p－multigrid的各種循環(huán)和五階格式相比收斂速度更快，脈動(dòng)范圍更小。圖11給出了壁面壓力最大偏差（瞬時(shí)值與5階精度方法的收斂值之差）迭代的演化情況，Pre＿V和FMG方法的收斂速度要明顯的快于V、S、W循環(huán)，主要是由于這兩種循環(huán)，可以看做是在V循環(huán)開始之前，使用三階格式和一階格式計(jì)算了一個(gè)較好的初場(chǎng)。

同時(shí)從圖10、圖11中可以看到在計(jì)算的初期pmultigrid方法可以有效地加速收斂，但是隨著流場(chǎng)逐漸收斂，在計(jì)算的后期p－multigrid方法加速收斂的效果不再明顯。

2.3 NACA0012翼型

采用圖12中所示的網(wǎng)格（C型網(wǎng)格：121×50，物面距離0.001）。計(jì)算條件為M＝0.3，Re＝30.0，攻角α＝0.0°。圖13給出了計(jì)算得到的壓力等值線。

從圖14殘差的比較中可知，S循環(huán)、W循環(huán)、V循環(huán)和FMG的收斂速度接近于三階格式的收斂速度，其中S循環(huán)和 W循環(huán)的收斂速度非常接近，F(xiàn)MG和V循環(huán)的收斂速度要稍快于三階格式。

圖15給出了計(jì)算過(guò)程前600s的V循環(huán)、W循環(huán)和五階格式壁面壓力的平均偏差的變化。兩種pmultigrid方法均明顯加快了收斂速度。

3 結(jié)束語(yǔ)

參考DG方法中廣泛使用的p－multigrid方法的基本概念，在結(jié)合有限差分法特點(diǎn)的基礎(chǔ)上，基于1階迎風(fēng)格式、3階加權(quán)格式和5階加權(quán)緊致格式WCNS－E－5，開展了以加快高精度格式收斂速度的pmultigrid方法研究，給出了三種基本循環(huán)（V、W、S形式）和兩種組合循環(huán)（Pre＿V、FMG形式）。通過(guò)典型算法的對(duì)比計(jì)算，得到以下階段性結(jié)論：

（1）p－multigrid方法的各種循環(huán)方式，主要影響收斂歷程，而最終計(jì)算結(jié)果同最高階算法精度一致；

（2）p－multigrid方法中所采用的1階精度算法收斂最快，5階精度算法收斂最慢，p－multigrid方法的收斂速度整體上接近三階精度的收斂速度，具有一定的收斂加速效果，與5階格式相比，能夠減少30%左右的CPU時(shí)間；

（3）p－multigrid方法明顯的加速收斂效果主要體現(xiàn)在計(jì)算過(guò)程的初期，當(dāng)計(jì)算結(jié)果接近收斂解時(shí)，p－multigrid方法加速收斂的效果不再明顯。

在下一步研究中，將在優(yōu)化p－multigrid方法循環(huán)過(guò)程和融合多重網(wǎng)格技術(shù)等方面開展研究工作，提高p－multigrid方法對(duì)高精度算法加速收斂的效果。

［1］MADAY Y，MUNOZ R.Spectral element multigrid，Part2：theoretical justification［R］.Tech Rep.88－73，ICASE，1988.

［2］BASSI F，REBAY S.Numerical solution of the Euler equations with a multiorder discontinuous finite element method［A］.Proceedings of the Second International Conference on Computational Fluid Dynamics［C］，Sydney，Australia，15－19July 2002.

［3］HONG LUO，JOSEPH D.BAUM，PAINALD L?HNER.Ap－multigrid discontinuous Galerkin method for the Euler equations on unstructured grids［J］.Journalof ComputationalPhysics，2006，211：767－783.

［4］KRYSZTOF J.FIDKOWSKI，TODD A.OLIVER，JAMES LU，DAVID L.Darmofal，P－multigrid solution of highorder discontinuous Galerkin discretizations of the compressible Navier－Stokes equations［J］.JournalofComputationalPhysics2005，207：92－113.

［5］KRIS VAN DEN ABEELE，TIM BROECKHOVEN，CHRIS LACOR.Dispersion and dissipation properties of the 1Dspectralvolume method and application to apmultigrid algorithm［J］.J.Comp.Phys.，2007，224：616－636.

［6］ALEKSANDAR JEMCOV，JOESPH P.MARUSZEWSKI.Acceleration and stabilization of algebraic multigrid solver applied to incompressible flow problems［R］.AIAA 2007－4330.

［7］XIAOGANG，DENG，HANXIN ZHANG，H.Develo－ping high－order weighted compact nonlinear schemes［J］.J.Comput.Phys.，2000，165：22－44.

［8］XIAOGANG DENG.A class of high－order dissipative compact schemes［R］.AIAA 96－1972，1996.

［9］劉昕，鄧小剛，毛枚良.高精度非線性格式WCNS的分析研究與其應(yīng)用［J］.計(jì)算力學(xué)學(xué)報(bào)，2006，24（3）：264－368.

［10］鄧小剛，劉昕，毛枚良，等.高精度加權(quán)緊致非線性格式的研究進(jìn)展［J］.力學(xué)進(jìn)展，2007，37（3）：417－427.