雙脈沖發(fā)動機中金屬膜片式隔艙設計方法①

劉偉凱，惠 博

(1.西北工業(yè)大學燃燒、熱結構和內流場重點實驗室，西安 710072;2.航天動力技術研究院，西安 710025)

0 引言

雙脈沖固體火箭發(fā)動機是固體火箭發(fā)動機的一項重大研究成果，其脈沖藥柱具有獨立的點火器，可以單獨點燃，II脈沖藥柱的點火時間可有不同的延遲［1］。目前，其推力向量控制和推力終止技術已得到較滿意的解決，但仍存在一些問題需要進一步研究。文獻［2］設計了一種軟質隔層結構，利用擴展有限元技術研究了其在脈沖發(fā)動機中的承壓和破壞過程，通過單項試驗驗證了數(shù)值仿真的準確性，說明其結構承壓、打開及密封性能均滿足設計要求，其結構可以應用于實際的脈沖發(fā)動機之中。文獻［3－6］主要研究了雙脈沖發(fā)動機的內流場特點:II脈沖工作時，由于級間通道的收縮導致燃氣在I脈沖燃燒室產(chǎn)生后臺階流動，從而使燃氣產(chǎn)生漩渦，強化了I脈沖的對流換熱及粒子沖刷。

金屬膜片式隔艙是利用輪輻式支撐件來減小重量，密封膜片與支撐件緊密貼實。為了減小打開壓強，在金屬膜片一側設有預制缺陷槽，同時在金屬膜片外側附著一層絕熱層進行絕熱。由于該類隔艙兼具有結構設計簡單、研制周期短、可靠性高等優(yōu)點，在國際上被廣泛應用于脈沖發(fā)動機領域［1］。

本文通過圓板大撓度理論和斷裂力學理論推導出金屬膜片預制缺陷處應力強度因子的計算式，從而得到金屬膜片在內壓作用下的設計方法。通過三維虛擬裂紋閉合法數(shù)值計算了預制缺陷處的應力強度因子并與公式計算結果對比，驗證公式的準確性;通過推導多孔圓板關心區(qū)域的應力大小，得到支撐件強度校核方法。建立三維有限元模型，計算關心點的應力大小，驗證校核方法的可靠性。最后通過實際的熱流試驗，驗證了文中推導方法的有效性。

1 隔艙設計方法

1.1 隔艙物理模型

本文的金屬膜片式隔艙與文獻［7］保持一致，具體如圖1所示，隔艙組件包括一個多孔支撐件和一個高強度易變形的金屬模片。為了得到可靠的打開形式，預制缺陷槽(V型槽)一般設計為“十字型”或“米字型”，本文選取“米字型”預制缺陷膜片作為研究對象。為減輕支撐件重量以及增加通氣面積，將支撐件設計為輪輻式結構。

圖1 金屬膜片式隔艙結構Fig.1 Schematic diagram of metal diaphragm PSD

1.2 金屬膜片設計式推導

為推導膜片預制缺陷處的應力強度因子，首先進行了如下假設:

(1)由于本文所設計金屬膜片的預制缺陷槽為均勻放射狀，各缺陷槽尺寸、受力狀態(tài)完全相同，忽略各缺陷槽之間的影響，取其中任一條缺陷作為研究對象;

(2)結合前期摸底試驗發(fā)現(xiàn)，本文的金屬膜片破壞形式均為I型裂紋擴展破壞，因此本文只求解預制缺陷處的I型應力強度因子KI;

(3)本文選取垂直于缺陷槽的任一截面作為研究對象，由于垂直于截面上的應力對應力強度因子不產(chǎn)生影響，故將該截面簡化為二維板條結構［8］。

(4)由于膜片實際破壞過程為瞬間動態(tài)過程，材料來不及發(fā)生塑性屈服，假設在破壞過程中材料都表現(xiàn)為線彈性屬性，因此近似認為膜片的動態(tài)破壞過程為線彈性斷裂問題。

為了量化膜片的破壞打開壓強與膜片結構尺寸之間的關系，其中包括預制缺陷深度a、缺陷V型槽開口角度α、膜片半徑R、膜片厚度h等。

本文將膜片承受燃燒室壓強的變形過程簡化為相同尺寸薄板(不含預制缺陷)的大撓度問題，圓板大撓度微分方程組如下［8］:

利用伽遼金法求解均布載荷下大撓度圓板應力場分布。設圓板的撓度符合式(3)，把代入基本微分方程組(2)中，并結合圓板周邊固定邊界條件，得圓板中心的最大撓度計算式(4)。

圓板任意點的切向應力包括兩部分，即薄膜應力和彎曲應力，具體見式(5):

式中各參數(shù)具體含義與文獻［7］相同。

通過式(5)即可計算得到預制缺陷相應位置的應力值，簡化后截面應力分布如圖2所示，利用疊加原理得到預制缺陷處應力強度因子的計算方法［7］。

通過疊加原理可得到中心處應力強度因子的計算式:

具體如下:

式中各參數(shù)具體含義與文獻［7］相同。由于板中心的應力最大，所以預制缺陷中心位置的應力強度因子也最大，應力強度因子的分布形式如圖3所示。

圖2 應力強度因子疊加原理圖Fig.2 Superimposed schematics of stress intensity factor

圖3 預制缺陷處應力強度因子的分布Fig.3 Stress intensity factor distribution of the prefab defect

1.3 支撐件強度校核方法

支撐件為“多孔+輪輻”形式，計算其應力分布時需要將其簡化為“當量無孔圓板”［9］，此時需要對材料的彈性常數(shù)(E、μ)進行修正［9］。其應力的計算方法是:先應用一般實心圓平板的應力和變形公式，計算當量實心圓平板的應力和變形。再將當量實心圓平板的應力，乘以應力乘數(shù)和面積削弱系數(shù)，就可以得到多孔板的應力強度實際值，并以此作為強度計算的依據(jù)。

由于支撐件厚度較厚，承壓變形過程中撓度較小，屬于圓板的小撓度問題。首先將支撐件簡化為周邊簡支的當量圓板，在均布壓力P作用下，其上下表面的應力為

應力計算時，沿最小管孔帶寬度方向取平均值，但并不沿板厚度方向取平均值，這時的有效應力為

式中 σave為由于外壓載荷作用而產(chǎn)生的徑向應力σr和環(huán)向應力σθ中最大者;κ=S0/S為實心圓板關心區(qū)域與有孔板關心區(qū)域面積與之比;K為應力乘數(shù)，無因次量［9］。

則支撐件強度校核方法:有效應力［σ］小于材料的強度極限σb。

2 數(shù)值計算驗證

2.1 應力強度因子計算方法

虛擬裂紋閉合技術(VCCT)［10］最先由 Rybicki和Kanninen提出，用來計算二維裂紋問題，它通過對有限元分析結果進行后處理得到所需要的裂紋擴展的能量釋放率;后被Shivakumar等推廣至三維裂紋，在應用三維虛擬裂紋閉合法時，為了保證計算二維裂紋擴展的方法能夠直接推廣應用至三維裂紋，裂紋前緣處的網(wǎng)格應采用六面體單元。以有限塊體中的裂紋(圖4)為例介紹三維虛擬裂紋閉合法的原理，有限塊體裂紋長度為a，裂紋擴展長度為Δa，網(wǎng)格采用8節(jié)點六面體單元，裂紋上表面節(jié)點和下表面節(jié)點一一對應，具有相同的坐標，則 I型、II型、III型應變能釋放率 GI、GII、GIII可通過式(10)計算得到。

式中 XLi、YLi、ZLi分別為節(jié)點 Li 3個方向的節(jié)點力;wLl、uLl、vLl為節(jié)點 Ll 3 個方向的位移;wLl*、uLl*、vLl*為節(jié)點Ll*3個方向的位移;ΔA為單元裂紋面的面積;節(jié)點Ll和Ll*初始坐標相同。

需要指出的是，以上結果都是基于裂尖的局部坐標系所得到的結果，當裂尖的局部坐標系和全局坐標系不一致時，要將所有應力和位移轉化到局部坐標系下。

圖4 三維裂紋示意圖Fig.4 Diagram for G calculation by 3D-VCCT

在線彈性情況下，G和SIF有如下關系:

式中 KI、KII、KIII分別為裂尖處的 I型、II型、III型應力強度因子;對于平面應力狀態(tài)E＇=E;對于平面應變狀態(tài)E＇=E/(1－ν2);E為材料的彈性模量;μ為材料的剪切模量。

ABAQUS已將VCCT功能集成到有限元程序中，本文直接利用ABAQUS的3D-VCCT計算功能，不需用戶子程序開發(fā)就可以得到滿意的計算結果。

2.2 計算模型

建立簡化后的金屬膜片的三維有限元模型，為簡化計算、有利于網(wǎng)格的劃分，簡化金屬膜片只有一條預制缺陷。合理簡化邊界條件，模擬膜片在發(fā)動機中的實際連接形式。整個膜片全部采用結構化網(wǎng)格，單元類型為C3D8R，單元總數(shù)為3.5萬。具體計算模型如圖5所示，膜片的材料選用LY12超強鋁合金。

將整個加載過程定義為靜態(tài)分析步:在隔膜片II脈沖一側施加壓力載荷，模擬II脈沖燃燒室初始工作壓強，分析步為靜態(tài)分析步，載荷的大小隨時間線性增加，計算預制缺陷在2.0 MPa下的應力強度因子。

建立支撐件及膜片的三維有限元模型及相應當量實心圓板三維有限元模型，在膜片一側分別施加壓力載荷，研究支撐件及實心圓板中心受拉一側應力變化情況。三維有限元模型如圖6所示，邊界條件為簡支，單元類型為C3D8R減縮積分單元，單元總數(shù)為20萬。

圖5 帶1條預制缺陷的圓平板模型Fig.5 Circular plate model with one prefab gap

圖6 有限元計算模型Fig.6 Finite element model

3 公式解與數(shù)值計算結果對比

3.1 膜片的應力強度因子對比

根據(jù)式(7)計算預制缺陷上的應力強度因子分布，并與有限元計算結果對比如圖7所示。

圖7 公式與數(shù)值計算結果對比Fig.7 Formula and numerical calculation results

對比圖7可知，2種方法計算的應力強度因子變化趨勢一致，在圓板中心位置達到最大，沿徑向方向逐漸減小，與圖3定性分析結果一致。同時，二者之間存在一定差異，在圓心位置公式計算結果為21.4 MPa·，數(shù)值計算結果為20.1 MPa·，誤差為 6.5%。在預制缺陷末端公式計算結果為3.4 MPa·，數(shù)值計算結果為2.5 MPa·，誤差達到34.2%。

分析認為，造成二者之間主要差別的原因是:

(1)公式的假設引入誤差。公式推導過程中作了大量假設，如假設所取截面為二維平面應變結構，必然會引入一定的誤差。

(2)數(shù)值計算結果對模型敏感。利用有限元數(shù)值計算裂紋應力強度因子時，計算結果對網(wǎng)格尺寸及裂紋前緣選取較為敏感，不同的網(wǎng)格劃分形式，及不同的裂紋前緣選取都會導致數(shù)值計算結果的差異(本文計算時選取裂紋前緣與裂尖一致)。

由于本文所設計的膜片結構破壞過程都是從圓心位置最先達到材料的斷裂韌性，從中心開始破壞，并擴展致整個預制缺陷。如果材料的斷裂韌性為20 MPa·，則根據(jù)公式及數(shù)值計算結果，在2 MPa內壓作用下，膜片已經(jīng)從中心位置破壞。

3.2 支撐件的應力大小對比

計算所得到各點的最大主應力，其中實心圓板的應力按式(8)計算，支撐件的應力按式(9)計算。由于本文主要關心圓板中心的應力狀態(tài)，在計算S0/S比值時，按圖8所示的2個面積進行計算，S0/S≈2，將公式計算結果與有限元計算結果對比如圖9所示。如圖9所示，公式計算結果與有限元計算結果保持一致，變化趨勢相同，但略高于有限元計算結果，說明利用該方法進行強度校核時，所得結果偏于“保守”，有利于提高隔艙整體的承載能力。

圖8 關心區(qū)域對應的面積Fig.8 The areas of interest region

圖9 不同壓力下計算結果對比Fig.9 Calculation results under different pressures

4 驗證試驗

4.1 膜片的打開壓強驗證試驗

考慮到某雙脈沖發(fā)動機的直徑，同時保證II脈沖藥柱能穩(wěn)定點燃(一般要求初始壓強控制在1.5～2.5 MPa)，本文設計II脈沖初始工作壓強為2.0 MPa，根據(jù)前文設計公式，設計膜片具體尺寸如下:制缺陷深度a=1 mm，缺陷V型槽開口角度α=π/2，膜片半徑R=142 mm，膜片厚度h=3 mm。

本文進行了5次膜片的破壞打開單項試驗［7］，膜片的打開壓強分別為:2.30、1.90、2.33、2.00、1.95 MPa。膜片的平均打開壓強為 2.1 MPa，比公式計算結果2.0 MPa高約5%。可見公式計算結果可較好預測膜片的真實動態(tài)打開壓強。

4.2 支撐件強度校核試驗

考慮到某雙脈沖發(fā)動機I脈沖的實際工作壓強為16 MPa，本文隔艙的承壓極限要求大于22 MPa，利用式(9)對前期的多種結構形式及材料的支撐件進行強度校核，最終確定支撐件的結構形式如圖6(a)所示，材料為LC9超強鋁合金，隔艙整體厚度為20 mm，根據(jù)式(9)隔艙的極限承載壓強為24 MPa。

本文對該隔艙進行了熱流承壓試驗，試驗裝置為模擬發(fā)動機，由于試驗成本較高，本文只進行了3發(fā)承壓試驗，發(fā)動機工作壓強分別為 21.5、22.3、20.5 MPa。試驗完成后，拆分模擬發(fā)動機，觀察支撐件一側，發(fā)現(xiàn)支撐件結構完整。說明本文的強度校核方法可以較好的預示隔艙的承載極限，為隔艙的設計提供一定依據(jù)。

5 結論

(1)根據(jù)本文推導的公式計算預制缺陷上的應力強度因子分布，與有限元計算結果一致性較好，在圓心位置公式計算誤差僅為6.5%。根據(jù)公式預估結果，膜片的破壞壓強為2.0 MPa，金屬膜片5次試驗的平均打開壓強為2.10 MPa，公式誤差僅為5%，進一步驗證了公式的準確性。

(2)根據(jù)支撐件強度校核公式得到的關心點的應力值略高于有限元計算值，說明校核公式更偏于保守。通過熱流承壓試驗，驗證了隔艙的承載性能，進一步說明本文所推導公式可以用來校核隔艙支撐件的強度，得到隔艙極限承載壓強。

(3)本文利用圓板大撓度理論及斷裂力學理論，得到了金屬膜片的設計公式，公式計算所得到的打開壓強可以用來預估膜片的真實打開壓強。利用多孔圓板理論，得到支撐件強度校核方法，可以用來對支撐件結構形式及材料進行前期優(yōu)選。本文得到的公式可以為雙脈沖發(fā)動機中隔艙設計提供參考依據(jù)。

［1］Wang Chun-guang，Liu Hong-chao，Yang De-min.Study of development status and future development proposals for pulse soild rocket motor［J］.Applied Mechanics and Materials，2012(110-116):2354-2358.

［2］王春光，田維平，任全彬.脈沖發(fā)動機中隔層工作過程的數(shù)值分析及試驗［J］.推進技術，2012，33(5):790-794.

［3］孫娜，婁永春，孫長宏.某雙脈沖發(fā)動機燃燒室兩相流場數(shù)值分析［J］.固體火箭技術，2012，35(3):335-338.

［4］劉亞冰，王長輝，劉宇.雙脈沖發(fā)動機燃燒室局部燒蝕特性分析［J］.固體火箭技術，2010，33(4):453-456.

［5］曹熙煒，劉宇，任軍學.內孔隔板脈沖固體火箭發(fā)動機流場分析［J］.航空動力學報，2011，26(2):448-452.

［6］朱衛(wèi)兵，張永飛，陳宏.雙脈沖發(fā)動機內流場研究［J］.彈箭與制導學報，2012，32(1):114-118.

［7］王春光，任全彬，田維平，等.脈沖發(fā)動機中金屬膜片式隔艙動態(tài)破壞過程研究［J］.固體火箭技術，2013，36(1):573-577.

［8］劉鴻文.板殼理論［M］.浙江大學出版社，1987.

［9］謝桂蘭，張平，龔曙光.基于均勻化理論的管板有效彈性常數(shù)的研究［J］.應用力學學報，2004，21(3).

［10］解德，錢勤，李長安.斷裂力學中的數(shù)值計算方法及工程應用［M］.北京:科學出版社，2009.