基于自適應(yīng)采樣的高超聲速飛行器氣動(dòng)熱全局快速預(yù)示

楊國(guó)濤，岳振江，2，，劉莉，2

1．北京理工大學(xué) 宇航學(xué)院，北京 100081

2．北京理工大學(xué) 飛行器動(dòng)力學(xué)與控制教育部重點(diǎn)實(shí)驗(yàn)室，北京 100081

隨著航空航天技術(shù)的發(fā)展，高超聲速飛行器由于速度快、突防能力強(qiáng)、覆蓋范圍廣、隱身效果顯著等優(yōu)點(diǎn)，具有廣泛的應(yīng)用前景[1］。由于高超聲速飛行器在大氣層內(nèi)高速飛行，氣動(dòng)熱防護(hù)是需要妥善解決的關(guān)鍵問(wèn)題之一[2］。傳統(tǒng)氣動(dòng)熱獲取主要包括試驗(yàn)研究、工程計(jì)算、數(shù)值模擬3 種方法[3］。飛行試驗(yàn)和地面風(fēng)洞試驗(yàn)?zāi)軌颢@取飛行器氣動(dòng)熱的真實(shí)狀況，但長(zhǎng)周期、高成本限制了其在設(shè)計(jì)階段初期的應(yīng)用[3］。氣動(dòng)熱工程計(jì)算基于試驗(yàn)數(shù)據(jù)或邊界層方程相似解，通過(guò)合理假設(shè)，采用經(jīng)驗(yàn)公式、半經(jīng)驗(yàn)公式獲得，具有快速高效的優(yōu)點(diǎn)，但只能處理簡(jiǎn)單外形的飛行器[4］。數(shù)值模擬是利用計(jì)算流體力學(xué)（Computational Fluid Dynamic，CFD）方法求解Navier-Strokes 方程及其近似，對(duì)飛行器外形適應(yīng)性強(qiáng)，能充分考慮氣流黏性、真實(shí)氣體效應(yīng)等，但較低的計(jì)算效率在設(shè)計(jì)迭代過(guò)程中應(yīng)用受限[3］。針對(duì)上述傳統(tǒng)方法存在的局限性，為有效解決高耗時(shí)高成本氣動(dòng)熱獲取的工程問(wèn)題，出現(xiàn)了利用相對(duì)簡(jiǎn)單的數(shù)學(xué)模型來(lái)表達(dá)復(fù)雜全階空氣動(dòng)力學(xué)模型輸入輸出關(guān)系的方法[5］，以達(dá)到快速預(yù)示氣動(dòng)熱的目的，并逐漸成為近年來(lái)氣動(dòng)研究的熱點(diǎn)[6］。陳海昕等[7］總結(jié)了機(jī)器學(xué)習(xí)在氣動(dòng)優(yōu)化設(shè)計(jì)中的應(yīng)用，并進(jìn)一步探討了利用深度代理模型和深度降階模型進(jìn)行氣動(dòng)數(shù)據(jù)的歸納分析和流場(chǎng)重建。Yan 等[8］基于POD 和修正的切比雪夫多項(xiàng)式提出一種POD-Chebyshev 方法，應(yīng)用于飛行器控制面將平均誤差降到了2%以內(nèi)。張智超等[9］基于徑向基神經(jīng)網(wǎng)絡(luò)提出一種氣動(dòng)熱預(yù)測(cè)代理模型方法，應(yīng)用在NASA 火星實(shí)驗(yàn)室的橢圓鈍化高超聲速飛行器中，關(guān)鍵區(qū)域預(yù)示誤差小于10%。

氣動(dòng)熱數(shù)學(xué)近似模型主要通過(guò)采樣獲取物理系統(tǒng)的主要特征，因此基于采樣技術(shù)獲取的樣本質(zhì)量很大程度決定樣本所包含的系統(tǒng)信息，進(jìn)而影響模型的預(yù)示精度?，F(xiàn)有的采樣方法包括一 次 性 采 樣（One-Shot）和 序 列 采 樣[6，10］。One-Shot 一次性獲取具有良好填充均勻性和投影均勻性的樣本，建模次數(shù)少、效率高，常用的方法有全析因設(shè)計(jì)、正交設(shè)計(jì)、均勻設(shè)計(jì)、拉丁超立方設(shè)計(jì)（Latin Hypercube Design， LHD）和最優(yōu)LHD（OLHD）等[6，10-11］，但是存在欠采樣或過(guò)采樣的局限[12］。序列采樣利用少量初始樣本構(gòu)建近似模型，并通過(guò)迭代增加樣本更新近似模型，以達(dá)到有效控制采樣規(guī)模和兼顧模型精度的目的[12］。自適應(yīng)采樣是一類重要的序列采樣方法，在利用已有樣本空間信息的同時(shí)，還充分考慮系統(tǒng)的物理特征在具有較大誤差或高度非線性等局部欠采樣區(qū)域重點(diǎn)新增樣本，構(gòu)建更加準(zhǔn)確高效的近似模型，并達(dá)到節(jié)約計(jì)算成本的目的，相較于一次性采樣，自適應(yīng)采樣具有更加廣泛的應(yīng)用前景[13］。

根據(jù)預(yù)示模型的用途，自適應(yīng)采樣方法主要包括兩類，一類是面向優(yōu)化以提高模型在其關(guān)注的局部區(qū)域的預(yù)示精度[14-16］；另一類是面向全局近似以提高整體預(yù)示精度[17］。針對(duì)全局近似任務(wù)，自適應(yīng)采樣通過(guò)順序迭代獲取新增樣本信息，以盡可能小的采樣規(guī)模獲取整個(gè)參數(shù)空間內(nèi)充分的系統(tǒng)信息，支持響應(yīng)的準(zhǔn)確預(yù)示。全局探索與局部開(kāi)發(fā)的平衡對(duì)全局近似效果至關(guān)重要，有效的平衡不僅可以在已識(shí)別的區(qū)域中進(jìn)行充分的局部開(kāi)發(fā)，還可以及時(shí)地提供全局探索以避免錯(cuò)過(guò)一些未檢測(cè)到的區(qū)域[18］。Xiong 等[19］利用基于非平穩(wěn)協(xié)方差的 Kriging 模型的預(yù)測(cè)方差輔助采樣以減小模型的不確定性，從而提高預(yù)測(cè)質(zhì)量。Yao 等[20］提出一種基于梯度的徑向基神經(jīng)網(wǎng)絡(luò)方法，利用現(xiàn)有模型的梯度信息依次擴(kuò)展樣本集并細(xì)化模型，從而有效提高模型近似精度。Farhang-Mehr 等[21］提出了最大熵自適應(yīng)采樣的方法，旨在將樣本放置在響應(yīng)的不規(guī)則區(qū)域。Turner 等[22］采用模擬退火來(lái)平衡局部開(kāi)發(fā)和全局探索，以達(dá)到多準(zhǔn)則自適應(yīng)采樣的目的。Li 等[23］提出一種基于預(yù)測(cè)誤差的Kriging 近似自適應(yīng)方法，以改善響應(yīng)中近似不佳的部分。Aute 等[24］提出了基于交叉驗(yàn)證誤差的序列采樣方法，根據(jù)預(yù)測(cè)誤差在高度非線性區(qū)域采集更多的點(diǎn)。Chen 等[25］提出了一種基于模糊聚類的自適應(yīng)加點(diǎn)策略（Adding-Point Strategy based on Fuzzy Clustering， APSFC），通過(guò)對(duì)具有較大誤差的樣本集進(jìn)行聚類，并在聚類內(nèi)部產(chǎn)生新增樣本，一定程度提高了構(gòu)建高超聲速升力面氣動(dòng)熱預(yù)示模型的效率和精度。Crombecq 等[12］提出了一種通用可靠的LOLA-Voronoi 混合序列采樣方法，通過(guò)局部線性近似衡量樣本周圍非線性程度和通過(guò)Voronoi 圖衡量樣本周圍密度相結(jié)合的方法獲取新增樣本，所提方法針對(duì)不同模型類型都能取得很好的效果，但是在估計(jì)樣本梯度時(shí)需要進(jìn)行復(fù)雜的鄰域計(jì)算。Xu 等[26］將LOO（Leave-One-Out）交叉驗(yàn)證和Voronoi 圖相結(jié)合提出一種通過(guò)誤差追蹤逐步提高預(yù)示模型準(zhǔn)確性的CV-Voronoi 方法，與LOLA-Voronoi 相比，該方法效果更優(yōu)。Kaminsky 等[27］還 結(jié) 合k-fold 交 叉 驗(yàn) 證 和 批 量 采樣技術(shù)提出了一種改善的KFCV-Voronoi 方法，降低了采樣過(guò)程中的計(jì)算成本。

為加快高超聲速飛行器氣動(dòng)熱環(huán)境分析，本文提出了一種基于模糊聚類超球體影響域和當(dāng)?shù)卣`差評(píng)分系數(shù)加權(quán)拒絕域的自適應(yīng)采樣方法（adaptive sampling method based on Fuzzy Clustering HyperSphere influence domain and local error scoring coefficient Weighted refused domain， FCHSW），在此基礎(chǔ)上建立氣動(dòng)熱全局快速預(yù)示模型（Rapid Predicting Model，RPM）。首先，對(duì)一次采樣獲取的少量高精度工況構(gòu)建氣動(dòng)熱RPM；然后，評(píng)估模型的誤差分布特征，建立基于模糊聚類的超球影響域和基于局部誤差評(píng)分系數(shù)加權(quán)的拒絕域，進(jìn)而結(jié)合maxmin 準(zhǔn)則在綜合確定的重點(diǎn)采樣域新增高精度樣本，提高采樣質(zhì)量，實(shí)現(xiàn)RPM 的自適應(yīng)更新，高效提升全局預(yù)示精度。

1 氣動(dòng)熱全局快速預(yù)示框架

高超聲速飛行器氣動(dòng)熱快速準(zhǔn)確預(yù)示是高超聲速熱防護(hù)和防隔熱設(shè)計(jì)分析的重要前提。面向高精度高效率氣動(dòng)熱設(shè)計(jì)的需求，本文采用基于數(shù)據(jù)驅(qū)動(dòng)的氣動(dòng)熱全局快速預(yù)示框架，根據(jù)有限的高精度氣動(dòng)熱數(shù)值計(jì)算結(jié)果構(gòu)造簡(jiǎn)單數(shù)學(xué)模型來(lái)替代復(fù)雜CFD 模型的評(píng)估，實(shí)現(xiàn)參數(shù)空間內(nèi)飛行器氣動(dòng)熱環(huán)境快速準(zhǔn)確預(yù)示。

氣動(dòng)熱全局快速預(yù)示框架按順序分為CFD氣動(dòng)熱計(jì)算與RPM 構(gòu)建、氣動(dòng)熱RPM 測(cè)試與更新2 個(gè)關(guān)鍵步驟，基本流程如圖1 所示。

圖1 氣動(dòng)熱RPM 框架示意Fig.1 Schematic of aerothermodynamics RPM framework

1.1 CFD 氣動(dòng)熱計(jì)算與RPM 構(gòu)建

為有效控制飛行器高精度數(shù)值計(jì)算成本，需要從眾多流場(chǎng)控制參數(shù)中確定對(duì)物理系統(tǒng)有重要影響的最小維度飛行狀態(tài)參數(shù)設(shè)計(jì)空間，并利用采樣技術(shù)獲取少量具有良好空間填充性和投影均勻性的初始訓(xùn)練樣本，以代表整個(gè)飛行狀態(tài)空間。

根據(jù)選定的飛行器氣動(dòng)外形，基于CFD 計(jì)算策略獲取所有初始飛行工況的高精度氣動(dòng)熱數(shù)值計(jì)算結(jié)果。然后基于給定狀態(tài)下流場(chǎng)的輸入輸出特性確定數(shù)學(xué)近似模型，并采用參數(shù)估計(jì)法對(duì)數(shù)學(xué)模型進(jìn)行訓(xùn)練，完成飛行器表面氣動(dòng)熱流RPM 的構(gòu)造。

1.2 氣動(dòng)熱RPM 測(cè)試與更新

初始構(gòu)建的訓(xùn)練模型往往精度較低，為了使RPM 能充分實(shí)現(xiàn)對(duì)大規(guī)模高耗時(shí)復(fù)雜CFD 系統(tǒng)的近似，根據(jù)訓(xùn)練樣本的空間分布信息和訓(xùn)練模型的誤差分布特性，基于自適應(yīng)采樣方法添加新的訓(xùn)練樣本，并調(diào)用CFD 模型獲取新增樣本響應(yīng)，據(jù)此迭代更新訓(xùn)練模型以逐步提高其預(yù)示精度。模型的誤差分布特性可以通過(guò)預(yù)示方差、交叉驗(yàn)證誤差以及梯度等信息來(lái)衡量，其中交叉驗(yàn)證誤差由于易于獲取且不依賴模型，更廣泛的用于自適應(yīng)指導(dǎo)模型更新[13，24］。

為了提供系統(tǒng)的主要特征，要求RPM 與CFD 模型誤差應(yīng)當(dāng)足夠小，因此獲取額外的測(cè)試樣本對(duì)每一輪更新的模型進(jìn)行精度校驗(yàn)，以驗(yàn)證模型的準(zhǔn)確性，并判斷迭代終止。

由于訓(xùn)練樣本規(guī)模或高精度計(jì)算成本的限制，為了提升高精度樣本采樣質(zhì)量，進(jìn)一步改善RPM 的全局預(yù)示精度，更好的適用于可能存在的多極值、強(qiáng)非線性問(wèn)題，自適應(yīng)采樣是指導(dǎo)RPM高效更新的關(guān)鍵技術(shù)。

2 FCHSW 基本方法介紹

針對(duì)有限計(jì)算成本下提升RPM 全局預(yù)示精度的問(wèn)題，本文提出的FCHSW 方法基本流程如圖2 所示，影響域和拒絕域分別為新增樣本的潛在區(qū)域和拒絕區(qū)域，同時(shí)結(jié)合maxmin 準(zhǔn)則，綜合權(quán)衡參數(shù)空間的全局探索和局部開(kāi)發(fā)，實(shí)現(xiàn)RPM精度的高效提升。

圖2 FCHSW 基本流程Fig.2 Basic flow of FCHSW

FCHSW 具體步驟如下：

步驟 1確定nv維設(shè)計(jì)變量、參數(shù)空間Ω、測(cè)試樣本數(shù)nt及方法相關(guān)參數(shù)。方法相關(guān)參數(shù)主要包括：初始訓(xùn)練樣本數(shù)ns0、每輪新增樣本數(shù)nadd、波 動(dòng) 系 數(shù)μ以 及 終 止 條 件（預(yù) 期 誤 差 閾 值e?est和終止樣本數(shù)閾值nest）。

步驟 2基于OLHD 在參數(shù)空間生成相互獨(dú)立的初始訓(xùn)練樣本集x和測(cè)試樣本集xt，計(jì)算真實(shí)響應(yīng)y和yt，并初始化樣本數(shù)據(jù)庫(kù)。

步驟 3利用初始訓(xùn)練樣本集構(gòu)建RPM，利用測(cè)試樣本集進(jìn)行模型精度校驗(yàn)。判斷是否滿足終止條件，即當(dāng)預(yù)期誤差e?＜e?est或訓(xùn)練樣本數(shù)ns＞nest，則迭代終止，否則進(jìn)行步驟4。本文將平均相對(duì)誤差作為模型精度評(píng)估的停機(jī)準(zhǔn)則：

其中：上標(biāo)U 和L 分別表示真實(shí)響應(yīng)yt的最大值和最小值；下標(biāo)Nt 表示響應(yīng)在區(qū)間[1，2]的歸一化值，表示為

步驟 4利用k-fold 或LOO 獲取所有樣本的交叉驗(yàn)證誤差CV，挑選出大于平均誤差的樣本做模糊聚類，進(jìn)而得到每個(gè)聚類的影響域DEi(1 ≤i≤nc)，其 中nc為聚 類個(gè)數(shù)。

步驟 5利用超球體對(duì)參數(shù)空間進(jìn)行分割，獲取每個(gè)樣本的拒絕域DRi(1 ≤i≤ns)。

步驟 6根據(jù)獲取的影響域和拒絕域，新的設(shè)計(jì)域?yàn)槌跏紖?shù)空間對(duì)拒絕域的補(bǔ)集與影響域的交集?？紤]到新的設(shè)計(jì)域?yàn)椴灰?guī)則區(qū)域，采用Monte Carlo 仿真在整個(gè)參數(shù)空間生成隨機(jī)分布的若干樣本xrand，并利用接受-拒絕采樣將設(shè)計(jì)域Dc內(nèi)的隨機(jī)樣本作為候選樣本xcandidate；進(jìn)而根據(jù)maxmin 準(zhǔn)則[28］從xcandidate中選出每輪迭代的新增樣本xadd并計(jì)算響應(yīng)yadd；最后，更新樣本數(shù)據(jù)庫(kù)并返回步驟3 以實(shí)現(xiàn)RPM 的更新。

隨機(jī)樣本個(gè)數(shù)定義為N(xrand)=ns0×nv×ω，ω為控制參數(shù)，為保證方法有效性的同時(shí)盡可能降低計(jì)算成本，參考文獻(xiàn)[26］的取值，本文取ω=100。

圖3 展示了二維情況下FCHSW 的一輪迭代加點(diǎn)示意，實(shí)心圓為訓(xùn)練樣本，圓半徑越大表示誤差越大，將具有較大誤差的樣本集分為2 個(gè)聚類，三角形為聚類中心，實(shí)線圓表示樣本的拒絕域邊界，填充部分為新的設(shè)計(jì)域，五角星代表新增樣本。

圖3 FCHSW 加點(diǎn)過(guò)程示意Fig.3 Adding-point process of FCHSW

2.1 模糊聚類

基于模糊聚類[29］按照相似度對(duì)數(shù)據(jù)集進(jìn)行分類，建立樣本與類的不確定性描述，以實(shí)現(xiàn)對(duì)參數(shù)空間的有效劃分和縮減。模糊c-均值聚類（Fuzzy C-Means clustering method，F(xiàn)CM）將問(wèn)題歸結(jié)為一個(gè)非線性規(guī)劃問(wèn)題，通過(guò)極小化各數(shù)據(jù)點(diǎn)到聚類中心的加權(quán)平均和得到每個(gè)數(shù)據(jù)點(diǎn)對(duì)所有類中心的隸屬度，實(shí)現(xiàn)對(duì)樣本自動(dòng)分類[29］。

考 慮 樣 本 集X=[xi](i=1，2，…，ns)，將樣本集中元素根據(jù)模糊聚類的方法劃分為nc(2 ≤nc≤ns)個(gè)模糊子集，每個(gè)數(shù)據(jù)點(diǎn)對(duì)所有類中心的整體差異可表示為

式中：U=[uij] (i=1，2，…，nc；j=1，2，…，ns)為隸屬矩陣，uij表示第j個(gè)樣本xj在第i個(gè)聚類的隸屬度，滿足，表達(dá)式為

其中：d為影響隸屬度矩陣模糊化程度的指數(shù)權(quán)重（通常取d=2）。求解式（3）中的極小值問(wèn)題，即可得到聚類中心和隸屬度矩陣U*。

通過(guò)最大輪廓系數(shù)法來(lái)預(yù)先確定FCM 的最佳聚類個(gè)數(shù)，以實(shí)現(xiàn)方法的無(wú)監(jiān)督學(xué)習(xí)。輪廓系數(shù)定義如下[30］：

式中：a(i)為第i個(gè)樣本xi與其所在簇內(nèi)其他樣本的平均距離，代表簇內(nèi)凝聚度；b(i)為樣本xi與其他每個(gè)簇樣本平均距離的最小值，代表簇間分離度。輪廓系數(shù)滿足-1 ≤S≤1，越接近1 說(shuō)明聚類效果越好，反之越接近-1 聚類效果越差。

2.2 影響域

影響域是由模糊聚類所屬較大誤差樣本集拓展得到的一個(gè)參數(shù)空間，定義為超球體?？紤]第i 個(gè)影響域DEi(1 ≤i ≤nc)包含的樣本集為xci，如果xci有多個(gè)樣本，定義影響半徑rci為聚類中心到該聚類所有樣本距離的最大值，超球體影響域表示為

如果xci只有一個(gè)樣本，定義超球體影響域?yàn)橐栽摌颖镜木芙^半徑r（拒絕半徑在2.3 節(jié)中給出了介紹）作為半?yún)^(qū)間寬度的超立方體的外接超球體，表示為

2.3 拒絕域與誤差評(píng)分系數(shù)

為充分考慮樣本的空間分布特性，為每個(gè)樣本定義一個(gè)拒絕域，拒絕域?yàn)樾略鰳颖拘枰?guī)避的區(qū)域，從而避免其與初始樣本空間分布過(guò)近造成樣本冗余。

誤差評(píng)分系數(shù)與拒絕域和maxmin 準(zhǔn)則聯(lián)合使用，一方面，誤差越大的樣本其拒絕域應(yīng)當(dāng)越小，以允許候選樣本向其適當(dāng)靠近；另一方面，利用maxmin 準(zhǔn)則從若干候選點(diǎn)中挑選新增樣本時(shí)，通過(guò)誤差評(píng)分系數(shù)對(duì)樣本間的幾何距離進(jìn)行乘法加權(quán)，更容易將靠近較大誤差樣本的候選點(diǎn)挑選出來(lái)，從而在其附近增加更多的樣本，以達(dá)到有效權(quán)衡全局探索和局部開(kāi)發(fā)的目的。

樣本所屬的拒絕域是以該樣本點(diǎn)為球心的超球體，拒絕半徑為該點(diǎn)與剩余樣本點(diǎn)最小距離的一半乘以誤差評(píng)分系數(shù)，樣本xj(1 ≤j ≤ns)的拒絕半徑和拒絕域分別表示為

式中：λj為樣本xj的誤差評(píng)分系數(shù)，根據(jù)當(dāng)?shù)卣`差對(duì)每個(gè)樣本進(jìn)行評(píng)分。λj表示為

式中：μ 為波動(dòng)系數(shù)，當(dāng)μ=1 時(shí)不考慮誤差評(píng)分，此時(shí)誤差評(píng)分系數(shù)λj=1。

根據(jù)每個(gè)聚類的影響域和所有樣本的拒絕域，第i 個(gè)聚類的設(shè)計(jì)域即為

添加樣本時(shí)，根據(jù)誤差評(píng)分系數(shù)分配各設(shè)計(jì)域的新增樣本數(shù)，實(shí)現(xiàn)在較大誤差區(qū)域增添更多的樣本。第i 個(gè)設(shè)計(jì)域的誤差評(píng)分系數(shù)表示為

其 中：e?CVij表示 第i 個(gè) 設(shè) 計(jì) 域 第j 個(gè) 樣 本 的 相 對(duì) 誤差。第i 個(gè)設(shè)計(jì)域新增樣本數(shù)為niadd=λic·nadd。

此外，基于maxmin 準(zhǔn)則能夠獲取較好空間均勻性的樣本，但是直觀上非線性程度較大的區(qū)域?qū)颖镜拿芗纫蟾?，因此?jì)算候選樣本xi∈xcandidate與初始樣本xj∈x 的距離時(shí)，引入誤差評(píng)分系數(shù)對(duì)距離進(jìn)行加權(quán)，表示為

3 算例驗(yàn)證及分析

本文所提自適應(yīng)采樣方法是一種非侵入方法，可以適配不同的數(shù)學(xué)近似模型，為驗(yàn)證方法的效率和可靠性，以Kriging[31-32］為例，應(yīng)用于數(shù)值測(cè)試算例，并討論相關(guān)參數(shù)的影響。

3.1 數(shù)值測(cè)試算例

附錄A 給出了測(cè)試函數(shù)用例，采用LOO 獲取交叉驗(yàn)證誤差。設(shè)置每輪新增樣本數(shù)nadd=nv和波動(dòng)系數(shù)μ=2，收斂準(zhǔn)則為歸一化預(yù)期誤差閾值e?est=1% 或終止樣本數(shù)閾值nest=1 000。選取1 000 個(gè)測(cè)試樣本進(jìn)行精度校驗(yàn)，將所提方法與One-Shot、CV-Voronoi 和APSFC 方 法 進(jìn) 行 對(duì)比。為減小隨機(jī)性影響，更客觀地評(píng)價(jià)方法的效果，每個(gè)算例運(yùn)行50 次，分析其結(jié)果的統(tǒng)計(jì)特性。文獻(xiàn)[13］總結(jié)了初始樣本數(shù)的合適取值，考慮到初始樣本數(shù)不是本文方法的主要關(guān)注點(diǎn)，故不再討論其對(duì)采樣結(jié)果的影響，后續(xù)不做特殊說(shuō)明處，均選取一個(gè)較小的初始樣本數(shù)ns0=5nv計(jì)算。

圖4 展示了初始樣本為6 終止樣本為16 的具有大范圍平坦區(qū)域以及多個(gè)局部極值的一維算例結(jié)果，其中圓為One-Shot 采樣，叉為初始樣本，下三角、上三角、五角星分別表示CV-Voronoi、APSFC 和FCHSW 對(duì) 應(yīng) 新 增 樣 本。One-Shot 在整個(gè)參數(shù)空間均勻采樣，APSFC 有多個(gè)新增樣本落在平坦區(qū)域且新增樣本間存在空間分布過(guò)近造成一定程度的資源浪費(fèi)，F(xiàn)CHSW 和CV-Voronoi重點(diǎn)探索非線性區(qū)域，在平坦區(qū)域上的點(diǎn)較少，且FCHSW 預(yù)示響應(yīng)與真實(shí)響應(yīng)曲線基本重合，相比One-Shot、CV-Voronoi 和APSFC，F(xiàn)CHSW具有更好的全局近似效果。

圖4 一維測(cè)試算例結(jié)果對(duì)比Fig.4 Results comparison of one-dimensional numerical examples

表1 展示了9 個(gè)二維～六維測(cè)試算例在One-Shot、CV-Voronoi、APSFC 和FCHSW 這4 種方法下終止樣本數(shù)的收斂均值和均方差對(duì)比結(jié)果，對(duì)應(yīng)的缺口箱線圖如圖5 所示。本文所提FCHSW 方法在9 個(gè)測(cè)試算例中均取得最優(yōu)效果。FCHSW 與基于One-Shot 的采樣相比，采樣規(guī)模平均減少14.2%，其中最大減少幅度相比ES 達(dá)到35.7%；與基于CV-Voronoi 的采樣相比，采樣規(guī)模平均減少6.8%，其中最大減少幅度相比MZ 達(dá)到21.5%；與基于APSFC 的采樣相比，采樣規(guī)模平均減少26.1%，其中最大減少幅度相比MZ 達(dá)到45.5%。此外，根據(jù)均方差對(duì)比結(jié)果，APSFC魯棒性最差，F(xiàn)CHSW 則具有很好的魯棒性。

表1 FCHSW 與One-Shot、CV-Voronoi、APSFC 測(cè) 試結(jié)果對(duì)比Table 1 Results comparison of One-Shot， CV-Voronoi， APSFC and FCHSW

圖5 FCHSW 與One-Shot、CV-Voronoi、APSFC 測(cè)試結(jié)果缺口箱線圖對(duì)比Fig.5 Notched boxplots of calculations by One-Shot， CV-Voronoi， APSFC and FCHSW

為進(jìn)一步驗(yàn)證FCHSW 方法的有效性，選取具有局部高度非線性和大面積平坦分布特征的ES、PK 和MZ 測(cè)試函數(shù)為例展示了自適應(yīng)采樣詳細(xì)結(jié)果，其中ES 和PK 分別為單峰和多峰函數(shù)，MZ 除了擁有多個(gè)局部極值還具有狹長(zhǎng)波谷帶。圖6 展示了ES、PK、MZ 計(jì)算收斂迭代圖，與One-Shot、CV-Voronoi 和APSFC 相 比，F(xiàn)CHSW能夠迅速收斂到期望誤差閾值。圖7 展示了FCHSW 分 別 針 對(duì)ES、PK 和MZ 算 例 的 新 增 樣本結(jié)果，其中正方形表示初始樣本，圓表示新增樣本。圖7（a）和圖7（b）展示了樣本在中部非線性區(qū)域采樣密集在邊緣平坦區(qū)域采樣稀疏，但是在整個(gè)平坦區(qū)域樣本分布較為均勻。圖7（c）展示了FCHSW 在非線性波谷帶采樣密集，尤其是在2 個(gè)局部極值附近更為明顯，在平坦波峰帶區(qū)域采樣相對(duì)均勻稀疏。結(jié)果表明，F(xiàn)CHSW 方法對(duì)形態(tài)良好的區(qū)域?qū)崿F(xiàn)稀疏、均勻且全覆蓋的采樣；對(duì)多極值、強(qiáng)非線性區(qū)域?qū)崿F(xiàn)重點(diǎn)密集采樣，實(shí)現(xiàn)了全局探索和局部開(kāi)發(fā)的高效平衡。此外，針對(duì)不同的模型特征以及高維問(wèn)題，F(xiàn)CHSW 均表現(xiàn)出良好的性能，從而驗(yàn)證了方法的有效性和廣泛適用性。

圖6 典型算例下FCHSW 與One-Shot、CV-Voronoi、APSFC 的迭代收斂圖Fig.6 Convergence histories of One-Shot， CV-Voronoi， APSFC and FCHSW for typical numerical cases

圖7 FCHSW 所得典型算例的樣本分布圖Fig.7 Sample distribution of typical test cases using FCHSW

針對(duì)計(jì)算效率，通過(guò)代理模型建模次數(shù)與計(jì)算耗時(shí)2 個(gè)方面對(duì)方法進(jìn)行評(píng)價(jià)。One-Shot 方法作為一次性采樣方法理論上只需要進(jìn)行一次模型訓(xùn)練，但是面對(duì)未知問(wèn)題，通常事先難以估計(jì)達(dá)到期望精度所需樣本數(shù)，故實(shí)際操作中需要多次嘗試，且由于此類采樣方法通常追求分布均勻性與投影正交性，不能針對(duì)相關(guān)物理問(wèn)題進(jìn)行重點(diǎn)區(qū)域的針對(duì)性采樣，樣本利用效率無(wú)法保證，下文討論計(jì)算效率時(shí)不再考慮。

表2 展示了3 種自適應(yīng)采樣方法在測(cè)試算例所需的建模次數(shù)及訓(xùn)練時(shí)間的統(tǒng)計(jì)結(jié)果。在建模次數(shù)方面，本文FCHSW 方法在所有測(cè)試算例中得到最小的均值與方差。在計(jì)算時(shí)間方面，在多數(shù)算例中取得了最佳或相當(dāng)?shù)谋憩F(xiàn)，特別是在存在多平臺(tái)死區(qū)的復(fù)雜算例MZ 中，CV-Voronoi 和APSFC 平均耗時(shí)分別是本文所提FCHSW 方法的4.7 倍和11.7 倍。以上結(jié)果表明本文所提方法在獲得指定精度的代理模型過(guò)程中，不僅可以減少對(duì)樣本數(shù)量的需求，同時(shí)也具有良好的計(jì)算效率。

表2 計(jì)算過(guò)程中產(chǎn)生的建模次數(shù)及所需時(shí)間Table 2 Number of modeling times and time required in calculation process

3.2 進(jìn)一步討論

為了驗(yàn)證FCHSW 方法的有效性，本節(jié)對(duì)影響域的選擇及每輪新增樣本數(shù)nadd和波動(dòng)系數(shù)μ對(duì)結(jié)果的影響進(jìn)行進(jìn)一步討論。

作為對(duì)比，在其他參數(shù)設(shè)置相同的情況下，圖8 展示了FCHSW 方法的2 種簡(jiǎn)化形式。其中圖8（a）展示了不考慮誤差評(píng)分系數(shù)的拒絕域，此時(shí)各樣本的誤差評(píng)分系數(shù)λi=1(i=1，2，…，ns)，相應(yīng)的波動(dòng)系數(shù)μ=1。圖8（b）展示了影響域?yàn)槌⒎襟w的情況，超立方體影響域?yàn)榫垲愃鶎贅颖驹谒芯S度空間的跨度，定義為

圖8 FCHSW 方法的2 種簡(jiǎn)化形式Fig.8 Two simplified forms of FCHSW method

3.2.1 關(guān)于影響域的討論

針對(duì)超立方體采樣區(qū)域使得新增樣本趨于在已有樣本內(nèi)空間，不利于全局探索的局限性，所提FCHSW 方法采用超球域，可在給定誤差分布特性下拓展不同維度的采樣范圍，提高新增樣本的潛在可探索性，進(jìn)而改善RPM 的全局預(yù)示精度。

分別針對(duì)超球體和超立方體兩種影響域計(jì)算方式，對(duì)比測(cè)試算例集的自適應(yīng)采樣結(jié)果，如表3 所示。超球體影響域在9 個(gè)算例中均表現(xiàn)出不同程度的優(yōu)異性，其采樣規(guī)模相比超立方體影響域平均減少6.7%，其中最大性能提升在PK 算例中取得，達(dá)到了13.6%。

表3 影響域?qū)Ρ冉Y(jié)果Table 3 Results comparison of influential domain

3.2.2 每輪新增樣本數(shù)對(duì)結(jié)果的影響

為考察每輪新增樣本數(shù)對(duì)方法有效性影響，選取nadd=nv，2nv，3nv，4nv，5nv，結(jié)果對(duì)比如 圖9所示，具體計(jì)算結(jié)果如表4 所示。

表4 每輪新增樣本數(shù)的影響Table 4 Summary of test results about number of new samples per round

圖9 每輪新增樣本數(shù)對(duì)比結(jié)果Fig.9 Results comparison of number of new samples per round

FCHSW 整體上呈現(xiàn)隨著每輪新增樣本數(shù)增加，終止樣本數(shù)緩慢增加的趨勢(shì)，nadd=nv在9個(gè)測(cè)試算例中均取得最佳性能，且當(dāng)nadd分別取2nv、3nv、4nv、5nv與nadd=nv相比，采樣規(guī)模平均增加幅度分別達(dá)到5.0%、6.9%、10.1%、13.6%，F(xiàn)CHSW 針對(duì)每輪新增樣本數(shù)具有較強(qiáng)魯棒性。另外，隨著nadd增大，所需迭代次數(shù)和交叉驗(yàn)證計(jì)算成本成倍降低。為有效兼顧計(jì)算效率和采樣質(zhì)量，本文建議nadd取nv～3nv進(jìn)行計(jì)算。

3.2.3 波動(dòng)系數(shù)對(duì)結(jié)果的影響

為分析波動(dòng)系數(shù)對(duì)方法有效性的影響，依次選取μ=1，1.2，1.5，2，2.5，3，5，計(jì)算結(jié)果如表5所示，結(jié)果對(duì)比如圖10 所示。

圖10 波動(dòng)系數(shù)對(duì)比結(jié)果Fig.10 Results comparison of fluctuation coefficient

表5 波動(dòng)系數(shù)的影響Table 5 Summary of test results about fluctuation coefficient

FCHSW 整體上隨著波動(dòng)系數(shù)增大，終止樣本數(shù)先減小后增大。當(dāng)μ=1.2 時(shí)，F(xiàn)CHSW 針對(duì)TG、AK、GP 算例具有最佳性能，相較于不考慮誤差評(píng)分系數(shù)分別有1.1%、0.6%、1.4%的性 能 提 升；當(dāng)μ=2 時(shí)，F(xiàn)CHSW 針 對(duì)ES、PK、MZ、HT、DP 算例具有最佳性能，相較于不考慮誤差評(píng)分系數(shù)分別有31.7%、6.7%、2.2%、10.2%、2.8% 的 性 能 提 升；當(dāng)μ=2.5 時(shí)，F(xiàn)CHSW 針對(duì)CH 算例具有最佳性能，相較于不考慮誤差評(píng)分系數(shù)有8.2%的性能提升。事實(shí)上，增大波動(dòng)系數(shù)，意味著允許更大的樣本聚集度，進(jìn)而減小全局探索同時(shí)增大了局部開(kāi)發(fā)能力，應(yīng)選取一個(gè)合適的波動(dòng)系數(shù)值來(lái)有效權(quán)衡二者以達(dá)到縮減采樣規(guī)模的目的。綜合考慮，本文建議μ取1.2～2.5 進(jìn)行計(jì)算。

4 高超聲速飛行器氣動(dòng)熱全局快速預(yù)示

本節(jié)基于FCHSW 方法進(jìn)行高超聲速飛行器氣動(dòng)熱RPM 的高效構(gòu)建，實(shí)現(xiàn)高精度高效率低成本獲取氣動(dòng)熱的目的，對(duì)所提方法的工程有效性與實(shí)用性進(jìn)行驗(yàn)證。

4.1 氣動(dòng)熱物理建模與高精度計(jì)算

選取類HTV-2 高超聲速飛行器作為氣動(dòng)外形。HTV-2 氣動(dòng)布局具有乘波體與升力體相結(jié)合的特征，外形扁平，具有尖前緣、大后掠[33］，幾何尺寸參考文獻(xiàn)[34］，如圖11 所示。

圖11 類HTV-2 飛行器幾何模型Fig.11 Geometry of HTV-2 type shape

飛行狀態(tài)參數(shù)為馬赫數(shù)、攻角和飛行高度，參數(shù)空間如表6 所示。給定側(cè)滑角0°，設(shè)置壁面為等溫壁面300 K 計(jì)算冷壁熱流。遠(yuǎn)場(chǎng)邊界由飛行參數(shù)控制，為節(jié)省計(jì)算資源，根據(jù)流場(chǎng)對(duì)稱特性，采用對(duì)稱邊界進(jìn)行求解，計(jì)算域及邊界條件如圖12 所示。采用344.6 萬(wàn)網(wǎng)格作為最終計(jì)算網(wǎng)格，網(wǎng)格劃分如圖13 所示。

表6 飛行狀態(tài)參數(shù)空間Table 6 Design space of flight parameters

圖12 計(jì)算域及邊界條件Fig.12 Computational domain and boundary conditions

圖13 飛行器網(wǎng)格劃分示意Fig.13 Schematic of vehicle meshing

選擇定常求解，黏性控制為層流，本文不考慮真實(shí)氣體效應(yīng)等。利用有限體積法空間離散求解質(zhì)量守恒方程、Navier-Stokes 方程和能量方程，采用Ausmdv 空間離散格式和minmod 限制器格式，MUSCL 插值精度為2 階迎風(fēng)。時(shí)間離散由庫(kù)朗數(shù)控制，起始庫(kù)朗數(shù)、結(jié)束庫(kù)朗數(shù)、變庫(kù)朗數(shù)步數(shù)分別設(shè)為0.01、10、500。

4.2 氣動(dòng)熱RPM 構(gòu)建與更新

利用OLHD 選取15 個(gè)樣本作為初始訓(xùn)練集，選取10 個(gè)測(cè)試樣本用于模型精度校驗(yàn)。在飛行器表面每個(gè)離散化的網(wǎng)格節(jié)點(diǎn)單獨(dú)構(gòu)造一個(gè)RPM，然后所有網(wǎng)格節(jié)點(diǎn)的氣動(dòng)熱RPM 協(xié)同預(yù)示飛行器表面熱流。

根據(jù)預(yù)示誤差分布特性，利用所提FCHSW方法獲取新增訓(xùn)練樣本并更新RPM，直至模型預(yù)示精度滿足要求，其中FCHSW 相關(guān)參數(shù)設(shè)置為：nadd=2nv和μ=2，迭代終止條件為測(cè)試集平均相對(duì)誤差e?≤5%。

4.3 結(jié)果分析

圖14 為相對(duì)誤差迭代收斂過(guò)程，由圖可知每輪加點(diǎn)后模型精度得到顯著改善。7 次迭代后，F(xiàn)CHSW 方法通過(guò)51 個(gè)訓(xùn)練樣本將RPM 的平均誤差從103.3%降低到4.8%，作為對(duì)比，相同樣本數(shù)下APSFC 方法的誤差為8.3%；而APSFC需要通過(guò)10 次迭代后的69 個(gè)訓(xùn)練樣本達(dá)到4.4% 的精度。 為滿足RPM 的精度要求，F(xiàn)CHSW 所需樣本是APSFC 的73.9%，使得高精度模型計(jì)算成本明顯減小，提高了構(gòu)建RPM的效率。

此外，使用最終得到的RPM 對(duì)飛行器表面熱流進(jìn)行預(yù)示，單工況耗時(shí)小于3 s，而采用CFD進(jìn)行高精度計(jì)算，單個(gè)工況約需240 核時(shí)，表明本文基于FCHSW 建立的氣動(dòng)熱RPM 與高精度CFD 計(jì)算相比體現(xiàn)出明顯的計(jì)算效率優(yōu)勢(shì)。

表7 給出了RPM 在10 個(gè)測(cè)試工況下精度對(duì)比，駐點(diǎn)熱流相對(duì)誤差在5%以內(nèi)，整體平均相對(duì)誤差在10%以內(nèi)，結(jié)果表明基于FCHSW 方法構(gòu)建的氣動(dòng)熱RPM 在有限高精度計(jì)算成本下通過(guò)縮減采樣規(guī)模達(dá)到了高精高效預(yù)示氣動(dòng)熱的目的。

表7 RPM 在10 個(gè)測(cè)試集上誤差對(duì)比Table 7 Results comparison of 10 test cases using RPM

考慮到駐點(diǎn)預(yù)示誤差和整體平均誤差的最大值均在工況Test3 取得，不失一般性，對(duì)工況Test3 進(jìn)行分析。圖15 給出了工況Test3 的高精度CFD 計(jì)算熱流與預(yù)示熱流qw及誤差分布對(duì)比，可以看到在飛行器的頭部出現(xiàn)高熱流區(qū)，在側(cè)翼（前緣）出現(xiàn)狹長(zhǎng)熱流分布帶。RPM 能準(zhǔn)確預(yù)示熱流分布，且與高精度計(jì)算結(jié)果保持一致。圖16（a）和圖16（b）分別為RPM 與CFD在工況Test3 的迎風(fēng)中心線和前緣中心線熱流對(duì)比，RPM 在迎風(fēng)中心線和前緣中心線平均預(yù)示誤差分別為12.6% 和7.8%，能準(zhǔn)確實(shí)現(xiàn)高熱流區(qū)重建，提供系統(tǒng)的主要物理特征。相對(duì)誤差較大的部分則主要集中在飛行器后半段靠近對(duì)稱軸附近的低熱流區(qū)。事實(shí)上，對(duì)于工況Test3，在所有相對(duì)誤差大于15%的區(qū)域，最大熱流為75.3 kW/m2；進(jìn)一步對(duì)于相對(duì)誤差大于20%的區(qū)域，最大熱流僅為8.4 kW/m2，比駐點(diǎn)熱流小兩個(gè)數(shù)量級(jí)，因而可對(duì)認(rèn)識(shí)飛行器氣動(dòng)熱變化規(guī)律、支持熱防護(hù)系統(tǒng)設(shè)計(jì)提供有意義的結(jié)果。

圖15 Test3 工況下RPM 和CFD 表面熱流分布及誤差對(duì)比Fig.15 Comparison of surface heat flux and error contours between RPM and CFD （Test3）

圖16 Test3 工況下RPM 和CFD 迎風(fēng)中心線和前緣中心線熱流對(duì)比Fig.16 Comparison of heat flux along center 1ine on windward wall and leading edge between RPM and CFD （Test3）

5 結(jié) 論

1）針對(duì)現(xiàn)有自適應(yīng)采樣方法更新RPM 難以兼顧計(jì)算效率和計(jì)算精度問(wèn)題，提出一種基于模糊聚類超球體影響域和誤差評(píng)分系數(shù)加權(quán)拒絕域的批量自適應(yīng)采樣方法。數(shù)值測(cè)試算例表明，所提方法與One-Shot、APSFC、CV-Voronoi相比能有效減小RPM 構(gòu)建所需樣本規(guī)模，高效提升預(yù)示精度。

2）基于模糊聚類的超球體影響域可有效改善全局探索能力，提升采樣質(zhì)量；對(duì)FCHSW 的2個(gè)主要參數(shù)（每輪新增樣本數(shù)和波動(dòng)系數(shù)）進(jìn)行了對(duì)比分析，結(jié)果表明采用本文建議值可兼顧計(jì)算效率和采樣質(zhì)量，有效應(yīng)對(duì)不同維度、不同類型問(wèn)題，具有較強(qiáng)的適應(yīng)性。

3）類HTV-2 高超聲速飛行器氣動(dòng)熱快速預(yù)示實(shí)例表明，基于FCHSW 構(gòu)建的氣動(dòng)熱RPM 在有限采樣規(guī)模下能夠獲得良好的全局預(yù)示精度，其中駐點(diǎn)預(yù)示誤差和整體平均預(yù)示誤差分別在5%和10%以內(nèi)。且所提FCHSW 方法相比APSFC 方法節(jié)約了26.1%的高精度數(shù)值計(jì)算成本，有效提升了氣動(dòng)熱RPM 的構(gòu)建效率。

附錄A：數(shù)值測(cè)試函數(shù)

1）一維測(cè)試算例（多峰，高度非線性），如圖 A1 所示。

圖 A1 一維測(cè)試算例Fig.A1 1D test function

2）Trigonometric function（TG，二維，多峰，均勻分布），如圖 A2 所示。

3）Easom function （ES，二維，單峰，非線性），如圖 A3 所示。

圖 A2 TG 測(cè)試算例Fig.A2 Test function TG

圖 A3 ES 測(cè)試算例Fig.A3 Test function ES

4）Peaks function（PK，二維，多峰，高度非線性）， 如圖 A4 所示。

圖 A4 PK 測(cè)試算例Fig.A4 Test function PK

5）Ackley’s path function（AK，二維，多峰，高度非線性），如圖 A5 所示。

圖 A5 AK 測(cè)試算例Fig.A5 Test function AK

式中：a=20，b=0.2，c=2π。

6）Goldstein-Price’s function（GP，二維，多峰，尖峰）， 如圖 A6 所示。

圖 A6 GP 測(cè)試算例Fig.A6 Test function GP

7）Michalewicz’s function（MZ，二維，狹長(zhǎng)峰谷帶），如圖 A7 所示。

圖 A7 MZ 測(cè)試算例Fig.A7 Test function MZ

8）Hartman function（HT，三維，非線性）

式中：

9） Colville-Himmelblau function （CH，四維）

10）Dixon-Price function（DP，六 維，均 勻分布）