侵徹條件下巖石類材料RHT模型參數(shù)敏感性分析

聶錚玥， 彭 永， 陳 榮， 蔡鎮(zhèn)清， 李付剛

(國防科技大學(xué) 文理學(xué)院，長沙 410073)

RHT(Riedel-Hiermaier-Thoma)模型廣泛應(yīng)用于巖石、混凝土類材料在高應(yīng)變率、大變形問題中的動態(tài)力學(xué)響應(yīng)模擬[1]，國內(nèi)外學(xué)者基于RHT模型開展的爆炸沖擊、彈體侵徹等研究都取得了比較理想的結(jié)果[2-4]。目前多數(shù)研究的重點在于，針對混凝土材料，基于實驗數(shù)據(jù)對部分參數(shù)或模型進行修正與改進[5-8]。然而，當RHT模型應(yīng)用于巖石材料時，由于不同種類巖石間的力學(xué)性能差異明顯，一般需通過理論分析與力學(xué)實驗相結(jié)合的方法來確定模型參數(shù)，不宜直接采用某種特定巖石或混凝土材料的參考值[9]。此外，RHT模型參數(shù)多達34個，獲得參數(shù)所需力學(xué)實驗包括物理性質(zhì)測定、靜態(tài)力學(xué)及動態(tài)力學(xué)等多種類型，使得參數(shù)確定比較復(fù)雜。因此，對模型中參數(shù)進行敏感性分析，甄別敏感性較強的參數(shù)，可以優(yōu)化模型參數(shù)確定方案，提高數(shù)值模擬的效率。

對于RHT模型參數(shù)敏感性的研究相對較少，李洪超等[10-11]采用LS-DYNA數(shù)值軟件并通過正交試驗的方法，模擬分析了RHT模型中的部分復(fù)雜參數(shù)對大理巖及紅砂巖霍普金森桿(split Hopkinson pressure bar, SHPB)壓縮實驗應(yīng)力-應(yīng)變曲線結(jié)果的敏感性，通過極差分析，給出了影響明顯的材料模型參數(shù)排序。該方法實驗設(shè)計合理，且能夠在全局上獲得擬合SHPB沖擊曲線的最佳解。但由于SHPB沖擊曲線結(jié)果將受材料物理性質(zhì)影響而具有一定離散性，且不同應(yīng)變率條件下獲得的SHPB沖擊曲線也將不同，因此該方法有一定局限性。劉殿柱等[12]同樣利用正交試驗法，以侵徹深度為評價指標，利用Autodyn軟件對彈丸侵徹混凝土問題進行了部分參數(shù)的敏感性分析。該方法簡明直觀，但根據(jù)侵徹深度變化趨勢圖對敏感性進行排序時，由于部分參數(shù)影響相對不顯著，因此結(jié)果具有一定主觀性。辛健[13]通過模擬TNT在土中的爆炸過程，采用單因素變量法，得到了以損傷面積為評價指標的參數(shù)敏感度結(jié)果。該方法簡潔直觀，但參數(shù)水平設(shè)置情況通常有限，且RHT模型中參數(shù)并非完全獨立，因此分析結(jié)果可能存在一定誤差。為客觀評價模型參數(shù)的敏感性差異，需要對計算參數(shù)矩陣進行合理設(shè)計與分析。同時，上述研究結(jié)果表明，針對TNT爆炸、SHPB試驗與侵徹問題，RHT模型材料參數(shù)的敏感性不同。因此，針對特定問題，需要進行相應(yīng)的參數(shù)敏感性分析，不能直接采用其他物理問題的敏感性分析結(jié)果。

針對實心彈體中低速正侵徹單層半無限厚巖石靶體問題，本研究基于ISIGHT軟件中的優(yōu)化拉丁超立方設(shè)計 (optimal Latin hypercube design，Opt LHD)算法設(shè)計了參數(shù)試驗矩陣，保證了設(shè)計矩陣具有較優(yōu)的空間填充性。采用LS-DYNA軟件對侵徹過程進行了大量數(shù)值模擬，獲得了不同參數(shù)組合條件下的仿真結(jié)果。然后以彈體侵徹深度、彈體質(zhì)量損失率及靶體損傷度為評價指標，分析了模型參數(shù)的交互效應(yīng)及主效應(yīng)。最后通過建立侵徹深度、彈體質(zhì)量損失率及靶體損傷度與參數(shù)間的響應(yīng)面模型，定量計算了RHT模型中參數(shù)的貢獻率。

1 RHT模型簡介及參數(shù)分類

1.1 RHT模型介紹

RHT模型由Riedel等[14]提出，可以反映材料的應(yīng)變率效應(yīng)、失效面及損傷演化等特點。RHT模型包括狀態(tài)方程部分及本構(gòu)方程部分，如圖1所示。

圖1 RHT模型原理圖

狀態(tài)方程采用Herrmann[15]提出的P-α狀態(tài)方程，以描述多孔及疏松介質(zhì)中孔隙從壓碎到壓實的過程。本構(gòu)模型由3個發(fā)展階段組成，隨著應(yīng)力增加，材料首先經(jīng)過彈性階段到達彈性屈服面，此后材料發(fā)生塑性變形并進入線性強化階段，材料表現(xiàn)出應(yīng)變硬化特征，直至到達失效面。當?shù)刃?yīng)力強度超出失效應(yīng)力強度后，材料開始出現(xiàn)累積損傷量，進入損傷軟化階段，最后到達殘余強度面[14-16]。

根據(jù)RHT模型原理，基于參數(shù)物理意義，可對34個參數(shù)進行初步歸類，如表1所示。

表1 RHT模型參數(shù)分類

1.2 參數(shù)分類及范圍選取

根據(jù)Rankine-Hugoniot方程及Mie-Grüneisen狀態(tài)方程可以求得部分參數(shù)間關(guān)系為

(1)

式中：c0為材料聲速；s為材料參數(shù)，其值可以通過平板撞擊試驗中測定的材料沖擊波速度和波后粒子速度結(jié)果線性擬合確定。其余參數(shù)間亦存在下列關(guān)系[17-18]

pel=fc/3
βc=4/(20+3fc)
βt=2/(20+fc)

(2)

因此，可將c0，s，fc作為主定參數(shù)，而A1，A2，A3，B0，B1，T1，pel，βc，βt作為其從屬參數(shù)處理。綜合各參數(shù)定義及式(1)、式(2)，主定參數(shù)、從屬參數(shù)及定值參數(shù)的具體分類如表2所示，其中主定參數(shù)共計21個。

表2 主定參數(shù)、從屬參數(shù)、定值參數(shù)分類

本研究針對巖石類材料，參考花崗巖、玄武巖、砂巖、千枚巖、大理巖、片麻巖共六種巖石的相關(guān)指標值[19]，確定了部分主定參數(shù)的取值范圍。對于材料參數(shù)s的取值范圍則根據(jù)巖石材料的平板撞擊試驗結(jié)果確定[20-23]。對于無法直接確定的其余主定參數(shù)，則利用C30/37混凝土RHT模型參數(shù)值作為基準值，并設(shè)置參數(shù)變化范圍為基準值的-80%～80%。所有主定參數(shù)的取值范圍如表3所示。

表3 主定參數(shù)參考取值范圍

2 敏感性分析數(shù)值實驗設(shè)計

2.1 參數(shù)矩陣設(shè)計

ISIGHT是一款可以靈活設(shè)計并調(diào)用各項軟件進行計算，進而對結(jié)果進行后處理的計算機輔助優(yōu)化軟件[24]。本研究利用ISIGHT中的試驗設(shè)計(design of experiments, DOE)功能，采用優(yōu)化拉丁超立方設(shè)計算法(Opt LHD)進行參數(shù)矩陣設(shè)計，從而保證試驗點具有較好的空間填充性及均勻度。

根據(jù)RHT模型中21個主定參數(shù)設(shè)計參數(shù)矩陣，為滿足后續(xù)的交互效應(yīng)及主效應(yīng)分析，以及建立響應(yīng)面模型的精度需求，參數(shù)水平個數(shù)的設(shè)計需滿足最小采樣點數(shù)。對于一階、二階響應(yīng)面模型，其最小采樣點數(shù)分別通過式(3)、式(4)計算

Smin,1=M+1

(3)

(4)

式中：Smin為最小采樣點數(shù)；M為試驗設(shè)計的主定參數(shù)個數(shù)。

因此，考慮到最小采樣點數(shù)要求及計算成本，基于上節(jié)中的主定參數(shù)范圍，利用Opt LHD算法設(shè)置了400個采樣點，獲得了圖2所示的參數(shù)設(shè)計矩陣。

圖2 主定參數(shù)設(shè)計矩陣

2.2 數(shù)值實驗設(shè)計

利用LS-DYNA有限元軟件對實心彈體中低速正侵徹半無限厚單層巖石靶體問題進行仿真。為節(jié)省計算時間，采用1/4對稱模型進行仿真實驗，尺寸及網(wǎng)格設(shè)計如圖3所示，彈體在靶體中心處以600 m/s的初速度垂直向下侵徹。靶體的上表面設(shè)置為自由邊界，下表面設(shè)置為無反射邊界，以近似模擬半無限厚條件。

圖3 模型尺寸及網(wǎng)格設(shè)計(mm)

彈體的材料選擇為硬度較大的 30CrMnSiNi2A 鋼，使用Johnson-Cook模型描述其材料力學(xué)性能，模型參數(shù)如表4所示。靶體材料采用RHT模型，其中主定參數(shù)為參數(shù)設(shè)計矩陣值，其余參數(shù)確定如1.2節(jié)所述。

表4 彈體Johnson-cook模型參數(shù)[25]

3 RHT模型參數(shù)敏感性分析方法

為了定量描述參數(shù)敏感性大小以及參數(shù)變化對評價指標(即響應(yīng)量)的影響規(guī)律，需對參數(shù)的交互效應(yīng)、主效應(yīng)以及參數(shù)對響應(yīng)量的貢獻率進行分析。參數(shù)的交互效應(yīng)可以反映某一參數(shù)變化對響應(yīng)量的貢獻隨其他參數(shù)的取值不同而發(fā)生變化的情況。參數(shù)交互效應(yīng)顯著則此時參數(shù)作用非獨立。當交互效應(yīng)不顯著時，則通過參數(shù)主效應(yīng)分析，可得到參數(shù)對響應(yīng)量的影響規(guī)律。通過建立參數(shù)與響應(yīng)量之間的近似模型，進一步計算各參數(shù)貢獻率，評價參數(shù)的敏感性大小。

3.1 交互效應(yīng)、主效應(yīng)分析原理

參數(shù)交互效應(yīng)及主效應(yīng)的計算是先通過式(5)得到主定參數(shù)xi與響應(yīng)量y之間多項式(以二次多項式為例)

(5)

式中：xixj為參數(shù)交互項；εi為誤差項。

對式(5)微分可得

dy=c1dx1+…+2cn+1x1dx1+…+
c2n+1d(x1x2)+…+εi

(6)

當式(6)中交互項為零或值很小時，可得到各參數(shù)的主效應(yīng)。線性項的主效應(yīng)為

Mxi=cidxi

(7)

二階項的主效應(yīng)為

(8)

當交互項不可忽略時，則參數(shù)間的交互效應(yīng)如式(9)所示

Mxixj=cmd(xixj)

(9)

3.2 參數(shù)貢獻率計算原理

為了計算各個參數(shù)的貢獻率，需建立樣本點和評價指標之間的近似模型。首先通過式(10)將各主定參數(shù)值歸一化至[-1,+1]，以消除21個主定參數(shù)自身數(shù)值大小的影響。

(10)

主定參數(shù)的貢獻率采用響應(yīng)面模型方法來計算，模型可分為多元一次回歸模型、多元二次回歸模型等，表達式分別為

(11)

(12)

(13)

3.3 評價指標選取

如圖4所示，不同參數(shù)條件下，彈體侵徹深度、彈體磨蝕程度及靶體損傷程度均不同。因此，分別選擇侵徹深度H、彈體質(zhì)量損失率及靶體損傷度D作為評價指標，進行參數(shù)敏感性分析。

圖4 不同工況下的仿真結(jié)果

其中，侵徹深度根據(jù)計算結(jié)果直接獲得。彈體的質(zhì)量損失率通過提取彈體單元體質(zhì)量變化曲線計算得到。而對于靶體損傷度的表征，則是通過將侵徹過程中已經(jīng)刪除的單元損傷度定義為1，并提取剩余所有單元的損傷度進行累加獲得。典型靶體的全模型損傷分布如圖5所示。

圖5 工況41全模型的損傷情況

4 RHT模型參數(shù)敏感性結(jié)果

數(shù)值模擬獲得的有效工況共計355個。剩余45個工況則因為參數(shù)取值不合理，導(dǎo)致計算時靶體單元體中出現(xiàn)復(fù)雜聲速，或由于靶體強度過低，彈體將貫穿靶體而不采用。提取各工況的侵徹深度、彈體質(zhì)量損失率及靶體損傷度作為敏感性分析的評價指標進行后續(xù)分析。

4.1 參數(shù)交互效應(yīng)分析

圖6 對于不同評價指標的部分參數(shù)間的交互效應(yīng)圖

因此，對于同一評價指標，不同參數(shù)之間的交互效應(yīng)不同。同時，對于不同評價指標，某兩個參數(shù)之間的交互效應(yīng)也可能存在差異。為了進一步描述參數(shù)之間的交互效應(yīng)，取21個主定參數(shù)進行總體參數(shù)間交互效應(yīng)分析。

SPSS(statistical product and service solutions)軟件是一個功能強大的統(tǒng)計分析軟件，通過SPSS軟件可以對參數(shù)交互效應(yīng)進行計算，并獲得相應(yīng)的顯著性結(jié)果(significance，Sig)，當Sig值小于0.05時，說明此時參數(shù)交互效應(yīng)顯著，參數(shù)間將互相影響，反之則不顯著。結(jié)果表明：當評價指標為侵徹深度時，21個主定參數(shù)組合的交叉效應(yīng)Sig值為0.296>0.05；當評價指標為彈體質(zhì)量損失率時，交叉效應(yīng)Sig值為0.397>0.05；當評價指標為靶體損傷度時，交叉效應(yīng)Sig值為0.995>0.05。因此，對于侵徹深度、彈體質(zhì)量損失率與靶體損傷度而言，總體主定參數(shù)組合間的交互效應(yīng)均不顯著。

4.2 參數(shù)主效應(yīng)分析

為了研究參數(shù)變化對評價指標的影響規(guī)律，針對各評價指標進行參數(shù)主效應(yīng)分析，結(jié)果如圖7～圖9所示。其中，橫坐標表示參數(shù)取值從低水平(對應(yīng)1)依次變化至高水平(對應(yīng)2)，縱坐標表示計算得到的評價指標的統(tǒng)計值。對于侵徹深度與靶體損傷度，各主定參數(shù)的主效應(yīng)均呈非線性及非單調(diào)性變化。當各主定參數(shù)從低水平值變化至高水平值時，侵徹深度及靶體損傷度均將先增大后減小。對于彈體質(zhì)量損失率，除了材料參數(shù)s外，其余參數(shù)對指標的主效應(yīng)亦呈非線性及非單調(diào)性變化。而當s值增加時，彈體的質(zhì)量統(tǒng)計損失率將降低。

圖7 侵徹深度主效應(yīng)圖

圖8 彈體質(zhì)量損失率主效應(yīng)圖

圖9 靶體損傷度主效應(yīng)圖

4.3 近似模型及參數(shù)貢獻率分析

分別建立一階、二階、三階響應(yīng)面模型，并利用355個侵徹深度、彈體質(zhì)量損失率及靶體損傷度結(jié)果進行誤差分析，結(jié)果如圖10所示。其中，平行于橫坐標的水平線表示用于計算的所有樣本點的評價指標的實際平均值，對角線表示實際值與模型預(yù)測值相等的情況，點越靠近對角線則表明模型預(yù)測越精確。同時，為了定量評價各響應(yīng)面模型對所有實際評價指標值的回歸擬合效果，計算各模型的可決系數(shù)R2(coefficient of determination)如表5所示，當R2越接近1，則回歸擬合效果越好。由圖10可知，對于侵徹深度、彈體質(zhì)量損失率及靶體損傷度，均為一階響應(yīng)面模型誤差顯著大于二階與三階響應(yīng)面模型，即主定參數(shù)與評價指標間的函數(shù)關(guān)系具有較強的非線性。同時，由于二階模型與三階模型的可決系數(shù)差異較小，考慮到三階響應(yīng)面模型較二階模型更為復(fù)雜，故選擇二階響應(yīng)面模型進行后續(xù)參數(shù)貢獻率計算。

圖10 各階響應(yīng)面模型誤差分析圖

表5 各階響應(yīng)面模型可決系數(shù)R2結(jié)果

基于二階響應(yīng)面模型，分別對侵徹深度、彈體質(zhì)量損失率及靶體損傷度計算各參數(shù)項的貢獻率，并提取絕對值較大的前21項，如圖11所示。結(jié)果表明，對于侵徹深度及靶體損傷度而言，參數(shù)二次項對評價指標的影響普遍較一次項及交互項大。對于彈體質(zhì)量損失率而言，貢獻率較大的參數(shù)項同樣以二次項居多。

為了較全面地衡量各主定參數(shù)的貢獻率情況，按下式計算各個主定參數(shù)總體貢獻率

(14)

表6 各主定參數(shù)總體貢獻率

5 結(jié) 論

本文針對彈體侵徹半無限厚巖石靶體問題，利用優(yōu)化拉丁超立方設(shè)計算法設(shè)計了參數(shù)試驗矩陣，并進行數(shù)值模擬。根據(jù)355個有效結(jié)果，以侵徹深度、彈體質(zhì)量損失率及靶體損傷度為評價指標，對21個主定參數(shù)進行了參數(shù)敏感性分析，并得到了以下結(jié)論：

(1) 主定參數(shù)間存在一定交互效應(yīng)，且不同參數(shù)之間的交互效應(yīng)不同。同時，對于不同評價指標，參數(shù)間的交互效應(yīng)亦存在差異。取21個主定參數(shù)對侵徹深度、彈體質(zhì)量損失率及靶體損傷度進行總體參數(shù)組合間交互效應(yīng)分析，其結(jié)果均為交互效應(yīng)不顯著。

(2) 對于侵徹深度與靶體損傷度，各主定參數(shù)的主效應(yīng)曲線均呈非線性及非單調(diào)性變化。隨著參數(shù)水平值增加，侵徹深度及靶體損傷度的統(tǒng)計值均將先增大后減小。而對于彈體質(zhì)量損失率，材料參數(shù)s的值增加時，彈體質(zhì)量統(tǒng)計損失率將降低。其余參數(shù)的主效應(yīng)則與侵徹深度和靶體損傷度相似。

(3) 對于不同評價指標，各參數(shù)的主效應(yīng)不一致。對于侵徹深度及彈體質(zhì)量損失率而言，殘余應(yīng)力強度指數(shù)的主效應(yīng)最顯著。此外，當參數(shù)水平變化時，各主定參數(shù)間的主效應(yīng)相對顯著情況將發(fā)生變化。