孟利平， 程遠征， 張倫平， 王海坤， 汪 俊， 劉建湖

（中國船舶科學研究中心， 江蘇 無錫214082）

0 引 言

隨著計算方法和計算程序的不斷完善，通過數(shù)值仿真進行艦船結(jié)構(gòu)抗爆抗沖擊性能研究已越來越重要，而材料模型和相關(guān)參數(shù)的正確與否直接關(guān)系到仿真結(jié)果的正確性，進而影響結(jié)構(gòu)抗爆抗沖擊性能評估結(jié)論的正確性[1-3]。 Q345B 鋼是工程結(jié)構(gòu)中廣泛使用的低合金鋼，是船舶領(lǐng)域中軍輔船、民船的主要結(jié)構(gòu)材料，還是艦艇設(shè)備的結(jié)構(gòu)用材料。所以針對Q345B 鋼動態(tài)力學性能的研究對于軍民用船結(jié)構(gòu)的均衡設(shè)計、準確評估和防護優(yōu)化等都具有很重要的意義。

國內(nèi)對Q345B 鋼動態(tài)力學性能的研究卻不多。 部分研究者[4-5]進行了Q345 鋼的動態(tài)試驗研究，但都局限在高應變率（＞103/s）對屈服強度的影響，缺乏中低應變率的試驗數(shù)據(jù)，且對于各應變率下的大應變及斷裂特性等缺乏較系統(tǒng)的試驗研究。 對于船用鋼材的模型及參數(shù)選取，國內(nèi)大部分研究者均采用低碳鋼的CS 模型參數(shù)，姚熊亮[6]在對船用917 鋼的研究中指出，在進行水下爆炸等仿真研究時，針對不同材料應采用試驗等方法得到合理的參數(shù)。 此外，由于鋼材在拉伸破壞時會發(fā)生頸縮，這導致斷裂判據(jù)在有限元計算中受網(wǎng)格的影響很大， 如何避免這一影響在實際計算中還缺乏一個定量化的指導。

本文針對Q345B 鋼進行了應變率0.002/s～3 040/s 范圍內(nèi)的單向拉伸試驗，并進行了壓縮、扭轉(zhuǎn)及不同缺口半徑試件拉伸試驗， 得到了Q345B 鋼在不同應變率和不同初始應力三軸度下的變形至斷裂的實驗結(jié)果，據(jù)此分析并修正了CS 本構(gòu)模型，擬合了JC 失效模型，結(jié)合數(shù)值仿真，給出了失效判據(jù)在有限元計算中針對不同尺寸網(wǎng)格的使用方法，為Q345B 鋼結(jié)構(gòu)的工程設(shè)計與應用提供了依據(jù)。

1 試驗研究

1.1 不同應變率拉伸試驗

采用MTS 準靜態(tài)材料試驗機、Instron 高速材料試驗機以及Hopkinson 拉桿進行了應變率在0.002/s～3 040/s 范圍內(nèi)的拉伸試驗，具體的應變率值及相應加載方式見表1，每組試驗進行3～5 次重復試驗保證了數(shù)據(jù)的有效性。 所有試驗均采用相同標距的試件以消除由于試件尺寸不同可能帶來的影響。 試件形狀為板狀，厚度為4 mm，試件尺寸示意圖及實物圖如圖1 所示，試件兩端夾持段的長度依試驗機夾持要求略有不同。

表1 Q345B 鋼動態(tài)力學性能Tab.1 Dynamic mechanical property of Q345B steel

圖1 試件形狀尺寸示意圖及實物圖Fig.1 Shapes and sizes of test specimens

1.2 不同應力三軸度試驗

采用電子萬能試驗機和電子扭轉(zhuǎn)試驗機進行了準靜態(tài)壓縮、扭轉(zhuǎn)試驗，以及不同缺口半徑的拉伸試驗。 壓縮和扭轉(zhuǎn)試驗均采用圓柱型試件，壓縮試驗試件的高度與直徑尺寸為15 mm×10 mm，高徑比3：2；扭轉(zhuǎn)試驗試件標距段尺寸為10 mm×70 mm；帶缺口試件的平行段直徑φ=10 mm，缺口中心半徑a=3 mm，缺口半徑R 分別為2 mm，6 mm，8 mm 和18 mm，以及光滑圓棒（R→∞），各種試件如圖2 所示，引伸計標距為50 mm。

圖2 壓縮、扭轉(zhuǎn)及帶缺口試件圖Fig.2 Test specimens of pressure, torsion and the notched

2 不同應變率試驗結(jié)果及對本構(gòu)模型修正

2.1 不同應變率拉伸試驗結(jié)果與分析

由Q345B 鋼不同應變率拉伸試驗， 得到0.002/s～3 040/s 之間共七個應變率下的工程應力應變曲線，如圖3 所示。 可見Q345B 鋼具有明顯的應變率效應，隨著應變率增大，屈服強度和抗拉強度均明顯增大。 與準靜態(tài)相比，應變率達到100/s 時，屈服強度增大將近50%，強化比（σd/σ0）達到1.42；應變率達到103/s 的量級時，動態(tài)屈服強度將近準靜態(tài)時的兩倍。當應變率大于103/s 時，應變率1 680/s 和3 040/s 兩條曲線幾乎重合在一起，說明在本文試驗范圍內(nèi)，Q345B 鋼存在一個應變率飽和值，達到該值后，應變率效應不再顯著。

圖3 各應變率下的工程應力應變曲線Fig.3 Engineering stress-strain at different rates

圖4 σd/σ0 及σd/σb 與對數(shù)應變率的關(guān)系Fig.4 σd/σ0 and σd/σb as function of strain rate

表1 給出了Q345B 鋼主要力學性能的試驗結(jié)果。 其中，σd（MPa）為鋼材各應變率下的屈服強度；σ0（MPa）為準靜態(tài)（應變率0.002 s-1）時的屈服強度，此時有σd=σ0；σb（MPa）為抗拉強度；強化比為σd/σ0（準靜態(tài)時為1）；屈強比為σd/σb；A%為斷裂延伸率；εb%為頸縮前的均勻延伸率，對應最大載荷時試件的伸長率。 將數(shù)據(jù)繪于圖中，如圖4 可清楚看出，Q345B 鋼的強化比和屈強比隨著應變率的對數(shù)呈近似線性關(guān)系增大。 屈強比增大， 直觀上反映出應變率增大時Q345B 鋼塑性變形的應變硬化能力減弱。圖5 中，Q345B 鋼的均勻伸長率隨著應變率對數(shù)的增大而逐漸減小，表現(xiàn)為塑性失穩(wěn)不斷提前。材料的塑性變形和應力應變曲線的形狀是應變硬化和動態(tài)回復競爭的結(jié)果[7]，應變率提高后，塑性變形過程逐漸趨于絕熱過程，變形功產(chǎn)生的熱量來不及傳出去，引起材料內(nèi)部溫度升高，溫升軟化與應變硬化相互競爭，導致高應變率變形時應變硬化作用的弱化；高應變率下材料塑性變形分布不均勻，往往會在局部區(qū)域滑移量累積過大，導致塑性失穩(wěn)的提前。 圖5 中斷裂伸長率在高應變率區(qū)有所發(fā)散，但總體隨應變率增大呈減小趨勢，說明應變率提高會導致Q345B 鋼延性降低。

Q345B 鋼在不同應變率拉伸過程中均出現(xiàn)明顯頸縮， 然后很快被拉斷。 典型斷口形狀如圖6，存在明顯的剪切唇和纖維區(qū)，放射區(qū)不明顯，根據(jù)宏觀斷口可以看出Q345B 鋼具有良好的塑性，其斷裂模式屬于典型的韌性斷裂[8]。 宏觀斷口形貌隨應變率的變化沒有明顯差異，可初步判斷應變率對Q345B鋼的斷裂方式?jīng)]有重要影響。 對斷裂機理的深入剖析，還需做進一步的微觀研究。

圖5 均勻/斷裂伸長率與對數(shù)應變率關(guān)系Fig.5 Uniform/fracture elongation as function of strain rate

圖6 典型斷口照片F(xiàn)ig.6 Typical fracture graph

2.2 動態(tài)本構(gòu)模型分析與修正

2.2.1 應變率強化模型

在水下爆炸等數(shù)值仿真中，常用于金屬動態(tài)變形的應變率強化模型有JC 模型[9]和CS 模型[10]。 JC模型用一個參數(shù)描述應變率對屈服強度的強化效應，其表達式為

式中:σd為動態(tài)屈服強度，MPa；σ0為準靜態(tài)屈服強度，MPa；C 為材料參數(shù)；為無量綱應變率為塑性應變率，為參考塑性應變率。 取參考應變率為0.002/s，有σ0=364 MPa，根據(jù)各個應變率下的屈服強度和應變率對數(shù)的關(guān)系，利用最小二乘法擬合，得到C=0.046，擬合曲線如圖7。

CS 模型不考慮參考應變率，采用兩個參數(shù)描述應變率強化效應，其關(guān)系式為

式中D 和q 為材料參數(shù)。 由各應變率下的屈服強度和應變率的關(guān)系，利用最小二乘法擬合，得到D=13 709.45，q=5.55，擬合曲線如圖8。

圖7 JC 應變率強化模型擬合結(jié)果Fig.7 JC model fitting curve

圖8 CS 應變率強化模型擬合結(jié)果Fig.8 CS model fitting curve

比較圖7 和圖8，可以很清楚地看到，JC 模型采用屈服強度和應變對數(shù)的線性關(guān)系，擬合相關(guān)度較差，并不能很好地描述Q345B 鋼的應變率強化效應；而CS 模型采用兩個參數(shù)，能很好地預測Q345B鋼的應變率強化效應。

2.2.2 模型修正

通過上節(jié)分析可知，為較好反映Q345B 鋼的應變率強化效應，應采用CS 強化模型。 而CS 模型是理想剛塑性模型，其沒有考慮應變硬化效應，對于材料動態(tài)屈服后塑性變形性能的預測能力不夠。 在實際應用中，有研究者采用線性硬化模型與之配合使用[6]，但這只是簡化了的模型，不能真實反映材料的實際力學性能。 JC 模型采用冪律形式的函數(shù)描述材料的應變硬化效應，效果通常較好。

已知準靜態(tài)塑性變形是等溫過程，而對于高應變率變形，由于變形速度較快，材料內(nèi)部因變形而產(chǎn)生的熱量來不及傳出去，導致了材料內(nèi)部溫度升高，材料高應變率變形過程中的溫升為

式中:η 為塑性功-熱轉(zhuǎn)化因子，一般取η=0.9；ρ 為材料密度，cV為材料比熱容。 由于溫升將導致Q345B鋼受到軟化，則還需考慮溫升軟化效應。 結(jié)合JC 模型的乘積形式耦合應變率強化、應變硬化以及溫度軟化效應，給出以下修正模型:

式中:T*=ΔT/Tm，Tm為鋼的熔點溫度；σeq，εeq，分別為等效塑性應力、 等效塑性應變和等效塑性應變率，為εeq的導數(shù)，σeq，εeq形式如下:

通過（5）式可將張量形式的應力應變等效為標量，以適用一維拉壓等試驗得出的本構(gòu)方程。 （4）式中，σ0，D，q 取值與CS 模型一致，將σ0，D，q 值代入，通過準靜態(tài)試驗數(shù)據(jù)擬合參數(shù)B 和n，并選應變率10/s 的曲線擬合m。最終得到（4）式中:σ0=364 MPa，B=558 MPa，n=0.515，D=13 709.45，q=5.55，m=0.45。

需要說明，拉伸試驗中試件達到抗拉強度時會發(fā)生頸縮失穩(wěn)，頸縮之后變形不再均勻，而是處于三向應力狀態(tài)，需進行Bridgman 修正[11]以得到等效真實應力，真實應變也需通過最小截面積得到。 而試驗未能獲得頸縮域的實時尺寸，本文在擬合上述參數(shù)時未考慮頸縮后的情況，只選取抗拉強度之前的數(shù)據(jù)進行擬合。注意到Q345B 鋼在抗拉強度點已達到一定變形量，因此選取該段數(shù)據(jù)擬合的參數(shù)對于描述其真實變形情況已足夠。 圖9 為將參數(shù)代入修正模型后對各應變率曲線的預測結(jié)果，可以看到預測結(jié)果與試驗曲線吻合較好。

圖9 修正模型預測效果的比較Fig.9 Modified model predicted results

2.3 應變率對斷裂應變的影響

單向拉伸試驗由于試件出現(xiàn)頸縮， 導致應變集中在較小范圍，受試件尺寸的影響，由試驗曲線得出的延伸率并不能反映材料的真實斷裂應變，不能直接用于計算?；谒苄泽w積不可壓縮理論，并假定斷口截面應變均勻[12]，可由（6）式計算試件的真實斷裂應變:

式中A0和A 分別為試件的初始截面積和斷口截面積。通過量取試件斷口及原始試件截面積，得到Q345B 鋼斷裂應變與應變率對數(shù)的關(guān)系如圖10，可見隨著應變率增大，Q345B 鋼的真實斷裂應變整體呈減小趨勢，由準靜態(tài)時的1.27 降到應變率3 040/s 時的1.02。

圖10 斷裂應變與應變率關(guān)系Fig.10 Relationship of fracture strain and strain rate

3 不同應力三軸度試驗結(jié)果及失效模型擬合

3.1 不同應力三軸度試驗結(jié)果與分析

應力三軸度是表征材料內(nèi)部受力狀態(tài)的參量，表達式為η=σm/σeq，其中σm是靜水壓力，MPa，σeq是Mises 等效應力，MPa。根據(jù)Bridgman 的研究[11]，帶缺口試件中心處的初始應力三軸度可根據(jù)（7）式計算:

式中:R 是缺口半徑；a 是缺口中心處圓棒半徑。 則可得到各缺口半徑試件中心處的初始應力三軸度，見表2，表2 中同時列出了壓縮和扭轉(zhuǎn)的應力三軸度，以及通過有限元計算得到的斷裂時刻中心點處的應力三軸度和平均應力三軸度。 其中，初始值是原始試件最小截面中心處按（7）式計算的應力三軸度，斷裂值是試件斷裂時刻最小截面中心處有限元計算的應力三軸度，平均值是按（8）式計算的試件最小截面中心處在整個變形過程中的平均應力三軸度。 有限元計算過程見第3.2 節(jié)。

表2 Q345B 鋼不同應力三軸度試驗與計算結(jié)果Tab.2 Test and FEM results of Q345B steel at different stress triaxiality

圖11 給出了各缺口試件拉伸的載荷位移曲線，從圖中可以看出，隨著缺口半徑的減?。磻θS度增大），試件的斷裂伸長率隨之下降，同時最大載荷顯著提高，表現(xiàn)出明顯的缺口強化效應。 這是因為當應力三軸度提高時，受到其他兩個方向的力的限制，體積應力和體積應變起的作用越來越大，使得材料的塑性變形比較困難，要達到相同的應變需要的應力就要提高；由于材料滑移變形受阻，只能產(chǎn)生較小的頸縮，應變進一步發(fā)展所受限制隨缺口半徑的減小而增大[13]。

壓縮試驗時，試件兩端面受到摩擦力的作用，試件最終被壓成圓鼓狀，如圖12 所示，壓縮前后的試件尺寸如圖12(a)，試件未發(fā)生斷裂。 扭轉(zhuǎn)試件被均勻地扭轉(zhuǎn)很多圈后斷裂，如圖13 所示，斷口位于箭頭所指處，與軸向垂直，較為平齊，屬于典型的延性斷裂。

圖11 各缺口及圓棒試件的載荷位移曲線Fig.11 Load-displacement of notched /smooth specimens

圖12 壓縮試件變形前后對比Fig.12 Comparison of specimens before and after pressure

圖13 扭轉(zhuǎn)斷裂試件Fig.13 Specimens after torsion fracture

3.2 應力三軸度和斷裂應變的確定

對于缺口拉伸試驗， 一般采用試件最小截面中心處的初始應力三軸度或者斷裂時的應力三軸度作為整個過程的確定量， 這兩種做法都忽略了試件變形過程中應力三軸度隨應變的累計變化。 實際上，應力三軸度隨著變形的進行是不斷變化的，應該考慮應力三軸度在時間上隨應變的累積效應，則定義平均應力三軸度，通過公式（8）確定各試件的應力三軸度值:

式中，εf是斷裂應變。 應力三軸度η 是等效塑性應變的函數(shù)，由于試驗中無法獲取這個函數(shù)關(guān)系，本文借助有限元模擬來確定該關(guān)系。 根據(jù)對稱條件，建立各缺口試件的二維軸對稱模型，對缺口處進行網(wǎng)格細化，試件一端固定，另一端施加位移載荷，采用ABAQUS/Standard 對試件拉伸過程進行模擬，當計算輸出的載荷位移曲線與試驗的曲線一致后， 提取試件最小截面關(guān)鍵點的應力三軸度η εeq（ ）和等效塑性應變。

圖14 給出了各缺口試件斷裂時刻最小截面的應力三軸度和等效塑性應變沿徑向的分布。 由于斷裂時各試件的半徑不一樣，為了便于比較，將半徑無量綱化，若r 是最小截面上某點的半徑，a*是各試件斷裂時最小截面的半徑，則γ=r/a*。 由圖14（a）可見各缺口試件斷裂時，最小截面中心點的應力三軸度最大，沿徑向逐漸減小，則選取中心點的應力三軸度進行分析；由圖14（b）可見對于缺口半徑R≥6 mm 的試件，最小截面斷裂時的應變從中心沿徑向逐漸減小，而R=2 mm 的試件相反，斷裂應變從中心沿徑向逐漸增大。觀察拉伸過程看到，R=2 mm 的試件缺口根部出現(xiàn)了嚴重應力集中現(xiàn)象，這是導致其邊緣斷裂時應變比中心大的原因，但仍然是中心點發(fā)生斷裂。 圖15 給出各缺口試件最小截面中心處的應力三軸度隨應變的變化，可見在拉伸過程中，試件的應力三軸度不斷變化，總體呈增大趨勢。 對于扭轉(zhuǎn)和壓縮試件，同樣采用有限元模擬確定應力三軸度η εeq（ ）關(guān)系，根據(jù)（8）式即可計算得到各試件對應的平均應力三軸度。

圖14 斷裂時刻最小截面應力三軸度和等效塑性應變的徑向分布Fig.14 Stress triaxiality and effective plastic strain along the minimum section ratio at fracture

圖15 最小截面中心應力三軸度隨應變的變化Fig.15 Stress triaxiality as strain varies at the section minimum center

圖16 JC 失效模型的參數(shù)擬合及外插Fig.16 JC failure model parameters fitting and extrapolation

各缺口試件和扭轉(zhuǎn)試件的斷裂應變均采用有限元計算對應試驗斷裂時刻的最大等效塑性應變。壓縮試件試驗時沒有失效，說明負應力三軸度時斷裂應變較其他大，則本文直接通過外插由剪切和缺口試件擬合的曲線得到壓縮及其他負應力三軸度對應的斷裂應變值，如圖16。表2 中給出了相應的試驗及計算結(jié)果。

3.3 J-C 失效模型參數(shù)擬合

材料失效具有較復雜的機理，延性斷裂一般有空穴成核聚集斷裂、剪切斷裂等模式，在動力學有限元計算中，主要的失效模型有最大塑性應變、最大主應力、最大主應變、Johnson-Cook 失效模型等，其中最常用的是Johnson-Cook 失效模型，表達式為

式中:D 為損傷變量，當D=1 時材料失效；εf為斷裂應變，不考慮溫度的表達式為

第一大項反映了應力三軸度對斷裂應變的影響，第二大項反映應變率的影響，兩項效應疊加。

根據(jù)平均應力三軸度與斷裂應變的關(guān)系， 結(jié)合斷裂應變與應變率的關(guān)系， 利用最小二乘法擬合（10）式。 擬合參數(shù)時未計入壓縮時的斷裂應變，得到:D1=1.734，D2=-0.281，D3=1.134，D4=-0.013，將參數(shù)代入并外插得到壓縮及其他負應力三軸度對應的斷裂應變，如圖16。

4 失效判據(jù)的網(wǎng)格依賴性分析

在有限元仿真中，是通過單元的失效來模擬結(jié)構(gòu)的破壞，當單元滿足失效判據(jù)時，將單元刪除。JC失效模型以斷裂應變作為失效判據(jù)，通過材料試驗得到的斷裂應變是材料的真實力學性能，材料斷裂時變形高度集中在較小的范圍，要真實模擬材料的斷裂行為，斷裂處的單元需劃分很細。 而在實際結(jié)構(gòu)計算中，考慮計算成本等因素，通常失效單元的尺寸比材料斷裂區(qū)域的尺寸要大很多，如簡單使用材料斷裂試驗得到的失效判據(jù)，將會帶來較大誤差，這時就涉及到如何正確使用失效判據(jù)的問題。 針對此問題，運用ABAQUS/Standard 建立試件的有限元模型，選擇不同尺寸的單元進行單向拉伸試驗模擬，輸出計算的載荷位移曲線與試驗曲線比較，確定不同單元尺寸所對應的斷裂應變。

選擇光滑圓棒的準靜態(tài)拉伸試驗進行模擬，根據(jù)軸對稱條件建立試件的二維軸對稱模型，共對10種不同單元尺寸的模型進行模擬，具體單元尺寸見表3，有限元模型如圖17。 固定試件為一端，在另一端施加位移載荷。 計算得到的斷裂應變隨單元尺寸的關(guān)系見表3， 圖18 給出了單元0.05 mm 和0.1 mm 的計算結(jié)果與試驗結(jié)果對比，可見兩種尺寸基本完全重合，說明此時的單元尺寸已能模擬出試件真實的變形情況，則可認為0.05 mm 單元的斷裂應變即為材料的真實斷裂應變。

圖17 不同單元尺寸的有限元模型Fig.17 FEM model of different size of grids

圖18 兩種尺寸單元的計算結(jié)果與試驗對比Fig.18 Comparison of FEM results of two meshes with the test

表3 不同尺寸單元對應的斷裂應變Tab.3 Fracture strain with different size of grids

圖19 為斷裂應變與單元尺寸的關(guān)系，可見開始時，隨著單元尺寸增大，斷裂應變逐漸減小；當單元尺寸大于1 mm 時，斷裂應變迅速減??；當單元尺寸增大到約15 mm 以后，斷裂應變逐漸趨于穩(wěn)定。 記0.05 mm 單元的斷裂應變?yōu)棣舊0， 令修正系數(shù)為α=εf0/εfc，εfc為其他尺寸單元的斷裂應變， 則修正系數(shù)α 與單元尺寸s（mm）存在（11）式所示關(guān)系，擬合結(jié)果見圖19。

圖19 修正系數(shù)α 隨單元尺寸的變化關(guān)系Fig.19 The parameter α as function of size of grids

綜合以上分析，可知失效判據(jù)應是應變率、應力三軸度以及單元尺寸函數(shù)， 則進行考慮單元尺寸的JC 失效判據(jù)修正，得到

代入?yún)?shù)得到最終表達式為

5 結(jié) 論

針對Q345B 鋼進行了準靜態(tài)、中低應變率和高應變率以及不同應力三軸度的試驗研究和數(shù)值計算，修正得到適用于Q345B 鋼的本構(gòu)模型和斷裂判據(jù)，并形成如下具有工程應用價值的結(jié)論:

（1） 應變率對Q345B 鋼的屈服強度和抗拉強度影響顯著，當應變率增大到3 040/s 時，屈服強度增大為準靜態(tài)的1.74 倍，抗拉強度為準靜態(tài)的1.24 倍；

（2） Q345B 鋼的斷裂呈典型的韌性斷裂特征，斷裂伸長率隨應變率增大總體呈下降趨勢，當應變率增大到3 040/s 時，斷裂伸長率由準靜態(tài)的49.6%下降到41.8%；

（3） 應力三軸度對Q345B 鋼的斷裂應變影響顯著，當應力三軸度由-0.217 增大到1.085 時，斷裂應變由1.517 下降到0.791；

（4） 有限元計算中，材料的斷裂應變受單元尺寸的影響，隨著單元尺寸s（mm）的增大，單元斷裂應變逐漸減小，對JC 失效模型進行考慮單元尺寸的修正后得到（13）式。