重力環境下液體大幅晃動運動脈動球模型及實驗研究1)

盧煜 岳寶增,2) 馬伯樂 郝柏隆 常園園 張 薇

* (北京理工大學宇航學院,北京 100081)

? (北京宇航系統工程研究所,北京 100076)

引言

現代航天器通常需要攜帶大量的液體燃料來完成長期、復雜的飛行任務.航天器在進行發射升空、姿態機動、交會對接等過程中,由于外部激勵的作用,可能會導致液體燃料發生劇烈晃動,由此產生附加的晃動力和力矩會對航天器造成重要影響,甚至導致航天器運行姿態失控或儲箱結構損壞.因此,針對航天器中液體燃料大幅晃動問題的研究是非常有意義的[1-5].

由于流體動力學計算的復雜性和星載計算機計算能力的局限性,在實際工程中幾乎不可能采用計算流體動力學(computational fluid dynamics,CFD)方法來模擬液體大幅晃動特性,需要通過等效建模的方法來簡化模型、減少計算量,使得通過等效模型得到的晃動特性在可接受的誤差范圍內逼近真實情況.傳統的液體晃動等效建模方法通常基于單擺或彈簧振子模型,只適用于液體燃料晃動幅度遠小于儲箱尺寸的線性小幅晃動,很難真實地反應出液體大幅晃動的非線性行為[6-11].因此,需要對液體大幅晃動等效力學模型展開深入研究.

在液體大幅晃動等效力學模型研究方面,國內外學者目前已取得一些進展.Berry 等[12-14]提出了一種名為LAMPS (large amplitude slosh)的模型,將Cassini 儲箱等效為橢球面,液體燃料等效為一個在面上大范圍自由運動的質點,并通過實驗驗證了該模型預測液體大幅晃動作用力的有效性.黃華等[15-17]對LAMPS 模型進行了修正,使得等效的質點可以離開約束面,并引入Hertz 接觸模型計算晃動力.Liu 等[18]提出了一類適用于研究常/低重力環境下,球形貯腔內液體非線性晃動的等效力學三自由度剛體擺復合模型,將液體的等效轉動慣量考慮在內,能夠較為全面地描述航天器貯腔內液體可能發生的多種運動模式,包括液體相對于貯腔的整體性的剛體運動、橫向晃動、旋轉晃動以及液體自旋運動.荷蘭科學家Vreeburg 等[19-22]首次提出了運動脈動球模型(moving pulsating ball model,MPBM),將航天器儲箱中液體等效成質量不變、半徑發生變化、與儲箱壁始終保持接觸的均勻球體,該模型可以模擬出液體大范圍運動過程中慣性張量的變化特性,并在Sloshsat FLEVO 衛星上驗證了其合理性與可靠性.Deng 等[23-26]對MPBM 模型進行了改進和修正,在保留運動方程的基本形式的前提下將靜態表面張力和離心力對等效質量的影響引入模型中,進一步提高了MPBM 模型在失重環境中的準確性.

以上所描述的液體大幅晃動等效力學模型均適用于失重或微重環境中,然而在實際工程應用中,很多情況下液體所處的重力環境不可忽略不計,例如:航天器在行星表面著陸前的懸停避障階段,液體處于行星引力形成的重力場中;航天器變軌機動的過程中,會通過加速產生等效重力,使得液體燃料聚集在“出液口”附近,提高燃料利用率.重力的存在使得等效模型的建模方法需要進行改進或修正,因此,本文基于運動脈動球模型,將重力和液體所收到的慣性力引入其中、對等效質量和液體質心的位置加以修正,并設計實驗,通過實驗、仿真的結果同時證明推廣后的模型的有效性.

1 對MPBM 模型的描述與改進

引入慣性坐標系Cixiyizi和航天器本體坐標系Cxyz,設 ? 為航天器本體坐標系Cxyz相對于慣性坐標系Cixiyizi和的角速度,Vc為航天器本體坐標系Cxyz相對于慣性坐標系Cixiyizi的速度.如圖1 所示,儲箱形狀為球形,半徑為R,航天器本體坐標系Cxyz的原點與球形儲箱幾何中心重合;假設液體總質量為m,參與晃動的液體質量為ms,被等效為一個半徑始終發生變化的脈動球,其余液體的質量為m0=m-ms.定義參數表示參與晃動的液體質量分數.脈動球與儲箱壁面始終保持接觸,接觸點為P,相對于儲箱的速度和角速度分別為Vs和 ωs,脈動球對儲箱的晃動作用力FL和力矩TL在接觸點P處傳遞,靜止液體對儲箱的作用力為F0.脈動球質心為S,在本體坐標系上的矢徑為rs,其單位向量和模長分別為e和r,顯然,脈動球的半徑為R-r;不參與晃動的液體同樣被等效為一個小球,質心在本體坐標系上的矢徑為r0.脈動球在慣性坐標系中的平動動力學方程為

圖1 改進的MPBM 模型示意圖Fig.1 Schematic diagram of improved MPBM model

其中,Gs代表脈動球所受到的重力作用,FL由法向力Ne和摩擦力Fb矢量相加得到,即FL=Ne+Fb.

為了求解法向力的大小N,將脈動球在慣性坐標系中的平動動力學方程向接觸點P的法線方向投影,方程兩端同時點乘單位向量e,即

化簡可得

式中 σ 為液體的表面張力系數.將能量關系式(4)與脈動球平動動力學方程的投影式(3)聯立,最終可以得到法向力大小N的通用表達式

當脈動球半徑達到最小值R-rmax時,法向速度反向,此時=0,法向力大小

脈動球和儲箱之間的相互作用力矩TL根據文獻[14]中的經驗公式給出,摩擦力Fb可以根據經驗公式寫成以下形式[26]

其中,fb為待定系數,可通過實驗或數值方法標定,μ為流體的運動黏性系數,L為脈動球半徑,Vu為脈動球和儲箱壁在接觸點P處的相對速度.

如圖2 所示,脈動球相對于儲箱壁面的速度和角速度分別為Vs和ωs;y方向為儲箱壁面的切向,Vy為Vs沿著y方向的分量,ωr為 ωs沿著法向e的分量,Rs為此刻脈動球的半徑.則

圖2 MPBM 模型運動分解示意圖Fig.2 Schematic diagram of motion decomposition of MPBM model

脈動球在儲箱內壁做“既滑動又滾動”的運動,其中 ωr不引起脈動球的平動運動,ωs-ωr是引起脈動球質心平動的角速度分量.因此,脈動球和儲箱壁在接觸點P處的相對速度Vu可以寫成以下形式

至此,可以寫出重力環境下在航天器本體坐標系Cxyz中,脈動球的動力學方程的等價形式

不參與晃動的液體在慣性坐標系中的平動動力學方程為

因此,這部分液體對儲箱的作用力

儲箱受到液體的總作用力Ftotal=FL+F0.

由于有質量為m0的液體并未參與晃動,因此,修正的等效模型中并不能將脈動球的質心位置當作所有液體的質心位置,應當對質心的求解加以修正.脈動球質心在本體坐標系上的矢徑為rs,不參與晃動的液體質心在本體坐標系的矢徑為r0,液體整體質心的位置矢量為rmass_center,有如下關系式

初始時刻液體靜止,rmass_center可根據球形儲箱的充液比求得,脈動球由于重力的作用沉在儲箱底部,rs也已知,由此可以確定不參與晃動的液體質心位置矢量r0.在后續計算中,修正后的液體質心坐標便可以根據關系式(9)求得.

2 仿真校驗

2.1 大幅晃動算例

為了檢驗改進后MPBM 模型的有效性,對充液比為50%的圓球腔中液體的晃動進行了模擬.為方便與文獻中已有結果進行比較,取儲箱半徑R為0.148 m.儲箱沿x方向進行簡諧運動,振幅為A,振動頻率為fE,初相位為0,則儲箱的位移為

在航天器本體坐標系中,運動脈動球所受到的慣性力

靜止質量所受到的慣性力

將MPBM 模型仿真結果與文獻[27-28]中實驗的大幅晃動算例進行對比驗證.設定儲箱位移的振幅A=0.003 04 m,頻率fE=1.5 Hz,重力加速度大小g=9.8 m/s,沿 -z方向.儲箱內液體為水,充液比50%,密度 ρ=1000 kg/m3,表面張力系數σ=0.072 N/m,動力學黏性系數μ=0.001 Pa·s.

初始時刻脈動球質心位置矢量rs=[0,0,-rmax]T,rmax的大小可以參考文獻[18]中剛體擺等效模型的擺長進行取值并進行微小修正.相比于液體的小幅晃動,大幅晃動的算例中黏性不可忽略,另外,液體晃動的一階固有頻率也會降低、參與晃動的質量分數相對變大,導致MPBM 模型中rmax和fem的取值需要修正得更大一些.本算例中摩擦力Fb表達式中的待定系數fb取為0.2,參與晃動的液體質量分數fem=0.57,與脈動球質心位置矢量相關的參數rmax=0.097 2 m.

圖3 給出了沿外激勵方向晃動力的對比結果.從圖中可以看出,液體晃動顯示出拍振的形式.推廣的MPBM 模型仍然可以很好地預測液體晃動力的幅度和變化趨勢,晃動力的響應有大概0.2 s 的提前量.

圖3 液體大幅晃動的晃動力Fig.3 Sloshing force of large amplitude liquid sloshing

圖4 展示了改進的MPBM 模型對液體質心位置預測的良好效果.藍色虛線表示改進之前MPBM模型的計算結果,紅色實線表示本文中對質心位置修正之后的結果,黑色虛線是通過CFD 軟件Flow3D預測得到的液體質心坐標.從圖中可以看到,MPBM模型對液體質心位置變化趨勢的預測是非常準確的.雖然 MPBM 模型質心的運動幅度明顯大于Flow3D的計算結果,但是本文對質心修正之后,沿外激勵方向質心坐標的變化幅度更加逼近CFD 軟件得到的結果.

圖4 液體大幅晃動質心坐標x 方向分量Fig.4 x-direction component of centroid coordinate of large amplitude liquid sloshing

圖5 分別展示了CFD 軟件計算得到的4.2 s 和7.0 s 時液體在x-z平面上的壓強分布云圖和速度場的矢量圖.根據MPBM 模型的計算結果,4.2 s 時液體質心x坐標達到一個極大值,此時液體質心水平速度為0.從圖5(a)中可以看到,此刻在x-z平面上的速度分布幾乎處處為0,也進一步印證了MPBM模型的計算結果.同理,圖4 的計算結果表明,7.0 s時液體質心剛剛從x軸正方向運動到x軸負方向,這也在CFD 軟件繪制的云圖和速度場中得到了驗證.

圖5 液體在x-z 平面上的壓強和流場Fig.5 Pressure and vectors of liquid on x-z plane

2.2 零動量機動算例

在工程實際中,如行星探測器著陸或者火箭回收過程中懸停避障,航天器主剛體可能會受到大小相等、方向相反的橫向加速度脈沖激勵,由于整個運動過程前后系統動量不變,因此稱為零動量機動.基于這個工程背景,研究了MPBM 模型在零動量機動過程中晃動力的等效情況.設定重力加速度大小g=9.8 m/s,沿-z方向.儲箱內液體為水,充液比50%,密度 ρ=1000 kg/m3,表面張力系數σ=0.072 N/m,動力學黏性系數 μ=0.001 Pa·s.本算例中摩擦力Fb表達式中的待定系數fb取為0.25,參與晃動的液體質量分數fem=0.56,與脈動球質心位置矢量相關的參數rmax=0.096 m.儲箱脈沖加速度的幅值aE=0.5 m/s2,持續時間TE=1 s,即

圖6 為零動量機動過程中晃動力和液體質心坐標的x方向分量的變化情況.通過圖6 中與CFD 軟件Flow3D 的計算結果對比可以看出,對于脈沖激勵下的液體晃動,MPBM模型晃動力的計算結果仍然是準確且可靠的.

圖6 零動量機動的晃動力Fig.6 Sloshing force of zero momentum maneuver

基于上述仿真與驗證結果,使用運動脈動球模型對脈沖激勵的時序優化方案進行仿真研究.對于脈沖加速度的幅值 0.5 m/s2,持續時間與間隔時間均為1 s 的兩對零沖量激勵,分別按照“正正負負”和“正負正負”兩種方式排列,并分別命名為“方案1”和“方案2”.兩種方案所對應的加速度圖像和液體響應的對比分別如圖7 和圖8 所示.

圖7 兩種方案的加速度時序Fig.7 Acceleration excitation of two cases

從圖8 中的晃動力和液體質心坐標的仿真結果可以清晰地看出,方案2 所造成液體晃動的響應明顯小于方案1.由此得出結論,零動量機動過程中,方向相反的脈沖激勵相鄰排列,可以有效地削弱液體晃動響應對航天器主剛體的影響.

圖8 兩種方案的液體晃動響應Fig.8 Sloshing responses of liquid of two cases

3 實驗研究與驗證

為了研究復雜形狀充液儲箱液體晃動行為,同時也為了對等效模型進行驗證,建立了一套用于測量充液儲箱液體晃動力的實驗裝置和實驗方法,并對矩形儲箱內液體橫向的受迫振蕩進行實驗研究.

3.1 實驗裝置和方法

實驗裝置包括實驗基座、靜壓導軌、滑塊、平臺、長孔角碼、力傳感器、加速度傳感器、充液儲箱、軸承滾珠、伺服電機、伺服電機控制系統、絲杠、信號采集系統,如圖9 和圖10 所示.

圖9 實驗裝置示意圖Fig.9 Schematic diagram of experimental device

圖10 實驗裝置實物圖Fig.10 Photo of experimental device

伺服電機用于為儲液箱提供水平方向可變頻率的簡諧振動外激勵,生產廠家為日本Yaskawa 公司,型號SGM7G-13A7C61,額定扭矩8.34 N·m,速度響應頻率1500 Hz.力傳感器采用美國PCB 公司208C02壓電式力傳感器,測量范圍: (壓縮)100 lb(0.444 8 kN),靈敏度: ± 15%(11 241 mV/kN).加速度傳感器和信號采集系統均為國產品牌東華測試,型號分別為1B104壓阻式加速度傳感器和DH5922D 動態信號測試分析系統.實驗過程中信號采集頻率為20 Hz.

兩條靜壓導軌由螺栓固定在實驗基座上,實驗平臺底部安裝滑塊,使得實驗平臺可以沿著導軌水平滑動;一個直角的長孔角碼的水平段由螺絲固定在實驗平臺上,其豎直段上固定一個力傳感器和一個加速度傳感器,充液儲箱與力傳感器用黏接的方式固連.充液儲箱底部均勻鋪設直徑9 mm 的軸承滾珠,用于減小充液儲箱和實驗平臺之間的摩擦力對實驗的干擾.伺服電機由螺栓固連在實驗基座上,并通過絲杠與實驗平臺進行連接,絲杠將伺服電機的繞軸線旋轉運動轉化為水平簡諧振動.實驗平臺在靜壓導軌上進行水平往復運動的振幅和頻率由伺服電機控制系統進行調節.信號采集系統用于采集力傳感器、加速度傳感器的數字信號,并通過濾波、四則運算等信號處理,輸出儲箱受迫振動時液體的晃動力.

晃動力的具體測量方法如下: 力傳感器測得的信號F(t) 為充液儲箱和其中的液體受到的總力,加速度傳感器測得的信號a(t) 為實驗平臺在靜壓導軌上進行水平往復運動的加速度,M為充液儲箱的質量,則液體晃動所產生的晃動力為Fsloshing=F(t)-M×a(t).

基于以上描述的實驗裝置和測量方法,采用矩形充液儲箱進行一組實驗,用于驗證該實驗裝置測量晃動力的準確性.儲箱由亞克力材料制成,底面為矩形,長寬L均為20 cm,內部裝載液體為水,高度h是10 cm.根據線性晃動理論,可知液體晃動各階固有頻率的計算公式為[29-30]

式中,n為晃動的模態階數.通過計算得到液體晃動的一階固有頻率f1=1.89 Hz.通過調整伺服電機的轉動幅度和頻率,使得振動平臺在導軌上往復運動,振幅A=0.0033 m,頻率fE=1.6 Hz,根據上述實驗方法測量了30 s 內液體的晃動力.同時,在Flow3D 軟件中建立相應的CFD 模型,計算同樣工況下液體晃動力的變化情況,對比結果如圖11.

圖11 實驗數據和Flow3D 對比圖Fig.11 Comparison between experimental data and Flow3D

從對比圖中可以看出,實驗方法和CFD 方法得到的結果可以非常好地吻合,既驗證了實驗平臺的設計和實驗方法的可靠性,也一定程度上說明了Flow3D 軟件計算的準確性,達到互相驗證的目的.

3.2 矩形儲箱的等效MPBM 模型

由前文所述,MPBM 模型的推導與分析是基于圓球形充液儲箱的,對于非球形儲箱,需要將其等效為圓球形.以3.1 節中實驗所用的矩形儲箱為例來說明問題.

如圖12 所示,將長寬L均為20 cm 的矩形充液儲箱等效為半徑R為0.3 m 的球腔,儲箱位移的振幅A=0.003 3 m,頻率fE=1.6 Hz,重力加速度大小g=9.8 m/s,沿 -z方向.液體密度 ρ=1000 kg/m3,表面張力系數 σ=0.072 N/m,動力學黏性系數μ=0.001 Pa·s.本算例中摩擦力Fb表達式中的待定系數fb取為0.3,參與晃動的液體質量分數fem=0.7,與脈動球質心位置矢量相關的參數rmax=0.07 m,即初始時刻脈動球質心位置矢量rs=[0,0,-0.07 m]T.

圖12 矩形儲箱的等效MPBM 模型Fig.12 Equivalent MPBM model of rectangular tank

將30 s 內MPBM 模型晃動力的計算結果與3.1 節中實驗所測得的數據進行對比,如圖13 所示.

圖13 實驗數據和MPBM 模型對比圖Fig.13 Comparison between experimental data and MPBM model

圖13 中紅色實線為MPBM 模型計算結果,黑色虛線為實驗數據,可以看出MPBM 模型用于等效矩形儲箱的液體晃動也同樣可以很好地反應出晃動力的變化趨勢.

4 結論

本文對失重環境下的液體大幅晃動運動脈動球模型進行了改進和推廣,將液體受到的重力引入動力學方程,并對質心位置加以修正,使得MPBM 模型的應用場景更加廣泛.文中通過大幅晃動、零動量機動的算例,驗證了推廣的MPBM 模型的有效性并在此基礎上研究了脈沖激勵時序的優化方案.另外,設計了一種測量充液儲箱晃動力的實驗裝置與實驗方法,驗證了MPBM 模型在等效非球形儲箱的液體晃動時也同樣可以很好地反應出晃動力的變化趨勢.