高興龍 高慶玉 張青斌 唐乾剛

（國防科技大學(xué)航天科學(xué)與工程學(xué)院，長沙 410073）

1 引言

降落傘空投過程可分為自由墜落、拉直、充氣、穩(wěn)定下降、著陸幾個階段[1]，其中充氣過程的力學(xué)行為最為復(fù)雜，涉及到柔性織物的彈性變形與氣動力的耦合作用，是典型的流固耦合問題。該問題的研究一直是航空航天領(lǐng)域的前沿性課題。由于傘衣結(jié)構(gòu)的透氣性和大變形等力學(xué)特性，使得作用于傘衣表面的氣動力分布較為復(fù)雜，內(nèi)外流場呈現(xiàn)極度不規(guī)則性，依靠傳統(tǒng)的經(jīng)驗理論方法難以對傘衣和周圍流場的動力學(xué)行為進行準(zhǔn)確的預(yù)測。

近年來，伴隨著計算機技術(shù)和計算力學(xué)方法的發(fā)展，數(shù)值模擬技術(shù)在降落傘研究領(lǐng)域得到了廣泛的應(yīng)用。俄羅斯早在20世紀(jì)70年代就將數(shù)值模擬技術(shù)引入到降落傘研究領(lǐng)域。美國的T.Tezduyar教授及其領(lǐng)導(dǎo)的“先進流體仿真建模團隊”[2-4]基于“變空間域/穩(wěn)定時間-空間”(Deforming-Spatial-Domain/Stabilized Space-Time,DSD/SST)格式對十字型傘、密實織物傘的流固耦合計算問題進行了深入研究，并開發(fā)了針對傘衣結(jié)構(gòu)幾何透氣性的均質(zhì)建模方法。法國的政府國防采購部門的Christine Espinosa等人曾利用任意拉格朗日-歐拉（Arbitrary Lagrange-Euler,ALE）有限元方法開發(fā)了SINAP軟件對降落傘流固耦合現(xiàn)象進行了數(shù)值仿真，并基于有限元軟件LS-DYNA對降落傘空投試驗進行了數(shù)值模擬[5]。

相比之下，國內(nèi)研究者早期對降落傘充氣過程的機理研究較多[6-8]，近年來則逐漸側(cè)重于對降落傘的動力學(xué)計算仿真及工程技術(shù)驗證，文獻[9-11]基于LS-DYNA軟件和ALE方法對降落傘開傘過程進行了三維數(shù)值仿真。整體來說國內(nèi)降落傘領(lǐng)域的研究水平近年來得到了顯著提高，但研究傘型對象較單一，傘衣結(jié)構(gòu)簡單，缺乏對復(fù)雜結(jié)構(gòu)降落傘充氣過程的三維可視化建模及仿真技術(shù)的研究。

本文緊跟國際前沿，基于ALE有限元方法對開縫救生傘進行了充氣過程的流固耦合問題研究，數(shù)值模擬了降落傘的充氣過程，對比分析了開縫救生傘的傘衣變形、開傘動載問題，同時給出傘衣結(jié)構(gòu)在變化流場作用下的三維動態(tài)響應(yīng)及傘衣周圍流場的衍變特性，仿真結(jié)果符合降落傘充氣規(guī)律，驗證了本方法的有效性。

2 耦合動力學(xué)模型

本文主要基于ALE描述下的有限元方法建立降落傘結(jié)構(gòu)和周圍流場的動力學(xué)計算模型，考慮降落傘開傘過程的復(fù)雜物理特性，所采用的計算模型做如下假設(shè)：

1）開始充氣前傘衣織物處于軸對稱的折疊狀態(tài)，且預(yù)應(yīng)力為0；

2）開傘過程為無限質(zhì)量情況，忽略重力影響；

3）傘衣流體域入口進氣速度不變；

4）低速流動，流體為不可壓流；

5）流場為來流速度恒定的準(zhǔn)定常流場。

2.1 結(jié)構(gòu)動力學(xué)模型

傘衣織物是柔性大變形體，具有典型的非線性動力學(xué)特性，且為多孔滲透性介質(zhì)的薄膜結(jié)構(gòu)。令為空間結(jié)構(gòu)域，固體邊界由 ?表示，可以寫出結(jié)構(gòu)的控制方程為：

式中 y 是位移矢量；ρs為結(jié)構(gòu)材料密度；f 是作用在結(jié)構(gòu)上的體積力；s為結(jié)構(gòu)體Cauchy應(yīng)力張量；t為積分時間。

2.2 流場動力學(xué)模型

空投過程降落傘開傘速度一般相對較小，可以認為充氣過程流體為不可壓流。在ALE描述下，有限元網(wǎng)格可以在空間域任意移動，則在參考坐標(biāo)系下不可壓縮流體的N-S方程可寫為：

式中 u為流體速度；w 為參考構(gòu)型下網(wǎng)格點的移動速度；ρf為流體密度。

式中 e為物質(zhì)內(nèi)能。

顯然，ALE描述是包含了歐拉描述和拉格朗日描述的[12]。

對于牛頓流體，應(yīng)力張量f定義為：

式中p為壓力；I為二階單位張量；μ為動力粘性系數(shù)。

在確定控制方程后還需要給定邊界條件，Dirichlet和Neumann邊界條件為：

2.3 耦合信息傳遞

在降落傘的整個充氣過程的耦合計算過程中，關(guān)鍵是耦合界面信息的傳遞，本文基于歐拉-拉格朗日描述下的罰函數(shù)法[13]進行流體與傘衣結(jié)構(gòu)間節(jié)點力信息的傳遞，該方法允許非匹配的流體和結(jié)構(gòu)網(wǎng)格。

假設(shè)流體為滲透性介質(zhì)，采用顯式動力學(xué)積分方法，當(dāng)積分時間t=tn（n為積分時間步）時，結(jié)構(gòu)節(jié)點的穿透深度為dn，則下一時刻迭代的結(jié)果為：

式中v為節(jié)點速度矢量；Δt為積分時間步長。

式中vr為主從節(jié)點的相對速度；從節(jié)點速度即結(jié)構(gòu)節(jié)點速度為vs，主節(jié)點速度可看作流體單元節(jié)點在等參坐標(biāo)變化所得的某一流體質(zhì)點速度vf。

3 仿真模型

3.1 開縫傘模型

本文所研究的開縫型救生傘傘衣平面為圓形，整個傘衣被沿輻射方向?qū)ΨQ分布的8條縫分成8個扇形傘衣幅，共28根傘繩?？p隙從距傘衣中心一定距離開始，其中4條縱向貫穿縫開至傘衣底邊，另4條縱向非貫穿縫開到離底邊一定距離。傘頂開有排氣孔，傘衣平面幾何模型如圖1所示，圖中 ，分別為傘頂孔直徑和傘衣結(jié)構(gòu)直徑。

圖1 開縫型傘衣結(jié)構(gòu)平面示意圖Fig.1 Canopy structure of slots-parachute

降落傘傘衣初始模型的好壞將直接影響到數(shù)值模擬的可靠性和精確性，因此需要準(zhǔn)確的復(fù)制出降落傘初始創(chuàng)建時的幾何外形?？紤]到柔性折疊織物建模的復(fù)雜性，本文暫不考慮織物經(jīng)向折疊，僅考慮緯向折疊，即只考慮傘系統(tǒng)拉直后充氣階段的仿真計算模型，圖2（a）為傘衣初始折疊模型。在建立傘衣初始折疊模型時，需要考慮傘縫的布置，這里在幾何建模過程中，在開縫區(qū)域用ANSYS軟件中常用的“生死單元”進行填充，當(dāng)作用于傘衣表面的氣動力達到所設(shè)置的閾值時，則開縫區(qū)域單元被“殺死”，即將其剛度矩陣乘以一個很小的因子，死單元載荷為0，氣流便可通過開縫區(qū)域。圖2（b）為傘衣充滿時傘縫的張開情況。

圖2 開縫型降落傘3D仿真模型Fig.2 3D simulation model of slots-parachute

3.2 流場模型

為避免壁面反射應(yīng)力波對計算結(jié)果造成影響，應(yīng)對照降落傘的尺寸和位置適當(dāng)加大流場尺寸，同時在流出截面設(shè)置無反射邊界條件。本文所建立的流場模型尺寸為30m×30m×80m。同時為了節(jié)省計算資源，并避免傘衣結(jié)構(gòu)對流體的擾動影響，可采用漸變網(wǎng)格劃分流場，即對靠近傘衣結(jié)構(gòu)的流體網(wǎng)格進行局部細化，以獲得較精確的計算結(jié)果。流場網(wǎng)格的截面圖如圖3所示。

圖3 流場網(wǎng)格模型Fig.3 Fluid mesh

4 計算結(jié)果

基于所建立的仿真模型，數(shù)值計算得到傘衣外形變化、充氣時間、開傘動載及流場衍變等結(jié)果，本文針對低速空投過程的開傘過程進行研究，來流速度取40 m/s，大氣密度為1.18kg/m3。自開始充氣至傘衣首次完全充滿時間內(nèi)，根據(jù)傘衣外形變化程度，任意選取其中的6個時間點，得到傘衣充氣過程外形的三維結(jié)構(gòu)動力學(xué)響應(yīng)結(jié)果如圖4所示。

圖4 傘衣充氣過程三維變形圖Fig.4 3D change of canopy inflation

充氣過程中，傘衣的形狀和有效透氣性發(fā)生改變，傘衣投影面積變化如圖5，傘衣最大投影面積為41.2m2。與文獻[14]的無限質(zhì)量情況下傘衣直徑和開傘力隨時間變化曲線相比，本文所計算獲得的開縫傘傘衣外形變化結(jié)果基本符合一般降落傘充氣過程的變形規(guī)律。從充氣過程傘衣投影面積的變化可看出，在氣流流入時出現(xiàn)較明顯的彈性現(xiàn)象，當(dāng)傘衣首次充滿后（即達到最大投影面積時），傘衣充氣便達到了穩(wěn)定狀態(tài)，并無明顯的呼吸現(xiàn)象，這主要是由于開縫傘衣的幾何排氣特性所致。

圖5 傘衣投影面積曲線Fig.5 Change of projected area

圖6 開傘力曲線Fig.6 Change of opening force

圖6為降落傘充氣過程開傘力曲線，開傘力是指在傘衣充氣過程中，降落傘作用在掛重上的力，它是隨著充氣時間變化的，通常在傘衣充滿前達到最大值。本文所計算的最大開傘力為36 092N，文獻[1]給出了無限質(zhì)量情況下最大開傘動載Fkmax的經(jīng)驗公式：

式中Kd稱為動載系數(shù)，對于一定結(jié)構(gòu)型式的傘衣，Kd是常值，這里參考亞聲速范圍的開縫傘取值Kd=1.05[15]；(C A )s為阻力特征，是其阻力系數(shù)與相應(yīng)的特征面積之乘積，這里阻力特征所對應(yīng)的特征面積A0取傘衣的名義面積60m2，換算得到對應(yīng)的阻力系數(shù)C為0.63；ρf為大氣密度；vL為氣流速度。

由圖5和圖6可看出傘衣自初始充氣至首次充滿所用時間為1.41s，結(jié)合經(jīng)驗公式可對降落傘的充氣時間估算，對開縫傘目前還沒有實用的解析方法來計算充滿時間，文獻[15]給出了計算開縫傘充氣時間tm的經(jīng)驗公式：

式中 0.65λg實質(zhì)等于以名義直徑表示的傘衣充滿距離D0；Ks為傘衣名義直徑。由于本文所研究傘型的試驗數(shù)據(jù)資料有限，故參考文獻[15]中無限質(zhì)量條件下降落傘充滿距離中平面有縫傘的試驗數(shù)據(jù)，計算得到傘的充滿時間為1.398s。

充氣過程傘衣周圍流場在流固耦合作用下同樣會產(chǎn)生較復(fù)雜的衍變。圖7為流場計算結(jié)果的流線變化云圖。從圖中可以看出在傘衣的阻流體作用下，沿傘衣外緣氣流流速明顯高于傘衣頂部區(qū)域，并在傘頂出現(xiàn)兩個渦流，且隨著傘衣逐漸充滿，渦流隨之增長，并由對稱分布逐漸衍變?yōu)榉菍ΨQ分布。

圖7 流場流線變化云圖Fig.7 Contour of streamlinechange and fluid velocity

計算結(jié)果表明：本文所采用的仿真技術(shù)有效的模擬了開縫救生傘充氣過程的變化特性，包括較大的開傘沖擊，傘衣的過度膨脹，尤其是傘衣開縫處的計算結(jié)果收斂，未出現(xiàn)不穩(wěn)定或發(fā)散的現(xiàn)象。開縫傘在較高的總透氣量下仍保持著可靠的開傘能力，這是由于開縫的流通特性與織物小孔的流通特性不同的緣故，類似具有顯著的射流效應(yīng)的銳緣小孔的功用，因此出現(xiàn)了較大的流通阻力[15]。傘衣周圍流場衍變特性復(fù)雜，隨著傘衣充滿，傘衣周圍氣流不對稱性加劇，會影響降落傘充氣穩(wěn)定性。

5 結(jié)束語

本文主要針對空投用開縫救生傘的開傘過程進行了研究。針對較為復(fù)雜的開縫傘衣結(jié)構(gòu)，基于ALE有限元方法對開縫救生傘充氣過程的流體和結(jié)構(gòu)的耦合動力學(xué)響應(yīng)進行了數(shù)值模擬，計算得到降落傘充氣過程的傘衣外形和開傘力的變化，驗證了ALE耦合方法求解降落傘流固耦合問題的有效性，尤其是有縫隙傘衣的耦合結(jié)果收斂，對降落傘技術(shù)的工程應(yīng)用具有借鑒意義。本文所采用的技術(shù)同樣可用于其它型降落傘，如回收傘、航彈傘等的動力學(xué)問題研究。

（References）

[1] 王利榮.降落傘理論及應(yīng)用[M].北京:宇航出版社,1997:102-104.WANGLirong.Parachute Theory and Application[M].Beijing:China Astronautic Publishing House,1997:102-104.(in Chinese)

[2] Vinay K,Tayfun E,Tezduyar.A Parallel 3D Computational Method for Fluid-structure Interactions in Parachute Systems[J].Computer Methodsin Applied Mechanicsand Engineering,2000,190:321-332.

[3] Sathe S,Benney A R,Charles R,et al.Fluid-structure Interaction Modeling of Complex Parachute Designs with the Space-time Finite Element Techniques[J].Computers&Fluids,2005,45(7930):1-9.

[4] Kenji Takizawa Tayfun E,Tezduyar.Computational Methodsfor Parachute Fluid-structure Interactions[J].Arch Comput Methods Eng,2012,19:125-169.

[5] Christine E,Yves DL D P,Pascal B,et al.Fluid-structure Interaction Simulation of Parachute Dynamic Behaviour[C].19th AIAA Aerodynamic Decelerator Systems Technology Conference and Seminar,Williamsburg,VA,21-24 May,2007.

[6] 彭勇，張青斌，秦子增.降落傘主充氣階段數(shù)值模擬[J].國防科技大學(xué)學(xué)報,2004,26(2):13-16.PENG Yong,ZHANG Qingbin,QIN Zizeng.Simulation of Parachute Final Inflation Phase[J].Journal of National University of Defense Technology,2004,26(2):13-16.(in Chinese)

[7] 余莉，史獻林，明曉.降落傘充氣過程的數(shù)值模擬[J].航空學(xué)報,2007,28(1):52-57.YU Li,SHIXianlin,MINGXiao.Numerical Simulation of Parachute During Opening Process[J].Acta Aeronautica et Astronautica Sinica,2007,28(1):52-57.(in Chinese)

[8] 潘星，胡利，曹義華.降落傘主充氣階段的動態(tài)仿真及流場分析[J].航空動力學(xué)報,2008,23(1):87-93.PAN Xing,HU Li,CAO Yihua.Analysis of Dynamic Simulation and Fluid Field of Parachute in Inflation Stage[J].Journal of Aerospace Power,2008,23(1):87-93.(in Chinese)

[9] 賈賀，榮偉，陳國良.基于LS-DYNA軟件的降落傘充氣過程仿真研究[J].航天器環(huán)境工程,2010,27(3):367-373.JIA He,RONG Wei,CHEN Guoliang.A Study on Simulation of Parachute Inflation Process Based on LS-DYNA[J].Spacecraft Environment Engineering,2010,27(3):367-373.(in Chinese)

[10] 程涵，余莉，李少騰，等.折疊降落傘展開過程研究[J].航天返回與遙感,2012,33(2):1-5.CHENG Han,YU Li,LI Shaoteng,et al.A Study on the Opening Process of Folded Parachute[J].Spacecraft Recovery and Remote Sensing,2012,33(2):1-5.(in Chinese)

[11] 程涵，余莉，李勝全.基于ALE的降落傘充氣過程數(shù)值仿真[J].南京航空航天大學(xué)學(xué)報,2012,44(3):290-293.CHENG Han,YU Li,LIShengquan.Numerical Simulation of Parachute Inflation Process Based on ALE[J].Journal of Nanjing University of Aeronautics&Astronautics,2012,44(3):290-293.(in Chinese)

[12] HUA Cheng,FANGChao.Simulation of Fluid-solid Interaction on Water Ditching of an AIrplane by ALEMethod[J].Journal of Hydrodynamics,2011,23(5):637-642.

[13] Jason W,Nicolas A,Benjamin T,et al.Porous Euler-lagrange Coupling Application to Parachute Dynamics[C].9th International LS-DYNA Users Conference,2005.

[14] Desabrais K J.Velocity Field Measurements in the Near Wake of a Parachute Canopy[D].Worcester:Worcester Polytechnic Institute,2002.

[15]Ewing EG,Knacke T W,Bixby H W.Recovery Systems Design Guide[M].Beijing:Aviation Industry Press,1988.