基于“各向異性”四面體網(wǎng)格聚合的復(fù)雜外形混合網(wǎng)格生成方法

趙 鐘，張來平，2，赫 新，2

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

0 引 言

在CFD工業(yè)化應(yīng)用中，生成計(jì)算網(wǎng)格一直占據(jù)整個(gè)計(jì)算流程的60%～70%時(shí)間，這一點(diǎn)對(duì)于復(fù)雜外形高質(zhì)量粘性流動(dòng)計(jì)算網(wǎng)格的生成更為突出［1］。

依照網(wǎng)格的拓?fù)浣Y(jié)構(gòu)，計(jì)算網(wǎng)格分為結(jié)構(gòu)網(wǎng)格和非結(jié)構(gòu)網(wǎng)格。結(jié)構(gòu)網(wǎng)格的優(yōu)點(diǎn)是數(shù)據(jù)結(jié)構(gòu)簡單，存儲(chǔ)方便，計(jì)算簡單快捷，計(jì)算結(jié)果精度高；其缺點(diǎn)是難以處理復(fù)雜外形，對(duì)一些外形上的奇點(diǎn)需要進(jìn)行簡化或特殊處理。非結(jié)構(gòu)網(wǎng)格最早來源于固體力學(xué)中的有限元計(jì)算，主要指二維的三角形和三維的四面體網(wǎng)格，其優(yōu)點(diǎn)是易于處理復(fù)雜外形，數(shù)據(jù)結(jié)構(gòu)的隨機(jī)性使得其易于進(jìn)行網(wǎng)格自適應(yīng)；缺點(diǎn)是存儲(chǔ)量大，計(jì)算效率低。非結(jié)構(gòu)網(wǎng)格的缺點(diǎn)在高雷諾數(shù)計(jì)算時(shí)表現(xiàn)得尤為突出：高雷諾數(shù)邊界層模擬要求在物面法向有足夠的網(wǎng)格分辨率，這一點(diǎn)使得非結(jié)構(gòu)網(wǎng)格數(shù)量急劇增長。為了在一定程度上減少網(wǎng)格量，可以使用法向高度壓縮的“各向異性”四面體單元。然而，法向相鄰兩個(gè)各向異性四面體單元的中心連線和物面法向的偏離使得計(jì)算在某些情況下不夠精確（如利用格心型有限體積法求解單元內(nèi)的物理量梯度時(shí)）。

基于結(jié)構(gòu)網(wǎng)格和非結(jié)構(gòu)網(wǎng)格各自的優(yōu)缺點(diǎn)發(fā)展起來的混合網(wǎng)格技術(shù)結(jié)合了二者的優(yōu)勢(shì)：邊界層內(nèi)的網(wǎng)格一般采用半結(jié)構(gòu)的三棱柱單元（在物面法向是結(jié)構(gòu)的，而在貼體方向是非結(jié)構(gòu)的），邊界層以外單元是各向同性的四面體網(wǎng)格或Cartesian網(wǎng)格。與結(jié)構(gòu)網(wǎng)格相比，生成復(fù)雜外形的混合網(wǎng)格更加容易；和純非結(jié)構(gòu)網(wǎng)格相比，邊界層法線上有很高的網(wǎng)格分辨率和很好的網(wǎng)格正交性，又極大地減少了網(wǎng)格量（邊界層內(nèi)比非結(jié)構(gòu)網(wǎng)格減少大約三分之二）。這無疑對(duì)提高復(fù)雜外形的數(shù)值模擬精度和效率都是非常有利的。

混合網(wǎng)格的出現(xiàn)，在提高邊界層網(wǎng)格質(zhì)量的同時(shí)，也在一定程度上緩解了生成復(fù)雜外形網(wǎng)格帶來的困難。目前常用的三棱柱網(wǎng)格生成方法主要有幾何層推進(jìn)方法［2－3］和 求 解 雙 曲 型 方 程 的 方 法［4－5］。 這 些方法在實(shí)際工程應(yīng)用中得到了成功的應(yīng)用。但是，對(duì)于工業(yè)應(yīng)用中很多極端復(fù)雜的實(shí)際外形而言，要生成邊界層的三棱柱網(wǎng)格并非易事，往往會(huì)在幾何曲率變化劇烈的區(qū)域出現(xiàn)網(wǎng)格相交的情況，為此需要付出很多人工勞動(dòng)調(diào)整局部網(wǎng)格分布，極大地制約了CFD的應(yīng)用。

鑒于復(fù)雜外形三棱柱／四面體混合網(wǎng)格生成的困難，本文發(fā)展了一種自動(dòng)化的基于“各向異性”四面體網(wǎng)格聚合的混合網(wǎng)格生成方法（為了敘述方便，下文稱“聚合法”），即首先用已經(jīng)發(fā)展成熟的非結(jié)構(gòu)網(wǎng)格生成方法生成粘性層內(nèi)的各向異性四面體單元和無粘區(qū)的各向同性四面體單元，然后將邊界層內(nèi)的各向異性四面體單元聚合為三棱柱，而外場(chǎng)的各向同性四面體單元保持不變。這種方法生成的網(wǎng)格在整個(gè)計(jì)算域是由四面體、三棱柱和少數(shù)多面體單元組成的混合網(wǎng)格。通過對(duì)DLR－F6翼身組合體和復(fù)雜戰(zhàn)斗機(jī)外形的混合網(wǎng)格生成及數(shù)值模擬表明，該方法生成的混合網(wǎng)格具有網(wǎng)格質(zhì)量好、自動(dòng)化程度高等特點(diǎn)，對(duì)提高復(fù)雜外形飛行器關(guān)鍵氣動(dòng)參數(shù)的模擬精度具有重要意義。

1 非結(jié)構(gòu)四面體網(wǎng)格生成方法簡介

目前，非結(jié)構(gòu)的四面體網(wǎng)格生成方法已經(jīng)發(fā)展的比較成熟［6－7］，多種商業(yè)軟件業(yè)已集成了這種方法，比如Gridgen、ICEM－CFD等。這里僅做簡要介紹。非結(jié)構(gòu)的四面體網(wǎng)格包括粘性層內(nèi)的各向異性四面體網(wǎng)格和無粘區(qū)域的各向同性四面體網(wǎng)格。

生成各向異性四面體網(wǎng)格采用的是由陣面推進(jìn)法（Advancing Front Method）改進(jìn)而來的層推進(jìn)法（Advancing Layer Method）。在生成好物面的三角形單元后，通過確定三角形頂點(diǎn)的法向和推進(jìn)距離生成出拉伸的四面體單元，生成步驟如下：

（1）計(jì)算物面三角形單元頂點(diǎn)的推進(jìn)方向。通過面積加權(quán)點(diǎn)p周圍的面法向得到點(diǎn)p的推進(jìn)法向，并做一定次數(shù)迭代進(jìn)行光滑：

這里，Np表示點(diǎn)p周圍的頂點(diǎn)數(shù)目，Nf表示和點(diǎn)p相連的三角形數(shù)目表示p點(diǎn)的法向，nf表示面的法向，n表示第n次迭代。

（2）在確定了推進(jìn)方向后，需要確定推進(jìn)距離，一般采用幾何指數(shù)增長模式：

這里，δ0和δn分別表示第一層的網(wǎng)格高度和第n層的網(wǎng)格高度，r是法向網(wǎng)格距離變化率。

（3）在確定了推進(jìn)方向和推進(jìn)距離后，從每個(gè)頂點(diǎn)向計(jì)算域推出一層網(wǎng)格點(diǎn)。對(duì)于每個(gè)三角形而言，三個(gè)頂點(diǎn)向外推出三個(gè)點(diǎn)。分別將推出的三點(diǎn)和前一層的三個(gè)點(diǎn)相連，形成三個(gè)四面體單元。

在推進(jìn)過程中，為了避免單元之間相交，保證交界面的相容性，需要兩兩判斷是否相交，這在三維情況下會(huì)帶來網(wǎng)格生成效率急劇下降。為此，可以采用“彈簧法”來判斷相交性以提高效率，具體過程可參考文獻(xiàn)［7］。

在推出給定層數(shù)的各向異性四面體網(wǎng)格后，停止各向異性四面體網(wǎng)格生成過程，形成各向異性四面體網(wǎng)格外邊界面。各向異性四面體網(wǎng)格推進(jìn)過程中，在局部會(huì)由于可能相交而停止推進(jìn)，這些局部停止推進(jìn)而形成的邊界面連同其他邊界面（如遠(yuǎn)場(chǎng)、對(duì)稱面）一起，形成了各向同性四面體網(wǎng)格生成的初始陣面，從而進(jìn)入各向同性四面體網(wǎng)格生成階段。生成空間各向同性四面體單元一般采用Delaunay法或陣面推進(jìn)法。在現(xiàn)有的商業(yè)軟件中，一般將二者結(jié)合起來使用。這樣既可以提高四面體網(wǎng)格生成效率，同時(shí)可以保證初始陣面的相容性。具體方法請(qǐng)參見文獻(xiàn)［6－7］，這里不再詳述。

2 基于各向異性四面體聚合的混合網(wǎng)格生成方法

在四面體網(wǎng)格生成好后，就可以用本文發(fā)展的聚合法來生成混合網(wǎng)格。聚合法的基本思想是，將各向異性四面體聚合為三棱柱，而各向同性四面體保持不變。這種方法一方面充分利用了非結(jié)構(gòu)網(wǎng)格填充復(fù)雜外形計(jì)算域的自動(dòng)化特點(diǎn)，另一方面又能聚合得到高質(zhì)量的混合網(wǎng)格以適應(yīng)粘性流動(dòng)計(jì)算的需要。聚合過程主要包括四面體單元聚合、棱柱側(cè)面聚合和邊界條件處理。

2．1 四面體單元聚合

聚合四面體網(wǎng)格是整個(gè)算法的核心部分，其原理是充分利用層推進(jìn)法生成各向異性四面體算法的特點(diǎn)，將三個(gè)各向異性四面體單元聚合為一個(gè)三棱柱。聚合方法直接影響到聚合后的網(wǎng)格質(zhì)量，具體包括幾個(gè)步驟：

首先，提取單元的面幾何特征。對(duì)初始四面體網(wǎng)格，計(jì)算每個(gè)四面體的所有表面（三角形）的面積，并找出每個(gè)單元的面積最大面、次大面和最小面。

然后，根據(jù)四面體單元面的幾何特征判斷是否聚合，具體算法如下：

（1）根據(jù)單元的最小面和最大面的面積之比確定單元性質(zhì)：

若r＜α，則該單元為各向異性，否則為各向同性。這里，Smin和Smax分別表示最小面面積和最大面面積。α是一經(jīng)驗(yàn)參數(shù)，實(shí)踐表明，α一般取值小于0．4，α取值越大則三棱柱層數(shù)越多，反之則邊界層網(wǎng)格和無粘區(qū)域網(wǎng)格過渡越光滑；

（2）第一次聚合。遍歷初始網(wǎng)格的所有面，設(shè)面Face兩側(cè)的單元分別為C1和C2，若Face是C1或者C2的最大面，并且C1或者C2是各向異性網(wǎng)格，則將這兩個(gè)網(wǎng)格聚合為一個(gè)粗網(wǎng)格CC。定義每個(gè)粗網(wǎng)格里含有的子單元個(gè)數(shù)為粗網(wǎng)格聚合率，則CC的聚合率為2；

（3）第二次聚合。遍歷初始網(wǎng)格的所有面，找出其兩側(cè)單元，若面Face兩側(cè)單元中有一個(gè)（C1）已被聚合而另一個(gè)（C2）尚未聚合，其中C1所在粗網(wǎng)格為CC，F(xiàn)ace是C1或C2的最大面，C1或C2是各向異性單元，并且CC的聚合率小于3，則將C2聚合到CC，CC聚合率加1。

（4）第三次聚合。遍歷每個(gè)尚未聚合的各向異性單元C1，找出其最大面以及最大面另一側(cè)單元C2，和C2所對(duì)應(yīng)的粗網(wǎng)格為CC。若CC的聚合率小于等于3，則將C1聚合到CC。

（5）遍歷所有初始四面體單元Cell，若Cell尚未被聚合，則將之單獨(dú)聚合為一個(gè)粗網(wǎng)格。

聚合過程如圖1所示。第一次聚合和第二次聚合的主要目的是將三個(gè)四面體單元聚合為一個(gè)三棱柱單元。但是對(duì)于外形變化劇烈的局部網(wǎng)格幾何特征不明顯，導(dǎo)致其不能被聚合而成為孤立單元。如果放任之，則該單元和周圍單元的體積比約為1∶3，網(wǎng)格過渡不光滑。為了避免網(wǎng)格過渡不光滑，第三次聚合將孤立單元聚合到周圍的粗網(wǎng)格，使之和周圍單元體積比約為4∶3，從而使網(wǎng)格光滑過渡。

圖1 聚合網(wǎng)格過程Fig．1 Procedure of cell agglomeration

2．2 側(cè)面的聚合

初始網(wǎng)格在經(jīng)過聚合后形成的“三棱柱”側(cè)面還是兩個(gè)三角形。為了減少存儲(chǔ)開銷和計(jì)算量，需要將每個(gè)側(cè)面的兩個(gè)三角形合并為一個(gè)四邊形，如圖2。聚合面的關(guān)鍵在于判斷哪些面需要聚合為一個(gè)面，可以通過面的左右單元屬性和面法向夾角判斷。值得注意的是，在合并邊界面時(shí)，只有邊界類型相同的兩個(gè)面才能合并，比如對(duì)稱面，否則會(huì)破壞幾何外形和邊界條件。

圖2 合并三棱柱側(cè)面Fig．2 Procedure of face agglomeration

邊界條件的處理比較簡單，將初始網(wǎng)格的邊界條件類型直接賦予聚合后的粗網(wǎng)格即可。在處理好邊界條件后就可以直接輸出網(wǎng)格并進(jìn)行計(jì)算。

2．3 復(fù)雜外形混合網(wǎng)格生成實(shí)例

為了驗(yàn)證本文的混合網(wǎng)格生成方法對(duì)復(fù)雜外形的適應(yīng)性，本文首先選取了第二屆AIAA CFD阻力預(yù)測(cè)會(huì)議（DPW－II）［8］中的 DLR－F6翼身組合體進(jìn)行驗(yàn)證。如圖3所示為該外形數(shù)模和局部表面網(wǎng)格。初始網(wǎng)格是全四面體網(wǎng)格，包括各向異性四面體網(wǎng)格和各向同性四面體網(wǎng)格，共1743萬個(gè)四面體單元，3496萬個(gè)面，295萬個(gè)點(diǎn)。在聚合過程中，α取為0．4，聚合后得到781萬個(gè)網(wǎng)格單元，2080萬個(gè)面，網(wǎng)格單元總量減少大約一半，面的總量大約減少三分之一。如圖4所示為聚合后的混合網(wǎng)格在x＝0．13截面和z＝0．3截面的網(wǎng)格視圖。聚合后邊界層有40層三棱柱網(wǎng)格，外場(chǎng)是四面體網(wǎng)格，從邊界層到外場(chǎng)網(wǎng)格均勻過渡。通過聚合過程不僅網(wǎng)格質(zhì)量得到提高，網(wǎng)格數(shù)量也大幅減少，計(jì)算效率也會(huì)有很大提高。

圖3 DLR－F6數(shù)模和局部表面網(wǎng)格Fig．3 Configuration and surface grids of DLR－F6

為了說明聚合法的魯棒性，本文生成了具有復(fù)雜外形的類F16戰(zhàn)斗機(jī)全機(jī)外形的混合網(wǎng)格。如圖5（a）所示為該外形的表面網(wǎng)格，其幾何構(gòu)型包括機(jī)身、機(jī)翼、平尾、立尾、腹鰭。我們采用聚合法生成該外形的混合網(wǎng)格量約為全場(chǎng)500萬。圖5（b）－（d）分別顯示了典型截面的混合網(wǎng)格。通過該實(shí)例表明，聚合法有很好的魯棒性，生成的混合網(wǎng)格可以保證邊界層內(nèi)絕大部分是三棱柱，而且從三棱柱網(wǎng)格光滑地過渡到外場(chǎng)四面體網(wǎng)格。

3 混合網(wǎng)格生成方法的應(yīng)用

為了驗(yàn)證本文的混合網(wǎng)格生成方法的實(shí)用性，我們對(duì)DLR－F6（帶發(fā)動(dòng)機(jī)和翼吊）外形進(jìn)行數(shù)值計(jì)算，并與風(fēng)洞試驗(yàn)結(jié)果和他人計(jì)算結(jié)果進(jìn)行了比較。本文采用的流場(chǎng)解算器是作者所在團(tuán)隊(duì)開發(fā)的CFD解算軟件 USTAR［9－10］。USTAR 是一款基于非結(jié)構(gòu)／混合網(wǎng)格的格心型有限體積雷諾平均N－S方程定常／非定常解算器。軟件采用迎風(fēng)型格式進(jìn)行空間無粘項(xiàng)離散，粘性項(xiàng)采用中心格式離散。時(shí)間推進(jìn)采用隱式Block LU－SGS算法［10］。為了提高計(jì)算效率，采用多重網(wǎng)格技術(shù)加速；對(duì)于大規(guī)模計(jì)算問題，采用了基于METIS分區(qū)的并行計(jì)算方法。軟件集成了典型湍流模型，包括：一方程SA模型，二方程k－ε模型、k－ω模型和SST模型。該軟件已在多種復(fù)雜外形的數(shù)值模擬中得到成功的應(yīng)用。其中的具體算法請(qǐng)參見文獻(xiàn)［9－10］。

計(jì)算的來流條件為DPW－II中的CASE2，即Ma＝0．75，Re＝3．0×106，攻角從－3．0°到1．5°。計(jì)算過程中采用SST湍流模型，沒有考慮轉(zhuǎn)捩影響。圖6所示為聚合后混合網(wǎng)格和聚合前非結(jié)構(gòu)網(wǎng)格的殘差、俯仰力矩收斂曲線對(duì)比，可以看到混合網(wǎng)格不管在收斂性還是計(jì)算效率上都有很大提高。

圖6 收斂歷程曲線Fig．6 Convergence history

DPW會(huì)議的重點(diǎn)是氣動(dòng)力精確預(yù)測(cè)。對(duì)于氣動(dòng)力預(yù)測(cè)來說，網(wǎng)格質(zhì)量尤其是邊界層內(nèi)網(wǎng)格質(zhì)量的好壞直接影響到計(jì)算精度。為了驗(yàn)證聚合法生成的混合網(wǎng)格計(jì)算結(jié)果，本文不僅對(duì)比了風(fēng)洞實(shí)驗(yàn)結(jié)果，還與DPW－II中一些參與者的計(jì)算結(jié)果進(jìn)行了對(duì)比［11－14］，包 括 NASA Langley 的 USM3D 和FUN3D，Cessna的非結(jié)構(gòu)解算器NSU3D和日本Kawasaki Heavy Industries的UG3。圖7是聚合后的混合網(wǎng)格（USTAR）、聚合前的非結(jié)構(gòu)四面體網(wǎng)格（Unstructured）、USM3D、FUN3D、NSU3D 計(jì)算結(jié)果對(duì)比。由圖可見，聚合前的非結(jié)構(gòu)網(wǎng)格計(jì)算結(jié)果和其他計(jì)算結(jié)果普遍偏差較大，證明了聚合后有效提高了網(wǎng)格質(zhì)量。

圖7（a）是計(jì)算得到的升力系數(shù)曲線。由圖可見，USTAR計(jì)算結(jié)果在攻角大于0．5°時(shí)和實(shí)驗(yàn)結(jié)果符合較好，在攻角小于0．5°時(shí)與實(shí)驗(yàn)值有一定偏差，但是升力曲線和FUN3D、USM3D的結(jié)果符合較好。DPW－II會(huì)議中，大部分提交的計(jì)算結(jié)果升力曲線也符合這種趨勢(shì)，即升力計(jì)算值比實(shí)驗(yàn)值偏大，斜率比實(shí)驗(yàn)值偏小［15］，這種現(xiàn)象在DPW－I中也普遍存在［15］。

圖7（b）是升阻力極曲線。第一、二屆阻力預(yù)測(cè)會(huì)議中，多數(shù)參與者提交的計(jì)算結(jié)果中阻力都比實(shí)驗(yàn)值偏大，如圖所示，NSU3D和FUN3D的阻力都比實(shí)驗(yàn)值偏大，但USTAR計(jì)算的極曲線和實(shí)驗(yàn)值符合的較為一致。

圖7（c）是俯仰力矩曲線。圖中計(jì)算結(jié)果都和實(shí)驗(yàn)值符合的不是很理想，USTAR的俯仰力矩曲線和USM3D符合較好，處于FUN3D和NSU3D之間。和阻力預(yù)測(cè)會(huì)議收到的大部分提交結(jié)果一樣，除FUN3D以外，俯仰力矩都比實(shí)驗(yàn)值偏大［15］。

除了氣動(dòng)力之外，壓力分布也是本文觀察的一個(gè)重點(diǎn)。圖8是機(jī)翼沿展向不同站位處的壓力分布，圖中為USTAR在混合網(wǎng)格上的計(jì)算結(jié)果與實(shí)驗(yàn)值、USM3D、FUN3D的計(jì)算結(jié)果對(duì)比。在z／b＝0．150處，由于機(jī)翼上表面后緣翼身結(jié)合處分離渦的影響，計(jì)算結(jié)果和實(shí)驗(yàn)值都有一定偏差；在z／b＝0．331站位處，由于機(jī)翼下表面的發(fā)動(dòng)機(jī)內(nèi)側(cè)有分離渦，計(jì)算結(jié)果也和實(shí)驗(yàn)有差異，但USTAR和USM3D明顯要比FUN3D符合的好；在z／b＝0．638站位，該處外形簡單流場(chǎng)結(jié)構(gòu)單一，因此三者都和實(shí)驗(yàn)值均符合較好。

圖7 DLR－F6翼身組合體氣動(dòng)特性計(jì)算Fig．7 Aerodynamic coefficients of DLR－F6

圖8 機(jī)翼不同站位處壓力系數(shù)分布Fig．8 Pressure coefficient distributions at z／b＝0．150，0．331，0．638

為了觀察流場(chǎng)結(jié)構(gòu)，圖9畫出了物面流線。圖9（a）為機(jī)翼上翼面流線，圖9（b）為翼身結(jié)合處流線，可看到有明顯的分離現(xiàn)象。圖9（c）為USTAR計(jì)算的下翼面發(fā)動(dòng)機(jī)短艙內(nèi)側(cè)的流線，圖9（d）為UG3計(jì)算結(jié)果，二者現(xiàn)象一致，分離區(qū)大小也差不多，但是和實(shí)驗(yàn)結(jié)果（圖9e）相比，分離區(qū)明顯偏大，在第二屆阻力預(yù)測(cè)會(huì)議的總結(jié)文獻(xiàn)里也指出大多計(jì)算結(jié)果具有分離區(qū)偏大的特點(diǎn)［15］。

圖9 USTAR、UG3計(jì)算的物面流線和實(shí)驗(yàn)結(jié)果Fig．9 Streamlines over DLR－F6

4 結(jié)束語

本文針對(duì)復(fù)雜外形粘性流動(dòng)網(wǎng)格生成費(fèi)時(shí)費(fèi)力的問題，發(fā)展了一種魯棒較好、自動(dòng)化程度較高的混合網(wǎng)格生成方法。聚合法生成的網(wǎng)格集成三棱柱、四面體和少數(shù)的多面體，在減少網(wǎng)格數(shù)量的同時(shí)提高了網(wǎng)格質(zhì)量以增強(qiáng)對(duì)流場(chǎng)結(jié)構(gòu)的捕捉能力。通過采用第二屆阻力預(yù)測(cè)會(huì)議的標(biāo)準(zhǔn)算例對(duì)該方法進(jìn)行了驗(yàn)證，結(jié)果表明該方法生成的混合網(wǎng)格質(zhì)量適應(yīng)計(jì)算的需要。通過對(duì)類F16復(fù)雜戰(zhàn)斗機(jī)全機(jī)構(gòu)型的混合網(wǎng)格生成，進(jìn)一步說明其魯棒性較好，表明其具有廣闊的應(yīng)用前景。

［1］ BAKER T J．Mesh generation：art or science［J］．Progress in Aerospace Sciences．2005，41：29－63．

［2］ PIRZADEH S．Three－dimensional unstructured viscous grids by the advancing－layers method［J］．AIAA Journal，1996，34：43－49．

［3］ KALLINDERIS Y，KHAWAJA A．Hybrid prismatic tetrahedral grid generation for complex geometries［J］．AIAA Journal，1996，34：291－298．

［4］ MATSUNO K．Hyperbolic upwind method for prismatic grid generation［R］．AIAA2000－1003，2000．

［5］ MATSUNO K．High order upwind method for hyperbolic grid generation［J］．Computers ＆ Fluids，1999，28：825－851．

［6］ LOHNER R，PARIKH P．Generation of three dimensional unstructured grids by the advancing front method［J］．Int．J．Numer．Meth．Fluids，1988，8：1135－49．

［7］ PIRZADEH S．Unstructured viscous grid generation by Advancing－Layers Method［R］．AIAA 93－3453，1993．

［8］ 2nd AIAA CFD Drag Prediction Workshop［EB］．http：／／aaac．larc．nasa．gov／tsab／cfdlarc／aiaa－dpw／，June 2003．

［9］ ZHANG L P，YANG Y J，ZHANG H X．Numerical simulations of 3Dinviscid／viscous flow fields on Cartesian／unstructured／prismatic hybrid grids［R］．In Proceedings of the 4th Asian CFD Conference［C］，Mianyang，Sichuan，China，September 2000．

［10］ZHANG L P，WANG Z J．A block LU－SGS implicit dual time－stepping algorithm for hybrid dynamic meshes［J］．Computers ＆ Fluids，2004，33：891－916．

［11］YAMAMOTO K，OCHI A．CFD sensitivity of drag prediction on DLR－F6configuration by structured method and unstructured method［R］．AIAA2004－398，2004．

［12］MAVRIPLIS D J．Drag prediction of DLR－F6using the turbulent Navier－Stokes calculations with multigrid［R］．AIAA2004－397，2004．

［13］LEE－RAUSCH E M，MAVRIPLIS D J．Transonic drag prediction on a DLR－F6transport configuration using unstructured solvers［R］．AIAA2004－554，2004．

［14］SCLAFANI A J，DEHAAN M A．OVERFLOW drag prediction for the DLR－F6transport configuration：a DPW－II case study［R］．AIAA 2004－393，2004．

［15］LAFLIN K R，BRODERSEN O．Summary of data from the second AIAA CFD Drag Prediction Workshop［R］．AIAA2004－0555，2004．