胡俊,于菲,于勇
(北京理工大學 宇航學院,北京 100081)
降落傘作為一種可展開式氣動減速器,具有質(zhì)量小、折疊體積小、價格便宜、阻力特性好等優(yōu)點,被廣泛應用于空投、航天器回收、新型武器等領域.
降落傘工作過程包括拉直、充氣、減速和穩(wěn)降,其中的拉直和充氣階段歷時短、變形大,涉及幾何非線性和材料非線性并存的流固耦合問題[1].降落傘的實驗研究包括風洞實驗、飛行實驗等方法.如 LEVIN 等[2]通過風洞試驗測量了5 種不同傘面形狀的十字型傘的氣動力以及俯仰力矩,并與空投實驗數(shù)據(jù)進行了結(jié)果對比,得到不同形狀、不同臂長比十字型傘的氣動力規(guī)律.RENDER 等[3]、LUDTKE[4]、BRAUN等[5]主要對小型十字傘傘群進行了一系列風洞試驗,研究了物傘連接繩長度和傘面之間系留繩長度對傘群系統(tǒng)氣動力的影響,對2 種不同尺寸的傘群系統(tǒng)進行不同的雷諾數(shù)范圍內(nèi)的阻力測量,得出了傘群阻力系數(shù)與單傘阻力系數(shù)之間的關系式.NASA 蘭利研究中心[6]多次采用飛行實驗,對“行星進入降落傘計劃”(planetary entry parachute program,PEPP)中不同構型的降落傘,進行有關穩(wěn)定性、阻力性能和降落傘尾跡的研究.
隨著計算機技術的發(fā)展,數(shù)值模擬也成為降落傘研究的一種手段.通過將穩(wěn)定的降落傘視為剛體,模擬研究降落傘的氣動特性和周圍流場.如XUE等[7-8]采用CFD 方法針對探測器超聲速降落傘系統(tǒng)進行了詳細研究,觀察了探測器尾流與傘前激波、探測器激波與傘前激波的氣動干擾以及傘繩激波、風洞連桿激波對整體流場的影響情況.周彤等[9]、李堯堯等[10]也采用了CFD 的方法對傘-彈系統(tǒng)中的相關氣動問題進行了研究.BARNHARDT 等[11]則采用分離渦模擬方法研究了NASA“火星科學實驗室任務”相關的超聲速降落傘系統(tǒng),模擬了超聲速來流下的探測器-降落傘系統(tǒng),研究了前體探測器尾流對降落傘氣動性能的影響,分析了尾流對傘衣壓力、阻力的影響以及傘衣帶底部有可能會出現(xiàn)部分坍塌的原因.
降落傘在穩(wěn)定減速和穩(wěn)定下降階段時,將傘衣視為剛體是可行的,但是在降落傘的拉直和充氣階段,傘衣會出現(xiàn)結(jié)構大變形,不能再將傘衣視為剛體,此時數(shù)值模擬主要采用任意拉格朗日-歐拉(arbitrary Lagrangian-Eulerian,ALE)法.2011 年TUTT 等[12]首次采用ALE 方法模擬降落傘有限質(zhì)量充氣過程,結(jié)果驗證了ALE 方法對開傘力、下降速度及降落傘呼吸頻率有較好的預測,驗證了ALE 流固耦合方法模擬有限質(zhì)量開傘過程的可行性.GAO 等[13]采用ALE算法模擬了開縫傘的有限質(zhì)量充氣過程,用顯式中心差分法求解了不可壓縮流動的Navier-Stokes (N-S)方程,得到了傘衣變形過程、減速規(guī)律及開傘過載.
降落傘作為氣動減速器,其主要作用是對載荷進行減速,避免載荷在到達落點時產(chǎn)生過大的沖擊,從而達到保護載荷的作用.但傳統(tǒng)的降落傘在運動過程中無法對物傘系統(tǒng)的運動狀態(tài)進行主動調(diào)整.為了提高物傘系統(tǒng)落點的精準度,近年來開始出現(xiàn)主動控制降落傘的應用研究.通過傘繩的收放,改變傘面構型,從而影響改變傘面所受氣動力,達到控制物傘系統(tǒng)運動狀態(tài)的目的.
降落傘的主動控制最早由美國提出,運用于精確控制空投系統(tǒng)的軌跡控制.BERGERON 等[14-15]在十字型傘的其中一片傘臂上添加了擾流片,通過擾流片的關閉和打開改變傘衣內(nèi)部的流場結(jié)構,使傘衣變形,達到改變降落傘氣動性能的目的.FIELDS 等[16]通過在傘面安裝驅(qū)動裝置改變傘繩長度,控制降落傘的旋轉(zhuǎn)或滑翔,根據(jù)空投實驗數(shù)據(jù)開發(fā)四自由度及PID 控制模型,將著陸位置圓誤差降低2~4 倍;李龍恩[17]通過動力學建模的方法,研究收放局部傘繩實現(xiàn)對圓形降落傘控制的可行性,發(fā)現(xiàn)收放局部傘繩能夠?qū)崿F(xiàn)對降落傘的主動控制.
有關學者對于降落傘主動控制的研究主要集中在實現(xiàn)傘繩收放的控制方法與控制裝置等方面,但對于降落傘收縮傘繩產(chǎn)生的氣動力規(guī)律沒有系統(tǒng)地認識.掌握降落傘收縮傘繩之后的運動規(guī)律和氣動力變化規(guī)律是實現(xiàn)降落傘主動控制應用的前提條件,同時也影響著發(fā)送調(diào)控指令的時間節(jié)點.不同傘繩收縮長度時降落傘的氣動特性,尤其是側(cè)向控制力的規(guī)律是降落傘主動控制成功應用的重要因素.本文以十字型降落傘為研究對象,建立收縮不同傘繩長度的降落傘無限質(zhì)量充氣的有限元模型,采用ALE流固耦合方法,對收縮不同傘繩長度時十字型傘的氣動特性進行研究,為用十字型降落傘進行主動控制提供參考.
柔性降落傘在運動過程中和空氣的相互作用具有強非線性時變的特點,兩者的耦合是同時涉及幾何非線性和材料非線性的瞬間大變形雙向耦合問題.ALE 算法是一種基于罰函數(shù)的接觸耦合算法,結(jié)合了拉格朗日算法(主要用于求解固體力學)和歐拉算法(主要用于求解流體力學)的優(yōu)點,可以實現(xiàn)結(jié)構和流體的雙向耦合.
ALE 描述下的流體模型遵循Navier-Stokes 方程,其中包括連續(xù)性方程、動量守恒方程、能量守恒方程.
ALE 描述下的連續(xù)性(質(zhì)量守恒)方程為
ALE 描述下的動量守恒方程為
ALE 描述下的能量守恒方程為
式中: ρ為流體密度;t為時間;vi為流體速度;wi為參考構型下網(wǎng)格點移動速度;bi為單位質(zhì)量體力;e為單位質(zhì)量內(nèi)能; σij為應力張量,表達式為:
其中:p為流體壓力; δij為克羅尼克函數(shù); μ為流體動力黏度系數(shù).
ALE 方法對固體域采用Lagrange 描述.降落傘的結(jié)構控制方程為
式中: ρs為降落傘織物密度;u為單元位移; σs為柯西應力張量;f為單元所受體積力;g為重力加速度.
降落傘傘衣采用不考慮厚度的二維薄殼結(jié)構,材料選擇織物,傘衣織物的透氣量與透過傘衣的流體速度呈二次多項式關系[18]:
式中: Δp為織物兩側(cè)壓力差;e為織物厚度;a、b分別為黏性系數(shù)和慣性系數(shù),均是關于多孔介質(zhì)相對透氣量 ε和空氣動力黏度 μ的多項式.
建立傘衣結(jié)構網(wǎng)格和流體網(wǎng)格重疊的有限元模型,使用罰函數(shù)接觸耦合算法實現(xiàn)力學參量的傳遞,在任意tn時間步內(nèi),追蹤拉格朗日節(jié)點(從節(jié)點)穿透歐拉流體物質(zhì)(主節(jié)點)的相對位移dn,則下一時刻迭代結(jié)果為:
式中:vrel為主從節(jié)點相對速度;主從節(jié)點速度為等參坐標變化后得到的某一流體質(zhì)點速度vf;從節(jié)點速度為結(jié)構節(jié)點速度vs; Δt為時間步長.
十字型降落傘主要包括傘衣和傘繩兩部分,其中,傘衣部分是由5 片尺寸相同、材質(zhì)相同的織物縫制而成,十字型傘的幾何結(jié)構如圖1 所示.L為傘長臂,W為傘短臂,Ls為傘繩長度.十字型傘的氣動特性主要受其幾何參數(shù)的影響,其中主要幾何參數(shù)包括臂長比A=L/W和傘繩比Rs=Ls/L.
如圖2 所示建立坐標系,傘繩節(jié)點作為坐標系原點O,來流方向為坐標系的+Z軸,與Z軸垂直的平面為XY平面,建立風軸坐標系OXYZ;以壓力中心O′作為坐標系原點,傘軸所在直線OO′作為Z′軸,與Z′軸垂直的平面為X′Y′平面,建立體軸坐標系O′X′Y′Z′.
風軸坐標系Z軸與體軸坐標系Z′軸之間的夾角為 α.在風洞實驗中,通過固定降落傘傘頂進行降落傘的靜態(tài)氣動特性研究,該夾角就是降落傘靜態(tài)實驗中的攻角.
圖3 所示為降落傘的受力示意圖.風軸坐標系中降落傘的氣動力主要包括阻力FD和升力FL.體軸坐標系中,降落傘的主要受力為軸向力FT和法向力FN.
計算氣動力系數(shù)時所用到的參考面積為十字型傘的名義面積S0,表達式[19]為:
為驗證ALE 方法在十字型降落傘研究中的適用性,以文獻[20]中風洞試驗的十字型降落傘為例,采用ALE 算法對對稱十字型降落傘進行數(shù)值模擬,十字型降落傘的具體尺寸和計算條件如表1 所示.
表1 十字型降落傘的尺寸和來流條件[20]Tab.1 Parameters of cross parachute and flow condition [20]
采用上述降落傘尺寸,模擬了攻角為0°和15°時十字型傘由平面狀態(tài)充氣至穩(wěn)定運動的過程,取計算穩(wěn)定后的1~4 s 的軸向力數(shù)據(jù),計算時間平均值,結(jié)果如圖4 所示.
圖4 風洞試驗數(shù)據(jù)[20]和模擬數(shù)據(jù)對比Fig.4 Comparison of wind tunnel data [20] and simulation data
攻角為0°時采用數(shù)值方法計算出十字型傘的軸向力系數(shù)大小為0.822,風洞試驗所得軸向力系數(shù)為0.8,兩者誤差約為2.75%;攻角為15°時數(shù)值計算出軸向力系數(shù)為0.765,風洞試驗所得軸向力系數(shù)為0.712,兩者誤差約為7.4%.
綜上所述,數(shù)值模擬計算與風洞試驗數(shù)據(jù)誤差小于10%,認為本文所使用的數(shù)值模擬方法適用于十字型降落傘模擬計算.
本文基于某火箭彈小角度迅捷轉(zhuǎn)彎制定了一種十字型降落傘,具體幾何尺寸如表2 所示.
表2 十字型降落傘的幾何尺寸Tab.2 Geometric parameters of cross parachute
每個傘臂連接2 根傘繩,整個降落傘系統(tǒng)共有8 根傘繩.本文目的是研究十字型降落傘在穩(wěn)定狀態(tài)下的氣動特性,對開傘過程不做精確要求,所以將降落傘的初始狀態(tài)設置為十字形平面狀態(tài),見圖1.
傘面結(jié)構為薄殼單元,對傘面進行正四邊形網(wǎng)格劃分,保證每個傘面網(wǎng)格的質(zhì)量均等.傘繩結(jié)構采取離散梁單元,流場域采取固體單元.十字型傘的流場域采用直徑為5L、高度為6L的圓柱形計算域,采用結(jié)構性網(wǎng)格劃分,相關的有限元模型參數(shù)如表3所示.
表3 降落傘的有限元模型參數(shù)Tab.3 Finite element model parameters of parachute
流場域和傘面的網(wǎng)格劃分如圖5 所示.
圖5 網(wǎng)格示意圖Fig.5 Diagram of mesh
如圖6 所示,收縮傘臂A所連接的2 根傘繩,收縮長度定義為Δl=Ls-,單位為m,收縮比定義為δ=Δl/Ls.
圖6 傘繩收縮示意圖Fig.6 Diagram of suspension line retracting
來流速度為68 m/s,來流密度為1.205 kg/m3,當收縮長度 Δl為0、0.16、0.20、0.24、0.32、0.40 m 時,所對應的收縮比δ及其他計算條件如表4 所示.
表4 數(shù)值模擬計算條件Tab.4 Numerical simulation conditions
通過傘繩的收縮,十字型降落傘由對稱性幾何構型變成非對稱性幾何構型,同時降落傘的氣動力參數(shù)發(fā)生改變.
圖7 所示為0°、5°、10°攻角時,十字型傘的瞬時阻力系數(shù)隨時間的變化.
圖7 不同收縮比十字型傘的阻力系數(shù)Fig.7 Drag coefficients of cross parachutes with different retraction ratios
選擇計算穩(wěn)定后的1~4 s 的數(shù)據(jù),取時間平均值,計算得到平均阻力系數(shù),如圖8 所示.
圖8 平均阻力系數(shù)隨收縮比 δ變化圖Fig.8 Retraction ratio vs mean drag coefficient
從圖中可以看出,十字型傘的平均阻力系數(shù)隨傘繩收縮比 δ的增大而減小.降落傘的阻力系數(shù)主要受阻力面積的影響,而十字型傘的傘面構型相較于其他傘型,受同等數(shù)量的傘繩長度縮短的影響程度較大,所以傘繩收縮時,平均阻力系數(shù)的減小較為顯著.
從圖7(b)、7(c)可以看出,在5°、10°攻角時,傘繩收縮的十字型傘的瞬時阻力系數(shù)隨時間有較為明顯的波動.這是因為當十字傘處于非零攻角且傘繩同時收縮時,由于傘面的柔性,迎風面積會隨著十字傘自旋發(fā)生變化,圖9 中展示了當攻角為10°、收縮比δ=20.0%時,十字型傘發(fā)生自旋運動時,t1、t2時刻(如圖7(c)中收縮比20%、攻角10°工況,選取某個周期中阻力系數(shù)最大和最小的2 個典型時刻)下,十字型傘的傘面狀態(tài).由圖中可以明顯看出,t1時刻傘面的迎風面積明顯大于t2時刻的迎風面積.
圖9 t1 時刻和t2 時刻的傘面狀態(tài) (收縮比δ=20.0%)δ=20.0%Fig.9 State of canopy at t1 and t2 (retracting ratio )
圖10 所示為0°、5°、10°攻角時,十字型傘的瞬時升力系數(shù)隨時間的變化.
圖10 不同收縮比 δ時十字型傘的升力系數(shù)Fig.10 Lift coefficient of cross parachute vs retraction ratio δ
選擇計算穩(wěn)定后的1~4 s,取時間平均值,計算平均升力系數(shù),如圖11 所示.傘的平均升力系數(shù)基本沒有影響.在5°和10°攻角時,傘繩收縮會引起平均升力系數(shù)的降低.
圖11 平均升力系數(shù)隨收縮比 δ變化圖Fig.11 Retraction ratio δ vs mean lift coefficient
根據(jù)圖10(a)可以看出,0°攻角時,瞬時升力隨著自旋運動呈現(xiàn)出明顯的周期性,升力的方向隨著下拉傘臂A旋轉(zhuǎn)到不同位置而發(fā)生變化;同時從不同收縮比的工況可以看出,十字型傘瞬時升力系數(shù)的峰值與傘繩收縮比 δ密切相關.
由圖10(b)、10(c)所示,攻角為5°、10°時,升力系數(shù)受到攻角和收縮比 δ這2 個因素的共同影響.當傘繩收縮比 δ為0 時,升力主要受攻角影響,方向始終指向+X方向,并且升力為正值;當 δ不為0 時,升力除了受攻角影響之外,還包括來自傘繩收縮的影響.隨著下拉傘臂A的旋轉(zhuǎn),如圖12(b)、12(c)所示,傘臂A的位置從傘軸左側(cè)運動至傘軸右側(cè),導致升力時刻發(fā)生變化.隨著傘繩收縮比 δ的增加,瞬時升力系數(shù)的峰值增大,如t1時刻(在收縮比為20%、攻角
圖12 收縮比δ=0 與δ=20.0%的傘面壓力云圖 (10°攻角)Fig.12 Pressure contour of canopy with δ=0 and δ=20.0% (angle of attack was 10°)
從圖中可以看出,在0°攻角時,十字型傘的平均升力系數(shù)基本為0,傘繩收縮比 δ對0°攻角的十字型為10°工況,選取某個周期升力最大的典型時刻,如圖10(c)所示).隨著傘繩收縮比 δ的增加,瞬時升力系數(shù)的最小值減小.當 δ較大時,收縮傘繩對升力的影響大于攻角對升力的影響,瞬時升力系數(shù)最小值逐漸出現(xiàn)負值,如t2時刻(在收縮比為20%、攻角為10°工況,選取某個周期升力最小的典型時刻,如圖10(c)所示).
當攻角為10°、收縮比δ=0 和δ=20.0%的t1、t2時刻,Y=0 平面的壓力云圖如圖12 所示,速度矢量圖如圖13 所示.
圖13 收縮比δ=0 與δ=20.0%的速度矢量圖 (10°攻角)Fig.13 Velocity vector of canopy with δ=0 and δ=20.0% (angle of attack was 10°)
由圖12 可以看出,在10°攻角、收縮比δ=0工況時,傘衣內(nèi)部的壓強關于傘軸非對稱分布.在收縮比δ=20.0%工況時,t1時刻,收縮傘繩在傘軸左側(cè),加劇了傘面幾何構型的非對稱性,導致高壓區(qū)向傘軸右側(cè)產(chǎn)生了進一步的偏移,此時傘面受到的升力最大;t2時刻,收縮傘繩在傘軸的右側(cè),下拉傘臂A在來流的沖擊下,出現(xiàn)傘面底邊向內(nèi)凹陷的現(xiàn)象,由于內(nèi)凹傘面的阻擋,來流在A處發(fā)生繞流,如圖13(c)所示,傘臂A內(nèi)側(cè)出現(xiàn)漩渦區(qū),出現(xiàn)較為明顯的低壓區(qū),而高壓區(qū)出現(xiàn)在傘臂C內(nèi)側(cè),如圖12(c)所示,此時,傘面受到的升力最小,為負值.
本文對不同傘繩收縮比 δ的十字型降落傘進行了數(shù)值模擬研究,結(jié)果表明:
① 十字型降落傘的阻力系數(shù)隨收縮比 δ的增大而減小.
② 由于十字傘存在自旋運動,下拉傘臂A的位置隨著自旋發(fā)生改變,導致升力具有較強的時變特性.十字型傘瞬時升力峰值隨收縮比 δ增大而增大.
本文的研究內(nèi)容是十字型傘作為主動控制裝置的基礎性研究工作,后續(xù)將對收縮多根(超過2 根)傘繩以及不同傘繩收縮組合下引起的其他非對稱構型對側(cè)向氣動力的影響進行研究,制定出合理的傘繩收縮策略,實現(xiàn)十字型傘對彈箭系統(tǒng)運動姿態(tài)的調(diào)控.