唐國(guó)金，申志彬，田四朋，楊東

(1.國(guó)防科學(xué)技術(shù)大學(xué) 航天與材料工程學(xué)院，湖南 長(zhǎng)沙410073;2.成都飛機(jī)設(shè)計(jì)研究所，四川 成都610041;3.中國(guó)航天科技集團(tuán)公司 四院43 所，陜西 西安710025)

0 引言

固體火箭發(fā)動(dòng)機(jī)(SRM)藥柱貯存壽命預(yù)估具有重要意義。常規(guī)的藥柱貯存壽命預(yù)估一般采用“確定性模型”，不考慮隨機(jī)參數(shù)的影響，取其均值進(jìn)行分析。這降低了預(yù)估的準(zhǔn)確度和可信度，是預(yù)估壽命偏離真實(shí)壽命的重要原因之一。

20世紀(jì)60年代初，文獻(xiàn)［1］采用數(shù)理統(tǒng)計(jì)方法和應(yīng)力－強(qiáng)度干涉理論處理隨機(jī)因素，初步形成了概率壽命預(yù)估的思想。文獻(xiàn)［2－3］基于高溫加速老化試驗(yàn)，分析了固體推進(jìn)劑的可靠貯存壽命。上述研究?jī)H以固體推進(jìn)劑力學(xué)性能參數(shù)作為發(fā)動(dòng)機(jī)貯存壽命的標(biāo)志參數(shù)，脫離了藥柱真實(shí)幾何形狀和載荷歷程，不可避免存在較大的誤差。文獻(xiàn)［4－5］將發(fā)動(dòng)機(jī)藥柱簡(jiǎn)化為一個(gè)彈性約束的無(wú)限長(zhǎng)中空?qǐng)A筒，利用經(jīng)驗(yàn)公式和Monte Carlo 方法進(jìn)行隨機(jī)結(jié)構(gòu)分析，得到了藥柱結(jié)構(gòu)可靠度隨貯存時(shí)間的變化。文獻(xiàn)［6］則采用粘彈性Monte Carlo 隨機(jī)有限元法對(duì)隨機(jī)溫度載荷下導(dǎo)彈結(jié)構(gòu)進(jìn)行了概率貯存壽命預(yù)估。但由于計(jì)算手段和分析方法的限制，這些研究?jī)H能處理一些簡(jiǎn)化的藥柱模型，難以應(yīng)用于實(shí)際工程問(wèn)題。預(yù)估具體發(fā)動(dòng)機(jī)藥柱的概率貯存壽命，需借助三維粘彈性隨機(jī)有限元法，目前尚無(wú)采用該方法開(kāi)展SRM 概率貯存壽命預(yù)估的研究報(bào)道。

本文基于三維粘彈性隨機(jī)有限元法和推進(jìn)劑高溫加速老化試驗(yàn)，對(duì)某SRM 藥柱的概率貯存壽命進(jìn)行預(yù)估。

1 概率壽命預(yù)估模型

高溫加速老化試驗(yàn)數(shù)據(jù)符合統(tǒng)計(jì)原則，假設(shè)總體服從正態(tài)分布，各觀察值之間相互獨(dú)立［7］。加速老化過(guò)程中，固體推進(jìn)劑某些性能參數(shù)的變化與時(shí)間、溫度存在一定的相關(guān)關(guān)系，通過(guò)對(duì)老化數(shù)據(jù)進(jìn)行擬合，可找到一個(gè)與其相適應(yīng)的老化模型方程。航天行業(yè)標(biāo)準(zhǔn)QJ 2328A—2005［7］給出了用于數(shù)據(jù)處理的3 個(gè)老化模型:

式中:I 為固體推進(jìn)劑某一時(shí)刻的性能;I0為性能初始值(常數(shù));K 為與溫度有關(guān)的性能變化速度常數(shù);t 為老化時(shí)間，單位:周。

在老化過(guò)程中，假定固體推進(jìn)劑性能變化服從Arrhenius 化學(xué)動(dòng)力學(xué)方程:

式中:E 為活化能;H 為頻率因子;R 為氣體常數(shù)。

考慮預(yù)測(cè)過(guò)程中存在的隨機(jī)誤差，需求K 的置信界限。先求lnK 的標(biāo)準(zhǔn)差［7］

于是，(4)式的置信界限為

式中tα為置信概率為1－α 時(shí)t 分布值。

從而可得貯存溫度為T(mén)s時(shí)，K 的上限值

根據(jù)文獻(xiàn)［7］，可取各老化試驗(yàn)溫度下初始性能的平均值為(1)～(3)式的性能初始值，即

固體發(fā)動(dòng)機(jī)結(jié)構(gòu)可靠度

式中:Zp(t)為某時(shí)刻藥柱的功能函數(shù);R(t)和S(t)分別為該時(shí)刻與強(qiáng)度和結(jié)構(gòu)響應(yīng)有關(guān)的隨機(jī)參數(shù)。假定R(t)和S(t)均為簡(jiǎn)單隨機(jī)變量，服從正態(tài)分布且相互獨(dú)立，則Zp亦服從正態(tài)分布，其均值和方差分別為

式中:μR(t)、μS(t)和σ2R(t)、σ2S(t)分別為R(t)、S(t)的均值和方差。對(duì)(9)式進(jìn)行標(biāo)準(zhǔn)化變換，則［8］

式中Φ(·)為標(biāo)準(zhǔn)正態(tài)分布的累積分布函數(shù)。

失效概率

結(jié)構(gòu)可靠度分析中，通常令

β(t)即可靠指標(biāo)。不難看出，在R(t)和S(t)均為正態(tài)分布的條件下，β(t)唯一確定了結(jié)構(gòu)可靠度。

本文藥柱的破壞準(zhǔn)則采用與八面體剪應(yīng)變等效Von Mises 應(yīng)變準(zhǔn)則［9］，因此，R(t)和S(t)分別表示某貯存時(shí)間推進(jìn)劑的最大伸長(zhǎng)率和藥柱的最大Von Mises 應(yīng)變。按照(11)～(13)式給出的可靠度分析模型，可得藥柱結(jié)構(gòu)可靠度隨貯存時(shí)間的變化規(guī)律，進(jìn)一步可求取滿足一定可靠度要求的概率貯存壽命。

以可靠指標(biāo)β 為例，給出基于可靠指標(biāo)的藥柱概率壽命預(yù)估模型。對(duì)各個(gè)貯存期的藥柱進(jìn)行隨機(jī)結(jié)構(gòu)完整性分析得到藥柱最大Von Mises 應(yīng)變均值μS和標(biāo)準(zhǔn)差σS，結(jié)合推進(jìn)劑最大伸長(zhǎng)率各貯存時(shí)間的均值μR和標(biāo)準(zhǔn)差σR.利用(13)式可得藥柱各貯存期的可靠指標(biāo)β(t)，對(duì)其進(jìn)行參數(shù)擬合得

β(t)與推進(jìn)劑材料參數(shù)M，藥柱載荷參數(shù)L，幾何參數(shù)G，貯存溫度T 以及貯存時(shí)間t 有關(guān)，其中M、L、G和T 為參數(shù)，t 為變量。f(M，L，T，t)可以為多項(xiàng)式模型、指數(shù)模型或?qū)?shù)模型。給定滿足要求的可靠指標(biāo)β0，利用(14)式便可得到一定可靠指標(biāo)下藥柱的貯存壽命t.(14)式即為基于可靠指標(biāo)的藥柱概率壽命預(yù)估模型。同理，利用(11)～(12)式可得基于可靠度和失效概率的藥柱概率壽命預(yù)估模型。

2 概率壽命預(yù)估流程

采用基于高溫加速老化試驗(yàn)和結(jié)構(gòu)完整性分析相結(jié)合的藥柱壽命預(yù)估方法［9］，利用三維粘彈性隨機(jī)有限元方法，對(duì)SRM 藥柱進(jìn)行概率壽命預(yù)估，流程如圖1所示。本文采用三維粘彈性響應(yīng)面隨機(jī)有限元法對(duì)藥柱進(jìn)行隨機(jī)結(jié)構(gòu)分析［10］，基于Von Mises 應(yīng)變準(zhǔn)則計(jì)算藥柱的結(jié)構(gòu)可靠度。

圖1 藥柱概率貯存壽命預(yù)估流程Fig.1 Flow chart of grain probabilistic storage life prediction

3 算例

利用提出的概率壽命預(yù)估模型對(duì)某固體發(fā)動(dòng)機(jī)藥柱進(jìn)行概率壽命預(yù)估。

3.1 試驗(yàn)方法

依照文獻(xiàn)［7］對(duì)某HTPB 固體推進(jìn)劑進(jìn)行了高溫加速老化試驗(yàn)，選取了40 ℃、50 ℃、60 ℃、70 ℃四個(gè)老化試驗(yàn)溫度，每個(gè)老化溫度下取樣8 次。按照GB 4889—85 的要求，為了獲取固體推進(jìn)劑性能參數(shù)的統(tǒng)計(jì)特征，每次取樣數(shù)目為25 個(gè)，并根據(jù)GJB 770B—2005 規(guī)定的方法對(duì)試件進(jìn)行單向拉伸試驗(yàn)，異常數(shù)據(jù)按GB 4883—85 中規(guī)定的方法處理。由試驗(yàn)數(shù)據(jù)知，老化導(dǎo)致該推進(jìn)劑最大伸長(zhǎng)率降低，模量升高，且老化速度隨老化溫度升高加快。

3.2 老化性能統(tǒng)計(jì)特征隨時(shí)間變化規(guī)律

按照文獻(xiàn)［7］給出的原則，為了選擇合適的老化數(shù)學(xué)模型，如圖2～3 所示，給出了各老化溫度下固體推進(jìn)劑最大伸長(zhǎng)率εmax及初始模量E0均值的散點(diǎn)圖，并給出了一次多項(xiàng)式的最小二乘擬合曲線。

圖2 最大伸長(zhǎng)率均值與老化時(shí)間散點(diǎn)圖Fig.2 Scatter diagram for specific elongation and aging time

圖3 初始模量均值與老化時(shí)間散點(diǎn)圖Fig.3 Scatter diagram for initial modulus and aging time

由散點(diǎn)圖可看出，在各老化試驗(yàn)溫度下，固體推進(jìn)劑最大伸長(zhǎng)率和初始模量的均值與老化時(shí)間具有明顯的線性關(guān)系。為了進(jìn)一步對(duì)線性相關(guān)性進(jìn)行檢驗(yàn)，如表1所示，列出了固體推進(jìn)劑最大伸長(zhǎng)率和初始模量均值隨老化時(shí)間變化的回歸方程，同時(shí)給出了線性相關(guān)系數(shù)和相應(yīng)的置信概率。計(jì)算表明，在置信概率大于99%的條件下，該型號(hào)固體推進(jìn)劑老化性能參數(shù)與老化時(shí)間的線性關(guān)系成立，因此可采用(1)式的老化模型方程求解各老化試驗(yàn)溫度下的性能變化速度常數(shù)K，然后求解Arrhenius 方程中的其他系數(shù)，進(jìn)而外推出常溫貯存條件下固體推進(jìn)劑的性能變化速度常數(shù)。

表1 固體推進(jìn)劑老化性能參數(shù)隨時(shí)間變化的回歸方程Tab.1 Regression equation for propellant aging performance parameter variation with time

基于表1中給出的不同溫度下性能變化速度常數(shù)，采用文獻(xiàn)［7］給出的步驟對(duì)Arrhenius 方程(4)式中的參數(shù)進(jìn)行計(jì)算及相關(guān)性檢驗(yàn)，如表2所示。取置信概率為99%，有tα=4.540 7.從而在貯存溫度Ts=25 ℃的條件下，可得固體推進(jìn)劑最大伸長(zhǎng)率εmax和初始模量E0的變化速度常數(shù)Ks的上限值分別為8.538 9 ×10－5和0.002 273.如表3所示，給出了固體推進(jìn)劑性能變化的回歸方程。

表2 老化溫度與性能變化速度常數(shù)相關(guān)的回歸方程Tab.2 Regression equation for performance velocity constant variation with aging temperature

這里從固體推進(jìn)劑高溫加速老化性能測(cè)試結(jié)果的均值出發(fā)給出了常溫貯存條件下固體推進(jìn)劑性能變化回歸方程的計(jì)算過(guò)程，未考慮老化性能參數(shù)的隨機(jī)性。實(shí)際上，由于不同老化溫度和老化時(shí)間的性能參數(shù)具有隨機(jī)性，在常溫貯存條件下，固體推進(jìn)劑性能參數(shù)必然隨貯存時(shí)間隨機(jī)變化，本文采用Monte Carlo 方法解決此問(wèn)題，主要步驟為:1)對(duì)不同老化溫度和老化試驗(yàn)時(shí)間的性能測(cè)試數(shù)據(jù)進(jìn)行n次隨機(jī)抽樣;2)重復(fù)進(jìn)行以上計(jì)算，可得到n 條常溫貯存下固體推進(jìn)劑性能變化的回歸方程(樣本函數(shù));3)統(tǒng)計(jì)分析各回歸方程在不同貯存時(shí)間的函數(shù)值，即可得到固體推進(jìn)劑性能參數(shù)統(tǒng)計(jì)特征隨貯存時(shí)間的變化規(guī)律。

表3 常溫貯存條件下固體推進(jìn)劑性能變化的回歸方程Tab.3 Regression equation for propellant performance parameter variation under normal temperature

取隨機(jī)抽樣次數(shù)n=10 000，經(jīng)計(jì)算可得常溫貯存條件下固體推進(jìn)劑最大伸長(zhǎng)率和初始模量標(biāo)準(zhǔn)差隨貯存時(shí)間的變化曲線，如圖4～5 所示。

圖4 最大伸長(zhǎng)率標(biāo)準(zhǔn)差隨時(shí)間變化曲線Fig.4 Curve of specific elongation standard deviation variation with time

圖5 初始模量標(biāo)準(zhǔn)差隨時(shí)間變化曲線Fig.5 Curve of initial modulus standard deviation variation with time

計(jì)算表明，老化試驗(yàn)數(shù)據(jù)的隨機(jī)性對(duì)推進(jìn)劑性能變化回歸方程的初值I0影響不大，但對(duì)性能變化速度常數(shù)K 影響非常顯著，從而導(dǎo)致推進(jìn)劑最大伸長(zhǎng)率和初始模量標(biāo)準(zhǔn)差隨貯存時(shí)間的增加而增大。

3.3 藥柱概率壽命預(yù)估

以文獻(xiàn)［10］中給出的SRM 為例，取貯存溫度Ts=25 ℃，載荷工況為常溫點(diǎn)火發(fā)射，如表4所示，給出了藥柱未老化時(shí)的部分結(jié)構(gòu)參數(shù)，其中，推進(jìn)劑松弛模量Prony 級(jí)數(shù)的松弛時(shí)間τj= 4.303 5 ×10j－1，j=1，2，3.

為簡(jiǎn)化計(jì)算，本文僅考慮推進(jìn)劑松弛模量和最大伸長(zhǎng)率的老化特性，不考慮其他隨機(jī)參數(shù)隨貯存時(shí)間的變化。由于推進(jìn)劑松弛模量Prony 級(jí)數(shù)中的各參數(shù)是由一組(多個(gè)試樣)應(yīng)力松弛試驗(yàn)確定的，不能通過(guò)簡(jiǎn)單的統(tǒng)計(jì)分析得到。本文采用張海聯(lián)［11］給出的近似解法得到。固體推進(jìn)劑初始模量的老化統(tǒng)計(jì)特征已在3.2 節(jié)給出，其他松弛模量老化統(tǒng)計(jì)特征不一一列出，僅給出未老化時(shí)推進(jìn)劑藥柱的松弛模量統(tǒng)計(jì)特征。

表4 藥柱結(jié)構(gòu)分析參數(shù)Tab.4 Structure analysis parameters of grain

采用響應(yīng)面隨機(jī)有限元方法對(duì)藥柱進(jìn)行隨機(jī)結(jié)構(gòu)完整性分析，假定貯存時(shí)間為20年，藥柱危險(xiǎn)點(diǎn)Von Mises 應(yīng)變均值和標(biāo)準(zhǔn)差隨貯存時(shí)間的變化曲線如圖6～7 所示?？煽闯觯S著貯存時(shí)間的增加，由于固體推進(jìn)劑模量上升，SRM 藥柱危險(xiǎn)點(diǎn)Von Mises 應(yīng)變的均值和標(biāo)準(zhǔn)差均有所下降。

利用(13)式計(jì)算各貯存時(shí)間藥柱的可靠指標(biāo)β，可得到藥柱結(jié)構(gòu)可靠指標(biāo)隨貯存時(shí)間的變化曲線，如圖8所示，由圖可知，在貯存期間SRM 藥柱的結(jié)構(gòu)可靠度不斷降低，失效概率增加。如圖9所示，進(jìn)一步給出了藥柱結(jié)構(gòu)失效概率隨貯存時(shí)間的變化曲線。

可見(jiàn)，根據(jù)藥柱結(jié)構(gòu)可靠度(或失效概率)隨貯存時(shí)間的變化趨勢(shì)，在不同結(jié)構(gòu)可靠度水平下，可得到藥柱不同的壽命預(yù)估結(jié)果。

圖6 應(yīng)變均值隨貯存時(shí)間的變化曲線Fig.6 Curve of strain equally variation with storage time

圖7 應(yīng)變標(biāo)準(zhǔn)差隨貯存時(shí)間的變化曲線Fig.7 Curve of strain standard deviation variation with storage time

圖8 藥柱結(jié)構(gòu)可靠指標(biāo)隨貯存時(shí)間的變化曲線Fig.8 Curve of grain structural reliable index variation with storage time

利用三次多項(xiàng)式模型對(duì)圖8中數(shù)據(jù)點(diǎn)進(jìn)行參數(shù)擬合，(14)式可表示為

式中:A0、A1、A2、A3分 別 為4.417、－ 7.032 × 10－4、－4.617×10－6、2.639 ×10－9，相關(guān)系數(shù)為0.999 9.給定β 或可靠度，就可從(15)式中求出藥柱的概率壽命。如表5所示，給出了幾種結(jié)構(gòu)可靠度水平下該發(fā)動(dòng)機(jī)藥柱的貯存壽命。

圖9 藥柱失效概率隨貯存時(shí)間的變化曲線Fig.9 Curve of grain failure probability variation with storage time

顯然，在較高的結(jié)構(gòu)可靠度水平條件下預(yù)估的貯存壽命較短，反之亦然。因此，確定合理的結(jié)構(gòu)可靠度水平非常重要。

表5 SRM 藥柱的概率貯存壽命Tab.5 Probabilistic storage life of solid rocket motor grain

4 結(jié)論

以三維粘彈性隨機(jī)有限元法和固體推進(jìn)劑高溫加速老化試驗(yàn)為基礎(chǔ)，提出了SRM 藥柱概率貯存壽命預(yù)估模型，并實(shí)現(xiàn)了具體SRM 藥柱的概率貯存壽命預(yù)估。主要研究結(jié)論如下:

1)基于粘彈性隨機(jī)有限元法的SRM 藥柱概率貯存壽命預(yù)估模型，能夠綜合考慮環(huán)境/載荷譜、幾何構(gòu)型及參數(shù)不確定性等因素對(duì)壽命預(yù)估的影響，結(jié)果更為合理、可信。

2)SRM 藥柱概率貯存壽命預(yù)估能夠在給定結(jié)構(gòu)可靠度要求的條件下，給出藥柱的貯存壽命。

3)對(duì)具體SRM 藥柱進(jìn)行概率貯存壽命預(yù)估，需借助三維粘彈性隨機(jī)有限元法和結(jié)構(gòu)可靠度分析理論。

本文建立和實(shí)現(xiàn)了適于工程應(yīng)用的SRM 藥柱概率貯存壽命預(yù)估模型，具有重要的理論意義和應(yīng)用價(jià)值。

References)

［1］ Chappell R N，Jensen F R，Burton R W.Statistical service life prediction:minuteman third－stage propellant grain［J］.Journal of Spacecraft and Rockets，1968，5(1):42－46.

［2］ 劉德輝，賀南昌.復(fù)合推進(jìn)劑貯存壽命及其可靠性研究［J］.推進(jìn)技術(shù)，1993，14(6):63－67.LIU De－h(huán)ui，HE Nan－chang.An investigation on storage life of composite propellant and its reliability［J］.Journal of Propulsion Technology，1993，14(6):63－67.(in Chinese)

［3］ 劉子如，郝仲璋，謝俊杰，等.FH－94 固體推進(jìn)劑的可靠度評(píng)估和可靠貯存壽命預(yù)測(cè)研究［J］.兵工學(xué)報(bào)，1994，15(2):24－27.LIU Zi－ru，HAO Zhong－zhang，XIE Jun－jie，et al.Research on reliability evaluation and reliable life prediction of FH－94 solid propellant［J］.Acta Armamentarii，1994，15(2):24－27.(in Chinese)

［4］ Heller R A，Kamat M P.Service life prediction for solid propellant motors，AIAA－1978－486 ［R］.New York:American Institute of Aeronautics and Astronautics，1978.

［5］ Heller R A，Thangjitham S，Janajreh I M.Probabilistic service life prediction for solid propellant motors subjected to environmental thermal loads，AD－330303［R］.Virginia:Defence Technical Information Center，1997.

［6］ Cost T L.Probabilistic service life prediction of missile structures subjected to random thermal loads，AIAA－1980—406［R］.New York:American Institute of Aeronautics and Astronautics，1980.

［7］ 李彥麗，趙海泉，杜秋麗.QJ 2328A—2005 復(fù)合固體推進(jìn)劑高溫加速老化試驗(yàn)方法［S］.北京:中國(guó)航天標(biāo)準(zhǔn)化研究所，2005.LI Yan－li，ZHAO Hai－quan，DU Qiu－li.QJ 2328—2005 test method for temperature accelerated aging of composite solid propellant［S］.Beijing:China Astronautics Standards Institute，2005.(in Chinese)

［8］ 何水清，王善.結(jié)構(gòu)可靠性分析與設(shè)計(jì)［M］.北京:國(guó)防工業(yè)出版社，1993.HE Shui－qing，WANG Shan.Structural reliability analysis and design［M］.Beijing:National Defense Industry Press，1993.(in Chinese)

［9］ 袁端才，唐國(guó)金，李海陽(yáng)，等.基于加速老化與三維粘彈性有限元分析的固體導(dǎo)彈發(fā)動(dòng)機(jī)壽命預(yù)估［J］.兵工學(xué)報(bào)，2006，27(4):685－689.YUAN Duan－cai，TANG Guo－jin，LI Hai－yang，et al.Three－dimension viscoelastic analysis and accelerated aging approach to predict the service life of solid missile motor［J］.Acta Armamentarii，2006，27(4):685－689.(in Chinese)

［10］ 田四朋.固體火箭發(fā)動(dòng)機(jī)藥柱三維粘彈性響應(yīng)面隨機(jī)有限元分析［J］.固體火箭技術(shù)，2010，33(1):17－20.TIAN Si－peng.Three－dimensional viscoelastic response surface stochastic finite element analysis on solid rocket motor grain［J］.Journal of Solid Rocket Technology，2010，33(1):17－20.(in Chinese)

［11］ 張海聯(lián).固體火箭發(fā)動(dòng)機(jī)藥柱的粘彈性不確定結(jié)構(gòu)分析［D］.長(zhǎng)沙:國(guó)防科學(xué)技術(shù)大學(xué)，2002.ZHANG Hai－lian.Viscoelastic uncertain structural analysis for solid rocket motor grain［D］.Changsha:National University of Defense Technology，2002.(in Chinese)