鄭 健，韓 波，周長省

（南京理工大學(xué) 機械工程學(xué)院，南京210094）

端羥基聚丁二烯（hydroxyl terminated polybutadiene，HTPB）推進劑是一種常用的復(fù)合推進劑，其力學(xué)行為十分復(fù)雜，其性能依賴于溫度和應(yīng)變率，為了能夠準確地描述出該種推進劑裝藥在受力狀態(tài)下的力學(xué)行為，需要使用粘彈性本構(gòu)模型。粘彈性本構(gòu)參數(shù)獲取的準確性影響仿真分析的準確性，而常規(guī)粘彈性實驗參數(shù)的獲取方法均存在較大的實驗誤差。為了解決這一問題，Zapas－Phillips和Sorvaris等人［1－2］在研究了由試驗機加載產(chǎn)生的誤差問題之后，分別使用近似積分處理的方法來獲取較為準確的松弛模量數(shù)據(jù)。Zapas－Phillips法使用積分中值定理來獲得修正的松弛模量數(shù)據(jù)，但是僅能得到t≥t1／2時刻的松弛模量數(shù)據(jù)，其中，t為松弛時間，t1為達到恒應(yīng)變的時間。Sorvaris法在此基礎(chǔ)上進一步研究、推導(dǎo)出了t≥0時刻的松弛模量數(shù)據(jù)。孟紅磊等人［3］在Sorvari方法的基礎(chǔ)上提出了改進型的Sorvari方法，即引入修正因子α來更加精確地近似積分過程。許進升等人［4］基于直接積分方法來獲取復(fù)合HTPB推進劑的松弛模量數(shù)據(jù)。

本文提出了一種基于遺傳算法的HTPB復(fù)合推進劑粘彈性本構(gòu)參數(shù)獲取方法，為固體火箭發(fā)動機裝藥結(jié)構(gòu)完整性分析提供準確的本構(gòu)參數(shù)。

1 本構(gòu)參數(shù)獲取

1.1 問題描述

HTPB推進劑是一種典型的粘彈性材料，其力學(xué)特性呈現(xiàn)出彈性固體和粘性流體的性質(zhì)，在固體火箭發(fā)動機裝藥完整性分析中常使用線粘彈性本構(gòu)模型來表征推進劑的力學(xué)行為。三維積分型的各向同性線粘彈性材料應(yīng)力－應(yīng)變關(guān)系可以表示為

應(yīng)力－應(yīng)變關(guān)系的應(yīng)力偏量部分sij和體變部分σkk分別為

式中：λ（t），G（t），K（t）分別為拉梅常數(shù)對應(yīng)的松弛函數(shù)、剪切松弛函數(shù)和體積松弛函數(shù)；δij為位移分量，eij為偏應(yīng)變張量，εkk為應(yīng)變偏量。對于各向同性線粘彈性材料，必須具備2種材料函數(shù)才能準確描述其應(yīng)力－應(yīng)變特性。但是通過實驗直接獲得粘彈性材料的材料函數(shù)G（t）和K（t）比較困難，而材料的松弛模量E（t）較容易獲得。因此，如果獲得了材料函數(shù)之間的轉(zhuǎn)化關(guān)系，就可以通過較容易測量的材料函數(shù)來表征線粘彈性材料的力學(xué)行為。

將式（2）進行Laplace變換，得到：

式中：s為Laplace變換的復(fù)變量。彈性力學(xué)中材料常數(shù)存在以下關(guān)系：

式中：E為松弛模量，ν為泊松比。

將式（4）中的參數(shù)進行Laplace變換，得到：

對式（5）、式（6）進行Laplace逆變換可以得到線粘彈性材料中E（t），ν（t）和G（t），K（t）之間的關(guān)系。實際工程應(yīng)用中，可以根據(jù)實際情況簡化材料函數(shù)之間的變化關(guān)系。在固體推進劑有限元計算中經(jīng)常采用下列假設(shè)：由于泊松比ν（t）隨時間變化很小，故假設(shè)泊松比為常數(shù)，則根據(jù)松弛模量E（t）可以得到剪切松弛模量G（t）和體積松弛模量K（t）：

因此，如何準確地獲得推進劑的松弛模量對結(jié)構(gòu)完整性分析十分重要。

1.2 傳統(tǒng)實驗方法

現(xiàn)行的復(fù)合推進劑松弛模量測定方法是使用拉伸試驗機快速地將推進劑試樣拉至某一恒定應(yīng)變，同時測量出推進劑的松弛應(yīng)力，從而擬合出該推進劑材料的靜態(tài)松弛模量［1－2］。

圖1給出了階躍位移和加載示意圖。

圖1 階躍位移加載示意圖

圖1 （a）表示試驗機系統(tǒng)的理想加載過程，即不考慮加載歷程。但真實情況下試驗機通過調(diào)節(jié)控制系統(tǒng)的PID參數(shù)，實際的上升階段可能存在超調(diào)、緩慢爬坡等現(xiàn)象，即存在加載歷程，如圖1（b）所示。且HTPB推進劑的松弛響應(yīng)較快，因此采用圖1（a）這種數(shù)據(jù)處理方法會造成實驗獲得的靜態(tài)松弛模量數(shù)據(jù)偏小。為了解決這一問題，很多研究學(xué)者提出了各種數(shù)據(jù)處理方法，但這些方法均假設(shè)試驗機的加載階段近似為一個直線上升的階段，即如圖1（c）中所示的加載歷程。國外Zapas－Phillips采用積分中值定理獲得修正的t≥t1／2時刻的松弛模量數(shù)據(jù)，其計算式如下：

Sorvaris［1－2］在 Zapas－Phillips的研究基礎(chǔ)上，推導(dǎo)出了t≥0時刻的松弛模量數(shù)據(jù)表達式：

式中 ：σ為應(yīng)力為應(yīng)力變化率 ，ε為應(yīng)變?yōu)閼?yīng)變變化率。

國內(nèi)孟紅磊［3］在Sorvari方法的基礎(chǔ)上提出了改進型的Sorvari方法，即引入修正因子α來更加精確地近似積分過程：

許進升等人［4］基于直接積分方法來獲取復(fù)合HTPB推進劑的松弛模量數(shù)據(jù)，即通過理論推導(dǎo)出Prony級數(shù)形式的松弛函數(shù)在圖1（c）加載條件下松弛應(yīng)力的顯式表達式，再通過最小二乘法來獲得靜態(tài)松弛模量。

雖然這些方法在一定程度上提高了實驗數(shù)據(jù)的準確性，但是與實驗過程中真實的加載過程仍存在一定的差異。

1.3 基于遺傳算法的松弛模量獲取

使用基于遺傳算法的全階段松弛模量數(shù)據(jù)擬合方法解決實驗系統(tǒng)和數(shù)據(jù)處理帶來的松弛數(shù)據(jù)不準確問題。遺傳算法（genetic algorithm，GA）是由美國Michgan大學(xué)的Holland J教授于1975年最早提出，它是一種借鑒了生物學(xué)進化規(guī)律的算法，該算法模擬了生物學(xué)的遺傳、突變、自然選擇和雜交等現(xiàn)象。對于一個優(yōu)化問題，首先從解空間中隨機抽取一定數(shù)量的候選解形成一個初始種群，使用一種編碼方法對初始種群進行編碼。上一代種群內(nèi)部通過遺傳、變異、雜交等過程形成新一代種群，通過算法的適應(yīng)度函數(shù)來判斷新種群的優(yōu)劣，通過這一過程來模擬自然界的優(yōu)勝劣汰。經(jīng)過若干代的優(yōu)勝劣汰獲得最優(yōu)的種群，即最優(yōu)解。

遺傳算法的基本步驟如下：

①編碼。使用某種編碼方式將解空間映射到遺傳算法的編碼空間中。

②產(chǎn)生初始種群。在編碼空間中隨機產(chǎn)生一定數(shù)量的碼形成初始種群，這些初始種群對應(yīng)著一組初始解。

③適應(yīng)度評估。使用適應(yīng)度函數(shù)來判斷種群中個體的優(yōu)劣。

④遺傳算子的確定。遺傳算子中確定了選擇、交叉、變異的整個過程。

⑤算法的終止條件。確定遺傳算法的終止條件。

使用MATLAB編制了擬合HTPB推進劑松弛曲線的計算程序，以擬合函數(shù)和實驗曲線的重合度作為遺傳算法中的目標函數(shù)，目標函數(shù)的表達式為

式中：f為適應(yīng)度，n為計算過程中的迭代次數(shù)，即種群繁衍的代數(shù)；σl（X）代表優(yōu)化過程中解空間內(nèi)某一個體在給定載荷作用下得到的應(yīng)力值；σe為實驗獲得的應(yīng)力值。使目標函數(shù)表達式獲得最小值的X即為最優(yōu)解，即為所需要的松弛模量擬合參數(shù)。為避免目標函數(shù)在尋優(yōu)過程中數(shù)值過小影響計算精度，將式（11）乘以一個較大的系數(shù)k可以將目標函數(shù)適當(dāng)放大（此處k取10 000）。

2 實驗方法和結(jié)果處理

為了獲得準確的推進劑松弛模量數(shù)據(jù)，使用美國BOSS動態(tài)力學(xué)實驗機進行實驗，該實驗機具有加載速度快、采集系統(tǒng)精度高等特點。為了能夠滿足動態(tài)力學(xué)實驗機的材料安裝要求，參照文獻［5］中的單軸拉伸實驗方法設(shè)計了單軸應(yīng)力松弛試樣。

圖2為松弛實驗的試樣示意圖。使用刀片將HTPB推進劑切割成長條形試樣，試樣尺寸B＝5mm，W＝5mm，H＝30mm。使用改性丙烯酸膠粘劑將HTPB推進劑試樣上下兩端和木質(zhì)接頭粘接，試樣制作完成后放入密封干燥箱中在室溫下靜置24h進行固化。實驗設(shè)備如圖3所示，整個實驗系統(tǒng)由控制柜、動態(tài)力學(xué)分析儀DMA、保溫箱、液氮罐和非接觸式應(yīng)變測量系統(tǒng)組成。HTPB推進劑模量較低且容易造成損傷，不宜直接使用引伸計來測量拉伸應(yīng)變。直接使用實驗機的加載位移來反推試樣應(yīng)變存在較大的誤差，因此實驗中采用基于數(shù)字圖像相關(guān)性的非接觸式應(yīng)變測量系統(tǒng)來測量拉伸過程中推進劑的應(yīng)變。

圖2 松弛試樣示意圖

圖3 DMA實驗系統(tǒng)

以松弛實驗過程中的實際位移變化曲線作為輸入條件開始遺傳算法的計算，圖4給出了優(yōu)化計算過程中的優(yōu)化歷程，圖中，n為計算過程中種群繁衍的代數(shù)，f為種群在優(yōu)化過程中的適應(yīng)度。圖中給出了適應(yīng)度的最優(yōu)值和平均值。從圖中可以看出，在種群的繁衍過程中種群的適應(yīng)度不斷下降，這表明種群的質(zhì)量在逐步提高，經(jīng)過292代之后種群的最佳適應(yīng)度之間達到一個穩(wěn)定值，并認為第292代中的最優(yōu)個體即為最優(yōu)解。計算結(jié)果如表1所示。表中，E0＝8MPa（瞬態(tài)松弛模量），N為Prony級數(shù)的階數(shù)，EN，tN分別為各階松弛模量、松弛時間。

圖4 遺傳算法計算過程中適應(yīng)度最優(yōu)值和平均值變化情況

表1 松弛實驗數(shù)據(jù)

分別使用傳統(tǒng)計算方法和遺傳算法獲得了HTPB推進劑的松弛模量數(shù)據(jù)，將該結(jié)果應(yīng)用于真實松弛實驗的加載過程中，所得的應(yīng)力σ隨時間變化結(jié)果如圖5所示。從圖中可以看出，遺傳算法結(jié)果和實驗結(jié)果吻合良好，而傳統(tǒng)方法所得到的結(jié)果明顯小于實驗結(jié)果；隨著時間的增加，傳統(tǒng)方法和遺傳算法的差距逐漸減小。由于傳統(tǒng)方法忽略了實驗過程中加載階段對實驗結(jié)果的影響，因此其獲得的結(jié)果偏低。但是加載階段對于HTPB推進劑的平衡模量的影響隨著時間逐漸減小，因此傳統(tǒng)方法和遺傳算法的差距會隨著時間的推移逐漸減小，直至重合。

圖5 遺傳算法和傳統(tǒng)算法結(jié)果對比

3 結(jié)論

本文建立了基于遺傳算法的復(fù)合推進劑松弛模量獲取方法，為該推進劑的粘彈性本構(gòu)模型提供了較精確的本構(gòu)參數(shù)，獲得如下結(jié)論：①由于拉伸試驗機系統(tǒng)存在加載歷程，傳統(tǒng)松弛模量獲取方法獲取的數(shù)據(jù)較真實值偏??；②基于遺傳算法的松弛模量擬合方法考慮了真實實驗過程中的加載歷程，因此可以準確地獲得復(fù)合推進劑的松弛模量。

［1］SORVARI J，MALINEN M.On the direct estimation of creep and relaxation functions［J］.Mechanics of Time－dependent Materials，2007（11）：143－157.

［2］SORVARI J，MALINEN M.Numerical interconversion between linear viscoelastic material function with regularization［J］.Solids and Structures，2007，44：1 291－1 303.

［3］孟紅磊.改性雙基推進劑裝藥結(jié)構(gòu)完整性數(shù)值仿真方法研究［D］.南京：南京理工大學(xué)，2011.MENG Hong－lei.Research on numerical simulation method of structure integrity analysis for modified double base propellant grain［D］.Nanjing：Nanjing University of Science and Technology，2011.（in Chinese）

［4］許進升，鞠玉濤，鄭健.復(fù)合固體推進劑松弛模量的獲取方法［J］.火炸藥學(xué)報，2011，34（5）：58－62.XU Jin－sheng，JU Yu－tao，ZHENG Jian.Acquisition of the relaxation modulus of composite propellant［J］.Chinese Journal of Explosive ＆ Propellants，2011，34（5）：58－62.（in Chinese）［5］JUNG G D，YOUN S K.A nonlinear viscoelastic constitutive model of solid propellant［J］.International Journal of Solids and Structures，1999，36（25）：3 755－3 777.