王文， 閻超,*， 袁武， 黃宇， 席柯

1.北京航空航天大學(xué) 航空科學(xué)與工程學(xué)院， 北京 100083 2.中國(guó)科學(xué)院 計(jì)算機(jī)網(wǎng)絡(luò)信息中心 超級(jí)計(jì)算中心, 北京 100190 3.中國(guó)兵器工業(yè)導(dǎo)航與控制技術(shù)研究所, 北京 100089

魯棒的結(jié)構(gòu)網(wǎng)格自動(dòng)化重疊方法

王文1， 閻超1,*， 袁武2， 黃宇1， 席柯3

1.北京航空航天大學(xué) 航空科學(xué)與工程學(xué)院， 北京 100083 2.中國(guó)科學(xué)院 計(jì)算機(jī)網(wǎng)絡(luò)信息中心 超級(jí)計(jì)算中心, 北京 100190 3.中國(guó)兵器工業(yè)導(dǎo)航與控制技術(shù)研究所, 北京 100089

網(wǎng)格技術(shù)是目前數(shù)值模擬中的關(guān)鍵技術(shù)之一。重疊網(wǎng)格是一種放寬拓?fù)湟?、減小網(wǎng)格生成難度的網(wǎng)格技術(shù)。本文以結(jié)構(gòu)重疊網(wǎng)格為基礎(chǔ)，分別針對(duì)挖洞、尋點(diǎn)以及洞面優(yōu)化方法進(jìn)行了研究和改進(jìn)，同時(shí)完成物面網(wǎng)格重疊，形成了一套魯棒的、自動(dòng)化的網(wǎng)格重疊系統(tǒng)。在挖洞方面，結(jié)合“最小洞映射”方法，提出“復(fù)合式挖洞”方法，節(jié)省內(nèi)存開(kāi)銷(xiāo)；在尋點(diǎn)方面，通過(guò)構(gòu)建格心虛網(wǎng)格，保證搜索空間的連續(xù)性，同時(shí)結(jié)合“有效搜索”思想，排除部分對(duì)尋點(diǎn)無(wú)貢獻(xiàn)的網(wǎng)格點(diǎn)，進(jìn)而減少ADT叉樹(shù)節(jié)點(diǎn)；在洞面優(yōu)化上，改變填補(bǔ)判別法則并引入兩類(lèi)受保護(hù)洞內(nèi)點(diǎn)，確保兩層插值邊界建立，提高魯棒性；在物面網(wǎng)格重疊上，利用物面投影法完成坐標(biāo)修正，實(shí)現(xiàn)物面附近網(wǎng)格流動(dòng)變量的準(zhǔn)確傳遞。為驗(yàn)證本文方法，分別對(duì)定常翼身組合體DLR-F6繞流和非定常機(jī)翼掛載分離過(guò)程進(jìn)行了數(shù)值模擬，計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果吻合良好，表明該結(jié)構(gòu)重疊網(wǎng)格系統(tǒng)對(duì)多物體間定常、非定常擾流具有較好的數(shù)值模擬能力和較高的模擬精度，具有較高的工程應(yīng)用價(jià)值。

計(jì)算流體力學(xué)； 結(jié)構(gòu)重疊網(wǎng)格； 洞映射； 尋點(diǎn)； 洞面優(yōu)化方法； 翼身組合體； 機(jī)翼掛載分離

隨著人類(lèi)對(duì)航空航天飛行器性能的不斷追求，近幾十年來(lái)航空航天技術(shù)有了突飛猛進(jìn)的發(fā)展，與此同時(shí)，飛行器的構(gòu)型越來(lái)越復(fù)雜，并且經(jīng)常伴隨著相對(duì)運(yùn)動(dòng)，典型的例子有：高升力體外形(High-lift Configuration)升阻力預(yù)測(cè)[1]、外掛物與載機(jī)分離過(guò)程[2]、固體火箭助推器與航天飛機(jī)分離[3]、子母彈拋撒過(guò)程[4]、直升機(jī)旋翼運(yùn)動(dòng)[5]等。在單個(gè)或者多個(gè)復(fù)雜構(gòu)型數(shù)值模擬過(guò)程中，網(wǎng)格生成技術(shù)是關(guān)鍵技術(shù)之一，網(wǎng)格質(zhì)量的好壞直接影響CFD計(jì)算結(jié)果的精度。因此，如何快速地、高質(zhì)量地生成網(wǎng)格成為了目前CFD技術(shù)中亟待解決的問(wèn)題之一。

目前處理復(fù)雜構(gòu)型的網(wǎng)格主要包括：結(jié)構(gòu)網(wǎng)格、非結(jié)構(gòu)網(wǎng)格和混合網(wǎng)格。結(jié)構(gòu)網(wǎng)格由于數(shù)據(jù)結(jié)構(gòu)簡(jiǎn)單、計(jì)算技術(shù)成熟、計(jì)算效率和計(jì)算精度較高以及黏性模擬能力強(qiáng)等優(yōu)點(diǎn)而被廣泛應(yīng)用，但復(fù)雜構(gòu)型飛行器生成一體化結(jié)構(gòu)網(wǎng)格仍然顯得十分困難，且不易保證網(wǎng)格質(zhì)量，有時(shí)簡(jiǎn)單構(gòu)型物體組合在一起，空間拓?fù)湟矔?huì)變得異常復(fù)雜。與此同時(shí)，剛性的結(jié)構(gòu)網(wǎng)格也無(wú)法處理有較大相對(duì)運(yùn)動(dòng)的情況。為降低結(jié)構(gòu)網(wǎng)格生成難度，Steger等[6]在20世紀(jì)80年代提出了重疊網(wǎng)格(Chimera Grids)方法。重疊網(wǎng)格也稱(chēng)為覆蓋網(wǎng)格(Overlapping Grids)或嵌套網(wǎng)格(Overset Grids)，它是將計(jì)算區(qū)域劃分成多個(gè)獨(dú)立的子區(qū)域，各子區(qū)域網(wǎng)格隨部件作剛性運(yùn)動(dòng)，流場(chǎng)信息在網(wǎng)格重疊部分通過(guò)插值進(jìn)行傳遞。重疊網(wǎng)格放寬了網(wǎng)格拓?fù)湟?，極大降低了網(wǎng)格生成難度，提高了網(wǎng)格質(zhì)量，因此得到了大力的發(fā)展和廣泛的應(yīng)用，同時(shí)涌現(xiàn)出一大批著名的軟件，如Pegasus5[7]、Overgrid[8]、Suggar[9]等。國(guó)內(nèi)對(duì)重疊網(wǎng)格的研究起步相對(duì)較晚，但發(fā)展十分迅猛，在工程領(lǐng)域也得到逐步應(yīng)用。如陳十一等[10]采用結(jié)構(gòu)重疊網(wǎng)格技術(shù)，利用壁面約束的大渦模擬(CLES)、分離渦模擬(DES)以及雷諾平均Navier-Stokes(RANS)方法對(duì)商用飛機(jī)進(jìn)行全機(jī)數(shù)值模擬；曹義華等[11]采用結(jié)構(gòu)重疊網(wǎng)格對(duì)帶旋翼的直升機(jī)進(jìn)行了數(shù)值模擬并對(duì)空間數(shù)值格式進(jìn)行評(píng)價(jià)；蔡晉生等[12]發(fā)展了一種多層多塊隱式嵌套重疊網(wǎng)格技術(shù)并對(duì)翼身平尾組合體進(jìn)行了數(shù)值模擬與分析；田書(shū)玲等[13-15]開(kāi)展了基于非結(jié)構(gòu)網(wǎng)格方法的動(dòng)態(tài)重疊網(wǎng)格算法研究。

本文主要是對(duì)現(xiàn)有的結(jié)構(gòu)重疊網(wǎng)格技術(shù)進(jìn)行研究，提出了魯棒的重疊方法，主要包括4個(gè)方面：① 完善洞映射方法，提出“復(fù)合式挖洞”方法，節(jié)省內(nèi)存開(kāi)銷(xiāo)；② 提出“有效搜索”思想，提高尋點(diǎn)效率；③ 將混合洞面優(yōu)化方法推廣至非定常多體運(yùn)動(dòng)數(shù)值模擬中，并引入兩類(lèi)受保護(hù)洞內(nèi)點(diǎn)，提高魯棒性；④ 采用物面投影法對(duì)物面附近重疊網(wǎng)格進(jìn)行坐標(biāo)修正，保證流動(dòng)變量在物面附近正確傳遞。最后，選取翼身組合體DLR-F6擾流和機(jī)翼掛載分離算例對(duì)本重疊網(wǎng)格技術(shù)的性能進(jìn)行了檢驗(yàn)和驗(yàn)證。

1 挖洞方法

挖洞是重疊網(wǎng)格需要進(jìn)行的第1步操作，目的是在流場(chǎng)計(jì)算前屏蔽掉一些不必要或者無(wú)實(shí)際意義的點(diǎn)，如落入物面、對(duì)稱(chēng)面或人為指定曲面內(nèi)的點(diǎn)，并將其標(biāo)記為“洞內(nèi)點(diǎn)”，以區(qū)別于參與流場(chǎng)計(jì)算的“洞外點(diǎn)”，也稱(chēng)為“正常點(diǎn)”。網(wǎng)格域內(nèi)，區(qū)分“洞內(nèi)點(diǎn)”和“洞外點(diǎn)”的界線(xiàn)稱(chēng)為“洞邊界”。常見(jiàn)的挖洞方法有點(diǎn)矢法[16]、射線(xiàn)求交法[17]、洞映射法[18]、X射線(xiàn)(Object X-Rays)法[19]、叉樹(shù)法[20]、直接法[21]、隱式切割(Implicit Hole Cutting， IHC)法[22]等，本文選取發(fā)展較為成熟、使用靈活的洞映射法為挖洞方法。

洞映射方法將曲物面投射到輔助的直角笛卡兒網(wǎng)格中，得到由笛卡兒網(wǎng)格構(gòu)成的近似挖洞面。根據(jù)笛卡兒網(wǎng)格與物面的相對(duì)位置關(guān)系，又將其劃分為“外部單元”、“內(nèi)部單元”和“邊緣單元”，邊緣單元的集合即是近似挖洞面。在洞映射方法的應(yīng)用中發(fā)現(xiàn)，包圍物體所有物面單元建立洞映射模型包含了較多無(wú)效的空間位置。以泰坦四號(hào)(Titan IV)運(yùn)載火箭[23]為例，主火箭與助推級(jí)軸向長(zhǎng)度存在較大差異，助推級(jí)網(wǎng)格僅與主火箭中后部區(qū)域重疊，因此主火箭頭部附近物面所構(gòu)建的洞映射模型無(wú)實(shí)際意義。為削減洞映射模型內(nèi)存消耗量，在模型建立時(shí)先判斷各物體網(wǎng)格包圍盒相交的情況，僅在其他物體網(wǎng)格涵蓋的范圍內(nèi)建立洞映射模型。這樣，新的洞映射邊界就由物面、網(wǎng)格邊界和截?cái)嗝鏄?gòu)成，顯然這是一個(gè)可供其他網(wǎng)格挖洞的最“緊湊”區(qū)域，稱(chēng)為“最小洞映射[24]”。泰坦火箭所建立的最小洞映射模型如圖1所示，方框?yàn)楦鞫从成淠Ｐ停渲兄骰鸺从成淠Ｐ洼^之前有顯著減小，網(wǎng)格量削減了75%。

圖1 Titan IV運(yùn)載火箭洞映射模型Fig.1 Hole mapping model of Titan IV launch vehicle

當(dāng)物面網(wǎng)格發(fā)生重疊時(shí)，空間網(wǎng)格離散方式發(fā)生了一些變化，物體初始體網(wǎng)格也可能是洞內(nèi)點(diǎn)。以圖 2為例，圖中機(jī)身網(wǎng)格為去掉機(jī)翼單獨(dú)生成，拓?fù)浣Y(jié)構(gòu)簡(jiǎn)單，網(wǎng)格質(zhì)量較好，機(jī)翼網(wǎng)格也單獨(dú)生成，然后將機(jī)翼網(wǎng)格嵌套入機(jī)身網(wǎng)格中，這樣，落在機(jī)翼物面內(nèi)的機(jī)身網(wǎng)格實(shí)質(zhì)為洞內(nèi)點(diǎn)，“Inverse Mark”[24]標(biāo)識(shí)方法不再適用。為此，本文將已建好的洞映射單元標(biāo)識(shí)為臨時(shí)單元，將有物面網(wǎng)格落入的單元標(biāo)識(shí)為過(guò)渡單元；對(duì)洞映射模型的8個(gè)角單元進(jìn)行判斷，若當(dāng)前角單元不是過(guò)渡單元，則將該單元標(biāo)識(shí)為正常單元；感染所有正常單元的相鄰單元，若該單元為臨時(shí)單元，則將其標(biāo)識(shí)為正常單元，直至再?zèng)]有臨時(shí)單元轉(zhuǎn)換為正常單元為止；最后將剩余的臨時(shí)單元標(biāo)識(shí)為內(nèi)部單元，進(jìn)而完成了洞映射模型的標(biāo)識(shí)。

圖2 機(jī)身網(wǎng)格洞內(nèi)點(diǎn)Fig.2 Hole points of wing grid

由于洞映射方法的固有特點(diǎn)，對(duì)模型表面模擬得越精細(xì)，所需的笛卡兒網(wǎng)格量就越大，內(nèi)存開(kāi)銷(xiāo)就越大。PEGASUS5.2[25]就通過(guò)細(xì)分洞映射模型，區(qū)別各洞映射模型笛卡兒網(wǎng)格尺度，針對(duì)精細(xì)結(jié)構(gòu)洞映射模型采用較小尺寸的笛卡兒網(wǎng)格，如窄縫等，從而在保證對(duì)模型細(xì)節(jié)精細(xì)化描述的前提下盡可能節(jié)省網(wǎng)格開(kāi)銷(xiāo)。本文對(duì)落入洞映射模型過(guò)渡單元的網(wǎng)格點(diǎn)進(jìn)行一次再判斷，在不增加原有笛卡兒網(wǎng)格量的基礎(chǔ)上解決物面附近網(wǎng)格點(diǎn)的屬性判斷問(wèn)題，本文稱(chēng)之為“復(fù)合式挖洞方法”。圖 3為復(fù)合式挖洞方法示意圖，圖中黑色方框?yàn)槎从成淠Ｐ?，左下角白色單元為?nèi)部單元，為洞內(nèi)區(qū)域，灰色單元為過(guò)渡單元，右上角為外部單元，網(wǎng)格為該物體體網(wǎng)格。圖中包含四種點(diǎn)，落在外部單元的正常點(diǎn)P1，落在內(nèi)部單元的洞內(nèi)點(diǎn)P4，還有一種則是落在過(guò)渡單元的點(diǎn)，可能是洞內(nèi)點(diǎn)如P3，也可能是正常點(diǎn)如P2。加密洞映射也能夠?qū)2和P3進(jìn)行更精確的判斷，但增加了網(wǎng)格量。在此，本文通過(guò)對(duì)落入過(guò)渡單元的點(diǎn)進(jìn)行尋找貢獻(xiàn)單元的操作，進(jìn)而確定該點(diǎn)的屬性，具體的實(shí)現(xiàn)方法如下。

1) 建立最小洞映射模型，通過(guò)上述方法標(biāo)識(shí)笛卡兒網(wǎng)格的外部單元、內(nèi)部單元和過(guò)渡單元。

2) 通過(guò)洞映射模型對(duì)其他網(wǎng)格進(jìn)行挖洞，將落入外部單元的點(diǎn)標(biāo)識(shí)為正常點(diǎn)Fix，將落入內(nèi)部單元的點(diǎn)標(biāo)識(shí)成為受保護(hù)洞內(nèi)點(diǎn)Pthole。

3) 針對(duì)有物面網(wǎng)格重疊的物體，進(jìn)行非重疊的物面網(wǎng)格單元尋點(diǎn)操作，若貢獻(xiàn)單元與該網(wǎng)格同屬一個(gè)物體且貢獻(xiàn)單元網(wǎng)格層級(jí)較低，則將該貢獻(xiàn)單元各頂點(diǎn)標(biāo)識(shí)為受保護(hù)洞內(nèi)點(diǎn)Pthole；再對(duì)體網(wǎng)格單元進(jìn)行尋點(diǎn)操作，若貢獻(xiàn)單元與該網(wǎng)格同屬一個(gè)物體且貢獻(xiàn)單元網(wǎng)格層級(jí)較低，則將貢獻(xiàn)單元的正常點(diǎn)Fix標(biāo)識(shí)為普通洞內(nèi)點(diǎn)Hole。

4) 針對(duì)物面網(wǎng)格重疊的物體體網(wǎng)格，感染受保護(hù)洞內(nèi)點(diǎn)Pthole毗鄰的正常點(diǎn)Fix，將其標(biāo)識(shí)為受保護(hù)洞內(nèi)點(diǎn)Pthole，直至再?zèng)]有正常點(diǎn)Fix發(fā)生轉(zhuǎn)化為止。

5) 對(duì)落入其他物體洞映射過(guò)渡單元的網(wǎng)格點(diǎn)，若能找到有效的貢獻(xiàn)單元(存在貢獻(xiàn)單元且貢獻(xiàn)單元不為受保護(hù)洞內(nèi)點(diǎn)Pthole或者普通洞內(nèi)點(diǎn)Hole)，則該點(diǎn)為正常點(diǎn)Fix，否則為普通洞內(nèi)點(diǎn)Hole。

6) 針對(duì)拉伸背景網(wǎng)格，進(jìn)行網(wǎng)格單元尋點(diǎn)操作，若存在貢獻(xiàn)單元，則將該網(wǎng)格標(biāo)識(shí)為普通洞內(nèi)點(diǎn)Hole。

對(duì)于有物面網(wǎng)格重疊的情況，需要對(duì)同屬一個(gè)物體的網(wǎng)格進(jìn)行分層處理，以方便挖洞過(guò)程。如翼身組合體，將模型拆成機(jī)身部分和機(jī)翼部分，兩部分單獨(dú)生成拓?fù)浜?jiǎn)單的網(wǎng)格，其中機(jī)翼網(wǎng)格為高層網(wǎng)格，機(jī)身網(wǎng)格為低層網(wǎng)格。

圖3 洞映射與四類(lèi)點(diǎn)示意圖Fig.3 Sketch map of hole-mapping cells and four kinds of nodes

以較為常見(jiàn)的導(dǎo)彈為例，導(dǎo)彈由彈身及氣動(dòng)舵組成，舵與彈身之間一般存在著細(xì)微的縫隙，若洞映射單元尺寸過(guò)太，縫隙將被邊鄰單元填補(bǔ)，無(wú)法實(shí)現(xiàn)挖洞，因此在這種情況下洞映射單元尺寸得小于縫隙。傳統(tǒng)洞映射方法與最小洞映射方法生成的導(dǎo)彈洞映射截面如圖 4(a)所示，復(fù)合式挖洞方法生成的導(dǎo)彈洞映射截面如圖 4(b)所示。傳統(tǒng)洞映射方法、最小洞映射方法和復(fù)合式挖洞方法的洞映射笛卡兒?jiǎn)卧约巴诙磿r(shí)間對(duì)比見(jiàn)表1，從表可以看出復(fù)合式挖洞方法顯著減少了笛卡兒網(wǎng)格數(shù)目，同時(shí)相對(duì)其他兩種方法并沒(méi)有消耗過(guò)多時(shí)間。

通過(guò)復(fù)合式挖洞方法可以在不額外增加內(nèi)存開(kāi)銷(xiāo)的前提下對(duì)網(wǎng)格進(jìn)行挖洞處理，同時(shí)保證挖洞效率。

圖4 3種挖洞方法洞映射截面對(duì)比Fig.4 Comparison of hole-mapping section among three hole-cutting algorithms

表1 3種挖洞方法比較Table 1 Comparison among three hole-cutting algorithms

AlgorithmCartesiangridTimeofhole?cutting/sOriginal59790560．4744Minimum18440530．4712Composite1651720．5035

2 基于虛網(wǎng)格的格心ADT搜索方法

尋點(diǎn)是重疊網(wǎng)格的另一項(xiàng)關(guān)鍵技術(shù)，是確定網(wǎng)格點(diǎn)在其他網(wǎng)格中的位置、查找網(wǎng)格貢獻(xiàn)單元的統(tǒng)稱(chēng)。在整個(gè)重疊網(wǎng)格過(guò)程中，需要反復(fù)尋找插值點(diǎn)的位置及插值點(diǎn)的貢獻(xiàn)單元，因此準(zhǔn)確性和效率是尋點(diǎn)方法的關(guān)鍵。目前基于格心網(wǎng)格的有限體積法在CFD中應(yīng)用廣泛，而格心網(wǎng)格塊與塊之間存在間隙，空間上不連續(xù)，這就給ADT(Alternating Digital Tree)[26]叉樹(shù)結(jié)構(gòu)的建立造成了不便。目前一種常用的解決方法是通過(guò)虛網(wǎng)格將格心單元進(jìn)行延拓，填補(bǔ)塊與塊之間的空白部分，保證空間的連續(xù)性。本文先通過(guò)線(xiàn)性外插的方式，分別沿i、j、k3個(gè)方向?qū)Ω鱾€(gè)塊的格心網(wǎng)格進(jìn)行延拓，再利用對(duì)接邊界條件，通過(guò)查找方塊相應(yīng)點(diǎn)對(duì)虛網(wǎng)格坐標(biāo)進(jìn)行修正，使得分區(qū)邊界上的虛網(wǎng)格坐標(biāo)嚴(yán)格對(duì)應(yīng)于對(duì)應(yīng)塊實(shí)網(wǎng)格坐標(biāo)，實(shí)現(xiàn)格心虛網(wǎng)格的延拓。通過(guò)這種方式，既可以實(shí)現(xiàn)內(nèi)部網(wǎng)格點(diǎn)的準(zhǔn)確延拓，又能實(shí)現(xiàn)遠(yuǎn)場(chǎng)、對(duì)稱(chēng)面等邊界網(wǎng)格的嚴(yán)格對(duì)接。

在完成格心單元虛網(wǎng)格延拓之后，首先需要組織ADT叉樹(shù)，對(duì)空間網(wǎng)格單元進(jìn)行規(guī)整。以二維情況為例，圖 5中包含了7個(gè)網(wǎng)格單元，節(jié)點(diǎn)1包含整個(gè)區(qū)域，為ADT叉樹(shù)的根節(jié)點(diǎn)。在根節(jié)點(diǎn)下插入子節(jié)點(diǎn)——節(jié)點(diǎn)2，由于該節(jié)點(diǎn)落在節(jié)點(diǎn)1所指代區(qū)域的右半?yún)^(qū)域，因此將節(jié)點(diǎn)1所指代的空間二分，節(jié)點(diǎn)2指代右半?yún)^(qū)域。同理，節(jié)點(diǎn)3指代左半?yún)^(qū)域，以此類(lèi)推，形成了如圖 5所示的ADT叉樹(shù)結(jié)構(gòu)。圖中各節(jié)點(diǎn)旁的外邊框?yàn)锳DT叉樹(shù)所包含的整個(gè)區(qū)域，灰色部分為該節(jié)點(diǎn)所指代的區(qū)域，如節(jié)點(diǎn)3指代左半?yún)^(qū)域。假設(shè)有一待查詢(xún)點(diǎn)P，如圖 5所示。將點(diǎn)P代入ADT叉樹(shù)中進(jìn)行搜索，點(diǎn)P落入根節(jié)點(diǎn)——節(jié)點(diǎn)1所指代的整個(gè)區(qū)域，因此節(jié)點(diǎn)1是備選單元，同理，節(jié)點(diǎn)2、5、6所指代的區(qū)域同樣包含點(diǎn)P，故也是備選單元。然后通過(guò)單元的空間包圍盒是否相交排除一部分備選單元，剩余的備選單元通過(guò)等參變換，計(jì)算出點(diǎn)P在該單元的權(quán)值，再結(jié)合模版漫游(Stencil-walking[27])法，確定最終的貢獻(xiàn)單元以及插值系數(shù)。

圖5 ADT搜索方法Fig.5 ADT searching method

對(duì)于有多個(gè)物體重疊的情況，應(yīng)盡可能減小搜索范圍，以提高搜索效率。圖6為3個(gè)物體重疊時(shí)搜索空間示意圖，圖中實(shí)線(xiàn)為各物體網(wǎng)格外邊界，虛線(xiàn)為包圍盒。以O(shè)bj3為例，Obj3中處于陰影區(qū)域內(nèi)的網(wǎng)格不與其他物體網(wǎng)格相交，故該區(qū)域內(nèi)網(wǎng)格不會(huì)成為其他物體網(wǎng)格的貢獻(xiàn)單元，為此在ADT叉樹(shù)建立時(shí)排除該部分網(wǎng)格，以減少ADT叉樹(shù)中無(wú)效結(jié)點(diǎn)的數(shù)目，本文稱(chēng)之為“有效搜索”。

有效搜索對(duì)各物體間網(wǎng)格包圍盒相交較少的情況效果顯著。以泰坦四號(hào)運(yùn)載火箭為例，主火箭網(wǎng)格為200萬(wàn)，助推火箭網(wǎng)格為30萬(wàn)，網(wǎng)格如圖 7所示，其中助推火箭網(wǎng)格僅延伸至主火箭中部，主火箭前部網(wǎng)格以及外部網(wǎng)格對(duì)尋點(diǎn)毫無(wú)貢獻(xiàn)，故通過(guò)有效搜索能將其剔除。改進(jìn)前后ADT節(jié)點(diǎn)數(shù)目、ADT建立時(shí)間以及貢獻(xiàn)單元搜索時(shí)間如表 2所示，其中ADT節(jié)點(diǎn)數(shù)目、ADT叉樹(shù)建立時(shí)間及貢獻(xiàn)單元搜索時(shí)間顯著減少。

圖6 3個(gè)物體重疊時(shí)的搜索空間Fig.6 Overlapping search space of three objects

圖7 Titan IV運(yùn)載火箭網(wǎng)格重疊Fig.7 Overlapping grids of Titan IV launch vehicle

表2 原始的和改進(jìn)的尋點(diǎn)方法比較

Table2Comparisonbetweenoriginalandimproveddonorsearchalgorithms

AlgorithmADTnodeTimeofATDbuilding/sTimeofdonorsearch/sOriginal23929804．37586．1862Improved7023191．10933．3455

3 改進(jìn)的混合洞面優(yōu)化方法

洞面優(yōu)化方法目的是改善重疊區(qū)域，將“挖洞”過(guò)程后建立的初始插值邊界推離物面，從而避免插值點(diǎn)落入邊界層內(nèi)部降低流場(chǎng)求解精度。插值邊界的位置以及它與貢獻(xiàn)單元的相互匹配關(guān)系影響整個(gè)流場(chǎng)的計(jì)算，因此插值邊界的確定是整個(gè)洞面優(yōu)化過(guò)程的關(guān)鍵。

文獻(xiàn)[28]對(duì)典型的洞面優(yōu)化方法進(jìn)行了概述，同時(shí)通過(guò)對(duì)割補(bǔ)法的分析，針對(duì)其洞面切割后洞邊界位置不確定性及不可預(yù)知性問(wèn)題，提出將其與物面距離優(yōu)化準(zhǔn)則相結(jié)合的新型混合洞面優(yōu)化方法。在后續(xù)的算例測(cè)試中發(fā)現(xiàn)，由于填補(bǔ)階段回填終止判據(jù)過(guò)于嚴(yán)格，陣面會(huì)推進(jìn)異常。圖8 為RAE2822翼型測(cè)試算例，網(wǎng)格由三部分組成，翼型前部網(wǎng)格Grid 1、后部網(wǎng)格Grid 2以及拉伸的直角背景網(wǎng)格Grid 3。由于切割后網(wǎng)格縫隙正好也處于虛線(xiàn)框處，填補(bǔ)時(shí)該區(qū)域洞邊界單元較多，早期的判別準(zhǔn)則要求貢獻(xiàn)單元全部為正常點(diǎn)，而背景網(wǎng)格洞邊界插值點(diǎn)的貢獻(xiàn)單元大多為洞邊界單元，故陣面會(huì)持續(xù)推進(jìn)，進(jìn)而出現(xiàn)圖中網(wǎng)格異?，F(xiàn)象。而對(duì)于窄縫問(wèn)題，例如機(jī)翼掛載分離算例中掛載物與掛架之間的窄縫，填補(bǔ)階段陣面會(huì)推進(jìn)至洞內(nèi)點(diǎn)位置，正常點(diǎn)毗鄰洞內(nèi)點(diǎn)造成插值陣面部分缺失從而導(dǎo)致該區(qū)域重疊失敗。

圖8 RAE2822翼型網(wǎng)格重疊圖Fig.8 Overlapping grids of RAE2822 airfoil

本文對(duì)填補(bǔ)過(guò)程判別法則進(jìn)行修改，將以插值點(diǎn)貢獻(xiàn)單元全為正常點(diǎn)的判別條件改為直接設(shè)定陣面推進(jìn)步數(shù)，同時(shí)引入兩類(lèi)受保護(hù)洞內(nèi)點(diǎn)Ntpthole 1與Ntpthole 2，確保兩層插值邊界的建立，提高優(yōu)化過(guò)程的魯棒性。改進(jìn)型的混合洞面優(yōu)化方法具體過(guò)程如下。

1) 求解各網(wǎng)格點(diǎn)到本物體的物面距離d0，若該網(wǎng)格點(diǎn)隸屬于背景網(wǎng)格，則將d0記為一個(gè)較大的值，如1×106m。

2) 將所有網(wǎng)格點(diǎn)標(biāo)記成為正常點(diǎn)Fix。

3) 通過(guò)復(fù)合式挖洞方法對(duì)網(wǎng)格進(jìn)行挖洞。

4) 對(duì)所有的網(wǎng)格點(diǎn)進(jìn)行ADT尋點(diǎn)操作，找到其所有貢獻(xiàn)單元，且貢獻(xiàn)單元不與該點(diǎn)同屬于一物體，通過(guò)三線(xiàn)性插值方法計(jì)算出該點(diǎn)到其他物面的最小距離dmin，若d0

5) 對(duì)所有落在遠(yuǎn)場(chǎng)邊界條件上的正常點(diǎn)Fix進(jìn)行ADT尋點(diǎn)操作，若存在有效貢獻(xiàn)單元(貢獻(xiàn)單元不為洞內(nèi)點(diǎn)Pthole或Hole)，將其標(biāo)識(shí)為臨界點(diǎn)Crt，將其相鄰的非落在遠(yuǎn)場(chǎng)邊界上的正常點(diǎn)Fix標(biāo)記為洞邊界點(diǎn)Frg 1。

6) 將與洞內(nèi)點(diǎn)Pthole、Hole相鄰的正常點(diǎn)Fix標(biāo)識(shí)為洞邊界點(diǎn)Frg 1。

7) 對(duì)洞邊界進(jìn)行切割操作。針對(duì)所有洞邊界點(diǎn)Frg 1尋點(diǎn)，若存在有效的貢獻(xiàn)單元，則該點(diǎn)標(biāo)識(shí)為普通洞內(nèi)點(diǎn)Hole，與之相鄰的正常點(diǎn)Fix標(biāo)識(shí)為洞邊界點(diǎn)Frg 1，直至再?zèng)]有正常點(diǎn)Fix轉(zhuǎn)化為洞邊界點(diǎn)Frg 1為止。

8) 將與受保護(hù)洞內(nèi)點(diǎn)Pthole相鄰的非受保護(hù)洞內(nèi)點(diǎn)標(biāo)識(shí)為第1類(lèi)受保護(hù)洞內(nèi)點(diǎn)Ntpthole 1。

9) 將與第1類(lèi)受保護(hù)洞內(nèi)點(diǎn)Ntpthole 1相鄰的非受保護(hù)洞內(nèi)點(diǎn)且非第1類(lèi)受保護(hù)洞內(nèi)點(diǎn)標(biāo)識(shí)為第2類(lèi)受保護(hù)洞內(nèi)點(diǎn)Ntpthole 2。

10) 對(duì)洞邊界進(jìn)行回填操作。針對(duì)所有洞邊界點(diǎn)Frg 1，將該點(diǎn)轉(zhuǎn)化為正常點(diǎn)Fix，將其相鄰的普通洞內(nèi)點(diǎn)Hole轉(zhuǎn)化為洞邊界Frg 1。

11) 重復(fù)過(guò)程10)N次(經(jīng)驗(yàn)值N=4)。

12) 將毗鄰正常點(diǎn)Fix的第2類(lèi)受保護(hù)洞內(nèi)點(diǎn)Ntpthole 2標(biāo)記成為洞邊界點(diǎn)Frg 1。

13) 將毗鄰洞邊界點(diǎn)Frg 1的非正常點(diǎn)標(biāo)識(shí)為洞邊界點(diǎn)Frg 2。

14) 檢索所有的插值點(diǎn)Frg 1和Frg 2，將非正常插值點(diǎn)(不存在貢獻(xiàn)單元或貢獻(xiàn)單元不全是正常點(diǎn))用其相鄰正常插值點(diǎn)貢獻(xiàn)單元代替。

15) 將剩余的兩類(lèi)受保護(hù)洞內(nèi)點(diǎn)Ntpthole 1與Ntpthole 2轉(zhuǎn)化為普通洞內(nèi)點(diǎn)Hole。

16) 針對(duì)相對(duì)運(yùn)動(dòng)數(shù)值模擬，對(duì)所有的洞內(nèi)點(diǎn)Hole、臨界點(diǎn)Crt進(jìn)行尋點(diǎn)，若存在有效貢獻(xiàn)單元，則將其標(biāo)識(shí)為插值點(diǎn)Frg 3。

其中，正常點(diǎn)Fix為參與計(jì)算的網(wǎng)格點(diǎn)，洞內(nèi)點(diǎn)(包括受保護(hù)洞內(nèi)點(diǎn)Pthole和普通洞內(nèi)點(diǎn)Hole)、臨界點(diǎn)Crt以及洞邊界點(diǎn)Frg 1、Frg 2、Frg 3為不參與流場(chǎng)計(jì)算的點(diǎn)，但洞邊界點(diǎn) Frg 1、Frg 2、Frg 3的守恒變量通過(guò)向正常點(diǎn)Fix插值獲得更新。同時(shí)，網(wǎng)格點(diǎn)一旦被標(biāo)識(shí)為受保護(hù)洞內(nèi)點(diǎn)Pthole，其屬性在優(yōu)化的任一環(huán)節(jié)均不會(huì)被修改。洞面優(yōu)化完成后通過(guò)檢索普通洞內(nèi)點(diǎn)Hole、臨界點(diǎn)Crt，盡可能多地將其轉(zhuǎn)變?yōu)椴逯迭c(diǎn)Frg 3，使其在非定常相對(duì)運(yùn)動(dòng)數(shù)值模擬過(guò)程中守恒變量與當(dāng)?shù)刂灯ヅ?，盡可能減少網(wǎng)格運(yùn)動(dòng)至下一個(gè)時(shí)刻時(shí)普通洞內(nèi)點(diǎn)Hole、臨界點(diǎn)Crt直接變成正常點(diǎn)Fix對(duì)流場(chǎng)收斂帶來(lái)的影響。圖9為該洞面優(yōu)化過(guò)程，簡(jiǎn)要示意了該洞面優(yōu)化的主要過(guò)程，圖 9(a)為復(fù)合式挖洞后的網(wǎng)格示意圖，對(duì)應(yīng)洞面優(yōu)化階段1)～3)；圖9(b)為采用壁面距離準(zhǔn)

圖9 洞面優(yōu)化過(guò)程示意圖Fig.9 Sketch map of hole optimization process

則標(biāo)記普通洞內(nèi)點(diǎn)，對(duì)應(yīng)洞面優(yōu)化階段4)；圖 9(c)為建立第一層插值邊界及網(wǎng)格切割，二者對(duì)應(yīng)洞面優(yōu)化階段5)～7)；圖 9(d)為標(biāo)記兩類(lèi)受保護(hù)洞內(nèi)點(diǎn)以及網(wǎng)格回填，對(duì)應(yīng)洞面優(yōu)化階段8)～12)；圖 9(e)為建立其余的插值邊界，對(duì)應(yīng)洞面優(yōu)化階段13)～16)。圖 10為采用改進(jìn)后的混合洞面優(yōu)化方法優(yōu)化后的網(wǎng)格，圖 8中網(wǎng)格回填異常情況消失。

圖10 RAE2822網(wǎng)格改進(jìn)型混合洞面優(yōu)化結(jié)果Fig.10 RAE2822 grids optimized result by improved mixed overlapping optimization method

4 物面網(wǎng)格重疊技術(shù)

飛行器表面非常復(fù)雜，較難生成高質(zhì)量的多塊結(jié)構(gòu)網(wǎng)格，因此希望在物面網(wǎng)格重疊的基礎(chǔ)上取消空間網(wǎng)格的拓?fù)湎拗?，從而降低空間網(wǎng)格生成難度。在重疊網(wǎng)格生成過(guò)程中，若對(duì)物面網(wǎng)格進(jìn)行分塊后獨(dú)立生成，則有可能因?yàn)楦鞑糠值木W(wǎng)格幾何誤差、物面曲率分辨率和光滑程度等因素的不同而導(dǎo)致重疊區(qū)內(nèi)彼此描述的物面不唯一，即“物面失配”問(wèn)題。圖 11為網(wǎng)格物面重疊中比較常見(jiàn)的兩種情況。圖 11(a)中物面為凹曲面，二者物面網(wǎng)格點(diǎn)均不落入對(duì)方網(wǎng)格中，無(wú)法完成貢獻(xiàn)單元搜索；圖 11(b)中物面為凸曲面，盡管二者物面網(wǎng)格均落入對(duì)方網(wǎng)格中，但貢獻(xiàn)單元遠(yuǎn)離對(duì)方物面，由于邊界層內(nèi)流動(dòng)梯度大，因此容易產(chǎn)生較大的插值誤差進(jìn)而影響流場(chǎng)的求解。

圖11 網(wǎng)格物面與物理物面失配問(wèn)題Fig.11 Miss match of body and grid surface

本文采用PEGASUS5[7]中的物面投影方法，將物面網(wǎng)格向目標(biāo)物面投影計(jì)算偏移量，再對(duì)近物面網(wǎng)格進(jìn)行坐標(biāo)修正，使其能找到合理的貢獻(xiàn)單元，從而實(shí)現(xiàn)物面附近網(wǎng)格間流動(dòng)變量的正確傳遞。本文首先對(duì)各物體物面建立ADT叉樹(shù)，通過(guò)ADT查找物面網(wǎng)格潛在投影單元，然后再利用Stencil-walking法做進(jìn)一步的判斷。曲面上使用Stencil-walking法的計(jì)算公式為

(1)

(2)

式中：ε通常取值0.01；Δs為網(wǎng)格距離；f為修正協(xié)調(diào)函數(shù)，當(dāng)插值點(diǎn)在物面附近時(shí)取1，距離物面較遠(yuǎn)時(shí)取0，實(shí)現(xiàn)坐標(biāo)修正光滑過(guò)渡。

5 計(jì)算方法

流動(dòng)控制方程為三維可壓縮Navier-Stokes方程。一般曲線(xiàn)坐標(biāo)系下，無(wú)量綱化的方程守恒形式為[29]

(3)

式中：Q為守恒變量；F、G和H為對(duì)流通量；Fv、Gv和Hv為黏性通量；t為時(shí)間；Re為自由來(lái)流雷諾數(shù)。

采用有限體積法對(duì)控制方程進(jìn)行離散求解，空間離散采用Roe的FDS(Flux Difference Splitting)格式、MUSCL(Monotone Upstream-centred Schemes for Conservation Laws)插值方法和Van Albada以及minmod限制器用于獲得二階空間離散精度，湍流模型采用S-A(Spalart-Allmaras)模型，時(shí)間離散方法為L(zhǎng)U-SGS(Lower-Upper Symmetric Gauss-Seidel)方法，非定常時(shí)間推進(jìn)采用雙時(shí)間步方法。

6 算例及計(jì)算結(jié)果

6.1 翼身組合體DLR-F6

為驗(yàn)證該物面重疊網(wǎng)格修正方法的正確性以及靜態(tài)復(fù)雜算例的模擬能力，本文選取德國(guó)宇航公司的DLR-F6翼身組合體(簡(jiǎn)稱(chēng)DLR-F6 WB[30])跨聲速擾流作為驗(yàn)證算例。計(jì)算來(lái)流馬赫數(shù)Ma=0.75，雷諾數(shù)Re=3×106，升力系數(shù)CL=0.500。該計(jì)算模型為半模，總共分為粗、中、細(xì)以及極細(xì)四套網(wǎng)格，網(wǎng)格量分別為300萬(wàn)、500萬(wàn)、750萬(wàn)以及1 000萬(wàn)，圖 12為DLR-F6 WB機(jī)翼、機(jī)身網(wǎng)格重疊示意圖。

圖13為該工況下DLR-F6 WB模型表面及對(duì)稱(chēng)面壓力系數(shù)Cp分布圖，圖中機(jī)翼與機(jī)身網(wǎng)格重疊區(qū)域等壓力系數(shù)線(xiàn)過(guò)渡光滑，分布合理，機(jī)翼表面激波結(jié)構(gòu)清晰，說(shuō)明本重疊網(wǎng)格方法在物面附近的處理是正確合理的。圖 14為該工況下機(jī)翼η=33.1%和η=51.4%站位不同密度的壓力系數(shù)分布對(duì)比。計(jì)算結(jié)果表明，各種不同密度網(wǎng)格計(jì)算得到的壓力分布基本相同，均與實(shí)驗(yàn)值吻合良好，不同的是激波位置與激波強(qiáng)度。

圖 12 DLR-F6 WB網(wǎng)格重疊示意圖Fig.12 Sketch map of overlapping grids for DLR-F6 WB

圖13 DLR-F6 WB模型對(duì)稱(chēng)面及表面壓力系數(shù)云圖Fig.13 Pressure coefficient contours on symmetric plane and surface of DLR-F6 WB model

圖 14 DLR-F6 WB機(jī)翼表面壓力系數(shù)分布Fig.14 Surface pressure coefficient distribution for DLR-F6 WB wing

6.2 機(jī)翼掛載分離

機(jī)翼掛載分離模型是測(cè)試非定常多體相對(duì)運(yùn)動(dòng)數(shù)值模擬能力的一個(gè)標(biāo)準(zhǔn)算例，可用來(lái)驗(yàn)證動(dòng)網(wǎng)格、六自由度方程模型以及非定常數(shù)值算法的正確性和準(zhǔn)確性。該模型由3部分組成，機(jī)翼、掛架以及外掛物，其幾何外形參數(shù)及其他參數(shù)詳見(jiàn)文獻(xiàn)[31]。計(jì)算條件為：Ma=0.95、H=7.92 km、α=0°，H為飛行高度，α為迎角。其中網(wǎng)格由四部分組成，機(jī)翼網(wǎng)格、掛架網(wǎng)格、外掛物網(wǎng)格以及拉伸的直角背景網(wǎng)格，總網(wǎng)格量為565萬(wàn)，其中機(jī)翼網(wǎng)格為400萬(wàn)，掛架網(wǎng)格為15萬(wàn)，掛載物網(wǎng)格為50萬(wàn)，背景網(wǎng)格為100萬(wàn)，其網(wǎng)格重疊情況如圖15所示。機(jī)翼和掛架間物面和近物面網(wǎng)格相互重疊，如圖16所示。

圖17為分離過(guò)程0.1 s、0.2 s和0.3 s表面以及對(duì)稱(chēng)面壓力分布圖，由圖可以看出計(jì)算得到的激波結(jié)構(gòu)清晰，掛架附近物面處等值線(xiàn)過(guò)渡光滑。圖 18為掛載物三方向力系數(shù)Cx、Cy和Cz、力矩系數(shù)Cl、Cm和Cn、角速度ωx、ωy和ωz、角位移、速度u、v和w以及位移的計(jì)算值與實(shí)驗(yàn)結(jié)果對(duì)比圖，圖中CFD代表計(jì)算結(jié)果，Exp.代表實(shí)驗(yàn)結(jié)果，二者吻合良好，驗(yàn)證了本文非定常多體相對(duì)運(yùn)動(dòng)數(shù)值模擬的精確性。由于本算例彈射力作用方式與實(shí)驗(yàn)有差別，且彈射力量值較大，z方向較小的位移也會(huì)造成較大的滾轉(zhuǎn)力矩變化，從而出現(xiàn)圖18(c)中滾轉(zhuǎn)通道角速率與實(shí)驗(yàn)值有一定偏差的現(xiàn)象，這在文獻(xiàn)[32]中也有相應(yīng)的論述。

圖15 機(jī)翼/掛架/掛載物網(wǎng)格示意圖Fig.15 Sketch map of grids of wing/pylon/store

圖16 機(jī)翼/掛架網(wǎng)格重疊示意圖Fig.16 Sketch map of overlapping grids of wing and pylon

圖17 各時(shí)刻表面及對(duì)稱(chēng)面壓力分布Fig.17 Surface and symmetric plane pressure contours at each moment

圖18 計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果對(duì)比Fig.18 Comparison of calculation and experiment results

7 結(jié) 論

為降低CFD計(jì)算過(guò)程中網(wǎng)格生成難度，本文基于結(jié)構(gòu)重疊網(wǎng)格技術(shù)，提出了新的挖洞、尋點(diǎn)和洞面優(yōu)化方法，并實(shí)現(xiàn)了帶有物面網(wǎng)格重疊的靜態(tài)定常、動(dòng)態(tài)非定常數(shù)值模擬，重疊過(guò)程無(wú)需人工干預(yù)，自動(dòng)化程度和魯棒性高。

1) 在最小洞映射的基礎(chǔ)上修改洞映射模型標(biāo)識(shí)方法，確保洞映射模型準(zhǔn)確標(biāo)識(shí)，通過(guò)對(duì)落入洞映射模型過(guò)渡單元的點(diǎn)進(jìn)行再判斷，提出復(fù)合式挖洞方法，在不增加額外內(nèi)存消耗的前提下準(zhǔn)確標(biāo)識(shí)網(wǎng)格點(diǎn)屬性，實(shí)現(xiàn)網(wǎng)格挖洞過(guò)程。

2) 通過(guò)對(duì)格心坐標(biāo)虛網(wǎng)格進(jìn)行延拓，確保格心系統(tǒng)空間的連續(xù)性，實(shí)現(xiàn)格心系統(tǒng)下尋點(diǎn)，提出有效搜索思想，提高了尋點(diǎn)效率。

3) 改變混合洞面優(yōu)化方法填補(bǔ)階段判別法則，引入兩類(lèi)受保護(hù)洞內(nèi)點(diǎn)，確保存在窄縫等嚴(yán)苛條件下兩層插值邊界的建立，同時(shí)將合適的洞內(nèi)點(diǎn)Hole、臨界點(diǎn)Crt轉(zhuǎn)換為插值點(diǎn)Frg 3，以利于非定常多體相對(duì)運(yùn)動(dòng)數(shù)值模擬過(guò)程收斂。

4) 利用物面投影法對(duì)物面網(wǎng)格重疊中的近壁面網(wǎng)格坐標(biāo)進(jìn)行修正，使其能準(zhǔn)確找到貢獻(xiàn)單元，實(shí)現(xiàn)了物面附近流動(dòng)變量的正確傳遞。

[1] THOMAS H P, ANTHONY J S. High-lift OVERFLOW analysis of the DLR-F11 wind tunnel model： AIAA-2014-2697[R]. Reston: AIAA, 2014.

[2] MENDENHALL M R, LESIEUTRE D J, CARUSO S C, et al. Aerodynamic design of pegasustm, concept to flight with CFD： AGARD CP 43[R]. Paris: AGARD, 1980.

[3] LIEVER P A, HABCHI S D. Separation analysis of launch vehicle crew escape systems: AIAA-2004-4726[R]. Reston: AIAA, 2004.

[4] WOODEN P A, BRYAN W, JUBARAJ S B. Calibrating CFD predictions for use in multiple store separation analysis: AIAA-1998-0754[R]. Reston: AIAA, 1998.

[5] KIM J W, PARK S H, YU Y H. Euler and Navier-Stokes simulations of helicopter rotor blade in forward flight using an overlapped grid solver: AIAA-2009-4268[R]. Reston: AIAA, 2009.

[6] BENEK J A, STEGER J L, DOUGHERTY F C. A flexible grid embedding technique with applications to the Euler equations: AIAA-1983-1944[R]. Reston: AIAA, 1983.

[7] ROGERS S E, SUHS N E, DIETZ W E. PEGASUS 5: An automated preprocessor for overset-grid computational fluid dynamics[J]. AIAA Journal, 2003, 41(6): 1037-1045.

[8] WILLIAM M C. The overgrid interface for computational simulations on overset grids: AIAA-2002-3188[R]. Reston: AIAA, 2002.

[9] NOACK R W. SUGGAR: A general capability for moving body overset grid assembly: AIAA-2005-5117[R]. Reston: AIAA, 2005.

[10] CHEN S Y, CHEN Y C, XIA Z H, et al. Constrained large-eddy simulation and detached eddy simulation of flow past a commercial aircraft at 14 degrees angle of attack[J]. Sci China-Phys Mech Astron, 2013, 56(2): 270-276.

[11] ZHAO M, CAO Y H. Numerical simulation of rotor flow field based on overset grids and several spatial and temporal discretization schemes[J]. Chinese Journal of Aeronautics, 2012, 25(2):155-163.

[12] 徐嘉, 劉秋洪, 蔡晉生, 等. 基于隱式嵌套重疊網(wǎng)格技術(shù)的阻力預(yù)測(cè)[J]. 航空學(xué)報(bào), 2013, 34(2): 208-217.

XU J, LIU Q H,CAI J S, et al. Drag prediction based on overset grids with implicit hole cutting technique[J]. Acta Aeronautica et Astronautica Sinica, 2013, 34(2): 208-217 (in Chinese).

[13] 田書(shū)玲, 伍貽兆, 夏健. 基于非結(jié)構(gòu)重疊網(wǎng)格的二維外掛物投放模擬[J]. 空氣動(dòng)力學(xué)學(xué)報(bào), 2007, 25(2): 225-249.

TIAN S L, WU Y Z, XIA J. 2D store separation simulation using unstructured overlapping grid[J]. Acta Aerodynamica Sinica, 2007, 25(2):225-249 (in Chinese).

[14] 田書(shū)玲, 伍貽兆, 夏健. 基于動(dòng)態(tài)非結(jié)構(gòu)重疊網(wǎng)格法的直升機(jī)前飛非定常流場(chǎng)數(shù)值模擬研究[J]. 航空學(xué)報(bào), 2007, 28(5): 1047-1054.

TIAN S L, WU Y Z, XIA J. Numerical simulation research of unsteady flow field around helicopter in forward flight on dynamic overset unstructured grids[J]. Acta Aeronautica et Astronautica Sinica, 2007, 28(5):1047-1054 (in Chinese).

[15] 田書(shū)玲, 伍貽兆, 夏健. 用動(dòng)態(tài)非結(jié)構(gòu)重疊網(wǎng)格法模擬三維多體相對(duì)運(yùn)動(dòng)繞流[J]. 航空學(xué)報(bào), 2007, 28(1): 46-51.

TIAN S L, WU Y Z, XIA J. Simulation of flow past multri-body in relative motion with dynamic unstructured overset grid method[J]. Acta Aeronautica et Astronautica Sinica, 2007, 28(1): 46-51 (in Chinese).

[16] SLOTNICK J, KANDULA M, BUNING P. Navier-Stokes simulation of the space shuttle launch vehicle flight transonic flowfield using a large scale Chimera grid system: AIAA-1994-1860[R]. Reston: AIAA, 1994.

[17] LABOZZETTA W F, GATZKE T D. MACGS-towards the complete grid generation system: AIAA-1994-1923[R]. Reston: AIAA, 1994.

[18] CHIU I, MEAKIN R. On automating domain connectivity for overset grids: AIAA-1995-0854[R]. Reston: AIAA, 1995.

[19] MEAKIN R. Object x-ray for cutting holes in composite overset structured grids: AIAA-2001-2537[R]. Reston: AIAA, 2001.

[20] NOACK R W. Resolution appropriate overset grid assembly for structured and unstructured grids: AIAA-2003-4123[R]. Reston: AIAA, 2003.

[21] NOACK R W. A direct cut approach for overset hole cutting: AIAA-2007-3835[R]. Reston: AIAA, 2007.

[22] LEE Y L, BAEDER J D. Implicit hole cutting-a new approach to overset grid connectivity: AIAA-2003-4128[R]. Reston: AIAA, 2003.

[23] TAYLOR S, WANG J C T. Launch-vehicle simulations using a concurrent, implicit Navier-Stokes solver: AIAA-1995-0223[R]. Reston: AIAA, 1995.

[24] 袁武, 閻超, 席柯. 洞映射方法的研究和改進(jìn)[J]. 北京航空航天大學(xué)學(xué)報(bào), 2012, 38(4): 563-568.

YUAN W, YAN C, XI K. Investigation and enhancement of hole mapping method[J]. Journal of Beijing University of Aeronautics and Astronautics, 2012, 38(4): 563-568 (in Chinese).

[25] STUART E R. Improvements to the Pegasus5 overset CFD software[C]//12th Symposium on Overset Grids and Solution Technology, 2014.

[26] BONET J, PERAIRE J. An alternating digital tree (ADT) algorithm for 3D geometric searching and intersection problems[J]. International Journal for Numerical Methods in Engineering, 1991, 31(1): 1-17.

[27] BELK D M, MAPLE R C. Automated assembly of structured grids for moving body problems: AIAA-1995-1680[R]. Reston: AIAA, 1995.

[28] 王文, 閻超, 袁武, 等. 新型重疊網(wǎng)格洞面優(yōu)化方法及其應(yīng)用[J]. 航空學(xué)報(bào), 2016, 37(3): 826-835.

WANG W, YAN C, YUAN W, et al. Novel overlapping optimization algorithm and its applications[J]. Acta Aeronautica et Astronautica Sinica, 2016, 37(3): 826-835 (in Chinese).

[29] 閻超. 計(jì)算流體力學(xué)方法及應(yīng)用[M]. 北京: 北京航空航天大學(xué)出版社, 2006: 18-25.

YAN C. Computational fluid dynamic’s methods and applications[M]. Beijing: Beihang University Press, 2006: 18-25 (in Chinese).

[30] KELLY R, BRODERSEN O. Summary of data from the second AIAA CFD drag prediction workshop: AIAA-2004-0555[R]. Reston: AIAA, 2004.

[31] HALL L H, PARTHASARATHY V. Validation of an automated chimera/6-DOF methodology for multiple moving body problems: AIAA-1998-0767[R]. Reston: AIAA, 1998.

[32] LEE S, PARK M, CHO K W, et al. A new automated chimera method for the prediction of store trajectory: AIAA-1999-3131[R]. Reston: AIAA, 1999.

王文男， 博士研究生。主要研究方向： 計(jì)算流體力學(xué)， 重疊網(wǎng)格方法。

Tel: 010-82338071

E-mail: wangwenbuaa@126.com

閻超男， 博士， 教授， 博士生導(dǎo)師。主要研究方向： 空氣動(dòng)力學(xué)， 計(jì)算流體力學(xué)。

Tel: 010-82317019

E-mail: yanchao@buaa.edu.cn

URL：www.cnki.net/kcms/detail/11.1929.V.20160222.1027.004.html

Arobustandautomaticstructuredoverlappinggridapproach

WANGWen1,YANChao1,*,YUANWu2,HUANGYu1,XIKe3

1.SchoolofAeronauticScienceandEngineering,BeihangUniversity,Beijing100083,China2.SupercomputingCenter,ComputerNetworkInformationCenter,ChineseAcademyofSciences,Beijing100190,China3.InstituteofOrdnanceIndustryNavigationandControlTechnology,Beijing100089,China

Thegridtechnologyisacrucialaspectofnumericalsimulationuptillthepresentmoment.Overlappinggridtechnologycouldrelaxtherestrictionofstructuretopology,thusreducingthedifficultyofgridgeneration.Basedonthestructuredoverlappinggrid,hole-cutting,donorsearch,andoverlappingoptimizationareexploredandimproved,andsurfacegridassemblageiscompletedinthisthesis.Simultaneously,arobustandautomaticgridoverlappingsystemisconstructed.Inordertosavememory,basedontheminimumholemappingmethod,acompositehole-cuttingmethodisproposed.Regardingdonorsearch,thecontinuityofthesearchspaceispreservedwiththevirtualcell-centeredgrid.Moreover,thenodesofalternatingdigitaltree(ADT)aredecreasedbyanefficientsearchalgorithmwhichexcludessomenon-contributinggridpoints.Withrespecttooverlappingoptimization,twokindsofprotectedhole-pointsareinsertedandanewprincipleisdevelopedtoensuretheconstructionoftwointerpolatedlayersinpastephase,intensifyingtherobustnessofoverlappingoptimization.Thecoordinatesofnear-bodygridsarecorrectedbyprojectionmethodtomakesuretheaccuratetransmissionofflowvariables.Withtheenhancedalgorithm,asteadyDLR-F6wing-bodyflowandanunsteadywing/pylon/storeseparationflowareperformed,andexcellentagreementofcomputationalresultscomparedwithexperimentaldatahasbeenachieved.Thealgorithmshowsremarkablecapabilityandaccuracyinthesimulationofsteadyorunsteadymultiplebodiesflowandprovidesgreatapplicationvalueforengineering.

computationalfluiddynamics;structuredoverlappinggrid;hole-mapping;donorsearch;holeoptimizationmethod;wing-body;wing/pylon/storeseparation

2016-01-03；Revised2016-01-18；Accepted2016-01-29；Publishedonline2016-02-221027

.Tel.：010-82317019E-mailyanchao@buaa.edu.cn

2016-01-03；退修日期2016-01-18；錄用日期2016-01-29； < class="emphasis_bold">網(wǎng)絡(luò)出版時(shí)間

時(shí)間：2016-02-221027

www.cnki.net/kcms/detail/11.1929.V.20160222.1027.004.html

.Tel.：010-82317019E-mailyanchao@buaa.edu.cn

王文， 閻超， 袁武， 等． 魯棒的結(jié)構(gòu)網(wǎng)格自動(dòng)化重疊方法J． 航空學(xué)報(bào),2016,37(10):2980-2991.WANGW,YANC,YUANW,etal.ArobustandautomaticstructuredoverlappinggridapproachJ.ActaAeronauticaetAstronauticaSinica,2016,37(10):2980-2991.

http://hkxb.buaa.edu.cnhkxb@buaa.edu.cn

10.7527/S1000-6893.2016.0034

V211.3

A

1000-6893(2016)10-2980-12