關 宏 徐世杰
北京航空航天大學宇航學院,北京 100191
?
基于Fluent計算的火箭離軌姿態(tài)運動仿真與分析
關 宏 徐世杰
北京航空航天大學宇航學院,北京 100191
在航天器姿態(tài)運動分析中常遇到液體的流動、晃動及流/固耦合等問題,針對這類問題建立精確的數學模型是很困難的,通常需要利用數值計算和仿真手段進行分析。利用用戶互動功能(UDF),引入自定義的變量源函數,應用Fluent對運載火箭離軌姿態(tài)控制進行剩余燃料的流場數據分析。通過在流體運動方程中加入相關的牽連運動,得到在軌運行狀態(tài)下貯箱內剩余液體的運動參數和對貯箱的干擾力矩,為運載火箭系統(tǒng)的姿態(tài)運動分析和仿真提供運算參數,使得流場的變化與運載火箭的姿態(tài)運動相關聯,分析各個時刻流場運動狀態(tài)對箭體姿態(tài)的影響。
Fluent;離軌;剩余液體;用戶互動功能(UDF)
隨著運載火箭技術的進步,使用液體推進劑的火箭已經成為航天器的主要運載工具。一些一次性使用的運載火箭和重復使用的航天飛行器,在主發(fā)動機關機后,不可用推進劑、安全儲備、以及飛行混合比偏差所引起的剩余推進劑都將存留在貯箱里。這些殘存的液體燃料在貯箱里的運動將會影響到運載火箭的姿態(tài)運動,甚至引起火箭姿態(tài)失穩(wěn)。因此對貯箱里殘存液體的運動分析和計算是十分必要的。實際工程中,運載火箭在分離后的回收過程中常遇到剩余液體晃動產生的問題。一般是將剩余液體的運動作為對火箭的干擾來處理,根據剩余液體不同的分布規(guī)律,采用偏微分方程描述其運動變化規(guī)律[1-2],并作為姿態(tài)運動的流體分析干擾模型。但實際上,當殘留液體體積遠小于貯箱容積時,此方法難以反映貯箱-液體系統(tǒng)參數的瞬態(tài)分布特征,如果對貯箱和液體建立完整的動態(tài)數學模型,則存在建模困難、計算量大的問題,而且往往不能反映符合實際的參數分布特征。
目前,隨著計算流體力學、分布式仿真技術的發(fā)展[3],利用成熟的CFD仿真軟件Fluent,對工程中常遇到的具有非定常分布模型的流場能夠有效地進行計算,并且其建模簡單、計算精度和可靠性都很高,在工程預研、仿真驗證中得到了廣泛的應用。Fluent軟件具有很多優(yōu)點,但還存在一定的局限性:一方面,僅靠流場仿真無法揭示航天器姿態(tài)運動和貯箱內流體運動的相互作用;另一方面,需要將Fluent加載到其他系統(tǒng)仿真軟件中,若接口設計不當或網絡連接不穩(wěn)定將會影響仿真效率[4]。
本文根據Fluent軟件本身的特點,利用其用戶互動功能(UDF)[4],通過系統(tǒng)參數的實時同步傳輸,使姿態(tài)動力學與運動學模塊能夠實時地影響Fluent計算參數,從而建立了研究分布參數的閉環(huán)控制仿真平臺,為含有流體運動的姿態(tài)運動仿真提供新的途徑。通過利用這種仿真手段,本文對某運載火箭離軌前的姿態(tài)穩(wěn)定控制進行了數值計算與仿真,驗證了該方法的可用性和有效性,分析了含有剩余燃料的情況下箭體離軌姿態(tài)的變化規(guī)律。
本文選擇某運載火箭末級推進器為研究對象,在完成軌道任務后氧化劑與燃料貯箱內均剩余少量推進劑,體積遠小于貯箱容積。離軌過程中由于姿態(tài)調整以及箭體自旋運動,液體的分布難以確知,運動存在不可控因素,故采用Fluent流場模型對貯箱內液體進行仿真,同時通過UDF功能實現箭體姿態(tài)與液體運動的閉環(huán)仿真。
此處需要說明的是當箭體與液體作為整體系統(tǒng)看待時,不受外界干擾力矩作用,但是為了分析流體運動對箭體姿態(tài)的影響,需要將箭體與液體作為兩個子系統(tǒng),以分析其相互作用。
1.1 箭體系統(tǒng)動力學分析
首先以箭體系統(tǒng)作為研究對象,其轉動慣量沿自旋軸對稱,箭體本體坐標系Sb定義如圖1所示,坐標原點Ob在箭體質心,Xb軸沿本體縱向對稱軸方向,Yb軸沿箭體縱切面向上,Zb軸按右手法則定義。
圖1 箭體本體坐標系定義
定義系統(tǒng)慣量主軸沿坐標軸方向:繞橫軸的轉動慣量為Iy,Iz,且Iy=Iz,繞縱軸(自旋軸)的轉動慣量為Ix,且Ix (1) 其中ωx,ωy,ωz分別為繞Xb,Yb,Zb軸的角速度。 方程左側的Tx,Ty,Tz為箭體系統(tǒng)受到的外干擾力矩,本文即為剩余液體對箭體的作用力矩。 由式(1)以及軸對稱剛體運動的特性可知: 1)剛體自旋軸Xb繞角動量H的圓錐運動稱為空間章動,其空間章動速率Ω=H/It; 2)章動角γ為動量矩H與自旋軸之間的夾角,滿足關系式cosγ=Ixωx/H; 圖2 本體錐和空間錐 定義軌道參考系S系,與t=0時刻的本體系重合,x軸軌道切向,y軸沿t=0時刻箭體縱切面向上,z軸按右手法則定義。定義空間章動參考系S′系,坐標原點在航天器質心,x′軸沿本體縱向對稱軸方向,y′軸沿章動角速度Ω在本體系的投影方向,t=0時刻箭體縱切面向上,隨空間錐運動,z′軸按右手法則定義。在t=0時刻,此3坐標系重合。運動Δt時刻后,y軸在空間指向不變,y′軸運動α角,Yb軸運動到β角位置,α=∫ΔtΩdt為箭體空間章動運動角度,β=∫Δtωxdt為箭體自旋運動角度。即 1.2 液體子系統(tǒng)運動 得到箭體空間運動規(guī)律,再考慮液體系統(tǒng)的相應晃動影響。 液體子系統(tǒng)在微重力環(huán)境下運動,受擾動后不能斷定其確切位置,液氣界面穩(wěn)定性差。一旦航天器進行姿態(tài)控制(起旋至穩(wěn)定旋轉),剩余液體產生的擾動就會影響到姿態(tài)控制的精度,破壞航天器的運動穩(wěn)定性。又因為液體運動規(guī)律不易簡化,如不能準確仿真其運動狀態(tài),便無法找到適當的控制策略,會對回收造成惡劣的影響[6]。 本文流場模型選擇圖2中S′系為參考系,該參考系隨箭體在空間不斷運動,始終保證液體受空間錐離心加速度方向在y′軸上,網格坐標系原點為箭體質心位置,便于對Fluent定義相對加速度以及對質心求作用力矩。 基于質量守恒定律、動量定理、能量守恒定律、熱力學定律以及流體的本身物性,在流體力學中存在一組制約流體運動的基本方程組,對于黏性不可壓縮流體,滿足: (2) 其中,ρ為流體密度,v為流體相對運動速度,ρ=const,div(P)表示單位體積上應力張量的散度,f1為單位質量上的質量力分布函數。 從動量定理出發(fā),任取一體積τ的流體,它的邊界為A,其中流體動量的變化率等于作用于該體積上的質量力和面力之和。以f1表示作用在單位質量上的質量力分布函數,pn為作用在單位面積上的面力分布函數,則作用在τ上和A上的總質量力和面力為∫τρfδτ及∫SpnδA,而體積τ內的動量是∫τρvδτ,于是動量定理寫成 (3) 由于∫ApnδA=∫τdiv(P)δτ,故得到式(2)中右側第2項,又ρf1表示單位體積上的質量力。根據相對運動學原理,流體運動方程將寫成: (4) 在原有的流體運動方程中,流場網格模型可以計算出式(4)中的ar,即流體相對于S′系的相對加速度。但是要考慮與箭體運動的耦合運動,就需要補充由于旋轉相對運動產生的慣性牽連加速度ae,以及柯氏加速度2(ω×vr)。本文采用S′系作為流場仿真的參考系,其本身就具有相對本體系旋轉的相對運動,因此只需要考慮牽連運動,即此處ae。 如前節(jié)中指出箭體上各點相對于角動量H方向的向心加速度a=Ω2r,此處的ae=a。通過對箭體的姿態(tài)運動學求解,得到各個時刻的相對姿態(tài)角、姿態(tài)角速度,根據式(1)以及軸對稱剛體運動的特性,即可通過瞬時姿態(tài)信息,解得空間章動速率Ω、章動角γ,通過三角函數關系便可輕松的得到各點到角動量H的距離r。由此得到需要嵌入到網格計算模型中的自定義加速度a。 1.3 Fluent-UDF接口的設計 為了解決流場模型與箭體姿態(tài)運動模型數據交換問題,在保證流場模型和姿態(tài)運動模型實時同步運行的前提下,數據準確交換,提出使用Fluent-UDF數據交換接口實現交互仿真。UDF與Fluent數據交互流程如圖3所示。 圖3 使用UFD的流體模型與姿態(tài)運動模型交互示意圖 首先,根據流體網格中的初始化邊界條件(其中包括了對流場加速度a的初始化定義),進行1個步長的網格計算。得到這一步中網格內剩余液體對箭體的干擾力矩Tx,Ty和Tz。下面就進入到UDF程序計算中: 第1步:由Fluent流場模型解決流體運動的分布計算問題,根據流體對箭體(貯箱)的作用情況,定義UDF程序的入口:液體子系統(tǒng)輸出對模型壁面的作用力、作用力矩;通過流場模型得到剩余液體對于箭體系統(tǒng)的干擾Tx,Ty和Tz,即為UDF程序的調用入口參數; 第2步:當UDF從入口處讀出Tx,Ty和Tz后,再根據UDF程序中編寫的姿態(tài)控制方程,計算包括了仿真對象的姿態(tài)動力學與運動學方程的姿態(tài)運動狀況; 第3步:仍然是在UDF程序包中,調用加速度計算和轉換程序。根據相對作用原理,得到流場各點的下一個時刻的牽連加速度a1,通過UDF輸出接口賦值到Fluent網格計算環(huán)境中。 此時,完成一個步長的箭體與液體系統(tǒng)的交互仿真計算,根據此時刻的網格場內加速度a1,進行下一個步長的網格參數計算,回到第1步。其中主要涉及的UDF源文件定義如下: 對應第1步:源文件1: DEFINE_INIT() {定義箭體姿態(tài)運動初始條件} 對應第2步:源文件2: DEFINE_ADJUST() {加載所需參數,包括作用力矩,姿態(tài)運動信息,姿態(tài)運動變化率,中間參數; 判斷循環(huán)條件; 姿態(tài)動力學與運動學運動方程計算、疊加; 返回姿態(tài)運動變化參數,姿態(tài)運動信息;} 對應第3步,源文件3: DEFINE_SOURCE() {Real 加載所需參數; 調用加速度定義函數; Return 加速度;} 仿真開始后,當網格計算一步長時,源文件1運行,將姿態(tài)運動方程的初始化條件Tx,Ty和Tz載入;之后,源文件2運行,計算姿態(tài)動力學和運動學方程,得出這一步箭體的姿態(tài)運動狀態(tài);最后,源文件3調用上一步得到的姿態(tài)運動狀態(tài)參數,計算流場模型需要加載慣性加速度。每一步UDF計算開始時,都需要從網格載入上次循環(huán)時在流體系統(tǒng)計算中得到的作用力矩Tx,Ty和Tz,通過上述3個源文件再得到下一步的網格場加速度。將下一個時刻的網格場加速度a1輸出到網格邊界條件中,進行下一步網格計算。循環(huán)這一過程,就實現了Fluent計算和姿態(tài)動力學與運動學計算之間的實時同步數據傳輸過程。 圖4 液體晃動對航天器的作用力矩 由于液體的振動以及章動運動,S′系中y′軸方向上的作用力矩振蕩上升,始終在其正方向上。即,在Fluent仿真的網格計算模型中,通過自定義源函數加載的慣性加速度始終指向離心方向,符合客觀規(guī)律。通過Fluent流體系統(tǒng)仿真模型還可以看到運動過程中各個時間點的流場分布情況,如圖5。 圖5 不同時刻流場分布圖 隨著貯箱內的液體不斷地向一側壁面碰撞,產生對該側的作用力。圖7中左側的貯箱離質心較遠,這種碰撞作用力產生的力矩效果更加明顯。右側的貯箱內包含有坐標原點,液體因表面張力的作用環(huán)繞壁面流動,產生對壁面的剪應力,以及對質心的作用力矩。 圖5中是根據時間順序選擇的3幅典型的流場分布圖,其中左側貯箱內,如果對照圖5中各個狀態(tài)下的流場分布情況,便可輕易地分析圖4(a)中的作用力矩變化的原因。左側貯箱的液體越來越多地與貯箱的一側碰撞,使碰撞力的作用增大;右側貯箱的液體分布逐漸均勻,與貯箱壁面的接觸面積增大,摩擦力矩也隨之增大。因此,圖4(a)中的力矩變化剛好與圖5中的流場分布變化相吻合。 原來羅漠與任何人都是一樣的,他的獨特,他的深情,都只為一人而已。他的愛那么忠誠,經得起各種各樣的考驗,她其實是滿意的,可是她忍不住流了一路的眼淚。 在這種作用下,箭體的姿態(tài)運動受到的干擾應當是越來越大的,下面再對箭體系統(tǒng)進行分析。圖4(a)中將S′系中的作用力矩投影到本體系Sb中,其方向與本體系當時的相對位置有關。箭體系統(tǒng)受到的液體干擾力矩就是圖4(b)中的作用力矩。根據本文第1.1節(jié)中箭體動力學方程,得到在剩余液體干擾下箭體的姿態(tài)運動狀態(tài),如圖6所示。 圖6 液體干擾作用下的箭體自旋軸空間指向 從圖6分析整個系統(tǒng)的運動狀態(tài):在仿真時間內,箭體系統(tǒng)章動逐漸增加,因此旋轉離心作用增大,必然導致液體相對于箭體運動加強。根據液體運動對質心的力矩分析可知,液體干擾作用力矩也應當是逐漸增強的,與Fluent網格模型的仿真結果一致。并且隨著箭體運動穩(wěn)定性的惡化,繞最小慣量軸(自旋軸Xb軸)的自旋運動最終轉化為繞最大慣量軸(Yb軸或Zb軸)的不規(guī)則旋轉運動,因而箭體運動發(fā)散后在S′系中y′軸方向上的作用力矩出現下降并且幅值也產生不規(guī)則的變化。 本文選擇的仿真方法不僅可以測得流場與箭體的相互影響參數,通過編寫合理的UDF程序,還可以選擇模型中任意位置作為觀測點,測量其在流場中所受的作用力矩,或者對流場內的壓強、溫度、過載進行中間控制,通過適當的控制程序在相應時間修改這些參數。并且,也可以將包含控制判斷算法的程序以UDF程序的形式添加到閉環(huán)系統(tǒng)中,形成含有實時控制系統(tǒng)的姿態(tài)耦合運動仿真模型。 這種方法的優(yōu)點在于可以直觀形象地觀察到仿真各個時刻系統(tǒng)的變化和狀態(tài),便于分析影響姿態(tài)運動的各種因素。同時,這種選擇離心運動坐標系為仿真參考系的方法可以簡化箭體姿態(tài)運動對液體影響的描述,通過簡單的投影關系就可以得到清晰、解耦的作用關系。 由此可知,Fluent-UDF技術適用于更復雜的分布參數系統(tǒng),可以推廣到更復雜的包含流體計算的系統(tǒng)控制模型中。 [1] S.Dutta, M.K.Laha.Analysis of the Small Amplitude Sloshing of a Liquid in a Rigid Container of Arbitrary Shape Using a Low-order Boundary Element Method[J].International Journal for Numberical Methods in Engineering, 2000, 1(47): 1633-1648. [2] 王照林, 李磊.航天器內部液體晃動對交會對接動力學與控制的影響[J].航天控制, 1991,(2): 24-32.(WANG Zhaolin, LI Lei.Effect on Dynamics and Control of Rendezvous and Docking due to Liquid Sloshing Internal to Spacecraft[J].Aerospace Control, 1991, (2): 24-32.) [3] 耿浩, 蔣志文.基于網格的分布式仿真技術[J].航天控制, 2004, 22 (6): 46-48.(GENG Hao, JIANG Zhiwen.Distributed Simulation Technology Based on the Grid[J].Aerospace Control, 2004,22 (6): 46-48.) [4] 鮑文, 李偉鵬, 等.基于Fluent/Matlab接口的分布參數系統(tǒng)閉環(huán)控制仿真[J].北京:系統(tǒng)仿真學報, 2008, 20 (11):2851-2854.(BAO Wen, LEI Weipeng, et al.Close-loop Control Simulation Technology of Distributed Parameter System Based on FLUENT/MATLAB Interface[J].Journalof System Simulation, 2008,20 (11): 2851-2854.) [5] 肖業(yè)倫.航天器飛行動力學原理[M].北京: 宇航出版社, 1995: 212-217.(XIAO YeLun.Theory of Spacecraft Dynamics[M].Chinese Astronautics Press, 1995: 212-217.) [6] 廖少英.低重力或失重環(huán)境下剩余推進劑的管理和排放技術研究[J].上海航天, 1993, 5: 13-18,40.(LIAO Shaoying.Research of the Residual Propellant Dumping System under Micro-gravity[J].Aerospace Shanghai, 1993, 5: 13-18,40.) 《航天控制》雜志 □歡迎訂閱 □歡迎刊登廣告 □歡迎投稿 開本:大16開 雙月刊 頁碼:96頁 定價:20.00元(全年120元) 國內郵發(fā)代號:80-338 國外發(fā)行代號:BM4668 國內統(tǒng)一連續(xù)出版物號:CN11-1989/V 國際標準連續(xù)出版物號: ISSN 1006-3242 編輯部地址:北京142信箱402分箱 郵政編碼:100854 電話:(010) 68762264、68388585 E-mail: ht12@httx.com.cn 網址:http://htkz.chinajournal.net.cn The Analysis and Simulation of Attitude Motion for a De-orbit Rocket Based on Fluent Flow Field Calculation GUAN Hong XU Shijie School of Astronautics, Beihang University, Beijing 100191, China Duetotheattitudemotionanalysisofaspacecraft,thefluxionandsloshingofliquidarefrequentlyencountered.Sinceitisdifficulttoestablishaprecisemathematicalmodelforsuchissues,researchersprefertousenumericalcalculationandsimulationforanalysis.Inthispaper,theUDFisusedtointroducetheuserdefinedvariablefunctionandtherocketde-orbitcontrolwithflowfieldsimulationofresidualfuelbasedonfluentistokenforexampletoanalyze.Byaddingtransportmotiontothefluidmotionequations,theparametersanddisturbingtorqueofflueundertheon-orbit-statecanbeobtained.Theoperationparametersofattitudeanalysisandsimulationforlaunchvehiclesystemareprovidedsothatthechangesinflowandtheattitudemotionofthelaunchvehiclearecorrelated. Fluent;De-orbit;Residualfuel;User-definefunction(UDF) 2011-08-16 關 宏(1985-),女,北京人,博士研究生,主要研究方向為航天飛行器動力學與控制;徐世杰(1951-),男,吉林人,教授,博士生導師,主要研究方向為航天飛行器動力學與控制。 V525 A 1006-3242(2012)03-0093-062 仿真結果與分析
3 結論