基于LS-DYNA軟件的降落傘充氣過程仿真研究

賈 賀，榮 偉，陳國(guó)良

（北京空間機(jī)電研究所， 北京 100076）

0 引言

降落傘系統(tǒng)作為一種高效、可靠的回收裝置，在我國(guó)載人飛船、返回式衛(wèi)星的回收著陸以及美國(guó)火星探測(cè)器的著陸等過程中，都得到了廣泛的應(yīng)用。充氣是降落傘整個(gè)工作過程中最為復(fù)雜的一個(gè)階段，也是研究人員迫切地想通過理論分析來解釋的[1-2]。目前，對(duì)于降落傘充氣性能（除最大開傘動(dòng)載外）的了解或掌握主要是通過一些試驗(yàn)結(jié)果的統(tǒng)計(jì)分析來獲取的，充氣性能與開傘條件之間的關(guān)系也是通過試驗(yàn)結(jié)果統(tǒng)計(jì)、擬合出來的[3]。這使得充氣時(shí)間、傘衣投影面積變化的一些經(jīng)驗(yàn)公式存在著一定的局限性，只能在某些情況下適用。

隨著“嫦娥一號(hào)”的成功發(fā)射，我國(guó)對(duì)于火星探測(cè)的研究即將進(jìn)入實(shí)質(zhì)階段?；鹦谴髿鈱拥某煞?、物理性質(zhì)與地球的大氣層存在著較大的差別，因此，對(duì)于降落傘系統(tǒng)在火星大氣環(huán)境下的工作性能評(píng)定，不能像地球上回收航天器所用的降落傘系統(tǒng)那樣可以通過空投試驗(yàn)來直接驗(yàn)證，而需要通過理論分析或者數(shù)值仿真來進(jìn)行。目前回收系統(tǒng)工作過程的數(shù)值仿真系統(tǒng)所用到降落傘開傘過程的一些參數(shù)，如充氣時(shí)間、傘衣投影面積的變化一般是采用經(jīng)驗(yàn)公式來處理的。有關(guān)研究表明，對(duì)于火星這樣稀薄的大氣環(huán)境，這些常用的經(jīng)驗(yàn)關(guān)系式將不再適用[4]。因此，如何通過數(shù)值仿真來獲取研究降落傘充氣過程性能的主要參數(shù)，已經(jīng)被提上了降落傘理論研究和工程實(shí)踐的日程。

進(jìn)入21世紀(jì)以來，應(yīng)用LS-DYNA有限元分析軟件來分析降落傘的工作過程（充氣階段、穩(wěn)定階段等）成為一個(gè)熱點(diǎn)。國(guó)外在這方面做了大量的工作，但是對(duì)降落傘充氣過程的仿真還不能完全地解決[5-8]。國(guó)內(nèi)使用 LS-DYNA軟件對(duì)降落傘充氣過程仿真研究尚屬空白。本文基于LS-DYNA軟件對(duì)數(shù)值模擬降落傘充氣過程的方法進(jìn)行了初步研究。

1 降落傘的充氣過程

降落傘的充氣過程是指從傘系統(tǒng)全長(zhǎng)地拉直到傘衣第一次充滿（傘衣投影直徑第一次能達(dá)到穩(wěn)定下降所需的傘衣投影直徑，不考慮過度充氣）為止的整個(gè)工作過程[1-2]，如圖1所示。

圖1 降落傘的充氣過程Fig. 1 Parachute inflation process

降落傘的充氣性能主要由臨界開傘速度、傘衣充滿條件、充氣時(shí)間、充氣距離、傘衣投影面積變化和最大開傘動(dòng)載等參數(shù)來描述和表征[1,3]。本文基于LS-DYNA軟件對(duì)降落傘的充氣時(shí)間和傘衣投影面積變化的數(shù)值模擬分析，來研究降落傘的充氣過程。

2 LS-DYNA的流固耦合算法

LS-DYNA軟件是著名的顯式動(dòng)力分析程序，具有強(qiáng)大精確的有限元仿真性能，擁有大量不同種類的單元模型、材料模型和算法選擇，能夠很方便地處理各種高度非線性問題。到目前為止，它是對(duì)降落傘充氣過程進(jìn)行仿真的最優(yōu)工具之一。LS-DYNA 971是一個(gè)主要的程序版本。它增加了透氣性材料的流固耦合算法，這樣就可以通過控制流固耦合過程的各個(gè)參數(shù)，真實(shí)地模擬降落傘的充氣過程。

LS-DYNA與一般的CAE輔助分析程序操作過程相似，LS-DYNA的一個(gè)完整的顯式動(dòng)力分析過程包括前處理、求解以及后處理3個(gè)基本操作環(huán)節(jié)，如圖2所示。

圖2 LS-DYNA的一般計(jì)算流程Fig. 2 Computation process of LS-DYNA

在解決流固耦合問題時(shí)，需要一種將Lagrange和 Euler的優(yōu)點(diǎn)結(jié)合起來的方法，即Arbitrary Lagrange-Euler，簡(jiǎn)稱 ALE方法[9,10]。LS-DYNA程序中采用關(guān)鍵字*CONSTRAINED_ LAGRANGE_IN_SOLID實(shí)現(xiàn)固體與流體間的耦合。采用流固耦合算法模擬問題時(shí)，往往要對(duì)Lagrange算法中的固體結(jié)構(gòu)進(jìn)行約束，將固體結(jié)構(gòu)的相關(guān)參數(shù)傳遞給流體單元。在LS-DYNA中數(shù)值模擬降落傘的充氣過程時(shí)，固體結(jié)構(gòu)的主要約束方法是罰函數(shù)耦合約束。對(duì)于罰函數(shù)約束而言，其原理是耦合系數(shù)追蹤Lagrange節(jié)點(diǎn)（結(jié)構(gòu)，即從物質(zhì)）和Euler流體物質(zhì)（主物質(zhì)）位置間的相對(duì)位移d，如圖3所示[11]。

圖3 罰函數(shù)耦合算法Fig. 3 Penalty function coupled method

3 LS-DYNA仿真模型的建立

3.1 降落傘仿真模型的建立

3.1.1 基本假設(shè)

對(duì)降落傘仿真模型作如下假設(shè)：

1）仿真充氣過程前模型為軸對(duì)稱的結(jié)構(gòu)；

2）模型的初始形狀為有一定的進(jìn)氣孔尺寸，截面為梅花形狀的模型；

3）模型只考慮傘繩的連接，對(duì)于傘衣徑向帶、頂孔繩等不作具體處理。

3.1.2 基本參數(shù)設(shè)置

1）降落傘幾何模型

降落傘仿真模型選用美國(guó)的經(jīng)典傘型—C-9平面圓形傘，如圖4(a)所示，具體結(jié)構(gòu)參數(shù)如表1所示。

表1 C-9傘結(jié)構(gòu)參數(shù)表Table 1 Structure parameters of C-9 parachute

完整的傘衣和傘繩幾何模型如圖4(b)所示。

2）材料參數(shù)的選擇

降落傘仿真模型的傘衣材料選用美軍標(biāo)的MIL-C-7020 III材料，假設(shè)該材料為各向同性材料，密度為533.77 kg/m3，彈性模量為0.4309 GPa，泊松比為0.14，傘衣厚度為0.000 1 m。傘繩材料的密度為462.00 kg/m3，彈性模量為97.0 GPa。假設(shè)降落傘傘繩的橫截面為圓形，面積為4.91×10-6m2。

3）網(wǎng)格的劃分

LS-DYNA的前處理軟件用于對(duì)降落傘仿真模型進(jìn)行有限單元的網(wǎng)格劃分。傘衣單元選用薄殼單元，傘繩單元選用離散的梁/索單元，得到降落傘的有限元模型如圖4(c)所示。

圖4 降落傘模型Fig. 4 Parachute model

4）K文件的參數(shù)設(shè)置

降落傘仿真模型的傘衣單元為薄殼單元（*ELEMENT_SHELL_THICKNESS），單元算法的控制關(guān)鍵字為*SECTION_SHELL。傘衣材料屬性由*MAT_FABRIC關(guān)鍵字來定義。傘繩為離散的梁/索單元（*ELEMENT_BEAM），單元算法的控制關(guān)鍵字為*SECTION_BEAM。傘繩材料屬性由*MAT_CABLE_DISCRETE_BEAM來定義。

傘繩匯交點(diǎn)采用關(guān)鍵字*BOUNDARY_ SPC_NODE來進(jìn)行約束，即對(duì)匯交點(diǎn)的3個(gè)方向平動(dòng)自由度進(jìn)行約束。

5）求解文件的生成

降落傘仿真模型建立完成，生成LS-DYNA軟件可以求解的文件格式，記為K文件(parachute.k)。

3.2 流場(chǎng)仿真模型的建立

3.2.1 基本假設(shè)

對(duì)流場(chǎng)仿真模型作如下假設(shè)：

1）流場(chǎng)模型模擬的為不可壓縮流場(chǎng)，符合N-S方程；

2）流場(chǎng)模型中流體的速度方向恒定，且速度大小不變；

3）流場(chǎng)模型外表面的邊界約束條件與風(fēng)洞壁性質(zhì)一致。

3.2.2 基本參數(shù)設(shè)置

1）流場(chǎng)幾何模型

流場(chǎng)仿真模型形狀選用圓柱形流場(chǎng)，如圖5(a)所示。綜合文獻(xiàn)[12?14]，其基本的幾何尺寸確定如下：流場(chǎng)半徑 Df=12 m；流場(chǎng)長(zhǎng)度 Lf=31 m；流場(chǎng)前端面距傘繩匯交點(diǎn)的長(zhǎng)度Lq=6 m；未充氣流場(chǎng)后端面距傘衣頂端的長(zhǎng)度Lh=12 m。

運(yùn)用LS-DYNA的前處理軟件建立流場(chǎng)幾何模型如圖5(b)所示，流體速度為80 m/s。

2）材料參數(shù)的選擇

流場(chǎng)仿真模型的流體定義為理想氣體，滿足理想氣體狀態(tài)方程。氣體的密度為1.18 kg/m3，動(dòng)態(tài)粘度系數(shù)為1.745 6×10-5，氣體的比熱容比為γ=1.4。

3）網(wǎng)格的劃分

同樣，在LS-DYNA前處理軟件中對(duì)流場(chǎng)的幾何模型進(jìn)行有限單元的網(wǎng)格劃分。根據(jù)有限元網(wǎng)格劃分的基本原則，流場(chǎng)的有限元模型如圖 5(c)所示。

圖5 流場(chǎng)模型Fig. 5 Fluid model

4）K文件的參數(shù)設(shè)置

流場(chǎng)仿真模型的單元選擇體單元（*ELEMENT_ SOLID），流場(chǎng)仿真模型單元算法的關(guān)鍵字為*SECTION_SOLID。LS-DYNA要想真實(shí)地模擬降落傘的充氣過程，必須在流場(chǎng)前端增加一層流體壓力入口單元，主要作用是向流場(chǎng)中源源不斷地輸入恒壓、恒速的流體。流場(chǎng)單元算法參數(shù)的選擇對(duì)于解決流固耦合問題很關(guān)鍵，本文為流體單元選擇的是單點(diǎn)ALE方法的多物質(zhì)單元方程。對(duì)于已定義關(guān)鍵字*MAT_NULL的材料而言，必須包含狀態(tài)方程的定義，用來提供壓力行為應(yīng)力組件，表示材料的體積變形與壓力之間的關(guān)系。這樣兩者合起來提供材料的整個(gè)應(yīng)力張量。流體仿真模型的狀態(tài)方程用*EOS_LINEAR_POLYNOMIAL來定義。

流場(chǎng)仿真模型的流場(chǎng)入口單元區(qū)域用*BOUNDARY_PRESCRIBED_MOTION_SET來施加對(duì)單元的強(qiáng)制運(yùn)動(dòng)，使得流體壓力入口單元獲得恒定的速度。

5）求解文件的生成

流場(chǎng)仿真模型建立完成，最后生成LS-DYNA軟件可以求解的文件格式，記為K文件（fluid.k）。

4 LS-DYNA流固耦合控制參數(shù)

LS-DYNA軟件是一種比較理想的顯式求解工具，通過它的前處理器，生成一個(gè)可以遞交給LS-DYNA的文件，稱之為關(guān)鍵字文件或K文件，該文件是由一系列的關(guān)鍵字組成的[15]。LS-DYNA求解流固耦合問題，實(shí)質(zhì)上就是對(duì)流固耦合控制關(guān)鍵字進(jìn)行參數(shù)設(shè)置的過程。

但是，降落傘充氣過程與一般的流固耦合問題不完全相同，如果按照一般的流固耦合問題處理，僅僅設(shè)置一般問題的基本控制關(guān)鍵字來仿充氣過程，可以發(fā)現(xiàn)傘衣沒有完全充開，最終的充氣形狀也很不理想，如圖 6所示。可見，使用LS-DYNA對(duì)充氣過程進(jìn)行仿真時(shí)，不能僅僅采用一般問題的解決方法和思路，必須選擇適合充氣過程的控制關(guān)鍵字，并對(duì)關(guān)鍵字進(jìn)行參數(shù)值的優(yōu)化處理。

圖6 基本控制參數(shù)下的傘衣最終充氣形狀Fig. 6 Final canopy shape with basic control parameters

通過大量的研究分析，得到在仿真降落傘充氣過程中，主要的流固耦合過程控制關(guān)鍵字及最優(yōu)化的控制參數(shù)值，其中主要的控制關(guān)鍵字如表2所示[16]?？梢姡x擇適合降落傘充氣過程的控制關(guān)鍵字及設(shè)置最優(yōu)化的參數(shù)值是仿真充氣過程的重點(diǎn)和難點(diǎn)，對(duì)仿真過程的實(shí)現(xiàn)起著決定性作用。

表2 流固耦合過程的主要控制關(guān)鍵字Table 2 Major keywords of Fluid-Structure Interaction(FSI) process

5 降落傘充氣過程的仿真結(jié)果及分析

5.1 降落傘充氣過程的數(shù)值模擬

圖7中(a) (i)表示的是在LS-DYNA中數(shù)值模擬降落傘充氣的全過程。

可以很清楚地看到，在充氣開始的時(shí)候，降落傘的充氣形狀變化比較慢。在0.02 s時(shí)，氣流使得降落傘的傘衣頂部開始發(fā)生形狀的變化；從0.02 0.10 s，傘衣頂部變化顯著，氣流將整個(gè)降落傘的頂部完全頂開。從0.10 0.16 s，氣流從傘衣頂部向傘衣周圍流動(dòng)，使得傘衣的充氣形狀不斷發(fā)生變化；在0.16 s時(shí)，降落傘充氣達(dá)到最大外形形狀，即達(dá)到充滿狀態(tài)。從0.16 0.18 s降落傘的充氣形狀發(fā)生一次較為明顯的收縮。從0.18 0.26 s降落傘的充氣形狀會(huì)不斷發(fā)生微小的變化，傘衣產(chǎn)生較為明顯的“呼吸現(xiàn)象”。通過對(duì)傘衣充氣形狀的測(cè)量，發(fā)現(xiàn)在0.26 s時(shí)降落傘的充氣過程基本完成。

圖7 LS-DYNA數(shù)值模擬的降落傘充氣過程Fig. 7 Parachute inflation simulation with LS-DYNA

5.2 降落傘的充氣時(shí)間

5.2.1 充氣時(shí)間的理論求解

基于 LS-DYNA對(duì)降落傘充氣過程的數(shù)值模擬，其實(shí)就是對(duì)降落傘做類似于風(fēng)洞試驗(yàn)的研究。因此，數(shù)值模擬降落傘的充氣過程是一個(gè)無限質(zhì)量條件下的充氣過程?？赏ㄟ^理論方法求解無限質(zhì)量條件下C-9平面圓形傘的充氣時(shí)間tf。

[17，18]可知，對(duì)于無限質(zhì)量條件下的致密織物的降落傘充氣過程而言，降落傘的充氣時(shí)間指的是傘衣第一次達(dá)到最大外形狀態(tài)的時(shí)間。

設(shè)降落傘的名義直徑為D0，有效透氣量為C，開始充氣時(shí)的速度為 v0，無因次充氣時(shí)間為 T0，對(duì)于基本方程

而言，當(dāng)考慮無限質(zhì)量條件充氣時(shí)，可以把速度 v看作一個(gè)常數(shù)，并且等于充氣開始時(shí)的速度v0，則基本方程可以直接積分。把上述方程寫成積分形式為

如果只考慮名義直徑為D0的平面圓形傘，并假設(shè)降落傘傘衣在完全充滿時(shí)呈半球狀，則有

式中Vmax為降落傘傘衣充滿時(shí)所包容的空氣體積。

對(duì)式(2)進(jìn)行積分，并把式(3)代入，得到充氣時(shí)間的表達(dá)式

已知降落傘的名義直徑D0、有效透氣量C和降落傘開始充氣時(shí)的速度v0，就可以直接計(jì)算出在無限質(zhì)量條件下降落傘的充氣時(shí)間 tf。數(shù)值模擬充氣過程的C-9傘的名義直徑D0為8.5 m，開始充氣時(shí)的速度v0為80 m/s，故求解降落傘的充氣時(shí)間轉(zhuǎn)化為求解其有效透氣量的問題。由于C-9傘傘衣材料選用的MIL-C-7020 III型材料，由圖8所示的曲線，可以得到其壓差和有效透氣量的關(guān)系。

圖8 MIL-C-7020 III材料壓差與有效透氣量曲線Fig. 8 MIL-C-7020 III material pressure vs effective permeability

由于模擬無限質(zhì)量條件下的充氣過程可看作為風(fēng)洞試驗(yàn)，大氣密度ρ為1.18 kg/m3，則壓差為

由圖8可知，C-9傘傘衣的有效透氣量C約為0.09。由式(4)計(jì)算出其在無限質(zhì)量條件下充氣時(shí)間為tf=0.21 s。

5.2.2 充氣時(shí)間的數(shù)值模擬值與理論求解值的比較分析

數(shù)值模擬 C-9傘傘衣第一次達(dá)到最大外形狀態(tài)的充氣時(shí)間為0.16 s，這與理論求解的充氣時(shí)間0.21 s的偏差為23%，是在可接受的范圍之內(nèi)。并且數(shù)值模擬充氣過程時(shí)，可以求解得到傘衣最終達(dá)到穩(wěn)定狀態(tài)時(shí)的充氣時(shí)間為0.26 s。

5.3 降落傘的投影面積變化

5.3.1 投影面積變化的數(shù)值模擬

數(shù)值模擬的 C-9傘傘衣底邊投影面積的形狀變化如圖9所示。

圖9 傘衣底邊投影形狀的變化Fig. 9 Canopy project shape variations

5.3.2 數(shù)值模擬與試驗(yàn)結(jié)果的比較分析

從充氣過程開始到傘衣最終達(dá)到穩(wěn)定狀態(tài)，數(shù)值模擬 C-9傘充氣過程的傘衣投影面積比隨充氣時(shí)間的變化情況。如圖10(a)所示，C-9傘充氣過程的傘衣投影面積比隨無因次充氣時(shí)間的變化如圖10(b)所示。

圖10 傘衣投影面積比的變化Fig. 10 Canopy project area ratio variations

由圖 10(b)中可以看出，傘衣投影面積隨著充氣時(shí)間的延長(zhǎng)而不斷增加，在T=1時(shí)達(dá)到最大值，投影面積與傘衣名義面積之比為 0.42。而在 T=1之后傘衣有一個(gè)明顯的收縮過程。而從圖 10(a)可以更加清楚地看到，在tf=0.21 s時(shí)，傘衣的投影面積變化趨勢(shì)基本穩(wěn)定，雖然由于傘衣的“呼吸現(xiàn)象”會(huì)使得傘衣的投影面積產(chǎn)生一些微小的變化，但是在tf=0.21 0.26 s之間，可以認(rèn)為傘衣的投影面積不再發(fā)生變化，即降落傘達(dá)到充滿穩(wěn)定狀態(tài)。

與 D0=8.5 m的平面圓形傘試驗(yàn)結(jié)果[1]對(duì)比可以發(fā)現(xiàn)，數(shù)值模擬C-9平面圓形傘充氣過程的傘衣投影面積變化規(guī)律與試驗(yàn)得到的變化規(guī)律基本一致。在達(dá)到充氣時(shí)間時(shí)，試驗(yàn)的傘衣投影面積與名義面積比為0.4，數(shù)值模擬的傘衣底邊投影面積與名義面積比為0.42，面積比基本相等?？梢?，數(shù)值模擬降落傘充氣過程，得到的傘衣投影面積的變化情況與試驗(yàn)得到的結(jié)果基本一致。

6 結(jié)束語(yǔ)

本文基于 LS-DYNA軟件對(duì)降落傘的充氣過程進(jìn)行了數(shù)值模擬，得到了降落傘的充氣時(shí)間和投影面積變化的規(guī)律，并對(duì)這些主要參數(shù)進(jìn)行了分析，驗(yàn)證了使用 LS-DYNA軟件模擬充氣過程的可行性，為降落傘的設(shè)計(jì)及其性能分析提供了一種新的分析手段。由于流固耦合過程控制關(guān)鍵字的參數(shù)值還沒有完全達(dá)到最優(yōu)化，導(dǎo)致降落傘充氣時(shí)間的誤差較大，這也是下一步工作需解決的主要問題。

參考文獻(xiàn)（References）

[1] 王利榮. 降落傘理論與應(yīng)用[M]. 北京宇航出版社, 1997

[2] 《降落傘技術(shù)導(dǎo)論》編寫組. 降落傘技術(shù)導(dǎo)論[M]. 北京: 國(guó)防工業(yè)出版社, 1977

[3] Ewing E G, Bixby H W, Knacke T W. 回收系統(tǒng)設(shè)計(jì)指南[M]. 北京: 航空工業(yè)出版社, 1988

[4] 榮偉. 火星探測(cè)器減速著陸技術(shù)研究[D]. 中國(guó)空間技術(shù)研究院博士學(xué)位論文. 北京: 中國(guó)空間技術(shù)研究院, 2008

[5] Taylor A P. An investigation of the apparent mass of parachutes under post inflation dynamic loading through the use of fluid structure interaction simulations, AIAA 2003-2104[R], 2003

[6] Benjamin A T, Anthony P T. The use of LS-DYNA to simulate the inflation of a parachute canopy, AIAA 2005-1608[R], 2005

[7] Benjamin A T, Anthony P T, Jean C B, et al. The use of LS-DYNA to assess the performance of airborne system North America candidate ATPS main parachutes, AIAA 2005-1609[R], 2005

[8] Benjamin A T. The application of a new material porosity algorithm for parachute analysis[C]//9thInternational LS-DYNA Users Conference, 2006

[9] 李裕春, 等. ANSYS 11.0/LS-DYNA 基礎(chǔ)理論與工程實(shí)踐[M]. 北京: 中國(guó)水利水電出版社, 2008

[10] 趙海鷗. LS-DYNA 動(dòng)力分析指南[M]. 北京: 兵器工業(yè)出版社, 2003

[11] John O H. Theoretical Manual[M]. Livermore Software Technology Corporation(LSTC), 1998

[12] Lingard J S, Darley M G. Simulation of parachute fluid structure interaction in supersonic flow, AIAA 2005-1607[R], 2005

[13] Lingard J S, Darley M G, Underwood J C. Simulation of Mars supersonic parachute performance and dynamics, AIAA 2007-2507[R], 2007

[14] Christine Espinosa, Yves de Lassat de Pressigny, Pascal Bordenave, et al. Fluid-structure interaction simulation of parachute dynamic behaviour, AIAA 2007-2510[R], 2007

[15] Livermore Software Technology Corporation(LSTC). LS-DYNA keyword user’s manual[G]. Version 971, 2007

[16] 賈賀. 基于LS-DYNA的降落傘充氣過程研究[D]. 中國(guó)空間技術(shù)研究院碩士學(xué)位論文. 北京: 中國(guó)空間技術(shù)研究院, 2009

[17] 彭勇, 宋旭民. 降落傘充氣時(shí)間的計(jì)算方法[J]. 航天返回與遙感, 2004, 25(1): 17-20

[18] 返回技術(shù)翻譯組. 氣動(dòng)力減速器原理及設(shè)計(jì)手冊(cè)(上冊(cè))[M]. 北京: 國(guó)防工業(yè)出版社, 1976