非均布瞬態(tài)內(nèi)壓作用下固體推進(jìn)劑藥柱泊松比隨機(jī)結(jié)構(gòu)分析①

王君祺，職世君，張建偉

(1.中國(guó)空空導(dǎo)彈研究院四所，洛陽(yáng) 471009;2.航空制導(dǎo)武器航空科技重點(diǎn)實(shí)驗(yàn)室，洛陽(yáng) 471009;3.北京航空航天大學(xué)，北京 100191)

0 引言

點(diǎn)火發(fā)射是固體火箭發(fā)動(dòng)機(jī)所經(jīng)歷的最惡劣的載荷工況。因此，在發(fā)動(dòng)機(jī)的設(shè)計(jì)階段及后期發(fā)動(dòng)機(jī)貯存后的再評(píng)估階段，對(duì)固體火箭發(fā)動(dòng)機(jī)點(diǎn)火工況下的結(jié)構(gòu)完整性分析就顯得格外重要。國(guó)內(nèi)外對(duì)固體火箭發(fā)動(dòng)機(jī)內(nèi)壓載荷下結(jié)構(gòu)完整性分析的研究較多，但大多采用均勻分布載荷［1－3］，而發(fā)動(dòng)機(jī)在實(shí)際點(diǎn)火過程中壓力波的傳遞，使得發(fā)動(dòng)機(jī)內(nèi)壓力載荷存在瞬態(tài)分布特性［4］。因此，若采用瞬態(tài)分布?jí)毫M(jìn)行分析，則更接近發(fā)動(dòng)機(jī)的實(shí)際工況。

由于固體發(fā)動(dòng)機(jī)裝藥受制造、工藝、測(cè)量等因素的影響，使得藥柱力學(xué)性能參數(shù)不可避免具有一定的分散性。如泊松比接近于0.5，通過實(shí)驗(yàn)很難測(cè)定其準(zhǔn)確數(shù)值，多次的測(cè)量結(jié)果呈概率分布，具有隨機(jī)性［5］。若采用確定性的分析方法，得出的結(jié)果可信性較差。近年來，隨著計(jì)算機(jī)計(jì)算能力的不斷提高，為了能描述工程結(jié)構(gòu)的實(shí)際情況，考慮不確定性因素已成為一種趨勢(shì)［6］。因此，有必要對(duì)推進(jìn)劑藥柱結(jié)構(gòu)進(jìn)行不確定性分析［7］。

本文首先采用Abaqus商用有限元軟件計(jì)算了固體火箭發(fā)動(dòng)機(jī)點(diǎn)火過程壓力瞬態(tài)分布時(shí)裝藥的結(jié)構(gòu)完整性，然后結(jié)合粘彈性有限元法和Monte-Carlo方法，即構(gòu)成所謂的粘彈性Monte-Carlo方法，考慮了泊松比隨機(jī)分布對(duì)固體火箭發(fā)動(dòng)機(jī)點(diǎn)火過程裝藥結(jié)構(gòu)完整性的影響。

1 物理模型和計(jì)算方法

1.1 發(fā)動(dòng)機(jī)模型及網(wǎng)格劃分

發(fā)動(dòng)機(jī)的構(gòu)成主要有殼體、絕熱層和藥柱，假設(shè)藥柱和絕熱層為線粘彈性材料，殼體為彈性材料。線粘彈性本構(gòu)關(guān)系為

式中 Yijkl(t)為四階張量，具體材料參數(shù)參見文獻(xiàn)［8］。

發(fā)動(dòng)藥柱為管狀結(jié)構(gòu)，考慮發(fā)動(dòng)機(jī)幾何構(gòu)型的對(duì)稱性，取其1/4建模，如圖1(a)所示。其中，面1為點(diǎn)火藥包瞬間壓力填充區(qū)域，面2的壓力符合式(3)的分布。采用六面體單元對(duì)發(fā)動(dòng)機(jī)進(jìn)行網(wǎng)格劃分，如圖1(b)所示。

圖1 固體火箭發(fā)動(dòng)機(jī)1/4模型及網(wǎng)格劃分Fig.1 SRM 1/4 model and mesh

1.2 發(fā)動(dòng)機(jī)瞬態(tài)載荷及藥柱泊松比分布

1.2.1 發(fā)動(dòng)機(jī)瞬態(tài)載荷分布

發(fā)動(dòng)機(jī)點(diǎn)火藥位于發(fā)動(dòng)機(jī)頭部，點(diǎn)火后，發(fā)動(dòng)機(jī)燃燒室藥柱頭部壓力急劇增加，并以很快的速度向燃燒室尾部傳播。在此期間，發(fā)動(dòng)機(jī)內(nèi)壓力載荷呈典型的非均布和非穩(wěn)態(tài)特征。文獻(xiàn)［4］中提出的固體火箭發(fā)動(dòng)機(jī)點(diǎn)火過程內(nèi)壓瞬態(tài)非均布函數(shù)，即

式中 p0為平衡壓強(qiáng);v為壓力軸向傳播速度;z為軸向坐標(biāo);t為發(fā)動(dòng)機(jī)點(diǎn)火時(shí)間;a為壓力增長(zhǎng)速因子。

當(dāng)點(diǎn)火藥被引燃后，壓力會(huì)瞬間充滿發(fā)動(dòng)機(jī)頭部區(qū)域，若直接采用上述公式，發(fā)動(dòng)機(jī)燃燒室壓力分布為從z=0位置開始遞減，與實(shí)際情況有所出入。因此，假定發(fā)動(dòng)機(jī)頭部點(diǎn)火藥部位在增壓過程中壓力保持一致，即在式(2)中引入發(fā)動(dòng)機(jī)頭部點(diǎn)火距離的參數(shù)b，則式(2)可改寫為

文中分別取 a=500，b=45 mm，v=50 m/s，發(fā)動(dòng)總長(zhǎng)為500 mm。

1.2.2 藥柱泊松比隨機(jī)分布

采用Monte-Carlo法進(jìn)行結(jié)構(gòu)隨機(jī)有限元分析時(shí)，首先要根據(jù)隨機(jī)變量的已知概率分布進(jìn)行抽樣，然后再對(duì)每個(gè)抽取的樣本，利用有限元法求取發(fā)動(dòng)機(jī)藥柱的結(jié)構(gòu)響應(yīng)量，最后根據(jù)所有樣本的計(jì)算結(jié)果進(jìn)行統(tǒng)計(jì)分析。對(duì)隨機(jī)變量進(jìn)行抽樣時(shí)，若采用直接Monte-Carlo方法，則在不改變抽樣的均值和方差的前提下，樣本數(shù)目巨大，由其對(duì)于三維粘彈性藥柱的有限計(jì)算，其計(jì)算量太大。因此，為提高粘彈性Monte-Carlo隨機(jī)有限元的計(jì)算效率，采用Latin超立方抽樣技術(shù)，該方法在不改變抽樣的均值和方差的前提下，可有效減少樣本數(shù)目。

在對(duì)固體火箭發(fā)動(dòng)機(jī)裝藥點(diǎn)火工況下的結(jié)構(gòu)完整性進(jìn)行分析時(shí)，泊松比的影響往往較大［9］。考慮藥柱泊松比為隨機(jī)參數(shù)，對(duì)泊松比進(jìn)行隨機(jī)抽樣。泊松比的分布采用截?cái)嘈透咚狗植迹嚓P(guān)參數(shù)如表1所示。

表1 藥柱泊松比隨機(jī)分布參數(shù)Table 1 Random distribution parameters of Poisson＇s ratio

根據(jù)表1數(shù)據(jù)，對(duì)泊松比進(jìn)行Latin超立方抽樣300次，其均值隨抽樣次數(shù)的變化，如圖2所示。從圖2中可看出，當(dāng)抽樣300次時(shí)，其均值已趨于平穩(wěn)，說明對(duì)泊松比進(jìn)行300次抽樣，滿足固體火箭發(fā)動(dòng)機(jī)裝藥結(jié)構(gòu)完整性隨機(jī)分析的要求。

圖3所示為泊松比隨機(jī)分布的統(tǒng)計(jì)結(jié)果，圖形與截?cái)嘈透咚狗植际纸咏@也說明抽樣次數(shù)是足夠的。

圖2 藥柱泊松比均值隨抽樣次數(shù)的變化曲線Fig.2 Curve of mean of Poisson＇s ratio vs number of sampling

圖3 藥柱泊松比隨機(jī)抽樣柱狀分布圖Fig.3 Random sample histogram of Poisson＇s ratio

2 計(jì)算結(jié)果與討論

2.1 藥柱受非均布內(nèi)壓的計(jì)算結(jié)果及分析

利用三維粘彈性有限元法，計(jì)算發(fā)動(dòng)機(jī)在瞬態(tài)非均布?jí)毫ψ饔孟碌慕Y(jié)構(gòu)完整性。取5個(gè)時(shí)間點(diǎn)為研究對(duì)象，其壓力分布如圖4所示。固體火箭發(fā)動(dòng)機(jī)的Mises等效應(yīng)變分布如圖5所示。

從圖5中可看出，點(diǎn)火過程中不同時(shí)刻發(fā)動(dòng)機(jī)的Mises等效應(yīng)變分布。當(dāng)點(diǎn)火藥被引燃后，發(fā)動(dòng)機(jī)頭部瞬間充滿氣體，壓力增大。如圖5(a)所示，在發(fā)動(dòng)機(jī)藥柱頭部產(chǎn)生了較大的應(yīng)變。

圖4 不同時(shí)刻發(fā)動(dòng)機(jī)壓力分布圖Fig.4 Distribution of pressure at different times

圖5 壓力傳播過程Mises等效應(yīng)變分布Fig.5 Distribution of Mises strain during pressure diffusion

由于本文忽略了發(fā)動(dòng)機(jī)的脫粘處理，因此發(fā)動(dòng)機(jī)頭部的計(jì)算結(jié)果存在一定誤差，本文主要考慮發(fā)動(dòng)機(jī)圓管段在壓力傳播過程的Mises等效應(yīng)變分布。隨時(shí)間增加，壓力逐漸向發(fā)動(dòng)機(jī)尾部傳播。從圖5中可明顯看出，發(fā)動(dòng)機(jī)藥柱上存在應(yīng)變峰的移動(dòng)。當(dāng)發(fā)動(dòng)機(jī)內(nèi)部壓力趨于穩(wěn)定時(shí)，發(fā)動(dòng)機(jī)藥柱的應(yīng)變分布亦趨于穩(wěn)定。對(duì)比圖5中的計(jì)算結(jié)果可知，在壓力傳播過程中，發(fā)動(dòng)機(jī)藥柱的最大Mises等效應(yīng)變值比整個(gè)發(fā)動(dòng)機(jī)內(nèi)部壓力分布穩(wěn)定時(shí)的最大Mises等效應(yīng)變值大許多?？梢?，若采用發(fā)動(dòng)機(jī)內(nèi)部均勻壓力分布計(jì)算發(fā)動(dòng)機(jī)點(diǎn)火過程的結(jié)構(gòu)完整性是不妥當(dāng)?shù)摹?/p>

在壓力傳遞的過程中，發(fā)動(dòng)機(jī)藥柱的受力狀態(tài)是復(fù)雜的。為研究發(fā)動(dòng)機(jī)在點(diǎn)火過程藥柱的復(fù)雜受力狀態(tài)，取5.84×10－3s時(shí)的發(fā)動(dòng)機(jī)為研究對(duì)象，如圖5(c)所示。建立路徑并統(tǒng)計(jì)路徑上節(jié)點(diǎn)的軸向應(yīng)變和Mises等效應(yīng)變分布，如圖6所示。從圖6中可看出，沿路徑分布的Mises等效應(yīng)變?cè)趬毫鞑サ降奈恢酶浇嬖诿黠@的波峰和波谷。軸向應(yīng)變分布是導(dǎo)致這一現(xiàn)象的一個(gè)主要原因，在壓力傳播的過程中，藥柱受到發(fā)動(dòng)機(jī)頭部和尾部殼體的約束，使得藥柱圓柱段受到內(nèi)壓載荷作用的部位軸向受拉，而未受到內(nèi)壓載荷的部位軸向受壓。藥柱頭部由于過渡圓弧的結(jié)構(gòu)，在受發(fā)動(dòng)機(jī)工作壓力作用時(shí)，其過渡圓弧處受壓。因此，可從圖6的軸向應(yīng)變分布曲線中看出存在兩個(gè)明顯的波谷，分別位于藥柱頭部和藥柱圓柱段。最大Mises等效應(yīng)變的位置基本上與軸向應(yīng)變負(fù)向最大值的位置一致，即出現(xiàn)在藥柱軸向受壓最大的位置。由于發(fā)動(dòng)機(jī)藥柱從受拉狀態(tài)到受壓狀態(tài)存在一個(gè)過渡階段，即在發(fā)動(dòng)機(jī)藥柱受拉段和受壓段存在一區(qū)域的軸向應(yīng)變接近零，該區(qū)域的Mises等效應(yīng)變很小。因此，在最大Mises等效應(yīng)變位置之前呈現(xiàn)波谷狀分布，這從圖6中可明顯看出。

圖6 5.84×10－3s時(shí)刻發(fā)動(dòng)機(jī)路徑上應(yīng)變分布Fig.6 Strain along SRM path at 5.84 ×10 －3s

以發(fā)動(dòng)機(jī)中間部位的節(jié)點(diǎn)14 639為對(duì)象，研究發(fā)動(dòng)機(jī)藥柱固定點(diǎn)在內(nèi)壓載荷傳遞過程中的受力狀況，如圖7所示。

圖7 節(jié)點(diǎn)14 639 Mises等效應(yīng)變隨時(shí)間變化曲線Fig.7 Mises strain-time curve of note 14 639

從圖7中可看出，在0.005 5 s時(shí)，節(jié)點(diǎn)的Mises等效應(yīng)變達(dá)到最大值。因?yàn)榇藭r(shí)發(fā)動(dòng)機(jī)藥柱的Mises應(yīng)變峰剛好移動(dòng)到節(jié)點(diǎn)14 639處，當(dāng)應(yīng)變峰繼續(xù)往發(fā)動(dòng)機(jī)后端移動(dòng)時(shí)，節(jié)點(diǎn)的Mises等效應(yīng)變迅速變小，在0.006 3 s時(shí)，節(jié)點(diǎn)位置剛好處于如圖6中所示的Mises應(yīng)變波谷，此時(shí)節(jié)點(diǎn)的等效應(yīng)變最小，隨應(yīng)變峰的繼續(xù)移動(dòng)，節(jié)點(diǎn)的Mises等效應(yīng)變逐漸增加，最后趨于穩(wěn)定。

2.2 藥柱隨機(jī)有限元計(jì)算結(jié)果及分析

為考慮泊松比對(duì)點(diǎn)火工況下發(fā)動(dòng)機(jī)藥柱結(jié)構(gòu)完整性的影響，采用如圖2所示的泊松比隨機(jī)分布數(shù)據(jù)，對(duì)發(fā)動(dòng)機(jī)在瞬態(tài)壓力分布作用下進(jìn)行有限元計(jì)算，共300個(gè)算例。仍以節(jié)點(diǎn)14 639為研究對(duì)象，統(tǒng)計(jì)應(yīng)變峰經(jīng)過該點(diǎn)時(shí)的Mises等效應(yīng)變，即該節(jié)點(diǎn)在點(diǎn)火過程中的最大Mises等效應(yīng)變，其柱狀分布圖如圖8所示。

圖8 藥柱節(jié)點(diǎn)14 639最大Mises應(yīng)變柱狀分布圖Fig.8 Mises strain maximum histogram of note 14 639

藥柱節(jié)點(diǎn)14 639最大Mises應(yīng)變隨機(jī)分布統(tǒng)計(jì)參數(shù)如表2所示。從表2中可看出，該節(jié)點(diǎn)統(tǒng)計(jì)的Mises應(yīng)變最大值和最小值相差僅0.571 16%。根據(jù)節(jié)點(diǎn)應(yīng)變均值和標(biāo)準(zhǔn)差可知，其變異系數(shù)為1.104 6%?？梢?，在瞬態(tài)非均布載荷的作用下，泊松比的隨機(jī)分布對(duì)發(fā)動(dòng)機(jī)Mises等效應(yīng)變的影響較小，這與采用均布載荷計(jì)算的結(jié)果有明顯差異。

表2 藥柱節(jié)點(diǎn)14 639最大Mises應(yīng)變分布參數(shù)Table 2 Distribution parameters of Mises strain maximum of note 14 639

節(jié)點(diǎn)14 639的最大Mises等效應(yīng)變隨藥柱泊松比的變化，如圖9所示。從圖9中可看出，泊松比越大，Mises等效應(yīng)變值越大，且呈線性分布。這與均布載荷下的所得到的結(jié)論完全相反。這主要與瞬態(tài)非均布載荷的作用下藥柱復(fù)雜的受力狀態(tài)有關(guān)。因此，不能直接認(rèn)定泊松比越大，Mises等效應(yīng)變值就會(huì)越小，應(yīng)充分考慮發(fā)動(dòng)機(jī)點(diǎn)火過程的壓力傳播對(duì)藥柱結(jié)構(gòu)完整性的影響。

圖9 藥柱節(jié)點(diǎn)14 639的最大Mises等效應(yīng)變隨泊松比的變化曲線Fig.9 Mises strain maximum-Poisson＇s ratio curve of note 14 639

3 結(jié)論

(1)固體火箭發(fā)動(dòng)機(jī)藥柱在瞬態(tài)非均布載荷的作用下，其受力狀態(tài)十分復(fù)雜。在點(diǎn)火過程中，發(fā)動(dòng)機(jī)藥柱上存在應(yīng)變峰的移動(dòng)，且計(jì)算得到的Mises等效應(yīng)變遠(yuǎn)大于均布載荷的計(jì)算結(jié)果。

(2)算例中藥柱泊松比的隨機(jī)分布對(duì)發(fā)動(dòng)機(jī)藥柱在瞬態(tài)非均布載荷作用下的結(jié)構(gòu)完整性影響不大，Mises等效應(yīng)變的變異系數(shù)僅為1.104 6%，且Mises等效應(yīng)變隨泊松比的增大而增大，呈線性分布。而在均布載荷下的泊松比的影響往往較大，且Mises等效應(yīng)變隨泊松比的增大而減小?？梢?，在瞬態(tài)非均布載荷與均布載荷下的結(jié)論差別較大，采用均布載荷計(jì)算藥柱應(yīng)變，并將其計(jì)算結(jié)果作為結(jié)構(gòu)完整性分析的依據(jù)有其局限性。因此，在對(duì)固體火箭發(fā)動(dòng)機(jī)結(jié)構(gòu)完整性分析時(shí)，應(yīng)充分考慮發(fā)動(dòng)機(jī)點(diǎn)火過程的壓力傳播過程，即壓力的瞬態(tài)非均勻分布。

［1］Finne S，F(xiàn)utsaether C，Botnan J.Three-dimensional analysis of solid propellant grains using a nonlinear viscoelastic model［R］.AIAA 1990.

［2］Shiang-Woei，Chyuan.Dynamic analysis of solid propellant grains subjected to ignition pressurization loading［J］.Journal o f Sound and Vibration，2003:465-483.

［3］鐘濤，張為華.點(diǎn)火瞬態(tài)過程對(duì)復(fù)合固體推進(jìn)劑力學(xué)響應(yīng)特性的影響［J］.國(guó)防科技大學(xué)學(xué)報(bào)，2004，26(3):7-10.

［4］孟紅磊，周長(zhǎng)省，鞠玉濤，等.非均布瞬態(tài)內(nèi)壓作用下星孔藥柱應(yīng)力分析［J］.固體火箭技術(shù)，2010，33(3):289-293.

［5］張海聯(lián)，周建平.固體火箭發(fā)動(dòng)機(jī)藥柱隨機(jī)結(jié)構(gòu)分析及其相關(guān)函數(shù)研究［J］.固體火箭技術(shù)，2001，24(4):25-28.

［6］張海聯(lián)，周建平.固體推進(jìn)劑藥柱泊松比隨機(jī)粘彈性有限元分析［J］.推進(jìn)技術(shù)，2001，22(3):245-249.

［7］田四朋，唐國(guó)金，李道奎，等.固體推進(jìn)劑三維粘彈性Monte-Carlo隨機(jī)有限元法［J］.強(qiáng)度與環(huán)境，2007，34(2):51-57.

［8］職世君，孫冰，張建偉.考慮泊松比的固體發(fā)動(dòng)機(jī)裝藥貯存壽命預(yù)估［J］.固體火箭技術(shù)，2011，34(5):593-597.

［9］Huang-Ta Chu，Jung-Hua Chou.Poisson ratio effect on stress behavior of propellant grains under ignition loading［J］.Journal of Propulsion and Power，2011，27(3):663-668.