點火沖擊作用下含裂紋固體發動機流固耦合仿真研究

張 澤，高 峰，郜 陽，呂 銳

（空軍工程大學，西安，710051）

0 引 言

固體火箭發動機在生產制造、服役及點火發射過程中受到溫度變化、振動、重力等多種載荷的作用，不可避免地會對推進劑藥柱造成一些損傷。這些損傷缺陷的存在很可能在發動機工作過程中形成“超燃燒表面”，即燃面面積大于設計值，導致發動機內彈道性能偏離設計指標。更加嚴重的情況下，當“超燃燒表面”產生的瞬時壓強峰值超過發動機殼體所能承受的最大壓強時，發動機會發生爆炸。

為了評判和預估含裂紋固體火箭發動機的工作性能，國內外學者對裂紋的擴展及裂紋內的流場分布情況做了相應的試驗及仿真。Kumar［1］通過試驗研究了點火增壓速率、裂紋尺寸及推進劑燃速等因素之間的相互影響關系。邢耀國［2-4］對燃燒條件下裂紋及脫黏面的擴展情況及危害性做了研究。何國強［5］對燃燒增壓過程脫黏擴展條件進行了試驗分析。Kuo等［6］的研究表明，點火升壓速率只有達到或超過1 GPa/s時，點火燃氣才可能對裂紋產生影響，使得裂紋產生失穩擴展。桂曉波［7］在冷流沖擊條件下對藥柱進行流固耦合數值分析，計算結果表明在壓力上升過程中靠近頭部位置的藥柱內孔壁的應力和應變更大。孫博［8］對發動機點火時藥柱不同深度的裂紋進行J積分計算，總結得出裂紋的擴展受裂紋深度的影響，并存在一個臨界值，超過這個臨界值裂紋會發生擴展。

文獻［9］～［12］中采用了不同的手段對推進劑藥柱進行了流固耦合計算，本文在上述學者研究的基礎上，針對點火沖擊過程中推進劑藥柱未被點燃的情況，進行含裂紋缺陷發動機流固耦合仿真研究，并對裂紋內流場壓強和速度的變化對裂紋變形的影響以及裂紋變形對流場的反作用進行了分析。

1 物理模型

圖1為BATES發動機模型，該發動機進行過上千次點火試驗，擁有較高的精度，目前已成為美國空軍火箭推進實驗室（Air Force Rocket Propulsion Laboratory，AFRPL）衡量固體推進劑性能的標準發動機［13］。根據文獻［7］和文獻［8］所得的結論建立本文的裂紋計算模型，其中裂紋寬為2 mm，深度為30 mm，距離燃燒室頭部120 mm。

圖1 含裂紋發動機計算模型Fig.1 Calculation model of cracked motor

2 數學模型

2.1 基本假設

a）將燃氣視為組分凍結的理想氣體，遵循理想氣體狀態方程；

b）點火燃氣從燃燒室頭部進入燃燒室，升壓速率為1 Gpa/s；

c）由于氣體在流動過程中溫度變化幅度較小，不考慮溫度變化對氣體各項參數的影響；

d）僅考慮點火過程生成燃氣對藥柱及裂紋的沖擊作用，不考慮推進劑的燃燒及化學反應作用；

e）將推進劑視為各向同性的線性粘彈性材料，泊松比為常數；

f）不考慮重力對計算過程的影響。

2.2 氣相控制方程

二維不可壓黏性流動N-S方程可表達為

質量守恒方程：

動量守恒方程：

能量守恒方程：

其中，

式中ρ為氣體密度；ui為氣體速度分量；p為氣體壓強；xi為空間中點坐標分量；E為能量；qi為熱流量；τij為切應力分量。

固體火箭發動機燃燒室內燃氣流動內通道流，流動過程屬于高雷諾數流動，因此湍流模型采用可實現型k-ε兩方程模式，其中湍流能量輸運方程和能量耗散輸運方程如下：

式中μt為渦黏性；Sij為平均速度應變率張量；k為湍動能；ε為湍流耗散率；τtij為雷諾應力。

在求解流場時，采用基于密度的有限體積法離散N-S 方程，并使用SⅠMPLE 算法進行求解，對求解區域采用分塊結構化網格進行劃分，并對邊界層進行網格加密處理?？臻g離散采用二階迎風格式，時間離散采用四階龍格-庫塔方法。

2.3 固相控制方程

依據有限元原理，結構動力學方程為

式中 [M]為總體結構質量矩陣；[C]為總體結構阻尼矩陣；[K]為總體結構剛度矩陣；{d?}為節點加速度列矢量；{d?}為節點速度列矢量；g0gggggg為節點位移列矢量；[F(t)]為節點等效載荷列矢量。

對于固體區域，采用任意拉格朗日-歐拉（Arbitrary Lagrangian-Eulerian，ALE）有限元法對結構區域進行離散。對計算區域采用四邊形結構網格劃分，為了提高裂紋附近計算精度，對裂紋附近網格進行加密處理。在計算中，將前端面和外壁面約束固定，后端面和內壁面為流固耦合面，采用Newmark時域連續方法對動態問題進行求解。圖2為裂紋區域網格示意。

選取HTPB推進劑為研究對象，將推進劑作為黏彈性材料處理，材料密度為1 600 kg/m3，泊松比為0.495，表1為松弛彈性常數。

表1 Prony級數參數Tab.1 Prony parameters

粘彈性材料剪切松弛核函數G(t)和體積松弛核函數K(t)采用Prony級數表示為

式中G∞，Gi為剪切模量；K∞，Ki為體積模量；為各Prony級數分量的松弛時間。

G(t)和K(t)可用松弛模量E(t)表示為

將松弛模量E(t)表示為Prony級數形式［14］：

2.4 流固耦合計算方法

流固耦合計算方法通常分為強耦合和弱耦合。強耦合方法在同一控制方程內同時對流體域與固體域進行求解。弱耦合方法分時間步對流體域和固體域分別進行求解，通過兩個區域之間的網格插值實現耦合計算。相比于強耦合方法，弱耦合更容易實現計算收斂，并且節省計算資源，因此本文選取弱耦合的方法進行計算。

由于采用弱耦合方法進行計算，流體域和固體域之間網格節點信息不匹配，在進行插值運算時需要尋找與兩區域耦合面上距離最近的節點進行匹配。節點之間建立匹配之后，將流體區域的壓強等信息傳遞到固體區域，然后求解固體區域的節點位移、速度并通過插值傳遞至流體區。接著采用動網格技術分別對流體域與固體域的網格進行更新，實現耦合面的移動，完成一個時間步的運算，按上述步驟進行迭代運算，直至收斂。

為了保證計算精度，需要對收斂進行判斷，應力收斂判斷表達式為

式中σk+1，σk分別為k+1和k迭代步耦合界面處的結構應力；εu為設定的應力收斂值；ε0為一個非常小的常數。速度和應變的收斂判定與應力形式相同，達到收斂標準后即認為計算結束。

3 計算結果與分析

點火藥被點燃后，生成的燃氣從燃燒室頭部向噴管出口方向欠膨脹噴出，在燃燒室內原有氣體的阻滯作用下，形成激波以及膨脹波，并在燃氣流過區域迅速形成高壓。隨著燃氣不斷向前運動，高壓區域不斷擴大，燃氣前峰到達噴管后，與噴管壁面形成反射波，沿與主通道燃氣相反方向運動，并與主通道燃氣相互作用，圖3為不同時刻燃燒室內壓強分布云圖。

圖3 不同時刻燃燒室內壓強分布云圖Fig.3 Pressure contours in the combustion chamber at different moments

燃氣在流動過程中，在激波的沖擊以及燃氣壓力作用下，推進劑出現變形并產生應力。在t=0.15 ms時刻，燃氣前峰到達裂紋處，與裂紋壁面碰撞形成壓縮波，如圖4所示。隨著燃氣的繼續運動，裂紋腔內開始注入燃氣，裂紋入口處壓強開始快速升高，因此在該時刻，裂紋內最大壓強出現在入口處，并且隨著裂紋深度的增加壓強逐漸降低。不同時刻裂紋內部沿裂紋方向壓強分布如圖5所示。在t=0.25 ms時到達燃氣末端與壁面碰撞后被反射形成反射波，在裂紋末端形成高壓區域，該時刻裂紋腔內最大壓強出現在裂紋末端處，沿裂紋向出口方向逐漸降低。在高壓的作用下，裂紋末端產生的形變大于裂紋內其他區域，形成一個鼓形區域，裂紋變形呈現內寬外窄的形狀。隨著壓強反射波繼續向裂紋入口處運動，在t=0.4 ms時刻裂紋中部壓強開始升高，此時壓強最大區域出現在裂紋中部，向兩邊依次遞減，并且由高壓產生的鼓形區域也隨著反射波向裂紋開口處運動，裂紋基本呈現內外同寬的形狀。在t=0.6 ms時壓強反射波到達裂紋入口處時，形變最大區域也出現在裂紋入口處，此時裂紋呈現內窄外寬的形狀，由于裂紋出口的擴張，裂紋腔內的反射波釋放至燃燒室主流區域，因此出口處壓強略低于裂紋內部壓強，裂紋腔內壓強分布逐漸趨于均勻。不同時刻裂紋位移、應力及應變云圖分別如圖6～7所示。

圖4 燃氣前峰與裂紋作用壓強云圖Fig.4 Pressure contour of gas front and crack action

圖5 不同時刻裂紋內壓強分布曲線Fig.5 Curves of pressure distribution in cracks at different moments

圖6 不同時刻裂紋位移云圖Fig.6 Crack displacement contours at different moments

圖7 不同時刻裂紋應力云圖Fig.7 Crack stress contours at different moments

圖8 為裂紋腔內不同時刻燃氣速度分布。在t=0.15 ms時刻，燃氣開始注入裂紋，裂紋入口處速度較大，隨著深度增加燃氣速度逐漸趨于0。此后裂紋內燃氣速度不斷提高，直至t=0.25 ms 時壓強波到達裂紋末端后形成反射波，由反射波產生負的壓強梯度，使裂紋內燃氣速度開始減小。當t=0.6 ms時，整個裂紋內部燃氣速度幾乎為0，這是因為裂紋內壓強升高，主通道內燃氣進入裂紋受阻。

圖8 不同時刻裂紋內燃氣速度分布曲線Fig.8 Distribution curve of gas velocity in crack at different moments

上述仿真計算所得流場壓強及速度變化規律，與Kumar［15］在試驗及仿真中所得結論基本相符，驗證了本文所采用的計算方法的合理性。在本文的研究基礎上可對工程中所遇到的含裂紋缺陷的固體火箭發動機進行仿真計算，從而對其工作性能進行初步判斷，以便為其延壽或判廢工作提供參考。

4 結 論

本文研究得出以下結論：

a）點火生成燃氣進入燃燒室后，推進劑藥柱受到激波以及高壓的沖擊作用發生變形，變形導致流場邊界發生變化，從而影響燃氣氣流的流動，因此在發動機流場計算中不能忽視固體區域的變形。

b）燃氣進入裂紋后膨脹加速，在裂紋前端達到最大速度，在向末端運動的過程中不斷減速，到達裂紋末端后發生反射形成反射波，造成裂紋內部壓強進一步升高，導致裂紋變形加劇，反射波在運動過程中使裂紋最大變形區域不斷向入口處擴展，最終形成內窄外寬的擴張型裂紋。

c）通過流固耦合方法計算了裂紋在點火燃氣作用下產生的形變，并且分析了裂紋附近產生的應力?？稍诖嘶A上進一步采用斷裂力學的方法對裂紋的擴展問題進行判斷，從而為含裂紋缺陷的固體火箭發動機的工作性能的評估工作提供數據參考。