靳晨暉，王剛，陳鑫，周豪

(西北工業(yè)大學(xué) 航空學(xué)院，西安 710072)

0 引 言

嵌套網(wǎng)格技術(shù)在航空航天領(lǐng)域有著廣泛的應(yīng)用，針對(duì)一些復(fù)雜的幾何構(gòu)型，使用嵌套網(wǎng)格技術(shù)可以分區(qū)生成網(wǎng)格，減少網(wǎng)格生成的難度。對(duì)于飛機(jī)副油箱分離[1]，座艙蓋分離/座椅彈射[2]，多級(jí)火箭助推器分離[3]等多體分離問題,不同物面間的相對(duì)位置會(huì)發(fā)生變化，同時(shí)伴有強(qiáng)烈的非定常效應(yīng)。常規(guī)研究方法如風(fēng)洞和飛行實(shí)驗(yàn)常帶有一定的局限性，使用嵌套網(wǎng)格技術(shù)進(jìn)行非定常計(jì)算可以對(duì)這些過程進(jìn)行精確數(shù)值模擬。

在嵌套網(wǎng)格算法中，一般通過插值實(shí)現(xiàn)嵌套區(qū)域的流場信息傳遞，尤其是在解決帶有動(dòng)邊界的問題，嵌套區(qū)域的網(wǎng)格在每個(gè)非定常步內(nèi)都需要更新，插值模板也需要重新選取。因此，使用正確的插值方法是保證計(jì)算結(jié)果準(zhǔn)確的關(guān)鍵。按照插值方法的不同可分為守恒型插值與非守恒型插值。對(duì)于守恒型插值而言，網(wǎng)格尺度對(duì)其插值精度的影響較小，但插值過程十分復(fù)雜；而非守恒型插值的實(shí)現(xiàn)則較為簡單，但需要控制嵌套區(qū)域網(wǎng)格的匹配程度，否則會(huì)損失計(jì)算精度。為了減少非守恒型插值過程中的精度損失，一般需要從兩方面入手：一方面要改善嵌套插值區(qū)域網(wǎng)格匹配程度，如在網(wǎng)格生成時(shí)，需要對(duì)背景網(wǎng)格中有相對(duì)運(yùn)動(dòng)的區(qū)域進(jìn)行局部加密，同時(shí)使用高效的嵌套裝配算法[4]劃分嵌套區(qū)域，以此來保證物面運(yùn)動(dòng)的整個(gè)過程中插值單元與貢獻(xiàn)單元的匹配程度。另一方面需要研究不同插值方法的插值效果，C.Mastin等[5]研究了雙線性、三線性插值方法在嵌套區(qū)域的插值效果；田書玲等[6]研究了一種基于Lagrange插值的非結(jié)構(gòu)嵌套網(wǎng)格插值方法；周乃春等[7]使用了逆向距離權(quán)方法對(duì)嵌套區(qū)域進(jìn)行插值；黃宇等[8]使用了具有二階精度的插值方法，有效降低了由于在嵌套區(qū)域進(jìn)行插值而引入的誤差。僅僅增加嵌套區(qū)域的網(wǎng)格量有時(shí)也難以達(dá)到較好的精度，網(wǎng)格量增加會(huì)導(dǎo)致仿真周期成倍增加，占用過多計(jì)算資源；同時(shí)，線性插值方法在處理復(fù)雜問題時(shí)的精度相對(duì)較低，選擇合適的插值方法對(duì)于處理不同的問題至關(guān)重要。

徑向基函數(shù)插值方法因其存儲(chǔ)十分方便和數(shù)據(jù)結(jié)構(gòu)相對(duì)簡單的特點(diǎn)，在近二十年內(nèi)迅猛發(fā)展。研究者們?cè)诓煌念I(lǐng)域?qū)较蚧瘮?shù)的應(yīng)用進(jìn)行大量研究，如R.Franke[9]發(fā)現(xiàn)對(duì)散亂的數(shù)據(jù)進(jìn)行擬合時(shí)，使用徑向基函數(shù)進(jìn)行插值的精度要高于其他插值方法；林言中等[10]將徑向基函數(shù)插值方法應(yīng)用于網(wǎng)格變形技術(shù)中的界面數(shù)據(jù)傳遞，對(duì)AGARD 445.6機(jī)翼顫振問題進(jìn)行非定常計(jì)算；劉溢浪等[11]提出一種增量形式的RBF插值方法，并將其應(yīng)用于有限體積方法的流場重構(gòu)中；Wang G等[12-13]在RBF網(wǎng)格變形技術(shù)的基礎(chǔ)上，提出將拉普拉斯光順與網(wǎng)格變形相結(jié)合，有效改善了變形后的網(wǎng)格質(zhì)量。但是，徑向基函數(shù)插值半徑的選取帶有一定的經(jīng)驗(yàn)性，傳統(tǒng)徑向基函數(shù)前述插值半徑選取不合適會(huì)嚴(yán)重影響插值的精度與效率。

為了解決嵌套區(qū)域的插值精度問題，本文對(duì)不同插值函數(shù)進(jìn)行研究，同時(shí)采用一種改進(jìn)的非結(jié)構(gòu)嵌套網(wǎng)格徑向基函數(shù)(RBF)插值方法。在嵌套網(wǎng)格流場信息傳遞的過程中，除了使用線性插值(LINE)和以距離倒數(shù)作為權(quán)重插值(WA)兩種傳統(tǒng)插值方法，還對(duì)逆二次徑向基函數(shù)(IQ)和Wendland’s C2徑向基函數(shù)(C2)的插值效果進(jìn)行研究。針對(duì)嵌套邊界區(qū)域出現(xiàn)的數(shù)值奇性現(xiàn)象，通過控制徑向基函數(shù)的作用半徑來調(diào)整插值矩陣條件數(shù)的方式進(jìn)行改善。

1 非結(jié)構(gòu)嵌套網(wǎng)格求解技術(shù)

1.1 控制方程

為便于對(duì)固壁邊界的六自由度運(yùn)動(dòng)進(jìn)行描述，采用ALE(Arbitrary Lagrangeian-Eulerian)方法對(duì)流動(dòng)控制方程進(jìn)行描述。ALE方法允許計(jì)算網(wǎng)格的進(jìn)行剛體運(yùn)動(dòng)和變形，通過在流動(dòng)控制方程中加入網(wǎng)格運(yùn)動(dòng)速度，將流體力學(xué)中的Lagrange方法和Euler方法進(jìn)行統(tǒng)一描述。使用ALE方法描述的三維非定常雷諾平均Navier-Stokes方程(URANS)的積分形式為

(1)

式中：Ω為控制體；?Ω為控制體單元邊界；Q為守恒變量，Q=[ρρuρvρwρE]T，其中ρ為流體密度,u、v、w分別為旋轉(zhuǎn)彈在體軸系下三個(gè)軸向的運(yùn)動(dòng)速度E為總內(nèi)能;n為控制體單元邊界上沿外法線方向的分量；F(Q)和G(Q)分別為N-S方程中的無粘通量項(xiàng)和有粘通量項(xiàng)。

對(duì)于式(1)，使用格心有限體積法(FVM)進(jìn)行空間離散后，得到半離散形式為

(2)

(3)

式中：R為計(jì)算殘差；S為網(wǎng)格單元的面積；n為時(shí)間步；Vgrid為網(wǎng)格進(jìn)行剛體運(yùn)動(dòng)的速度。為了保證足夠的時(shí)間精度，使用二階向后的歐拉隱式方法[14]求解式(2)，如式(4)所示

(4)

式中：Δt為時(shí)間步長；上標(biāo)n為真實(shí)時(shí)間迭代的步數(shù)。直接求解式(4)較為復(fù)雜，一般采用偽時(shí)間迭代的方法對(duì)其進(jìn)行求解，為此在方程的左端加入一個(gè)守恒變量對(duì)虛擬時(shí)間τ的導(dǎo)數(shù)，從而方程的解可以等價(jià)為求解一個(gè)關(guān)于虛擬時(shí)間τ的一階常微分方程組在虛擬時(shí)間τ趨近無窮時(shí)的漸進(jìn)解：

(5)

(6)

由于式(5)的精度與偽時(shí)間迭代方法無關(guān)，因此，偽時(shí)間迭代過程可以使用一些加速收斂方法來提高流場的計(jì)算效率。具體方法可以參照文獻(xiàn)[15-16]。采用改進(jìn)的LU-SGS[17]算法對(duì)式(5)進(jìn)行向后的歐拉隱式時(shí)間迭代，同時(shí)使用基于OpenMP并行算法來提高計(jì)算的效率。

1.2 嵌套區(qū)域劃分與貢獻(xiàn)單元查找

嵌套網(wǎng)格裝配技術(shù)的兩個(gè)核心問題是如何確定嵌套邊界和如何查找貢獻(xiàn)單元。確定嵌套邊界也叫 “挖洞”，當(dāng)計(jì)算網(wǎng)格量較大時(shí)，會(huì)顯著增加預(yù)處理的時(shí)間，因此高效并且合理地劃分出嵌套區(qū)域?qū)⑹乔短拙W(wǎng)格裝配技術(shù)的關(guān)鍵。本文在嵌套區(qū)域劃分時(shí)采取的策略是根據(jù)壁面距離來確定嵌套邊界區(qū)域，同時(shí)計(jì)算子網(wǎng)格到自身物面和其他物面的最短距離，通過比較到不同物面的最短距離來確定嵌套邊界，相比于傳統(tǒng)的枚舉法，有效減少了前處理時(shí)間。節(jié)點(diǎn)屬性判定如圖1所示，i點(diǎn)是屬于A網(wǎng)格的網(wǎng)格節(jié)點(diǎn)，在求解距離時(shí)，距離A網(wǎng)格更近，則i點(diǎn)的屬性標(biāo)記為激活的狀態(tài)；而j點(diǎn)距離B網(wǎng)格更近，則j點(diǎn)的屬性標(biāo)記為非激活的狀態(tài)。

圖1 節(jié)點(diǎn)屬性判定示意圖Fig.1 Schematic diagram of node attribute determination

當(dāng)確定了嵌套邊界后，需要在嵌套區(qū)域查找貢獻(xiàn)單元。對(duì)于復(fù)雜模型來講，網(wǎng)格量較大，查找貢獻(xiàn)單元往往要花費(fèi)較多時(shí)間，因此如何快速的查找到貢獻(xiàn)單元是提高嵌套效率的關(guān)鍵。本文根據(jù)文獻(xiàn)[4]采用了一種基于非結(jié)構(gòu)嵌套網(wǎng)格的“蛇形查找方法”(如圖2所示)，在空間中以近乎最佳路徑的方式查找到了貢獻(xiàn)單元，有效縮短了查找貢獻(xiàn)單元的時(shí)間與距離。假設(shè)要尋找i點(diǎn)的貢獻(xiàn)單元，該方法的具體實(shí)現(xiàn)過程如下：

(1) 在空間中賦予i點(diǎn)一個(gè)初始的出發(fā)單元，這里假設(shè)是1號(hào)網(wǎng)格單元；

(2) 連接當(dāng)前的出發(fā)單元的網(wǎng)格格心和i點(diǎn)，判斷其穿過當(dāng)前出發(fā)單元的哪一條邊(三維為面)；

(3) 若穿過某一邊(面)，則將該邊(面)的另一側(cè)單元當(dāng)做出發(fā)單元，循環(huán)第二步；

(4) 若沒有穿過當(dāng)前的出發(fā)單元所有的邊(面)，則當(dāng)前的出發(fā)單元就是i點(diǎn)的貢獻(xiàn)單元，結(jié)束查找。

圖2 蛇形查找示意圖Fig.2 Schematic diagram of snake search

通過以上步驟，嵌套區(qū)域的網(wǎng)格裝配完畢。在進(jìn)行流場計(jì)算時(shí)，網(wǎng)格實(shí)際被劃分為三個(gè)區(qū)域，分別是計(jì)算區(qū)域、關(guān)閉區(qū)域和嵌套插值區(qū)域。在計(jì)算區(qū)域，進(jìn)行正常的CFD求解，當(dāng)接觸到邊界時(shí)，面外側(cè)位于嵌套區(qū)域的網(wǎng)格信息由貢獻(xiàn)單元插值得到；在關(guān)閉區(qū)域，只計(jì)算網(wǎng)格單元體積和網(wǎng)格單元面積等預(yù)處理操作；在嵌套插值區(qū)域，在每個(gè)偽時(shí)間迭代結(jié)束后進(jìn)行插值，嵌套插值區(qū)域的網(wǎng)格不進(jìn)行計(jì)算，而是被當(dāng)作嵌套邊界來處理，這類邊界具有插值后獲得的流場信息。

1.3 插值方法簡介

流場信息的插值傳遞方法較多，本文主要考察了線性插值、距離倒數(shù)權(quán)重插值兩種傳統(tǒng)插值方法，以及逆二次徑向基函數(shù)、Wendland’s C2徑向基函數(shù)插值方法。由于線性插值結(jié)構(gòu)較為簡單，此處不再贅述，下面重點(diǎn)對(duì)其余三種插值方法進(jìn)行簡要介紹。

以距離倒數(shù)作為權(quán)重來進(jìn)行插值的原理是：對(duì)待插值點(diǎn)附近點(diǎn)進(jìn)行加權(quán)平均，其中權(quán)重的大小與需要插值的點(diǎn)與相鄰點(diǎn)之間的距離密切相關(guān)，與距離呈反比，因此也叫做距離反比插值方法。其具體形式為

(7)

本文所使用的徑向基函數(shù)的基本形式為

(8)

式中：下標(biāo)i為支撐點(diǎn)；r為該點(diǎn)的位置矢量；N為網(wǎng)格節(jié)點(diǎn)的總數(shù)；φ(‖rp-ri‖)為所采用的基函數(shù)的形式；ωi為與第i個(gè)支撐點(diǎn)相關(guān)的權(quán)重系數(shù)。徑向基函數(shù)大致可分為以下三類：全局型(Global)徑向基函數(shù)、局部型(Local)徑向基函數(shù)、緊支型(Compact)徑向基函數(shù)，這里簡單介紹一下本文所使用的全局型/局部型(IQ)與緊支型(C2)兩種徑向基函數(shù)。

(1) 逆二次徑向基函數(shù)方法

(9)

在對(duì)嵌套區(qū)域的P點(diǎn)進(jìn)行流場信息傳遞時(shí)，根據(jù)之前查找得到的貢獻(xiàn)單元，將貢獻(xiàn)單元以及它的相鄰單元一起作為徑向基函數(shù)的插值模版構(gòu)建徑向基函數(shù)的系數(shù)矩陣。插值模版示意圖如圖3所示，i網(wǎng)格是P的貢獻(xiàn)單元，則使用i、j、k、l四個(gè)網(wǎng)格構(gòu)建徑向基函數(shù)系數(shù)矩陣。

圖3 插值模版示意圖Fig.3 Interpolation template schematic diagram

為了選取比較光滑的模版值，在選取模版時(shí)，每個(gè)物理量減去他們的平均值，最后再加上被減去的值。將各個(gè)點(diǎn)的流場信息做差后的ΔZ帶入到RBF插值函數(shù)之中，得到的線性方程組如下：

ΦW=ΔZ

(10)

(11)

W=[ωi,ωj,ωk,ωl]T

(12)

式中：‖ri-rj‖為歐式距離，即每個(gè)基網(wǎng)格格心之間的距離，求解方程后得到每個(gè)基網(wǎng)格的權(quán)重系數(shù)ωi，將需要插值的節(jié)點(diǎn)P坐標(biāo)rP帶入到徑向基函數(shù)中就可以重構(gòu)出該點(diǎn)的值。

(2) Wendland’s C2型徑向基函數(shù)

基本形式如下所示：

(13)

插值半徑的選取對(duì)于徑向基函數(shù)插值方法的效率有著很重大的影響。較大的插值半徑會(huì)擴(kuò)大物面插值的影響域，與此同時(shí)，也會(huì)使得基矩陣的維度增加并引起計(jì)算量的增加。選擇較小的插值半徑可以減小矩陣維度，降低計(jì)算量。

1.4 改進(jìn)的徑向基函數(shù)插值策略

徑向基函數(shù)插值精度高，靈活性很強(qiáng)，使用不同類型的插值函數(shù)的插值效果也有所差異。對(duì)于同一個(gè)徑向基函數(shù)，選取不同的參數(shù)d也會(huì)得到不同的插值效果。本文所使用的C2徑向基函數(shù)，作用半徑d的初始值為插值模版中格心和P點(diǎn)(嵌套區(qū)域需要進(jìn)行插值的點(diǎn))的最大距離的50倍。矩陣的條件數(shù)不僅是判斷矩陣病態(tài)與否的一種度量，也可以表示矩陣計(jì)算對(duì)誤差的敏感度強(qiáng)弱，條件數(shù)越大，數(shù)值穩(wěn)定性越差。對(duì)于C2徑向基函數(shù)的系數(shù)矩陣Φ，取二范數(shù)后的條件數(shù)可以表示為：

(14)

為了改善在一些網(wǎng)格上的條件數(shù)過大的情況，本文采用了一種基于控制條件數(shù)的徑向基函數(shù)插值算法，通過人為的改變C2徑向基函數(shù)的作用半徑d來調(diào)整插值矩陣的條件數(shù)，從而在條件數(shù)過大的插值模版上進(jìn)行有效的控制，有效改善了由于條件數(shù)過大使得穩(wěn)定性較差以及出現(xiàn)數(shù)值奇性的現(xiàn)象。具體實(shí)施的步驟為：當(dāng)流場求解開始時(shí)，人為輸入一個(gè)期望達(dá)到的條件數(shù)，當(dāng)使用徑向基函數(shù)插值進(jìn)行流場信息傳遞時(shí)，每次迭代對(duì)貢獻(xiàn)單元插值矩陣的插值半徑自動(dòng)乘以0.9，計(jì)算所有貢獻(xiàn)單元插值矩陣的條件數(shù)，不滿足要求則返回上一步繼續(xù)縮小插值矩陣半徑，直到所有貢獻(xiàn)單元插值矩陣的條件數(shù)都滿足要求。

2 計(jì)算結(jié)果與分析

2.1 插值精度分析

為了分析本文所使用的插值函數(shù)進(jìn)行解析的精度，使用解析函數(shù)作為插值精度驗(yàn)證算例。該函數(shù)的數(shù)學(xué)形式為：

(15)

嵌套網(wǎng)格組裝過程示意圖如圖4所示，圓心在(0,0)處,半徑為1。

(a) 嵌套組裝前網(wǎng)格

(b) 嵌套完成后網(wǎng)格圖4 嵌套網(wǎng)格組裝過程示意圖Fig.4 Diagram of overset grid assembly process

使用線性插值，以距離倒數(shù)作為權(quán)重插值，逆二次徑向基函數(shù)插值、Wendland’s C2徑向基函數(shù)插值四種不同的方法進(jìn)行插值，四種插值方法的誤差統(tǒng)計(jì)如表1所示。

表1 四種插值方式在嵌套區(qū)域的誤差統(tǒng)計(jì)Table 1 Error statistics of four interpolation modes in nest area

從表1可以看出：使用C2徑向基函數(shù)進(jìn)行插值具有較小的誤差,得到的結(jié)果和解析解最為接近。在整個(gè)插值的區(qū)域內(nèi)，C2徑向基函數(shù)插值的平均絕對(duì)誤差僅為0.007 254 9%，而IQ徑向基函數(shù)插值為0.036 369%，距離反比函數(shù)插值(WA)為0.197 52%，線性插值(LINE)為2.361 2%，說明使用C2徑向基函數(shù)進(jìn)行插值具有更高的計(jì)算精度，本文將優(yōu)先選擇使用C2徑向基函數(shù)插值進(jìn)行流場信息傳遞。

2.2 MD30P/30N多段翼流場數(shù)值模擬

當(dāng)直接使用C2徑向基函數(shù)用于嵌套網(wǎng)格插值過程中時(shí)，在某些點(diǎn)可能存在數(shù)值原因引起的振蕩或插值誤差過大的點(diǎn)。以MD30P/30N多段翼為例，計(jì)算馬赫數(shù)Ma=0.2，雷諾數(shù)Re=9.0×106，迎角α=8.010 9°，湍流模型選擇Spalart-Allmaras(SA)。圍繞前緣縫翼、主翼、后緣襟翼分別生成計(jì)算網(wǎng)格，其中前緣縫翼網(wǎng)格量約為6萬，后緣襟翼網(wǎng)格量約為11萬，主翼網(wǎng)格量約為22萬。三套網(wǎng)格按照之前距離判定的準(zhǔn)則自動(dòng)進(jìn)行網(wǎng)格嵌套組裝。三套網(wǎng)格均按照預(yù)期完成組裝示意圖如圖5所示。

圖5 MD30P/30N嵌套網(wǎng)格組裝示意圖Fig.5 Assembly schematic of MD30P/30N overset grid

未采用控制條件數(shù)與采用控制條件數(shù)后的C2徑向基函數(shù)插值得到的流場圖如圖6所示。從流場計(jì)算結(jié)果來看，直接使用C2函數(shù)進(jìn)行插值得到了圖6(a)所示的結(jié)果，可以看出：在嵌套區(qū)域存在有數(shù)值奇性。根據(jù)之前提出的條件數(shù)控制策略，對(duì)于不同的插值模板，通過調(diào)整他們的插值半徑來控制條件數(shù)。經(jīng)過條件數(shù)控制為10，即cond=10時(shí)，得到了圖6(b)所示的計(jì)算結(jié)果，整個(gè)流場的求解更加光滑，消除了之前出現(xiàn)的數(shù)值奇性現(xiàn)象。

(a) 未采用控制條件數(shù)

(b) 采用控制條件數(shù)圖6 MD30P/30N流場計(jì)算壓力分布云圖Fig.6 Pressure distribution cloud diagram of MD30P/30N flow field

MD30P/30N表面Cp分布與實(shí)驗(yàn)對(duì)比示意圖如圖7所示，可以看出，當(dāng)直接使用C2徑向基函數(shù)插值進(jìn)行流場信息傳遞時(shí)，表面Cp分布在主翼上表面與后緣襟翼下表面與實(shí)驗(yàn)數(shù)據(jù)差別較大，當(dāng)采取條件數(shù)控制策略的C2徑向基函數(shù)插值策略將矩陣條件數(shù)控制在10后，所得到的翼型表面Cp分布與實(shí)驗(yàn)數(shù)據(jù)吻合度較高，說明該插值方法能夠很好地解決直接使用C2徑向基函數(shù)插值帶來的數(shù)值奇性問題，且計(jì)算精度較高。MD30P/30N殘差收斂歷程對(duì)比示意圖如圖8所示，當(dāng)條件數(shù)控制在10時(shí)，殘差下降的速率較快，說明采取條件數(shù)控制策略的C2徑向基函數(shù)插值策略能有效改善殘差的收斂效果。

圖7 MD30P/30N表面Cp分布示意圖Fig.7 Schematic diagram of MD30P/30N surface Cp distribution

圖8 MD30P/30N殘差收斂歷程對(duì)比示意圖Fig.8 Comparison diagram of MD30P/30N residual convergence process

2.3 AEDC外掛物投放標(biāo)模驗(yàn)證

AEDC標(biāo)準(zhǔn)外掛物分離模型的計(jì)算參數(shù)[18-19]如圖9所示。計(jì)算馬赫數(shù)Ma=0.95，雷諾數(shù)Re=7.874×106，迎角α=0°，湍流模型選擇Spalart-Allmaras(SA)。圍繞機(jī)翼和外掛物分別生成計(jì)算網(wǎng)格，網(wǎng)格總數(shù)約450萬，兩套網(wǎng)格按照之前的距離判定的準(zhǔn)則自動(dòng)進(jìn)行網(wǎng)格嵌套組裝。

圖9 標(biāo)準(zhǔn)外掛物分離模型的計(jì)算參數(shù)Fig.9 The calculation parameters of standard separation model

為了考察使用控制條件數(shù)策略的C2徑向基函數(shù)的插值效果，將插值矩陣的條件數(shù)控制在100，即cond=100時(shí)，如圖10所示。從流場計(jì)算結(jié)果來看，可以看出：未采用控制條件數(shù)時(shí)，在嵌套區(qū)域的邊界上，出現(xiàn)了數(shù)值奇性，這些異常導(dǎo)致在尾跡區(qū)出現(xiàn)了和周圍顯著不同的區(qū)域；當(dāng)采用控制條件數(shù)后，尾跡區(qū)流場過渡的更加光滑和連續(xù)，這說明采用條件數(shù)控制策略的插值方法可以有效消除嵌套邊界區(qū)域出現(xiàn)的數(shù)值奇性的問題。

(a) 無條件數(shù)控制的z-x截面

(b) 采用條件數(shù)控制的z-x截面

(c) 無條件數(shù)控制的x-y截面

(d) 采用條件數(shù)控制的x-y截面圖10 外掛物截面流場對(duì)比示意圖Fig.10 Flow field comparison diagram of the cross section

外掛物表面Cp分布與實(shí)驗(yàn)值對(duì)比示意圖如圖11所示，可以看出：直接使用C2徑向基函數(shù)和未使用控制條件數(shù)策略的C2徑向基函數(shù)進(jìn)行插值得到的外掛物周向Cp分布均與實(shí)驗(yàn)值較為吻合，當(dāng)加入控制插值矩陣條件數(shù)后，外掛物在彈尾處的Cp分布與實(shí)驗(yàn)值吻合度更高一些。

(a) 外掛物兩側(cè)表面Cp分布

(b) 外掛物上下表面Cp分布圖11 外掛物表面Cp分布與實(shí)驗(yàn)對(duì)比示意圖Fig.11 Schematic diagram of Cp distribution on the surface of the plug-in and comparison with the experiment

AEDC外掛物殘差收斂歷程對(duì)比曲線如圖12所示，可以看出：未采用控制條件數(shù)時(shí)，殘差收斂的歷程波動(dòng)較為劇烈；采用控制插值矩陣條件數(shù)后，殘差收斂效果較好，呈現(xiàn)不斷下降的趨勢(shì)，說明采用條件數(shù)控制策略的C2徑向基函數(shù)插值策略能提高殘差收斂的效果。

圖12 殘差收斂歷程對(duì)比Fig.12 Comparison of residual convergence process

非定常數(shù)值模擬得到的外掛物質(zhì)心處的位移Xg、Yg、Zg與實(shí)驗(yàn)值的對(duì)比結(jié)果如圖13所示，可以看出：采用控制條件數(shù)策略的C2徑向基函數(shù)插方法數(shù)值模擬得到的水平位移Xg整體與實(shí)驗(yàn)吻合度更高一些。兩種不同的插值方式得到外掛物沿豎直方向的位移與實(shí)驗(yàn)結(jié)果吻合較好。

圖13 質(zhì)心位移與實(shí)驗(yàn)值的對(duì)比圖Fig.13 Comparison diagram of center of mass displacement and experimental value

外掛物分離過程中的姿態(tài)角變化曲線如圖14所示。

圖14 姿態(tài)角與實(shí)驗(yàn)值的對(duì)比圖Fig.14 Comparison diagram of attitude angle and experimental value

從圖14可以看出：使用控制條件數(shù)策略的C2徑向基函數(shù)進(jìn)行插值得到的俯仰角θ變化曲線在分離中后期與實(shí)驗(yàn)吻合更好一些，而偏航角ψ與滾轉(zhuǎn)角φ使用兩種插值方法數(shù)值模擬得到的結(jié)果較為接近，均與實(shí)驗(yàn)結(jié)果吻合度較高。同時(shí)可以看出彈射力矩與氣動(dòng)力的共同作用，使得外掛物的俯仰角θ會(huì)呈現(xiàn)先增加后減小的趨勢(shì)；外掛物受到三角翼下洗與側(cè)洗的影響，偏航角也會(huì)逐漸增加。

3 結(jié) 論

(1) 使用解析函數(shù)作為插值模板，對(duì)線性插值、距離倒數(shù)作為權(quán)重插值，逆二次徑向基函數(shù)插值、Wendland’s C2徑向基函數(shù)插值四種不同的嵌套網(wǎng)格插值方法的精度進(jìn)行驗(yàn)證。計(jì)算結(jié)果表明C2徑向基函數(shù)的插值精度較高，同時(shí)收斂較快。

(2) 改進(jìn)后的C2徑向基函數(shù)插值方法能夠有效消除嵌套插值區(qū)域數(shù)值奇性，同時(shí)計(jì)算收斂的速率更快，魯棒性較好。