苗 楠， 劉戰(zhàn)合， 王曉璐， 王 菁， 田秋麗

(1.鄭州航空工業(yè)管理學(xué)院航空工程學(xué)院， 鄭州 450046； 2.鄭州航空工業(yè)管理學(xué)院經(jīng)濟學(xué)院， 鄭州 450046)

為提高飛機的續(xù)航性能，通常在機身或機翼下面加掛副油箱[1]。在非滿油狀態(tài)下，飛機副油箱內(nèi)的燃油將隨飛機的機動飛行而發(fā)生晃動。燃料的晃動將施加在飛機上一定作用力和力矩，同時影響飛機的質(zhì)心位置，從而影響其飛行品質(zhì)和飛行安全[2]。同時，隨著燃料的消耗，副油箱內(nèi)燃油的晃動動力學(xué)特性也將發(fā)生改變。因此，燃油晃動動力學(xué)行為的描述是航空工程中頗為關(guān)心的問題。

對液體晃動行為的分析通常從理論研究、數(shù)值模擬和實驗研究三方面展開[3-4]。相比于理論和實驗研究，數(shù)值方法在對不同形狀容器的適用性及研究成本等方面有其獨特優(yōu)勢。計算流體動力學(xué)(computational fluid dynamics, CFD)方法通過直接求解Navier-Stokes方程或其簡化方程，在時域內(nèi)對液體晃動進行數(shù)值模擬[5]。用于建模驗證的CFD軟件Flow-3D基于流體體積法(volume of fluid, VOF)編寫，在自由液面追蹤和流動仿真方面較為成熟，在液體晃動仿真中得到了廣泛應(yīng)用[6-7]。

此外，考慮到計算效率等方面的要求，航天工程中常用等效力學(xué)模型描述飛行器內(nèi)液體的晃動動力學(xué)行為[8]。常用的等效力學(xué)模型有單擺模型和彈簧-質(zhì)量模型兩種形式，如圖1[9]所示，二者在液體發(fā)生小幅晃動時等價，且等效力學(xué)模型參數(shù)之間存在簡單的轉(zhuǎn)換關(guān)系。

圖1 等效力學(xué)模型[9]Fig.1 Equivalent mechanical model[9]

然而，由于飛機油箱外部形狀及內(nèi)部結(jié)構(gòu)的復(fù)雜性，航空工程中的燃油晃動問題目前較多地采用數(shù)值方法進行研究[10-11]。在燃油晃動的等效力學(xué)模型建立方面，主要困難在于結(jié)構(gòu)復(fù)雜導(dǎo)致等效模型參數(shù)難以獲取。為提高等效建模方法對于不同油箱形狀和結(jié)構(gòu)的適用性，基于勢流理論和有限元法[5]進行液體晃動等效建模，以飛機副油箱為例建立了固定充液高度及考慮燃油消耗兩種情況下的單擺等效力學(xué)模型，并將其計算結(jié)果與CFD軟件Flow-3D的仿真結(jié)果進行對比驗證。

1 油箱及液體模型

飛機副油箱多采用流線型設(shè)計，將副油箱近似為中間圓柱段、兩端半橢球的軸對稱形狀，如圖2所示。圓柱段半徑和長度分別取為r=0.25 m、lc=2.4 m。兩端半橢球大小相同，均為繞長軸(x軸)旋轉(zhuǎn)而成的旋轉(zhuǎn)橢球體的一半，長半軸和短半軸分別取為a=1 m、b=r=0.25 m。

圖2 油箱形狀及坐標(biāo)系定義Fig.2 Tank shape and coordinate system definition

建立油箱本體坐標(biāo)系Obxyz，其中坐標(biāo)原點Ob點位于副油箱的幾何中心，x軸即縱軸為油箱對稱軸并以指向機頭方向為正，y軸即立軸以向上為正，z軸根據(jù)右手法則確定，以指向右翼為正。假設(shè)初始時刻，副油箱內(nèi)的燃油在沿y軸負向的重力作用下聚集在油箱底部，則充液高度范圍可表示為(-0.25, 0.25) m。為進行燃油晃動等效建模、且考慮到實際工程中的油箱使用情況，針對充液高度(相對于油箱中心Ob點)ys=-0.22：0.02：0.22 m的不同情形，在有限元軟件MSC.Patran中建立液體初始狀態(tài)的幾何模型，并進行網(wǎng)格劃分。建模及網(wǎng)格劃分結(jié)果以ys=±0.06 m和ys=±0.2 m共4種充液情形為例給出，如圖3所示。

圖3 液體初始形態(tài)及網(wǎng)格劃分Fig.3 Liquid initial configuration and mesh generation

同時，為對所提出的等效建模方法的有效性和準(zhǔn)確性進行驗證，在CFD軟件Flow-3D中同樣建立油箱模型。在代表流體流動區(qū)域的空腔外部建立長方體實體以包裹油箱，并以立方體單元對含油箱的長方體區(qū)域進行網(wǎng)格劃分。立方體單元各邊尺寸均設(shè)定為0.02 m，網(wǎng)格數(shù)量為240×28×28個。假設(shè)副油箱吊掛于右側(cè)機翼下，在飛機飛行過程中隨飛機發(fā)生平動且繞飛機重心轉(zhuǎn)動，因此將CFD仿真的坐標(biāo)原點平移至飛機重心，此時副油箱中心Ob點的坐標(biāo)表示為(-0.2,-0.5, 2) m。沿重力方向即-y軸方向設(shè)置充液條件，充液高度ys=0 m和ys=±0.06 m的結(jié)果如圖4所示。液體選為航空煤油RP-3[12]，密度為791.5 kg/m3，動力黏度為1.289×10-3(N·s)/m2。

圖4 不同充液情況下CFD油箱及液體模型Fig.4 CFD tank and liquid model under different fill situation

2 固定充液高度下的晃動等效建模

首先對固定充液高度或固定充液比下的副油箱內(nèi)燃油晃動進行等效建模。根據(jù)等效原則，等效系統(tǒng)與原系統(tǒng)的總質(zhì)量、質(zhì)心位置、轉(zhuǎn)動慣量對應(yīng)相等，對油箱的作用力、作用力矩和運動動能等也應(yīng)分別對應(yīng)相等。由此，可根據(jù)文獻[5]中給出的方法，基于勢流理論和有限元網(wǎng)格劃分結(jié)果計算出晃動質(zhì)量的相關(guān)參數(shù)，如x軸方向第一階晃動模態(tài)所對應(yīng)的晃動質(zhì)量即單擺擺球質(zhì)量mx，x軸方向單擺擺長lx及懸掛點相對于Ob點的高度hx，x軸晃動質(zhì)量關(guān)于Ob點的慣量張量Ix及晃動模態(tài)阻尼比ζx等。類似地，可計算得到z軸方向第一階晃動模態(tài)所對應(yīng)的晃動質(zhì)量mz，單擺擺長lz及懸掛點高度hz，關(guān)于Ob點的慣量張量Iz及阻尼比ζz等。

根據(jù)等效關(guān)系，有：

mliq=mx+mz+m0

(1)

mliqhc=mxhx+mzhz+m0h0

(2)

Jliq=Ix+Iz+I0+I′0

(3)

式中：mliq為液體總質(zhì)量；m0為靜止質(zhì)量；hc和h0分別為液體總質(zhì)心和靜止質(zhì)量塊質(zhì)心相對于Ob點的高度；Jliq代表液體等效剛體關(guān)于Ob點的慣量；I0和I′0分別為靜止質(zhì)量塊相對于其質(zhì)心的慣量張量和其質(zhì)心關(guān)于Ob點的慣量。根據(jù)副油箱的對稱性可知，x軸和z軸晃動質(zhì)量的平衡位置及靜止質(zhì)量塊的質(zhì)心均位于Oby軸上。由上述等效關(guān)系，即可求得晃動質(zhì)量和靜止質(zhì)量塊的所有參數(shù)，并計算等效系統(tǒng)對油箱的作用、力矩和質(zhì)心位置等。

針對不同充液高度，可由該方法計算得到三維副油箱內(nèi)的液體晃動的單擺等效力學(xué)模型參數(shù)。由于同樣激勵作用下，充液比在40%～60%時會施加給油箱更大的作用力和力矩，因此以充液高度ys=0 m和ys=±0.06 m為例，給出單擺等效力學(xué)模型的主要參數(shù)(表1)，并以俯仰激勵為例對等效模型進行驗證(參見4.1節(jié))。

表1 等效力學(xué)模型主要參數(shù)Table 1 Main parameters of equivalent mechanical model

3 考慮燃料消耗的晃動等效建模

由計算結(jié)果可知，隨著副油箱內(nèi)燃油的消耗，不同充液高度將對應(yīng)不同的等效力學(xué)模型參數(shù)。若考慮燃油從滿箱到接近空箱的全過程，則必須考慮燃油消耗過程中等效模型參數(shù)的變化。為研究等效模型參數(shù)的變化過程，首先計算不同充液高度上x軸和z軸的晃動質(zhì)量和擺長，然后由等效條件進一步計算得到靜止質(zhì)量等參數(shù)，主要參數(shù)隨充液高度的變化情況及其線性插值結(jié)果如圖5所示。

圖5 主要等效模型參數(shù)線性插值結(jié)果Fig.5 Linear interpolation results of main equivalent model parameters

由圖5可見，隨充液高度的不同，各個等效模型參數(shù)變化范圍均較大，但變化較為平緩，故可嘗試基于線性插值結(jié)果近似描述充液比連續(xù)變化時導(dǎo)致的等效模型參數(shù)變化[13]。具體來說，可根據(jù)燃油剩余質(zhì)量實時計算得到副油箱充液比和充液高度，并以充液高度ys為基準(zhǔn)進行線性插值，求得此時單擺等效力學(xué)模型中晃動質(zhì)量和靜止質(zhì)量的所有參數(shù)，并計算對油箱的作用力、作用力矩及液體總質(zhì)心位置。通過與CFD軟件的對應(yīng)計算結(jié)果進行對比分析，驗證考慮燃料消耗的變充液比燃油晃動等效建模方法的可行性和準(zhǔn)確性(參見4.2節(jié))。

4 算例驗證

分別以俯仰激勵和滾轉(zhuǎn)激勵為例，對固定充液比和考慮燃料消耗(變充液比)情形下飛機副油箱內(nèi)的燃油晃動等效力學(xué)模型進行算例分析。給出燃油對油箱的作用力(晃動力)、作用力矩(晃動力矩)、液體質(zhì)心位置等的等效模型預(yù)測結(jié)果，并將其與CFD軟件Flow-3D的仿真結(jié)果進行對比驗證。

4.1 固定充液比算例驗證

首先對第2節(jié)中建立的固定充液比單擺等效力學(xué)模型進行驗證。假設(shè)飛機繞其重心做俯仰簡諧運動，俯仰角幅值為15°，周期為10 s，因此俯仰角運動規(guī)律可表示為θ=15°sin(0.2πt)。

以充液比50%即充液高度ys=0 m為例，利用CFD軟件進行100 s的仿真，期間部分時刻副油箱xy對稱面上的液體形態(tài)及壓強分布如圖6所示。

圖6 俯仰激勵下部分時刻液體形態(tài)及壓強分布 (ys=0 m)Fig.6 Liquid configuration and pressure distribution at some moments under pitch excitation (ys=0 m)

為驗證所建立單擺等效力學(xué)模型的有效性和準(zhǔn)確性，給出晃動力、晃動力矩和液體質(zhì)心位置的預(yù)測結(jié)果，與CFD軟件仿真結(jié)果的對比情況分別如圖7～圖9所示。

圖7 俯仰激勵下晃動力對比 (ys=0 m)Fig.7 Slosh forces comparison under pitch excitation (ys=0 m)

圖8 俯仰激勵下晃動力矩對比 (ys=0 m)Fig.8 Slosh moments comparison under pitch excitation (ys=0 m)

圖9 俯仰激勵下液體質(zhì)心位置對比 (ys=0 m)Fig.9 Liquid center of mass comparison under pitch excitation (ys=0 m)

從圖7給出的晃動力對比結(jié)果可見，俯仰運動將引起x軸方向單擺的晃動，并產(chǎn)生x軸方向晃動力，等效力學(xué)模型的預(yù)測結(jié)果與CFD軟件的仿真結(jié)果有較好吻合。在y軸方向上也會產(chǎn)生一定程度上作用力的波動，初始階段等效力學(xué)模型與CFD軟件吻合較好，但由于CFD仿真的誤差累積，其后半段仿真結(jié)果與基準(zhǔn)值相比有了一定程度上的偏離。由于俯仰激勵理論上不會引起z軸方向單擺的運動，因此z軸方向上的晃動力理論值和等效力學(xué)模型的預(yù)測值均為零。

如圖8所示，俯仰運動將引起較大的俯仰力矩即z軸方向的力矩，同時x軸方向的晃動力將產(chǎn)生一定的y軸方向晃動力矩。從y軸和z軸晃動力矩的對比結(jié)果可見，等效力學(xué)模型與CFD軟件的預(yù)測結(jié)果吻合較好。y軸方向的晃動力將引起x軸晃動力矩，因此在仿真后半段出現(xiàn)了一定的偏差。

由圖9給出的液體質(zhì)心位置對比結(jié)果可見，燃油質(zhì)心主要沿x軸方向運動，在y軸和z軸方向上變化較小，且等效力學(xué)模型在質(zhì)心位置上同樣有較為準(zhǔn)確的預(yù)測?；谏鲜鼋Y(jié)果，可認為等效力學(xué)模型對該組俯仰激勵下的液體動力學(xué)響應(yīng)進行了較為準(zhǔn)確的預(yù)測。

為進一步驗證該模型在其他充液高度上的有效性，對同組激勵下充液高度ys=0.06 m的情況進行計算。僅給出較為關(guān)心的x軸方向晃動力、z軸方向晃動力矩(俯仰力矩)和x軸方向質(zhì)心變化的預(yù)測結(jié)果，對比情況如圖10所示。可以看出，雖然同樣存在一定的誤差累積，但整體而言，等效力學(xué)模型可以對晃動方向上的作用力、力矩和質(zhì)心位置進行準(zhǔn)確預(yù)測。

圖10 俯仰激勵下液體動力學(xué)響應(yīng)對比(ys=0.06 m)Fig.10 Liquid dynamic responses comparison under pitch excitation (ys=0.06 m)

4.2 考慮燃料消耗算例驗證

對第3節(jié)中建立的變充液比等效力學(xué)模型進行驗證，考慮燃料消耗過程中等效力學(xué)模型各個參數(shù)的變化。為快速驗證模型有效性，假設(shè)副油箱內(nèi)初始燃油體積為0.586 4 m3(充液比為80%)，單位時間內(nèi)流出燃油體積0.005 m3，經(jīng)過100 s仿真時間后剩余燃油體積為0.086 4 m3(充液比為11.8%)。假設(shè)飛機在燃油消耗過程中繞其重心做滾轉(zhuǎn)運動，以如下兩組簡諧激勵為例進行模型驗證。

4.2.1 滾轉(zhuǎn)激勵(1)

滾轉(zhuǎn)角γ=30°sin(0.4πt)。利用CFD軟件進行仿真，期間部分時刻副油箱yz對稱面上液體形態(tài)及壓強分布如圖11所示，可明顯看到燃油的消耗過程，且燃油沿z軸有較大幅度晃動。

圖11 滾轉(zhuǎn)激勵(1)下部分時刻液體形態(tài)及壓強分布Fig.11 Liquid configuration and pressure distribution at some moments under roll excitation (1)

類似地，為驗證等效力學(xué)模型的有效性和準(zhǔn)確性，給出晃動力、晃動力矩和液體質(zhì)心位置的CFD仿真對比結(jié)果。由圖12所示的晃動力對比結(jié)果可見，滾轉(zhuǎn)激勵理論上不引起x軸方向單擺的晃動；z軸方向單擺的運動將引起y軸和z軸方向作用力的改變，且等效力學(xué)模型的預(yù)測結(jié)果與CFD軟件的仿真結(jié)果吻合良好。此外，隨著燃料的消耗，y軸方向上的作用力大小呈下降趨勢，等效力學(xué)模型同樣對這一影響進行了準(zhǔn)確預(yù)測。

圖12 滾轉(zhuǎn)激勵(1)下晃動力對比Fig.12 Slosh forces comparison under roll excitation (1)

由圖13給出的晃動力矩對比結(jié)果和圖14給出的液體質(zhì)心位置對比結(jié)果可見，等效力學(xué)模型對燃油消耗和燃油晃動所引起的動力學(xué)響應(yīng)都有著非常準(zhǔn)確的預(yù)測。從力、力矩、質(zhì)心位置等參數(shù)的結(jié)果都可以看到，等效力學(xué)模型的預(yù)測結(jié)果隨時間即燃油消耗而平緩變化，從一定程度上驗證了等效模型參數(shù)插值處理的合理性和可行性。

圖13 滾轉(zhuǎn)激勵(1)下晃動力矩對比Fig.13 Slosh moments comparison under roll excitation (1)

圖14 滾轉(zhuǎn)激勵(1)下液體質(zhì)心位置對比Fig.14 Liquid center of mass comparison under roll excitation (1)

4.2.2 滾轉(zhuǎn)激勵(2)

為進一步驗證模型有效性，改變滾轉(zhuǎn)激勵幅值和周期，令滾轉(zhuǎn)角γ=20°sin(0.2πt)。僅給出與滾轉(zhuǎn)激勵更為相關(guān)的z軸方向晃動力、x軸方向晃動力矩(滾轉(zhuǎn)力矩)和z軸方向質(zhì)心變化的預(yù)測結(jié)果，與CFD仿真結(jié)果的對比情況如圖15所示?？梢钥闯觯m然CFD仿真結(jié)果存在一定的高頻振動現(xiàn)象，但等效力學(xué)模型可以對作用力、力矩和質(zhì)心位置等的變化趨勢及變化范圍進行準(zhǔn)確計算。

圖15 滾轉(zhuǎn)激勵(2)下液體動力學(xué)響應(yīng)對比Fig.15 Liquid dynamic responses comparison under roll excitation (2)

5 結(jié)論

以飛機副油箱為研究對象，分析燃油在油箱內(nèi)的晃動行為，并基于勢流理論和有限元法建立了高效準(zhǔn)確描述燃油晃動動力學(xué)行為的單擺等效力學(xué)模型。以俯仰激勵和滾轉(zhuǎn)激勵為例，分別以不同充液高度的固定充液比副油箱和考慮燃料消耗的變充液比副油箱為對象，對所提出的燃油晃動等效建模方法進行了驗證。

通過將燃油晃動中最主要的動力學(xué)響應(yīng)如對油箱的作用力和力矩、燃油質(zhì)心位置變化等與CFD軟件的仿真結(jié)果進行對比，驗證了單擺等效力學(xué)模型的有效性和準(zhǔn)確性，以及模型參數(shù)插值處理的合理性和可行性。該建模方法可推廣至簡單形狀整體油箱內(nèi)的燃油晃動描述，并為含防晃裝置的副油箱及復(fù)雜結(jié)構(gòu)機翼油箱內(nèi)的燃油晃動建模提供重要參考。