航天器火工作動裝置流固耦合過程的數值研究

水　龍，楊震春，杜永剛，楊　勇，劉軼鑫

（蘭州空間技術物理研究所 真空技術與物理重點實驗室，蘭州　730000）

航天器火工作動裝置流固耦合過程的數值研究

水龍，楊震春，杜永剛，楊勇，劉軼鑫

（蘭州空間技術物理研究所 真空技術與物理重點實驗室，蘭州730000）

航天器上配備的火工作動裝置用于完成關鍵程序動作與任務，具有很高的可靠性與安全性要求。針對復雜結構火工作動裝置的工作過程，建立非線性流固耦合動力學模型，在推力與拉力負載兩種工況下，采用有限體積法、有限元法進行數值計算，得到火工作動裝置工作過程中流場變化規律與輸出性能。數值模擬結果表明，航天器火工作動裝置流固耦合過程的數值分析能夠模擬其工作過程，火工作動裝置的負載對輸出性能有較大影響。

航天火工裝置；流固耦合；有限元；ALE方法

0　引言

航天器入軌后始終運行在真空環境中，通常配備有多個火工作動裝置，以完成關鍵程序動作與任務［1］。火工作動裝置以裝藥燃燒產生的高溫高壓氣體作為驅動源，將化學能轉化為機械能并輸出線性作動力。在實際工程應用中，對火工作動裝置具有很高的可靠性與安全性要求，且需要考慮到火工沖擊對航天器結構的影響。火工作動裝置的設計主要依賴工程經驗，一般研制程序為經驗/半經驗設計-大量試驗-改進設計。由于火工作動裝置關鍵參數的微調能對輸出性能產生巨大影響［2］，實際研制過程中總是需要進行大量試驗與反復修改設計，導致研制周期長、成本高。另一方面，由于火工作動裝置的體積較小、工作時間極短（毫秒級），對其工作過程性能的全面測試有較大困難。因此，為了提高火工作動裝置的設計水平、縮短研制周期、降低研制成本、全面而準確獲得其工作過程性能，對火工作動裝置工作過程進行數值模擬是十分必要的。

20世紀50年代末，美國、俄羅斯就開始在航天器上采用各種火工裝置，已積累了豐富的工程經驗，并進行了一定的相關理論研究。1993年，Kuo等［3］分析了由NASA標準電起爆器驅動拔銷器的動態特性，分別采用C語音、MESA-2D代碼程序建立了兩個理論分析模型，分析結果能與試驗數據更好地吻合。此后，Goldstein等［4］采用NASA-2D和DYNA 3D軟件對拔銷器和電爆閥門的工作過程進行了動力學仿真分析，為結構受力和變形的研究提供了依據。1994年，Gonthier等［5-6］以NASA標準電起爆器驅動的拔銷器為研究對象，采用LSODE標準程序對所建立的理論模型進行求解計算，對該拔銷器火藥（Zr/KClO4）燃燒過程、活塞運動過程進行了分析。美國的一些專業火工裝置生產廠也開發了自己的性能分析和模擬手段［7］，如Scot公司能夠對火藥燃燒過程、分離作動過程、溫度、壓力等進行計算機模擬仿真，并進行設計優化。

南京理工大學的王濤等［8］基于經典內彈道和氣體動力學理論，建立了二級活塞式拋放彈射機構的理論模型，采用Godnov差分格式對該彈射機構的工作過程進行了數值模擬計算，分析了不同參數對其彈射效果的影響。高濱［7，9］基于經典內彈道理論，建立了火工作動裝置的性能計算模型，利用性能仿真模型對一種彈射裝置進行了分析，計算結果與試驗結果基本吻合。北京理工大學的葉耀坤等［10］對一種用于高速導彈分離系統的楔塊式火工解鎖螺栓動作過程建立了內彈道模型，并利用MATLAB/Simulink進行了仿真計算，可以反映該火工解鎖螺栓的分離運動特性。綜上所述，火工作動裝置的仿真分析模型關注火藥的燃燒過程，采用牛頓第二定律描述活塞的運動過程，均沒有考慮火工作動裝置工作過程中的非線性流固耦合等本質特性。

針對火工作動裝置工作過程中的非定常、高速可壓縮高溫高壓氣體與活塞之間的非線性流固耦合問題，以伸長型火工作動裝置為研究對象，引入任意拉格朗日-歐拉（ALE）方法描述流場控制方程，建立火工作動裝置工作過程的流固耦合系統動力學模型。在不同負載工況條件下，結合有限體積法與有限元法進行求解計算，獲得火工作動裝置的流場變化規律、活塞運動位移和速度等輸出性能。

1　流固耦合系統動力學模型的建立

火工作動裝置燃燒室內火藥燃燒產生的高溫高壓氣體流經腔室后作用于活塞，活塞克服負載開始運動，輸出滿足要求的推力，如圖1所示。火工作動裝置的工作過程是一個復雜的物理、化學變化過程［7］，工作時間極短，涉及到高溫高壓氣體的超音速流動、幾何非線性（活塞的大位移運動）、狀態非線性（筒壁與活塞的接觸）、流固熱多場耦合等重要問題，這是一個強瞬時性、強非線性和強耦合的復雜系統。對火工作動裝置工作過程性能的數值模擬實質就是非線性流固熱耦合系統動力學建模及其求解問題。

圖1　火工作動裝置工作原理圖

1.1非線性偏微分方程組

給出三個假設條件：（1）火藥燃燒產生的高溫高壓氣體為三維非定常、高速可壓縮流動氣體；（2）關于任一單元體，高溫高壓燃氣滿足完全氣體狀態方程；（3）整個工作過程忽略熱量損失。

建立笛卡爾坐標系O-xyz，對于流場內任意一點（x，y，z），在t時刻的速度為v=［u v w］T，密度為ρ，壓力為p，溫度為T；流體的分子粘性系數、第二粘性系數分別為ρ、λ，熱傳導系數為k，單位質量的體積力為f，通過表面的熱通量為q，單位質量的體積加熱率為r˙，單位質量的內能為e。根據質量守恒、動量守恒以及能量守恒定律，可以分別得到連續性方程、動量方程與能量方程為：

這一組控制方程就是守恒形式的納維-斯托克斯（N-S）方程，寫成更簡潔的形式為：

火工作動裝置內火藥燃燒產生的高溫高壓氣體作用于活塞，高溫高壓燃氣流場的變化影響活塞的變形和運動，而活塞的變形和運動又會影響燃氣流場的變化［11］。由于火工作動裝置活塞的運動會導致流場網格的運動，這里采用任意拉格朗日-歐拉（ALE）方法來描述網格的幾何變形與運動。

流場中任一控制體V的控制面為S，控制面的外法向矢量為S，在t時刻動網格的運動速度為w。

將式（2）寫成在ALE坐標系中的積分形式：

式中：Δv=v-w。

為了得到封閉方程組，補充氣體狀態方程：

聯立式（3）、（4），得到一組非線性耦合偏微分方程組，該方程組的求解變量為守恒變量U。

1.2邊界條件和初始條件

火工作動裝置非線性流固耦合系統包括三個求解域：火藥燃燒產生高溫高壓氣體形成的流場，火工作動裝置的筒壁與活塞形成的結構場，流場與結構場的耦合界面。火工作動裝置在工作過程中，其流場各點的速度、壓力、溫度等變量之間是強耦合關系，每個邊界條件都需要對變量進行綜合考慮。注意到火工作動裝置的對稱性，取原模型的1/2作為分析模型，并在截斷面上施加對稱邊界條件。定義火工作動裝置的筒壁作為高溫高壓流場的固定壁面，并采用無滑移條件，即氣體在筒壁上的速度v=0，忽略熱量損失。定義氣體與活塞界面為流固耦合邊界，該邊界需要同時滿足運動學與動力學條件，實現流場與結構場的耦合。結構場內在火工作動裝置的筒壁與活塞之間定義接觸面，在活塞端部施加與實際情況符合的負載，并定義與流場相對應的流固耦合邊界條件。邊界方程離散化后，可以整合到系統動力學模型中進行求解計算。

根據火工作動裝置的裝藥結構，分別定義火藥組與空氣組的初始條件，可以模擬火藥瞬間燃燒產生高溫高壓燃氣的過程。

因此，采用ALE方法描述了火工作動裝置內流場、流場與結構場的耦合過程，并給出流固耦合過程的邊界條件和初始條件，得到火工作動裝置非線性流固耦合系統的動力學模型。

2　流固耦合系統的求解

對于前面得到的流固耦合系統動力學模型，是一個強非線性、強耦合復雜系統，無法直接得到其解析解，結合有限元法和有限體積法進行數值求解，得到變量的變化規律。

對火工作動裝置的結構域、流場域離散化。火工作動裝置的筒體與活塞結構復雜，均采用四面體單元進行自由網格劃分；流場域與結構域互補，也采用四面體單元進行網格劃分。如圖2所示，每個四面體單元包含4個控制體、4個節點。

圖2　四面體單元網格劃分圖

取時間步長為Δt，對時間項的積分采用歐拉方法，取積分參數α=1。流場內任一控制體V，對式兩邊進行積分，并將積分項用平均值表示，整理后可得：

式中：Fn是矢通量F在法矢量n方向的分量，即Fn= nF；同理，Gn=nG。

下面采用有限體積法求解包含對流項和壓力的非粘性項Fn，有限元法求解包含粘性和傳導項的粘性項Gn。

2.1計算非粘性項

對于流場中某一單元體的一個面，定義該三邊形的三個頂點分別為1、2、3，這三個頂點的位移為ri（i=1，2，3）、運動速度為wi（i=1，2，3）。于是有:

流場體積變量為：

式中：w0=（w1+w2+w3）/3。

在t～t+Δt時間段內，該單元面的平均運動速度為：

分別用L、R表示共用一個單元面的兩個控制體，Roe平均變量為：

矢通量增量ΔF用Roe平均變量表示為：

式中：P、P-1為特征矩陣。

其中：Δu=u-wn，c為聲速。

因此，平均非粘性矢通量為：

2.2計算粘性項

采用有限元法計算粘性項，關于控制體的面S積分有：

于是，熱通量可用式（14）計算：

式中：hi為單元形函數。

2.3求解步驟

基于ALE描述的火工作動裝置工作過程流固耦合系統有限元方程，結合有限體積法與有限元法進行數值計算。求解步驟為：

（1）對求解域內所有單元循環求解：

a.計算形函數及其導數；

b.對求解域內所有控制體循環求解：

c.結束控制體循環；

圖3　工作過程中活塞的位移變化圖

e.結束控制體面循環；

（2）結束單元循環。

3　結果與分析

建立火工作動裝置的流場、結構場有限元模型，在結構場分別施加推力與拉力兩種工況的負載，并進行直接耦合求解。值得注意的是，流場有限元模型的網格劃分對計算效率和解的精度有直接關系，對計算收斂性有重要影響。

在推力負載作用下，火工作動裝置工作過程中活塞的位移隨時間的變化如圖3所示，活塞的運動速度變化如圖4所示。

圖4　工作過程中活塞的速度變化圖

由圖4可知，在點火后經過0.02 s，活塞開始運動，活塞的位移隨著時間的增加而增加，且位移與時間近似呈線性關系；試驗結果與分析結果基本一致，驗證了仿真分析的正確性。活塞的速度先波動后趨于穩定，達到穩定狀態后隨著時間的增加逐漸減小。經過0.345 s后活塞運動到位，火工作動裝置的行程為0.106 m，活塞的到位速度為0.29 m/s，活塞到位時會對筒體與周圍結構產生較大沖擊。

在火工作動裝置腔室空間內取點1、點2、點3、點4以及點5，如圖5所示。

圖5　腔室空間內點的示意圖

考察在火工作動裝置工作過程中所選取各點壓力、溫度隨時間的變化關系，分別如圖6、7所示。

圖6　腔室空間內各點的壓力變化圖

圖7　腔室空間內各點的溫度變化圖

由圖6可知，腔室空間內5個點的壓力隨時間的變化關系完全一致，腔室內壓力先增大，經過較小波動后逐漸減小。在0.028 s時刻，腔室空間內壓力達到最大值13.12 MPa。這5條完全重合的曲線說明，在火工作動裝置工作過程的各個時刻，其腔室空間內各點的壓力處處相等，即可以認為腔室空間內壓力瞬時平衡。由圖7可知，腔室空間內5個點的溫度隨時間的變化趨勢是相同的，均先增大后逐漸減小，在活塞運動到位時刻，位于腔室空間后端的點2、點3的溫度明顯高于位于前端的點1、點4，腔室空間中間的點5的溫度介于兩者之間且靠近前端空間點。

在拉力負載作用下，火工作動裝置工作過程中活塞的位移變化如圖8所示，活塞的運動速度變化如圖9所示。

圖8　推力與拉力負載作用下活塞位移變化的對比圖

圖9　拉力負載作用下活塞的運動速度變化圖

由圖9可知，與在推力負載作用下相同，點火后經過0.02 s，活塞開始運動，活塞的位移隨著時間的增加而增加，且位移與時間近似呈線性關系；活塞的速度先波動后趨于穩定，達到穩定狀態后隨著時間的增加逐漸減小。經過0.196 s后活塞運動到位，火工作動裝置的行程為0.106 m，活塞的到位速度為0.56 m/s。與推力負載作用下工況相比，受到拉力負載作用時，火工作動裝置活塞的運動速度更快，到位時間更短。從圖9可以看出，活塞的運動速度在經歷波動后明顯出現突然增大后突然減小的現象，這是由于活塞在與其運動方向相同的拉力作用下開始快速運動，火工作動裝置通過調節給活塞一個反向作用力使得其運動速度減小，因此該火工作動裝置具有較好的負載自適應性。另外，在拉力負載作用下，活塞到位時會對筒體與周圍結構產生更大沖擊。

4　結論

以復雜結構火工作動裝置為研究對象，引入任意拉格朗日-歐拉（ALE）描述方法，建立了火工作動裝置工作過程的非線性流固耦合動力學模型，采用有限體積法、有限元法進行了數值計算。計算結果表明：

（1）采用有限體積法、有限元法對火工作動裝置的工作過程進行數值模擬是有效的，理論模型表征了火工作動裝置工作過程的非線性流固耦合本質特征，通過求解計算得到了其流場變化規律與輸出性能；

（2）火工作動裝置在推力、拉力兩種負載工況下，其輸出性能基本保持一致，具有較強的負載自適應性；

（3）由于火工作動裝置的體積較小，工作時間極短，工作過程中腔室空間內各點壓力變化規律一致，可以認為腔室空間內各點壓力是瞬時平衡的；

（4）火工作動裝置活塞運動到位時，會對筒體產生較大沖擊，進而對火工作動裝置及其周圍結構產生影響。因此，需要考慮采取相應的緩沖措施。

［1］王希季.航天器進入與返回技術（下冊）［M］.北京：中國宇航出版社，1991.

［2］Bement L J，Schimmel M L.A Manual for Pyrotechnic Design Development and Qualification［R］.NASA Technical Memorandum110172，1995.

［3］Kuo J H.Dynamic analysis of NASA Standard Initiator driven pin puller［C］//AIAA，Joint Propulsion Conference and Exhibit，29th，Monterey，CA，1993.

［4］Goldstein S，Lu Y M，Wong T E.Importance of enhanced test dataforcomputermodelingofexposivelyactuateddevices［C］// AIAA，Joint Propulsion Conference and Exhibit，31 st，San Diego，CA，1995.

［5］Gonthier K A，Powers J M.Formulations，predictions，and sensitivity analysis of a pyrotechnically actuated pin puller model［J］.Journal of propulsion and power，1994，10（4）: 501-507.

［6］Gonthier K A，Kane T J，Powers J M.Modeling Pyrotechnic Shock in a NASA Standard Initiator Driven Pin Puller［R］. AIAA94-3054，1994.

［7］高濱.火工分離裝置的性能研究［D］.長沙：國防科學技術大學，2005.

［8］王濤，廖振強，賈安年，等.二級活塞式彈射機構動態仿真與分析［J］.彈道學報，2002，14（4）:36-39.

［9］高濱.火工作動裝置設計參數的敏感性分析［J］.航天返回與遙感，2006，27（3）:57-59.

［10］葉耀坤，嚴楠.低沖擊火工解鎖螺栓的內彈道特性分析［J］.北京工業大學學報，2012，38（9）:1332-1336.

［11］邢景棠，周盛，崔爾杰.流固耦合力學概述［J］.力學進展，1997，27（1）:19-38.

NUMERICAL SIMULATION ON NONLINEAR FLUID-STRUCTURE INTERACTION PROCESS OF PYROTECHNICALLY ACTUATED MECHANISM

SHUI Long，YANG Zhen-chun，DU Yong-gang，YANG Yong，LIU Yi-xin
（Science and Technology on Vacuum Technology and Physics Laboratory，Lanzhou Institute of Space Technology and Physics，Lanzhou730000，China）

The pyrotechnically actuated mechanisms equipped on spacecraft are used for achieving key procedures and missions，and need to satisfy very high reliability and security requirements.This paper focused on the working process of pyrotechnically actuated mechanism with complex structure，the dynamics model considering nonlinear fluid-structure interaction has been established and solved by using finite volume method and finite element method under two different load working conditions.The change of variables in fluid field and output performance of pyrotechnically actuated mechanism has been obtained.The numerical simulation results showed the availability of numerical simulation on pyrotechnicallyactuatedmechanism，andtheoutputperformanceofpyrotechnicallyactuatedmechanismwas relatedtotheloads.

space pyrotechnic device；fluid-structure interaction；finite element analysis；ALE method

V423

A

1006-7086（2015）01-0042-06

10.3969/j.issn.1006-7086.2015.01.010

2014-10-28

水龍（1989-），男，碩士研究生，主要從事非線性多場耦合仿真研究工作。E-mail：shuilonghit@126.com。