王新光， 毛枚良，， 鄧小剛， 涂國華

（1.中國空氣動力研究與發(fā)展中心，四川 綿陽621000；2.空氣動力學國家重點實驗室，四川 綿陽621000）

0 引 言

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

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

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

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

1 基本理論

1.1 空間離散

考慮標量雙曲型線性方程：

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

式（2）中一階導數(shù)的離散采用3種精度的格式進行離散，一階迎風格式

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

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

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

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

1.2 P－multigrid

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

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

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

3）將u＊0→u＊；

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

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

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

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

2 算例考察

2.1 一維粘性Burgers方程

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

其精確解為：

其中通過下式通過迭代求解：

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

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

表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 二維無粘圓柱

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

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

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

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

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

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

2.3 NACA0012翼型

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

從圖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給出了計算過程前600s的V循環(huán)、W循環(huán)和五階格式壁面壓力的平均偏差的變化。兩種pmultigrid方法均明顯加快了收斂速度。

3 結(jié)束語

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

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

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

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

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

［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］.計算力學學報，2006，24（3）：264－368.

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