基于數(shù)據(jù)同化的氣動(dòng)壓力稀疏重構(gòu)方法

黃俊，郭雨欣，冀晶晶，黃永安

華中科技大學(xué) 機(jī)械科學(xué)與工程學(xué)院，武漢 430074

0 引言

風(fēng)洞作為空氣動(dòng)力學(xué)研究的重要地面實(shí)驗(yàn)設(shè)備，在先進(jìn)飛行器研制和基礎(chǔ)空氣動(dòng)力學(xué)問(wèn)題研究中發(fā)揮了不可替代的作用。實(shí)驗(yàn)流體力學(xué)（EFD）著眼于風(fēng)洞實(shí)驗(yàn)數(shù)據(jù)的采集、分析和處理，是風(fēng)洞系統(tǒng)工程中不可分割的部分。作為風(fēng)洞技術(shù)測(cè)量目標(biāo)之一，表面壓力的測(cè)量直接關(guān)系到飛行器的升力系數(shù)與姿態(tài)控制，同時(shí)還可以通過(guò)壓力分布判斷飛行器表面的分離、轉(zhuǎn)捩等信息?，F(xiàn)有的壓敏材料（PSP）利用非接觸方式獲得連續(xù)、大范圍的模型表面壓力分布[1]，是風(fēng)洞實(shí)驗(yàn)中有效的測(cè)壓方法，然而，基于光學(xué)原理的PSP 存在觀測(cè)死角（如飛機(jī)短艙內(nèi)表面）問(wèn)題。傳統(tǒng)的在模型表面布設(shè)測(cè)壓孔的方法技術(shù)成熟度高、測(cè)量準(zhǔn)確性好、可按需布置[2]，最新發(fā)展的柔性智能蒙皮可以粘貼于模型內(nèi)外表面，實(shí)現(xiàn)多物理量同步測(cè)量[3-4]。但測(cè)壓孔和柔性傳感陣列往往數(shù)量有限，需發(fā)展基于稀疏測(cè)點(diǎn)的場(chǎng)重構(gòu)技術(shù)，實(shí)現(xiàn)壁面繞流及空間流場(chǎng)的感知。

流場(chǎng)數(shù)據(jù)具有多維度且復(fù)雜的特點(diǎn)，而測(cè)量值來(lái)自壁面稀疏測(cè)點(diǎn)，如何構(gòu)建有限數(shù)據(jù)與全空間流場(chǎng)的映射關(guān)系是流場(chǎng)重構(gòu)的難點(diǎn)所在。國(guó)內(nèi)外學(xué)者從流場(chǎng)物理模型與數(shù)據(jù)驅(qū)動(dòng)的角度進(jìn)行了理論嘗試，Vamsi Krishna 等[5]基于快速畸變理論和Taylor假設(shè)建立雙向加權(quán)模型，基于先驗(yàn)速度場(chǎng)對(duì)湍流演化進(jìn)行重構(gòu)。Callaham 等[6]探索了基于二維流場(chǎng)稀疏表示和數(shù)據(jù)庫(kù)的流場(chǎng)重構(gòu)法，嘗試了數(shù)據(jù)驅(qū)動(dòng)的流場(chǎng)重構(gòu)技術(shù)，并討論了需要的隨機(jī)測(cè)點(diǎn)最小數(shù)量。Sun 等[7]結(jié)合Bayesian 神經(jīng)網(wǎng)絡(luò)與先驗(yàn)物理模型重構(gòu)了二維鐘形與T 形流場(chǎng)，并對(duì)結(jié)果不確定性進(jìn)行估計(jì)。李靜等[8]采用本征正交分解方法，以較少階模態(tài)高精度再現(xiàn)非定常圓柱繞流的完整流場(chǎng)。這些重構(gòu)方法大多局限在簡(jiǎn)單幾何和穩(wěn)態(tài)流場(chǎng)重構(gòu)，與實(shí)際工程應(yīng)用之間仍存在一定距離。得益于計(jì)算流體力學(xué)（CFD）的發(fā)展，通過(guò)數(shù)值方法求解Navier–Stokes 方程[9]，能夠得到豐富的流場(chǎng)先驗(yàn)知識(shí)，為全流場(chǎng)重構(gòu)提供良好的條件。然而，CFD 的計(jì)算過(guò)程通常不能充分考慮風(fēng)洞實(shí)際流動(dòng)過(guò)程中的不確定性[10]，導(dǎo)致確定性的計(jì)算機(jī)仿真難以復(fù)現(xiàn)風(fēng)洞實(shí)驗(yàn)對(duì)物理參數(shù)的影響，從而無(wú)法實(shí)現(xiàn)與風(fēng)洞實(shí)測(cè)數(shù)據(jù)的匹配和融合。因此，結(jié)合CFD 與實(shí)測(cè)稀疏數(shù)據(jù)對(duì)風(fēng)洞的真實(shí)流場(chǎng)進(jìn)行重構(gòu)極富挑戰(zhàn)。

數(shù)據(jù)同化（Data Assimilation，DA）方法能夠結(jié)合估計(jì)值和測(cè)量值進(jìn)行預(yù)測(cè)[11]，可以作為聯(lián)系CFD理論先驗(yàn)信息和EFD 風(fēng)洞實(shí)測(cè)數(shù)據(jù)的橋梁，其在流場(chǎng)重構(gòu)中已有發(fā)展。Chandramouli 等[12]提出一種變分?jǐn)?shù)據(jù)同化方法，通過(guò)圓柱繞流在過(guò)渡狀態(tài)下的三維湍流尾流案例進(jìn)行方案驗(yàn)證，與PIV 實(shí)驗(yàn)數(shù)據(jù)對(duì)比，驗(yàn)證了變分?jǐn)?shù)據(jù)同化在不可壓縮湍流重構(gòu)中的有效性。Belligoli 等[13]通過(guò)數(shù)據(jù)同化方法實(shí)現(xiàn)了二維翼型的迎角及馬赫數(shù)的修正，降低了實(shí)驗(yàn)測(cè)量值與RANS 模擬計(jì)算值之間的誤差。集合變換卡爾曼濾波（Ensemble Transform Kalman Filter，ETKF）是卡爾曼濾波（Kalman Filter，KF）針對(duì)非線性系統(tǒng)模型的擴(kuò)展方法[14]，在復(fù)雜系統(tǒng)中有良好的數(shù)據(jù)同化效果。ETKF 以蒙特卡洛方法為依托，以貝葉斯原理為核心，利用稀疏測(cè)量數(shù)據(jù)所包含的信息對(duì)物理場(chǎng)先驗(yàn)預(yù)報(bào)數(shù)據(jù)進(jìn)行濾波，借此給出對(duì)不確定參數(shù)的連續(xù)后驗(yàn)估計(jì)[15]。該方法源自地球物理領(lǐng)域，是基于海面離散探測(cè)數(shù)據(jù)對(duì)海洋與大氣的真實(shí)演化狀態(tài)進(jìn)行連續(xù)修正的重要工具[16]。近年來(lái)，國(guó)內(nèi)外學(xué)者將該方法引入空氣動(dòng)力學(xué)領(lǐng)域，美國(guó)加州大學(xué)通過(guò)無(wú)限長(zhǎng)薄板周圍隨機(jī)擾動(dòng)無(wú)黏渦流模擬，印證了ETKF 在流場(chǎng)數(shù)據(jù)同化中的有效性[17]。日本宇航研究開發(fā)機(jī)構(gòu)（JAXA）應(yīng)用ETKF 對(duì)風(fēng)洞中存在不確定性的馬赫數(shù)、迎角、湍流黏度進(jìn)行估計(jì)[18]，成功重構(gòu)剛性機(jī)翼繞流場(chǎng)，證實(shí)了該方法對(duì)風(fēng)洞中復(fù)雜流場(chǎng)的估計(jì)能力。國(guó)內(nèi)上海交通大學(xué)在預(yù)測(cè)直升機(jī)轉(zhuǎn)子三維流場(chǎng)特征時(shí)，采用了ETKF 優(yōu)化剪切應(yīng)力傳輸模型常數(shù)，為修正逆壓梯度下邊界層流動(dòng)分離提供了參考[19]。

本文面向風(fēng)洞實(shí)驗(yàn)的壓力測(cè)量應(yīng)用，使用集合變換卡爾曼濾波方法，以二維翼型RAE 2822 和二維對(duì)稱翼型NACA 0012 為研究對(duì)象，進(jìn)行機(jī)翼表面氣動(dòng)壓力重構(gòu)，以達(dá)到風(fēng)洞實(shí)驗(yàn)對(duì)全域數(shù)據(jù)重構(gòu)精度的要求。先通過(guò)CFD 計(jì)算得到全域先驗(yàn)分布，結(jié)合機(jī)翼表面有限數(shù)量的壓力測(cè)量值，對(duì)機(jī)翼的迎角及馬赫數(shù)進(jìn)行修正，重構(gòu)得到高精度的機(jī)翼周圍全域壓力場(chǎng)。再采用ETKF 進(jìn)行流場(chǎng)數(shù)據(jù)同化：一方面，可以充分利用稀疏測(cè)量數(shù)據(jù)，將這些高精度數(shù)據(jù)泛化到整個(gè)流域，使得最終展示的流場(chǎng)更加接近真實(shí)的流場(chǎng)，為空間全域感知提供可能；另一方面，ETKF 可以作為一種流場(chǎng)風(fēng)洞干擾的修正方法，基于真實(shí)測(cè)量數(shù)據(jù)，直接進(jìn)行流場(chǎng)分布的修正，獲取高精度的流場(chǎng)分布。

1 集合變換卡爾曼濾波方法

1.1 系統(tǒng)模型

在流場(chǎng)計(jì)算中，從邊界條件到物理量的分布不是簡(jiǎn)單的線性系統(tǒng)，而是離散的非線性系統(tǒng)，其中狀態(tài)變量估計(jì)可通過(guò)求解Navier–Stokes 方程得到。本文通過(guò)軟件FLUENT 對(duì)流場(chǎng)進(jìn)行數(shù)值計(jì)算，選用Spalart–Allmaras（S–A）湍流模型。采用Navier–Stokes 方程作為流動(dòng)控制方程，并將此方程作為卡爾曼濾波中的系統(tǒng)模型：

式中：W 為守恒狀態(tài)矢量，包含了密度、速度和能量；V 為控制體體積；Fc為對(duì)流通量；Fv為黏性通量；S 為控制體表面積；n 為控制體外法線方向的單位矢量；τ為時(shí)間。

1.2 狀態(tài)空間矩陣

系統(tǒng)模型的估計(jì)與實(shí)驗(yàn)測(cè)量均存在誤差，卡爾曼濾波的數(shù)據(jù)同化手段通過(guò)狀態(tài)空間模型將誤差代入系統(tǒng)當(dāng)中：

式中：下標(biāo)t 表示迭代次數(shù)，xt、yt向量分別表示系統(tǒng)模型和實(shí)驗(yàn)測(cè)量數(shù)據(jù)的狀態(tài)向量，vt、ωt向量分別表示系統(tǒng)模型的噪聲和觀測(cè)噪聲，假定其符合高斯分布。非線性算子F 是從第t-1次迭代到第t 次迭代的映射，在數(shù)據(jù)同化中，由系統(tǒng)模型計(jì)算得到。H矩陣是將系統(tǒng)估計(jì)的數(shù)據(jù)矩陣投影到實(shí)驗(yàn)觀測(cè)數(shù)據(jù)的投影矩陣。

向量xt包含了迎角α、馬赫數(shù)Ma 和所有（n 個(gè)）計(jì)算網(wǎng)格節(jié)點(diǎn)上的密度ρ、笛卡爾速度分量u、v 和壓力p，其維數(shù)l=4n+2，向量xt的表達(dá)式為：

向量yt由實(shí)驗(yàn)測(cè)點(diǎn)對(duì)應(yīng)位置的測(cè)量值構(gòu)成。在稀疏壓力重構(gòu)中，選取實(shí)測(cè)壓力p，向量yt維數(shù)與實(shí)驗(yàn)測(cè)點(diǎn)的數(shù)量（m 個(gè)）一致：

1.3 集合變換卡爾曼濾波方法及實(shí)現(xiàn)流程

式中：H為式（3）中的投影矩陣，在翼型網(wǎng)格和實(shí)驗(yàn)?zāi)Ｐ偷膶?duì)應(yīng)位置設(shè)置為1，其余位置設(shè)置為0。

圖1 ETKF 流程圖Fig.1 The flowchart of ETKF

在各物理量符合收斂條件判定準(zhǔn)則之后，將最后一次迭代的集合元素均值作為整個(gè)數(shù)據(jù)同化過(guò)程的來(lái)流條件修正值，代入FLUENT 進(jìn)行重新計(jì)算，得到修正后的流場(chǎng)分布。

2 實(shí)驗(yàn)過(guò)程和結(jié)果

RAE 2822 是一個(gè)二維跨聲速湍流流動(dòng)的經(jīng)典翼型，被許多國(guó)外的項(xiàng)目合作組（如EUROVAL）選作經(jīng)典算例[21]；NACA 0012 作為經(jīng)典的對(duì)稱翼型，有實(shí)驗(yàn)數(shù)據(jù)可用來(lái)驗(yàn)證其數(shù)值模擬的準(zhǔn)確性，故本文實(shí)驗(yàn)使用二維翼型RAE 2822 和NACA 0012 進(jìn)行驗(yàn)證。RAE 2822 翼型實(shí)驗(yàn)主要用來(lái)展現(xiàn)重構(gòu)的收斂過(guò)程，ETKF 對(duì)迎角、馬赫數(shù)的影響，以及與其他方法在精度上的對(duì)比；NACA 0012 翼型實(shí)驗(yàn)則面向風(fēng)洞測(cè)量中基于有限測(cè)點(diǎn)重構(gòu)的需求，用以探究影響ETKF 重構(gòu)精度的因素。

2.1 RAE 2822 翼型壓力重構(gòu)

2.1.1 網(wǎng)格和計(jì)算條件設(shè)置

使用ICEM 軟件對(duì)RAE 2822 進(jìn)行網(wǎng)格劃分（圖2），采用四邊形結(jié)構(gòu)化網(wǎng)格，機(jī)翼模型弦長(zhǎng)c=1 m，計(jì)算域在弦長(zhǎng)方向（即x 方向）為40 m，弦長(zhǎng)垂直方向（即y 方向）為40 m，網(wǎng)格質(zhì)量如表1所示。本文笛卡爾坐標(biāo)系的原點(diǎn)位于機(jī)翼最前緣點(diǎn)。將FLUENT 軟件的計(jì)算模型設(shè)置為S–A 湍流模型，此模型對(duì)跨聲速流動(dòng)的求解有較好的效果，且計(jì)算成本相對(duì)較低。自由流空氣設(shè)置為可壓縮空氣，計(jì)算方式為穩(wěn)態(tài)計(jì)算。實(shí)驗(yàn)數(shù)據(jù)來(lái)源為Cook 等[22]的案例6 中RAE 2822 跨聲速實(shí)驗(yàn)壓力數(shù)據(jù)，其實(shí)驗(yàn)工況為：Ma=0.725，α=2.92°，Re=6.2×106。共有103 個(gè)測(cè)量點(diǎn)，其中一部分測(cè)點(diǎn)與機(jī)翼網(wǎng)格點(diǎn)難以匹配，會(huì)給數(shù)據(jù)同化帶來(lái)較大誤差，將其去除后剩余的壓力測(cè)點(diǎn)個(gè)數(shù)為75（即m=75）。

表1 RAE 2822 網(wǎng)格質(zhì)量及節(jié)點(diǎn)數(shù)量Table 1 Quality and node number of RAE 2822 mesh

圖2 RAE 2822 翼型網(wǎng)格Fig.2 Mesh of RAE 2822 airfoil

2.1.2 RAE 2822 翼型ETKF 數(shù)據(jù)同化結(jié)果

2.1.2.1 迎角及馬赫數(shù)

圖3 集合成員迎角在ETKF 同化前后的分布Fig.3 Distribution of ensemble angel of attack before and after ETKF assimilation

圖4 集合成員馬赫數(shù)在ETKF 前后的分布Fig.4 Distribution of ensemble Mach number before and after ETKF assimilation

圖5 集合的迎角及馬赫數(shù)均值隨迭代過(guò)程的變化Fig.5 Mean of angle of attack and Mach number changes with iterations

圖6 為迭代過(guò)程中所有網(wǎng)格點(diǎn)的密度（ρ）、速度分量（u、v）和壓力系數(shù)（Cp）的均方誤差MSE 變化曲線。從圖中可以看出，在進(jìn)行了第3 次迭代之后，所有物理量的均方誤差已經(jīng)小于第一次迭代的1%，這表明在第3 次迭代時(shí)ETKF 就已經(jīng)滿足收斂要求，繼續(xù)迭代發(fā)現(xiàn)均方誤差可以進(jìn)一步降低，于是迭代繼續(xù)進(jìn)行。

圖6 ρ、u、v、Cp 的MSE 隨迭代過(guò)程的變化Fig.6 MSE of ρ,u,v and Cp changes with iterations

隨著迭代的進(jìn)行，集合中所有成員的邊界條件都在產(chǎn)生變化并且更加集中（圖3 和4）。計(jì)算集合成員的方差如表2所示，可以看到在數(shù)據(jù)同化過(guò)程中，隨機(jī)抽樣產(chǎn)生的所有集合成員在向更為準(zhǔn)確的修正值變化，且馬赫數(shù)的集中效果要比迎角的集中效果更好。最終，集合成員迎角及馬赫數(shù)均值在第10 次迭代之后分別收斂為2.434°和0.732 8（圖5）。

表2 ETKF 前后集合成員迎角、馬赫數(shù)的均值及方差Table 2 Mean and variance of ensemble angle of attack and Mach number before and after ETKF

2.1.2.2 壓力系數(shù)、升力系數(shù)和力矩系數(shù)

基于圖2 網(wǎng)格，將收斂后集合成員的迎角及馬赫數(shù)均值（表3）作為邊界條件代入FLUENT 進(jìn)行重新計(jì)算，得到ETKF 同化后的壓力系數(shù)Cp曲線分布，此結(jié)果即為ETKF 數(shù)據(jù)同化的最終結(jié)果。將ETKF 同化結(jié)果與線性理論修正值進(jìn)行對(duì)比，結(jié)果如圖7 和8所示。

表3 RAE 2822 case 6 邊界條件Table 3 Boundary condition of RAE 2822 case 6

圖7 ETKF 同化后的壓力系數(shù)曲線與線性理論修正曲線對(duì)比Fig.7 Comparison of the pressure coefficient curve after assimilation of ETKF with the linear theory correction

圖8 激波位置的壓力系數(shù)對(duì)比Fig.8 Comparison of pressure coefficient at shock position

從圖7 和8 中可以看出：與線性理論修正后的計(jì)算結(jié)果相比，ETKF 同化后的壓力系數(shù)曲線更加靠近實(shí)驗(yàn)數(shù)據(jù)；線性理論的迎角修正值過(guò)大，導(dǎo)致其計(jì)算誤差與ETKF 相比更大，尤其是在激波位置。從數(shù)值上分析（式（25））可得，ETKF 同化的壓力系數(shù)平均相對(duì)誤差e 比線性理論修正后降低了約3%。e 的計(jì)算方法如下：

式中：j 為實(shí)驗(yàn)測(cè)點(diǎn)的數(shù)量，Cp,exp為實(shí)驗(yàn)測(cè)點(diǎn)的壓力系數(shù)。

表4 ETKF 同化和線性理論修正后CL 和Cm、激波位置Cp 平均相對(duì)誤差與實(shí)驗(yàn)值的對(duì)比Table 4 Comparison of CL,Cm and Cp average relative error near the shock wave position between ETKF,linear theory and experiment

2.2 NACA 0012 翼型壓力重構(gòu)

2.2.1 網(wǎng)格和計(jì)算條件設(shè)置

如圖9 和表5所示，NACA 0012 翼型采用的計(jì)算網(wǎng)格為四邊形結(jié)構(gòu)化網(wǎng)格，機(jī)翼模型弦長(zhǎng)c=0.152 4 m，計(jì)算域在弦長(zhǎng)方向?yàn)?.5 m，在垂直弦長(zhǎng)方向?yàn)? m，網(wǎng)格質(zhì)量高，F(xiàn)LUENT 計(jì)算時(shí)選用S–A湍流模型。實(shí)驗(yàn)數(shù)據(jù)來(lái)源于NASA 進(jìn)行的NACA 0012 翼型風(fēng)洞實(shí)驗(yàn)[23]，選取TEST 119 的壓力數(shù)據(jù)進(jìn)行數(shù)據(jù)同化的測(cè)試，其實(shí)驗(yàn)工況為：Ma=0.402 2，α=4.020 8°，Re=6.070 6×106。共有46 個(gè)壓力測(cè)量點(diǎn)，由于靠近機(jī)翼前緣位置的數(shù)據(jù)在同化時(shí)難以與網(wǎng)格節(jié)點(diǎn)進(jìn)行匹配，容易帶來(lái)較大誤差，故在同化時(shí)將其去除，僅使用其中44 個(gè)數(shù)據(jù)（即m≤44）。

表5 NACA0012 網(wǎng)格質(zhì)量及節(jié)點(diǎn)數(shù)量Table 5 Quality and node number of NACA 0012 mesh

圖9 NACA 0012 翼型網(wǎng)格Fig.9 Mesh of NACA 0012 airfoil

NACA 0012 翼型實(shí)驗(yàn)主要面向風(fēng)洞測(cè)量中基于有限測(cè)點(diǎn)重構(gòu)的應(yīng)用需求，對(duì)比分析基于少數(shù)測(cè)點(diǎn)進(jìn)行壓力稀疏重構(gòu)的可行性。圖10 給出了實(shí)驗(yàn)壓力測(cè)點(diǎn)的位置及序號(hào)。實(shí)驗(yàn)共進(jìn)行4 組，其中3 組實(shí)驗(yàn)分別選用上翼面（序號(hào)1～22）的不同位置6 個(gè)實(shí)驗(yàn)測(cè)點(diǎn)進(jìn)行數(shù)據(jù)同化，第4 組實(shí)驗(yàn)作為同化結(jié)果的對(duì)標(biāo)，使用所有數(shù)據(jù)測(cè)點(diǎn)進(jìn)行同化處理。

圖10 NACA 0012 實(shí)驗(yàn)測(cè)點(diǎn)位置及序號(hào)Fig.10 Test point location and serial number of NACA 0012

2.2.2 NACA 0012 翼型ETKF 數(shù)據(jù)同化結(jié)果

在數(shù)據(jù)同化之前進(jìn)行初始集合的抽樣，因?yàn)榇藢?shí)驗(yàn)無(wú)具體參考修正值，故基于實(shí)驗(yàn)原始工況進(jìn)行初始值的抽樣，并設(shè)定較大的抽樣范圍。設(shè)置集合成員數(shù)量k=40，在3.52～4.52 范圍進(jìn)行迎角抽樣，在0.372 2～0.432 2 范圍進(jìn)行馬赫數(shù)抽樣，得到初始集合的迎角均值為4.020 8°，馬赫數(shù)均值為0.402 1。采用第1 節(jié)的理論與流程進(jìn)行數(shù)據(jù)同化，同化實(shí)驗(yàn)1～4 均從此初始值開始分別進(jìn)行迭代，以對(duì)比使用不同實(shí)驗(yàn)測(cè)量數(shù)據(jù)得到的迭代效果。4 組實(shí)驗(yàn)所用的壓力測(cè)點(diǎn)及迭代至收斂所需的迭代次數(shù)見表6。

表6 4 組實(shí)驗(yàn)的壓力測(cè)點(diǎn)序號(hào)及迭代次數(shù)Table 6 The number of pressure measuring points and iteration times of 4 groups of experiments

實(shí)驗(yàn)結(jié)果如表7 和圖11所示，表7 中壓力系數(shù)平均相對(duì)誤差由式（25）計(jì)算得到，4 組同化實(shí)驗(yàn)均使用上翼面實(shí)驗(yàn)測(cè)量值進(jìn)行誤差的計(jì)算（即j=22）。如圖11所示，相比于深綠色點(diǎn)劃線表示的未同化實(shí)驗(yàn)數(shù)據(jù)，經(jīng)過(guò)ETKF 同化的壓力系數(shù)均有一定的優(yōu)化效果。同化實(shí)驗(yàn)4 使用了所有實(shí)驗(yàn)測(cè)量點(diǎn)進(jìn)行同化，其結(jié)果的誤差最小，只有1.05%，說(shuō)明ETKF 方法的同化精度與使用的數(shù)據(jù)量有關(guān)，測(cè)點(diǎn)數(shù)據(jù)越多，得到的同化效果越好；只采用6 個(gè)測(cè)點(diǎn)進(jìn)行同化的實(shí)驗(yàn)1～3 也達(dá)到了小于6%的平均相對(duì)誤差，證實(shí)了基于有限測(cè)點(diǎn)進(jìn)行壓力稀疏重構(gòu)、實(shí)現(xiàn)全域感知的可行性。且對(duì)比同化實(shí)驗(yàn)1～3 的結(jié)果可以發(fā)現(xiàn)，重構(gòu)精度與初始測(cè)點(diǎn)位置的選取密切相關(guān)，初始測(cè)點(diǎn)位置不同，相對(duì)誤差從2.42%到5.74%不盡相同，這說(shuō)明在風(fēng)洞實(shí)驗(yàn)時(shí)可以基于重構(gòu)算法進(jìn)行傳感測(cè)點(diǎn)位置的優(yōu)化布置。

表7 4 組實(shí)驗(yàn)ETKF 后的迎角、馬赫數(shù)及壓力系數(shù)平均相對(duì)誤差Table 7 Angle of attack,Mach number and pressure coefficient average relative error after ETKF in 4 groups of experiments

圖11 NACA 0012 的ETKF 同化結(jié)果與實(shí)驗(yàn)值的對(duì)比Fig.11 Comparison of ETKF assimilation results with experimental values in NACA 0012

3 結(jié)論

本文針對(duì)風(fēng)洞測(cè)量中有限測(cè)點(diǎn)進(jìn)行氣動(dòng)壓力重構(gòu)開展研究，使用集合變換卡爾曼濾波方法對(duì)二維翼型RAE 2822 和NACA 0012 的湍流流場(chǎng)進(jìn)行數(shù)據(jù)同化，結(jié)合稀疏的實(shí)驗(yàn)測(cè)量壓力數(shù)據(jù)，對(duì)兩種翼型的實(shí)驗(yàn)迎角及馬赫數(shù)進(jìn)行修正，通過(guò)修正后的邊界條件進(jìn)行重新計(jì)算得到了全域流場(chǎng)。將ETKF 同化后的機(jī)翼表面壓力系數(shù)及機(jī)翼的升力系數(shù)、力矩系數(shù)與修正前數(shù)據(jù)和線性理論修正值作對(duì)比，得到以下結(jié)論：

1）ETKF 可以通過(guò)有限的測(cè)量信息進(jìn)行流場(chǎng)的全域重構(gòu)，其重構(gòu)的精度與測(cè)量點(diǎn)的選取有直接關(guān)系，在實(shí)驗(yàn)測(cè)量位置進(jìn)行壓力系數(shù)的修正，其精度比經(jīng)典線性理論更加靠近實(shí)驗(yàn)測(cè)量值。

2）使用ETKF 修正后的迎角及馬赫數(shù)計(jì)算得到的機(jī)翼升力系數(shù)及力矩系數(shù)與實(shí)驗(yàn)測(cè)量值之間的誤差很小，這表明了基于ETKF 數(shù)據(jù)同化的有效性。

3）由于ETKF 進(jìn)行同化的數(shù)據(jù)基礎(chǔ)是實(shí)驗(yàn)測(cè)量數(shù)據(jù)，其對(duì)于流場(chǎng)的估計(jì)精度取決于此測(cè)量值的準(zhǔn)確性。

本文研究表明，作為一種用于湍流流場(chǎng)的數(shù)據(jù)同化方法，集合變換卡爾曼濾波能夠有效利用CFD數(shù)值計(jì)算數(shù)據(jù)和EFD 高精度稀疏測(cè)量數(shù)據(jù)，重構(gòu)得到更高精度的全域流場(chǎng)，適用于湍流場(chǎng)穩(wěn)態(tài)計(jì)算的分布預(yù)測(cè)。本文在給定了實(shí)驗(yàn)迎角及馬赫數(shù)的條件下進(jìn)行了壓力場(chǎng)稀疏重構(gòu)，后續(xù)可以進(jìn)一步進(jìn)行其他分布條件下（如剪應(yīng)力分布）的數(shù)據(jù)同化實(shí)驗(yàn)，進(jìn)行多數(shù)據(jù)融合下的流場(chǎng)重構(gòu)，推進(jìn)風(fēng)洞實(shí)驗(yàn)測(cè)量技術(shù)向高效化發(fā)展。