趙璇 常思江 倪旖

摘 要：為提高炮彈氣動(dòng)力計(jì)算的精度， 常需要通過(guò)計(jì)算流體力學(xué)數(shù)值仿真反復(fù)迭代計(jì)算炮彈的氣動(dòng)參數(shù)， 優(yōu)化效率較低。 本文提出一種基于多可信度代理模型的氣動(dòng)外形快速優(yōu)化方法， 在保證計(jì)算精度的條件下大幅提高計(jì)算效率。 以某滑翔制導(dǎo)炮彈的鴨舵外形優(yōu)化為例， 采用兩種可信度樣本訓(xùn)練出多可信度代理模型代替耗時(shí)的計(jì)算流體力學(xué)仿真獲得氣動(dòng)參數(shù)， 依據(jù)滑翔升阻比最大、 穩(wěn)定性與機(jī)動(dòng)性相匹配等設(shè)計(jì)要求， 利用遺傳算法搜尋鴨舵的最優(yōu)外形參數(shù)。 與初始方案相比， 優(yōu)化方案在升阻比方面提升顯著。 通過(guò)與數(shù)值模擬結(jié)果對(duì)比， 該方法對(duì)升阻比的預(yù)測(cè)平均誤差為1.94%， 精度較高， 同時(shí)計(jì)算量大大降低， 驗(yàn)證了其可行性和有效性。

關(guān)鍵詞： 滑翔制導(dǎo)炮彈; 氣動(dòng)外形; Co-Kriging代理模型; 計(jì)算流體力學(xué); 遺傳算法; 優(yōu)化設(shè)計(jì)

中圖分類號(hào)： TJ760; TJ413+.6? 文獻(xiàn)標(biāo)識(shí)碼：??? A 文章編號(hào)： 1673-5048（2021）05-0099-07

0 引? 言

滑翔制導(dǎo)炮彈是通過(guò)在常規(guī)尾翼彈的基礎(chǔ)上增加制導(dǎo)控制系統(tǒng)而形成的遠(yuǎn)程精確打擊武器[1]。 在彈道升弧段， 鴨舵折插在彈體內(nèi)， 經(jīng)過(guò)彈道頂點(diǎn)某位置后， 鴨舵張開(kāi)進(jìn)行滑翔控制， 從而保證炮彈在正攻角下始終以最大升阻比狀態(tài)飛行。 可見(jiàn)， 鴨舵是滑翔制導(dǎo)炮彈在滑翔過(guò)程中極為重要的部件， 其氣動(dòng)特性影響了全彈的升阻比、 滑翔效率及穩(wěn)定性等。 因此， 有必要對(duì)滑翔制導(dǎo)炮彈的鴨舵外形開(kāi)展優(yōu)化設(shè)計(jì)， 以期提高滑翔制導(dǎo)炮彈的性能。

目前， 國(guó)內(nèi)外關(guān)于制導(dǎo)炮彈氣動(dòng)外形優(yōu)化開(kāi)展了一些工作。 文獻(xiàn)[2]利用DATCOM軟件和粒子群算法來(lái)優(yōu)化遠(yuǎn)程制導(dǎo)炮彈控制舵的外形參數(shù)， 并通過(guò)Cart3D軟件分析評(píng)估了其飛行性能; 文獻(xiàn)[3]提出了一種基于試驗(yàn)設(shè)計(jì)和響應(yīng)面方法的制導(dǎo)炮彈外形設(shè)計(jì)優(yōu)化框架， 得到了一組射程最遠(yuǎn)的帕累托最優(yōu)設(shè)計(jì)集; 文獻(xiàn)[4]研究了一種多學(xué)科優(yōu)化設(shè)計(jì)方法， 用于設(shè)計(jì)并評(píng)估無(wú)尾式遠(yuǎn)程制導(dǎo)炮彈的性能。 國(guó)內(nèi)方面， 文獻(xiàn)[5]對(duì)某滑翔增程制導(dǎo)炮彈進(jìn)行了氣動(dòng)數(shù)值計(jì)算， 并通過(guò)風(fēng)洞實(shí)驗(yàn)進(jìn)行驗(yàn)證; 文獻(xiàn)[6]建立了優(yōu)化設(shè)計(jì)模型， 并采用復(fù)形調(diào)優(yōu)法對(duì)制導(dǎo)炮彈的鴨舵外形進(jìn)行優(yōu)化; 文獻(xiàn)[7]針對(duì)常規(guī)氣動(dòng)布局的制導(dǎo)炮彈提出了適用于概念設(shè)計(jì)階段的氣動(dòng)外形優(yōu)化方法， 采用參數(shù)化建模方法和非線性面元法快速準(zhǔn)確地預(yù)估炮彈氣動(dòng)特性。

為提高優(yōu)化效率， 國(guó)內(nèi)外目前對(duì)滑翔制導(dǎo)炮彈的外形優(yōu)化設(shè)計(jì)常采用工程算法計(jì)算氣動(dòng)參數(shù)， 但獲得氣動(dòng)參數(shù)的方法還包括數(shù)值模擬、 風(fēng)洞試驗(yàn)和靶道自由飛行試驗(yàn)等。 在不具備開(kāi)展風(fēng)洞試驗(yàn)或自由飛行試驗(yàn)條件的初步設(shè)計(jì)階段， 工程算法和計(jì)算流體力學(xué)（Computational Fluid Dynamics，? CFD）數(shù)值模擬是重要的氣動(dòng)參數(shù)預(yù)測(cè)手段。 一般而言， CFD數(shù)值模擬的計(jì)算精度高于工程算法， 如文獻(xiàn)[8]所述， 采用工程算法數(shù)據(jù)計(jì)算某型炮彈的軌跡時(shí)， 射程誤差可達(dá)10%， 而采用CFD數(shù)據(jù)所得的誤差為2%左右。 但是， CFD數(shù)值計(jì)算耗時(shí)很長(zhǎng)， 難以實(shí)現(xiàn)快速優(yōu)化。 為此， 研究者們開(kāi)發(fā)了一些代理模型來(lái)代替耗時(shí)的CFD仿真。 其本質(zhì)是采用機(jī)器學(xué)習(xí)技術(shù)， 從設(shè)計(jì)空間內(nèi)的少量樣本信息中學(xué)習(xí)目標(biāo)函數(shù)隨設(shè)計(jì)變量變化的規(guī)律。 目前比較常用的代理模型有響應(yīng)面模型、 徑向基函數(shù)模型、 人工神經(jīng)網(wǎng)絡(luò)模型和Kriging模型等[9]。 文獻(xiàn)[10]采用人工神經(jīng)網(wǎng)絡(luò)模型實(shí)現(xiàn)了對(duì)導(dǎo)彈的氣動(dòng)參數(shù)預(yù)測(cè); 文獻(xiàn)[11]則發(fā)展了一套基于Kriging模型的機(jī)翼氣動(dòng)優(yōu)化設(shè)計(jì)方法。 構(gòu)造一個(gè)精確的氣動(dòng)力代理模型需要大量的樣本點(diǎn)， 若每個(gè)樣本點(diǎn)都調(diào)用CFD數(shù)值模擬， 其計(jì)算量同樣不容小覷。 考慮引入計(jì)算快速的工程算法數(shù)據(jù)輔助建立代理模型， 將CFD數(shù)值模擬和工程算法結(jié)合起來(lái)， 用較少的CFD數(shù)據(jù)實(shí)現(xiàn)高精度的模型， 從而提高計(jì)算效率。

本文以滑翔制導(dǎo)炮彈的鴨舵為對(duì)象開(kāi)展外形參數(shù)優(yōu)化， 充分考慮工程算法的快速性和CFD數(shù)值模擬的高精度性， 構(gòu)建代理模型近似代替對(duì)不同外形方案滑翔制導(dǎo)炮彈的氣動(dòng)力計(jì)算， 并以此為基礎(chǔ)建立了鴨舵外形優(yōu)化的數(shù)學(xué)模型， 通過(guò)遺傳算法求解該模型得到最優(yōu)外形方案， 驗(yàn)證所提方法的可行性和有效性。

1 氣動(dòng)力快速計(jì)算方法

對(duì)于復(fù)雜的多設(shè)計(jì)變量氣動(dòng)外形優(yōu)化問(wèn)題， 氣動(dòng)力計(jì)算結(jié)果的精準(zhǔn)度影響著炮彈外形優(yōu)化后的氣動(dòng)特性， 若采用CFD仿真， 那么氣動(dòng)力的反復(fù)計(jì)算在迭代優(yōu)化的過(guò)程中就會(huì)耗時(shí)很長(zhǎng); 若采用工程算法， 計(jì)算精度則較低。 如何在保證計(jì)算精度的前提下提高計(jì)算速度， 就顯得尤為重要。

本文采用Co-Kriging模型來(lái)代替在優(yōu)化過(guò)程中對(duì)滑翔制導(dǎo)炮彈氣動(dòng)參數(shù)的計(jì)算， Co-Kriging模型是在Kriging模型的基礎(chǔ)上發(fā)展起來(lái)的地質(zhì)統(tǒng)計(jì)學(xué)方法[12]， 通過(guò)該模型可以用少量的CFD仿真數(shù)據(jù)修正大量易獲取的工程算法數(shù)據(jù)， 從而建立高精度代理模型。 相較于直接使用CFD和工程算法計(jì)算， 該方法不但可以保證較高的精度而且計(jì)算時(shí)間大大降低， 從而實(shí)現(xiàn)對(duì)鴨舵氣動(dòng)外形的快速優(yōu)化。

1.1 Co-Kriging代理模型

Co-Kriging模型是利用少量高可信度數(shù)據(jù)來(lái)修正大量的低可信度數(shù)據(jù)， 從而建立模型來(lái)近似代替炮彈外形參數(shù)與氣動(dòng)參數(shù)之間的映射關(guān)系 y=f（x） 。 那么如何獲得該映射關(guān)系最佳的近似函數(shù) f^（x） 就是代理模型方法的關(guān)鍵所在。

假定在樣本集Xe處抽樣獲得精度較高的高可信度數(shù)據(jù) n e個(gè)， 其函數(shù)值記為Ye， 在樣本集Xc處抽樣獲得精度較低的低可信度數(shù)據(jù) n c個(gè)， 其函數(shù)值記為Yc。 將樣本集X處的函數(shù)值視為一個(gè)穩(wěn)態(tài)隨機(jī)過(guò)程， 于是， 對(duì)于Co-Kriging模型， 可得到如下的隨機(jī)變量場(chǎng)：

Y=Yc （Xc ）

Ye （Xe ）=

Yc （X（1）c ）Yc （X（nc ）c ）Ye （X（1）e ）Ye （X（ne ）e ） （1）

構(gòu)造的代理模型的近似估計(jì)可以表示為

Ze （x）=ρZc （x）+Zd （x） （2）

式中： ?Zc （x） 為對(duì)低可信度樣本數(shù)據(jù)構(gòu)造的Kriging模型; ?Zd （x） 為對(duì)高、 低可信度樣本數(shù)據(jù)差值構(gòu)造的Kriging模型。 這里實(shí)質(zhì)上是將高可信度數(shù)據(jù)近似看作低可信度數(shù)據(jù)乘以一個(gè)比例因子常量 ρ 后， 再加上一個(gè)表示 ρZc （x） 和 Ze （x） 之間差別的高斯隨機(jī)過(guò)程 Zd （x） 。 考慮兩組獨(dú)立的高、 低可信度數(shù)據(jù)， 可以構(gòu)建Co-Kriging完整的協(xié)方差矩陣：

C=σ2c Ψc（Xc， Xc） ρσ 2cΨc（Xc， Xe）

ρσ 2cΨc（Xe， Xc） ρσ 2cΨc（Xe， Xe）+ σ 2dΨd（Xe， Xe）（3）

式中： ?σ 2c， ?σ 2d為 Z c （x） 和 Z d （x） 的方差; ψc（Xe， Xc）為Xe和Xc之間的相關(guān)性矩陣， 這些相關(guān)性依賴于樣本點(diǎn)間的絕對(duì)距離以及參數(shù) θ 和 p 。 與Kriging模型不同的是， 這里有兩種不同的空間相關(guān)函數(shù)ψc和ψd， 就需要訓(xùn)練更多的模型參數(shù)： ?θ c， ?θ d， ?p c， ?p d以及比例因子 ρ 。 為了估計(jì)參數(shù) θ d， ?p d和 ρ 的值， 首先定義：

d=Ye- ρ Yc（Xe）（4）

其中， Yc（Xe）是在點(diǎn)Xe處低可信度模型的值Yc。 這樣就可以采用與標(biāo)準(zhǔn)Kriging類似的優(yōu)化過(guò)程來(lái)求得對(duì)θd， pd及 ρ 的估計(jì)[13]。 Co-Kriging模型預(yù)估值表達(dá)式如下：

e （x）=+c T C-1（y-1） （5）

考察一下采用該預(yù)估值表達(dá)式來(lái)預(yù)測(cè)某一高可信度樣本點(diǎn)處的函數(shù)值。 假設(shè)有 x（n+1）=x（i） e， 且c是C的第 nc+i 列， 那么cTC-1是一個(gè)第 nc+i 個(gè)元素為1、 其余元素為0的單位向量， 則 e （x（i）e ）=+y（nc +i）-=y（i） e。 可以認(rèn)為， 該式是對(duì)高可信度數(shù)據(jù)的一個(gè)插值模型。 從另一個(gè)角度， 一定程度上也可認(rèn)為是對(duì)低可信度數(shù)據(jù)的回歸， 并使得在高可信度樣本點(diǎn)處的預(yù)測(cè)值與已知的函數(shù)值Ye一致。

1.2 數(shù)值方法與驗(yàn)證

高可信度樣本點(diǎn)的響應(yīng)值由CFD數(shù)值計(jì)算得到。 為保證計(jì)算精度， 采用結(jié)構(gòu)化網(wǎng)格劃分流場(chǎng)并在彈箭頭部和翼面等氣流變化劇烈處加密網(wǎng)格， 以提高網(wǎng)格質(zhì)量和數(shù)量。 求解器中， 以N-S方程作為流場(chǎng)計(jì)算的主控方程， 采用realisable K-epsilon湍流模型、 二階迎風(fēng)差分格式以及密度基隱格式， 對(duì)滑翔制導(dǎo)炮彈的繞流場(chǎng)進(jìn)行數(shù)值模擬。

為了驗(yàn)證擬采用的計(jì)算方法， 選取120 mm迫擊炮彈為對(duì)象進(jìn)行氣動(dòng)力計(jì)算驗(yàn)證， 如圖1所示。

所得計(jì)算結(jié)果與文獻(xiàn)[14]中提供的氣動(dòng)數(shù)據(jù)規(guī)律基本一致， 且升力系數(shù)導(dǎo)數(shù)的誤差小于5%、 俯仰力矩系數(shù)導(dǎo)數(shù)的誤差小于9%。 這表明擬采用的CFD數(shù)值計(jì)算方法可行有效。

1.3 代理模型構(gòu)建過(guò)程

作為代替真實(shí)氣動(dòng)力計(jì)算方法的工具， 代理模型的建立必然離不開(kāi)先驗(yàn)數(shù)據(jù)的獲取與計(jì)算， 其構(gòu)建過(guò)程如圖2所示。

（1） 確定鴨舵各個(gè)外形參數(shù)的取值范圍， 形成一個(gè)多維設(shè)計(jì)空間。

（2） 在設(shè)計(jì)空間中采用全因子試驗(yàn)選取大量的低可信度樣本點(diǎn)構(gòu)成樣本集Xc， 并均勻選取少量的高可信度樣本點(diǎn)構(gòu)成樣本集Xe。 其中， 高可信度樣本點(diǎn)是低可信度樣本點(diǎn)的子集（XeXc）。

（3） 對(duì)低可信度樣本點(diǎn)采用文獻(xiàn)[6]中的工程算法計(jì)算其氣動(dòng)參數(shù)。 高可信度樣本點(diǎn)的氣動(dòng)參數(shù)則采用CFD數(shù)值模擬來(lái)求解。

（4） 基于以上高、 低可信度數(shù)據(jù)組成的樣本庫(kù)， 采用第1.1節(jié)方法構(gòu)建代理模型。

所建立的Co-Kriging代理模型， 其輸入為鴨舵的外形設(shè)計(jì)變量xi， 輸出為炮彈的氣動(dòng)參數(shù)。 利用該模型， 設(shè)計(jì)空間內(nèi)的其他外形方案， 就無(wú)需經(jīng)過(guò)復(fù)雜的CFD仿真而能快速預(yù)估其氣動(dòng)參數(shù)并保證足夠的精度。 相較于傳統(tǒng)的氣動(dòng)力計(jì)算方法， 這樣一個(gè)有兩種數(shù)據(jù)來(lái)源的多可信度代理模型不僅有著更快的計(jì)算速度， 而且可以過(guò)濾數(shù)值噪聲并顯著提高優(yōu)化效率。

2 鴨舵優(yōu)化設(shè)計(jì)模型

滑翔制導(dǎo)炮彈通常是在尾翼彈的基礎(chǔ)上通過(guò)加裝制導(dǎo)裝置和控制機(jī)構(gòu)而來(lái)， 其氣動(dòng)外形設(shè)計(jì)可在彈體確定的前提下， 在彈頭部添加具有適當(dāng)外形的舵面。 由于鴨舵對(duì)制導(dǎo)炮彈滑翔性能具有重要影響， 本文選取鴨舵作為優(yōu)化對(duì)象， 對(duì)其主要外形參數(shù)進(jìn)行優(yōu)化設(shè)計(jì)， 從而使炮彈具有較高的滑翔控制能力。

2.1 設(shè)計(jì)變量及目標(biāo)函數(shù)的選取

依據(jù)舵面參數(shù)確定的原則[15]， 為了兼顧舵面有較好的氣動(dòng)特性、 較小的舵面鉸鏈力矩和舵面結(jié)構(gòu)的強(qiáng)度和剛度等， 選用后掠的梯形舵面。 如圖3所示， 其外形參數(shù)主要包括舵面的展長(zhǎng) B c（本文以展長(zhǎng)的一半作為設(shè)計(jì)變量， 即 12B c）、 根弦長(zhǎng) C c、 舵面前緣至彈頂?shù)木嚯xXc和前緣后掠角 χ c， 以上4個(gè)主要幾何參數(shù)就是鴨舵外形優(yōu)化的設(shè)計(jì)變量， 即

X=[x1， x2， x3， x4] T（6）

對(duì)于滑翔制導(dǎo)炮彈， 為提高射程， 通常保持在最大升阻比狀態(tài)下滑翔飛行。 因此， 本文選取制導(dǎo)炮彈的升阻比K作為優(yōu)化目標(biāo)函數(shù)。 當(dāng)4個(gè)設(shè)計(jì)變量確定后， 就可以在一定的馬赫數(shù)和雷諾數(shù)下， 計(jì)算出制導(dǎo)炮彈的升阻比， 其函數(shù)關(guān)系為

K=RyRx=f（X）=f12Bc ， Cc， ?Xc ， χc（7）

式中： K為全彈升阻比; ?Ry， Rx 為飛行彈道上全彈的升力和阻力。

2.2 約束條件

約束條件是對(duì)優(yōu)化設(shè)計(jì)問(wèn)題提出的一些條件限制， 以使設(shè)計(jì)方案在滿足這些限制下達(dá)到較優(yōu)。 對(duì)于滑翔增程制導(dǎo)炮彈， 所參加調(diào)優(yōu)的設(shè)計(jì)方案應(yīng)保證如下條件。

（1） 穩(wěn)定性約束條件

鴨舵張開(kāi)不偏轉(zhuǎn)時(shí)， 全彈靜穩(wěn)定性儲(chǔ)備量 ε 應(yīng)滿足： ?c1≤ε≤c2 ， 其中 c1 ， ?c2 為預(yù)定常數(shù)。 靜穩(wěn)定儲(chǔ)備量可以表示為

ε=XCP -XCG l （8）

式中： ?X CP為全彈壓心至彈頂距離; ?XCG為全彈重心至彈頂距離; ?l 為彈長(zhǎng)。

（2） 機(jī)動(dòng)性約束條件

制導(dǎo)炮彈在飛行過(guò)程中， 假設(shè)在任意時(shí)刻都處于平衡狀態(tài)， 即舵面偏轉(zhuǎn)時(shí)作用的力矩在每個(gè)瞬時(shí)都處于平衡狀態(tài)。 由此可得

qSm Cαyα（XCP -XCG ）=qSm Cδzyδzlδ （9）

即

αδz=CδzylδCαyXCP -XCG? （10）

式中： ?q=ρv22 ; ?S m為彈體的特征面積; ?Cαy ， ?Cδzy 為彈身和鴨舵的升力系數(shù)導(dǎo)數(shù); ?αδz 為單位舵偏角所能引起的平衡攻角， 該值越大， 則穩(wěn)定性越小， 機(jī)動(dòng)性增加。 鴨舵張開(kāi)進(jìn)行控制時(shí)， 應(yīng)滿足機(jī)動(dòng)性要求： ?c3≤αδz≤c4 ， 其中 c3 ， ?c4 為預(yù)定常數(shù)。 對(duì)于制導(dǎo)炮彈在舵面實(shí)施控制時(shí)， 需要兼顧炮彈的穩(wěn)定性和機(jī)動(dòng)性。

（3） 尺寸約束條件

制導(dǎo)炮彈在發(fā)射前， 通常將尾翼和舵面折疊在彈體內(nèi)， 因此， 希望舵面的尺寸、 體積滿足火炮發(fā)射條件， 即

Bc ≤B*c Sc ≤S*c （11）

式中： ?B *c為1對(duì)鴨舵展長(zhǎng)的最大允許值; ?S *c為1個(gè)舵面的面積最大值。

（4）幾何約束條件

受鴨式氣動(dòng)布局限制及要求， 鴨舵位置受一定幾何關(guān)系約束， 即

x*c≤ x c≤ X CG- C c（12）

式中： ?x *c為舵面前緣至彈頂?shù)脑试S最近距離; ?X CG為全彈重心至彈頂距離。

2.3 優(yōu)化流程

綜上所述， 滑翔制導(dǎo)炮彈鴨舵的外形優(yōu)化數(shù)學(xué)模型可以表示為

min F（X）=-f（X）=-k? X∈R∩E4

gi（X）≥0，? ?i=1， 2， 3， 4 ??（13）

式中： ?E4 為四維歐氏空間; R為滿足約束條件的實(shí)數(shù)空間; ?gi（X） 為上述約束條件。

綜合考慮， 采用遺傳算法進(jìn)行目標(biāo)函數(shù)尋優(yōu)， 該算法具有通用性強(qiáng)、 穩(wěn)定性好等特性， 且有著良好的全局搜索能力。 滑翔制導(dǎo)炮彈鴨舵外形的優(yōu)化過(guò)程大致如下：

（1） 確定鴨舵外形優(yōu)化的設(shè)計(jì)變量、 約束條件及目標(biāo)函數(shù)等;

（2） 利用遺傳算法在鴨舵外形參數(shù)構(gòu)成的設(shè)計(jì)空間中生成初始種群;

（3） 通過(guò)訓(xùn)練好的Co-Kriging代理模型快速計(jì)算種群中對(duì)應(yīng)不同鴨舵外形的升阻比K;

（4） 在約束條件下對(duì)目標(biāo)函數(shù)， 即升阻比K進(jìn)行適應(yīng)度評(píng)估， 若不滿足收斂條件則使種群交叉變異生成新種群并返回步驟（3）， 直至收斂輸出最優(yōu)結(jié)果;

（5） 采用CFD數(shù)值模擬對(duì)遺傳算法的收斂結(jié)果進(jìn)行驗(yàn)證， 若誤差較大則重新構(gòu)建代理模型， 若誤差較小則確定為最優(yōu)外形方案。

3 算例仿真

3.1 算例分析

本文以某大口徑滑翔制導(dǎo)炮彈為例， 該炮彈尾翼3對(duì)、 舵面2對(duì)， 舵面與尾翼以“++”形式配置。 該炮彈鴨舵的各個(gè)設(shè)計(jì)變量及其取值范圍如表1所示。

在鴨舵外形參數(shù)構(gòu)成的4維設(shè)計(jì)空間中（表1）， 采用4×5×3×3的全因子試驗(yàn)設(shè)計(jì)均勻采樣。 全因子設(shè)計(jì)可以均勻填充樣本空間， 當(dāng)因子數(shù)不超過(guò)5時(shí)較為合適。 共采集低可信度樣本點(diǎn)180個(gè)， 并在設(shè)計(jì)狀態(tài)： 來(lái)流馬赫數(shù)Ma=0.7、 攻角 α=6 °的飛行條件下選用工程算法計(jì)算。

為方便分析， 利用控制變量法篩選部分樣本點(diǎn)， 得到關(guān)于各個(gè)外形參數(shù)對(duì)全彈升阻比K及靜穩(wěn)定儲(chǔ)備量 ε 的影響趨勢(shì)曲線， 如圖4所示。 由圖可知， 全彈升阻比隨著鴨舵展長(zhǎng)的增大而不斷增大， 隨著鴨舵前緣后掠角的增大而減小， 并且隨著鴨舵根弦長(zhǎng)的增大而緩慢升高。 鴨舵至彈頂?shù)木嚯x越短， 全彈的壓心不斷前移， 則炮彈的靜穩(wěn)定度就越小。 通過(guò)對(duì)工程算法的計(jì)算結(jié)果分析， 可以認(rèn)為在4維設(shè)計(jì)空間內(nèi)存在一個(gè)點(diǎn)， 可以使得升阻比達(dá)到一個(gè)最大值。

如表2所示， 在設(shè)計(jì)空間中以樣本點(diǎn)2為中心點(diǎn)均勻采集高可信度樣本點(diǎn)10個(gè)。 對(duì)表2中每一個(gè)高可信度樣本點(diǎn)所對(duì)應(yīng)的滑翔制導(dǎo)炮彈分別建立模型、 劃分流場(chǎng)網(wǎng)格并采用CFD計(jì)算氣動(dòng)參數(shù)， 計(jì)算條件與上述低可信度樣本點(diǎn)一致。 以某外形為例， 彈體表面網(wǎng)格如圖5所示， 結(jié)構(gòu)化網(wǎng)格數(shù)量900萬(wàn)左右， 網(wǎng)格質(zhì)量0.7以上， 以此保證足夠的精度。

利用上文10個(gè)高可信度數(shù)據(jù)和180個(gè)低可信度數(shù)據(jù)構(gòu)建Co-Kriging代理模型， 實(shí)現(xiàn)在優(yōu)化設(shè)計(jì)過(guò)程中對(duì)不同鴨舵外形的滑翔制導(dǎo)炮彈的氣動(dòng)力快速計(jì)算。 該模型輸入為鴨舵的外形參數(shù)， 輸出為氣動(dòng)特性參數(shù)（即全彈升阻比）。 采用遺傳算法作為優(yōu)化方法在約束條件下搜索最大升阻比， 經(jīng)反復(fù)試算， 設(shè)置種群規(guī)模為20、 變異概率為0.2、 超過(guò)進(jìn)化代數(shù)30代時(shí)優(yōu)化停止， 可取得較好優(yōu)化效果。

3.2 仿真結(jié)果

對(duì)構(gòu)建的代理模型進(jìn)行驗(yàn)證。 取9個(gè)點(diǎn)（與構(gòu)建代理模型的樣本點(diǎn)不重復(fù)）分別計(jì)算其高精度CFD仿真結(jié)果和Co-Kriging代理模型預(yù)測(cè)值， 結(jié)果如圖6所示。 觀察圖6可以發(fā)現(xiàn)， 代理模型的預(yù)測(cè)結(jié)果與CFD仿真計(jì)算結(jié)果差距不大， 樣本點(diǎn)5處的誤差最大， 為3.4%， 誤差最小處（樣本點(diǎn)2）只有0.08%， 平均誤差為1.94%。 同時(shí)， 代理模型對(duì)隨外形參數(shù)變化的樣本點(diǎn)的預(yù)測(cè)結(jié)果趨勢(shì)與CFD計(jì)算趨勢(shì)相同。 這表明， 本文構(gòu)建的多可信度代理模型有較高的精度， 可以代替CFD數(shù)值模擬用于對(duì)制導(dǎo)炮彈的氣動(dòng)參數(shù)預(yù)測(cè)。

在多可信度代理模型基礎(chǔ)上， 按照制導(dǎo)炮彈的鴨舵氣動(dòng)外形優(yōu)化模型編制計(jì)算程序。 遺傳算法對(duì)目標(biāo)函數(shù)的優(yōu)化收斂歷程如圖7所示， 當(dāng)種群更新至第七代時(shí)開(kāi)始逐漸收斂。

表3為優(yōu)化前后具體結(jié)果， 優(yōu)化后鴨舵外形如圖8所示， 實(shí)體部分為鴨舵外形優(yōu)化后的炮彈頭部， 虛影部分為優(yōu)化前鴨舵的初始外形方案。 與初始方案相比， 優(yōu)化方案的半展長(zhǎng)增加了24%， 根弦長(zhǎng)增加了10.76%， 至彈頂距離縮小至199.85 mm， 后掠角變?yōu)?2.78°。

據(jù)優(yōu)化結(jié)果， 其升阻比為2.34， 采用CFD數(shù)值仿真驗(yàn)證其阻力系數(shù)0.41、 升力系數(shù)0.95， 升阻比為2.317， 升阻比誤差在0.99%， 表明優(yōu)化結(jié)果有較高可信度。 優(yōu)化方案相較于初始方案， 升阻比提升了11.4%。 該外形下的制導(dǎo)炮彈， 單位舵偏角下產(chǎn)生的平衡攻角為 α/δ=0.73 ， 依據(jù)文獻(xiàn)[15-16]判斷其滿足操縱性約束條件。 鴨舵張開(kāi)后其穩(wěn)定儲(chǔ)備量為5.76%， 滿足在2%～6%之間的穩(wěn)定性設(shè)計(jì)要求[17]。

對(duì)較優(yōu)外形采用CFD數(shù)值仿真計(jì)算， 圖9為全彈升阻比在不同攻角下隨馬赫數(shù)的變化曲線。 由圖可知， 全彈升阻比在亞音速范圍內(nèi)較大， 在跨聲速時(shí)期隨馬赫數(shù)增加迅速減小， 這是因?yàn)樵诳缏曀贂r(shí)由于激波的出現(xiàn)阻力陡增。 全彈的最大升阻比可以達(dá)到2.5左右， 滿足鴨式布局制導(dǎo)炮彈在亞聲速段機(jī)動(dòng)飛行的總體要求。

經(jīng)以上計(jì)算獲得制導(dǎo)炮彈最優(yōu)外形下氣動(dòng)參數(shù)后， 采用有控質(zhì)點(diǎn)彈道模型[18]對(duì)該滑翔制導(dǎo)炮彈進(jìn)行彈道仿真： 假設(shè)炮彈初始質(zhì)量為45 kg， 火箭助推發(fā)動(dòng)機(jī)平均推力為1 200 N， 發(fā)動(dòng)機(jī)裝藥質(zhì)量耗散速率為0.535 7 kg/s， 采用標(biāo)準(zhǔn)氣象條件， 初速930 m/s， 射角61°， 出炮口后， 助推火箭點(diǎn)火工作14 s， 到達(dá)彈道最高點(diǎn)附近時(shí)張開(kāi)鴨舵保持全彈在6°正攻角下滑翔飛行。 飛行彈道軌跡如圖10， 無(wú)控彈射程僅42 km， 在鴨舵控制下初始方案射程為64 km， 優(yōu)化方案的射程可達(dá)到72 km。 相比于初始方案， 鴨舵外形優(yōu)化后的射程提升了12.5%， 在僅僅優(yōu)化鴨舵外形參數(shù)的條件下， 這一增程效果是較為可觀的。

值得注意的是， 對(duì)于上述算例， 對(duì)10個(gè)高可信度樣本點(diǎn)的仿真評(píng)估過(guò)程共花費(fèi)30 h， 對(duì)180個(gè)低可信度樣本點(diǎn)的評(píng)估、 訓(xùn)練代理模型以及程序優(yōu)化時(shí)間一共只需25 s左右（采用Intel Core四核處理器、 主頻2.6 GHz、 8 G內(nèi)存）。 若所有樣本點(diǎn)均采用CFD數(shù)值模擬來(lái)評(píng)估氣動(dòng)特性， 則計(jì)算量是巨大的， 需要耗時(shí)數(shù)天乃至數(shù)月[8]。 后續(xù)利用該快速優(yōu)化方法， 可以考慮更多的設(shè)計(jì)變量、 約束條件及目標(biāo)函數(shù)等， 以滿足復(fù)雜多樣的作戰(zhàn)需求。

4 結(jié)? 論

為了兼顧氣動(dòng)參數(shù)的計(jì)算精度和效率， 本文設(shè)計(jì)了一種基于多可信度代理模型的氣動(dòng)力高精度快速計(jì)算方法， 能夠快速預(yù)估不同外形參數(shù)對(duì)應(yīng)的炮彈氣動(dòng)力系數(shù)。 與基于CFD的氣動(dòng)外形優(yōu)化相比， 本文所建立的優(yōu)化過(guò)程大幅縮短了優(yōu)化時(shí)間， 升阻比的計(jì)算誤差約為1.94%。 以某大口徑滑翔制導(dǎo)炮彈為算例， 對(duì)其鴨舵外形的主要幾何參數(shù)進(jìn)行優(yōu)化， 在滿足機(jī)動(dòng)性和穩(wěn)定性的條件下， 炮彈的滑翔升阻比提升了11.4%， 射程增加了12.5%， 驗(yàn)證了本文方法的可行性和有效性。

參考文獻(xiàn)：

[1] 陳琦. 滑翔增程制導(dǎo)炮彈彈道優(yōu)化及制導(dǎo)控制方法研究[D]. 南京： 南京理工大學(xué)，? 2017.

Chen Qi. Investigation on Trajectory Optimization and Guidance & Control Scheme for a Type of Gliding Guided Projectiles[D]. Nanjing： Nanjing University of Science and Technology，? 2017. （in Chinese）

[2] Vasile J D，? Bryson J，? Fresconi F. Aerodynamic Design Optimization of Long Range Projectiles Using Missile DATCOM[C]∥AIAA Scitech 2020 Forum，??? 2020.

[3] Fowler L，? Rogers J. Airframe Performance Optimization of Guided Projectiles Using Design of Experiments[J]. Journal of Spacecraft and Rockets，? 2015，? 52（6）： 1603-1613.

[4] Vasile J D，? Bryson J，? Gruenwald B C，? et al. A Multi-Disciplinary Approach to Design Long Range Guided Projectiles[C]∥AIAA Scitech 2020 Forum，? 2020.

[5] 雷娟棉，? 吳甲生. 增程制導(dǎo)炮彈氣動(dòng)外形設(shè)計(jì)[J]. 航空學(xué)報(bào)，? 2005，? 26（3）： 294-297.

Lei Juanmian，? Wu Jiasheng. Aerodynamic Configuration Design for Extended Range Guided Munition[J]. Acta Aeronautica et Astronautica Sinica，? 2005，? 26（3）： 294-297.（in Chinese）

[6] 常思江，? 王中原，? 劉鐵錚. 鴨式布局防空炮彈的鴨舵設(shè)計(jì)優(yōu)化[J]. 兵工學(xué)報(bào)，? 2010，? 31（4）： 521-524.

Chang Sijiang，? Wang Zhongyuan，? Liu Tiezheng. Optimum Design of Canard Parameters of Air Defense Trajectory Correction Projectile Equipped with a Pair of Canard Surfaces[J]. Acta Armamentarii，? 2010，? 31（4）： 521-524.（in Chinese）

[7] 孫世巖，? 梁偉閣， ?陳俊丞. 適合概念設(shè)計(jì)的制導(dǎo)炮彈氣動(dòng)外形優(yōu)化方法[J]. 海軍工程大學(xué)學(xué)報(bào)，? 2020，? 32（2）： 81-86.

Sun Shiyan，? Liang Weige，? Chen Juncheng. An Aerodynamic Shape Optimization Method for Projectile Conceptual Design[J]. Journal of Naval University of Engineering，? 2020，? 32（2）： 81-86.（in Chinese）

[8] Arnoult G，? Zeidler M，? Garnier E. Control Surface Geometry Surrogate-Based Optimization for Spin-Stabilized Projectile Course Correction[J]. AIAA Journal，? 2019，? 58（2）： 550-560.

[9] 韓忠華，? 許晨舟，? 喬建領(lǐng)，? 等. 基于代理模型的高效全局氣動(dòng)優(yōu)化設(shè)計(jì)方法研究進(jìn)展[J]. 航空學(xué)報(bào)，? 2020，? 41（5）： 30-70.

Han Zhonghua，? Xu Chenzhou，? Qiao Jianling，? et al. Recent Progress of Efficient Global Aerodynamic Shape Optimization Using Surrogate-Based Approach[J]. Acta Aeronautica et Astronautica Sinica，? 2020，? 41（5）： 30-70.（in Chinese）

[10] 原智杰，? 張公平，? 崔茅，? 等. 基于神經(jīng)網(wǎng)絡(luò)的導(dǎo)彈氣動(dòng)參數(shù)預(yù)測(cè)[J]. 航空兵器，? 2020，? 27（5）： 28-32.

Yuan Zhijie，? Zhang Gongping，? Cui Mao，? et al. Prediction of Missiles Aerodynamic Parameters Based on Neural Network[J]. Aero Weaponry，? 2020，? 27（5）： 28-32.（in Chinese）

[11] 劉俊. 基于代理模型的高效氣動(dòng)優(yōu)化設(shè)計(jì)方法及應(yīng)用[D]. 西安： 西北工業(yè)大學(xué)，? 2015.

Liu Jun. Efficient Surrogate-Based Optimization Method and Its Application in Aerodynamic[D]. Xian： Northwestern Polytechnical University，? 2015. （in Chinese）

[12] Forrester A I J，? Sóbester A，? Keane A J. Engineering Design via Surrogate Modelling[M]. Chichester： John Wiley & Sons， 2008.

[13] Toal D J J，? Keane A J. Efficient Multipoint Aerodynamic Design Optimization via Cokriging[J]. Journal of Aircraft，? 2011，? 48（5）： 1685-1695.

[14] Mccoy R L. Modern Exterior Ballistics[M]. Atglen： Schiffer Pub Ltd，? 2012.

[15] 史金光，? 王中原，? 許厚謙，? 等. 滑翔增程彈鴨式舵的氣動(dòng)設(shè)計(jì)與分析[J]. 彈道學(xué)報(bào)，? 2006，? 18（4）： 33-37.

Shi Jinguang，? Wang Zhongyuan，? Xu Houqian，? et al. Aerodynamic Design and Analysis of Canard Rudder for Gliding Extended Range Projectile[J]. Journal of Ballistics，? 2006，? 18（4）： 33-37.（in Chinese）

[16] 雷娟棉，? 吳甲生. 遠(yuǎn)程制導(dǎo)炮彈氣動(dòng)設(shè)計(jì)[C]∥中國(guó)航空學(xué)會(huì)青年科技論壇論文集， 2005： 158-163.

Lei Juanmian，? Wu Jiasheng. Extended Range Guided Munition Aerodynamic Configuration Design[C]∥Proceedings of the China Aviation Society Youth Science and Technology Forum， 2005：158-163.（in Chinese）

[17] 雷娟棉. 制導(dǎo)兵器氣動(dòng)特性工程計(jì)算方法[M]. 北京： 北京理工大學(xué)出版社，? 2015.

Lei Juanmian. Engineering Prediction Methods of Aerodynamic Characteristics for Guided Weapon [M]. Beijing： Beijing Insititute of Technology Press，? 2015.（in Chinese）

[18] 錢杏芳， 林瑞雄， 趙亞男.導(dǎo)彈飛行力學(xué)[M].北京： 北京理工大學(xué)出版社， 2000.

Qian Xingfang，? Lin Ruixiong，? Zhao Yanan.Missile Flight Aerodynamics [M]. Beijing： Beijing Institute of Technology Press，? 2000. （in Chinese）

Research on the Rapid Aerodynamic Shape Optimization of

Canards of Gliding Guided Projectile

Zhao Xuan1， Chang Sijiang1*， Ni Yi2

（1. College of Energy and Power Engineering，? NUST，? Nanjing 210094，? China;

2. Shanghai Institute of Mechanical and Electrical Engineering，? Shanghai 201109，? China）

Abstract： ?In order to improve the accuracy of? calculation，? it is often necessary to calculate the aerodynamic parameters of the projectile iteratively through computational fluid dynamics simulation，? so the optimization efficiency is low. In this paper，? a rapid aerodynamic shape ptimization method based on multi-fidelity surrogate model is proposed，? which improves the calculation efficiency greatly under the condition of ensuring the calculation accuracy. Taking the canard shape optimization of gliding guided projectile as an example，? the multi-fidelity surrogate is trained by two kinds of reliability samples to replace the computational fluid dynamics simulation process to calculate the aerodynamic parameters. According to the design requirements of maximum gliding lift-to-drag ratio，? matching stability and maneuverability，? the genetic algorithm is used to search for the optimal design of canard aerodynamic shape. Compared with the initial scheme，? the lift-to-drag ratio of the optimized scheme are improved significantly. Compared with the results of numerical simulation，? the average error of the method for the lift-to-drag ratio prediction is 1.94%，? the accuracy is higher，? and the amount of calculation is greatly reduced，? which verifies the feasibility and effectiveness of the method.

Key words： gliding guided projectile; aerodynamic shape; Co-Kriging surrogate model; CFD; genetic algorithm; optimal design