基于微元法的侵徹體彈頭摩擦升溫計算方法

邵 偉，范錦彪，耿宇飛，王 瑋

(中北大學儀器科學與動態測試教育部重點實驗室，山西 太原 030051)

0 引言

侵徹過程中彈靶間發生的高速滑動摩擦主要集中在彈頭處，其中摩擦力做功的90%以上將轉化為熱能,導致彈頭表面溫度劇變[1-3]。

該部分摩擦溫度需在瞬態條件下獲取，難度較大。文獻[4]先將彈頭錐面近似分解為無數個的等邊三角形，通過彈體侵徹過程中摩擦力和法向力，計算出每個三角形上的熱流強度，再乘以熱流持續時間，得到彈頭單位面積上的摩擦升溫值；文獻[5]先在彈體侵徹中承受介質阻抗外力公式的基礎上分析得到彈體所受摩擦力的表達式，再分析了摩擦力與溫度的轉化公式，將摩擦力表達式帶入到轉化公式中，得到彈頭處摩擦升溫值；文獻[6]對侵徹過程中彈體所受的等效摩擦應力、滑動率、摩擦功轉熱系數和摩擦功耗散功率四項求積得到摩擦升溫值。

由上述介紹可知，這些方法都是將彈體侵徹混凝土靶摩擦升溫瞬態過程簡化為穩態過程，彈頭處受到摩擦力、法向力和侵徹速度等按均值處理,其分析過程不嚴密必然導致計算精度較低。本文針對該問題，根據Forrestal侵徹阻力分析法和摩擦熱力學，提出基于微元法的侵徹體彈頭摩擦升溫計算方法。

1 Forrestal侵徹阻力分析法及摩擦 熱力學

1.1 Forrestal侵徹阻力分析法

文獻[7—9]基于對侵徹體侵徹土壤類介質的摩擦阻力分析，提出侵徹混凝土力學方程。當侵徹目標靶為混凝土時，侵徹體侵徹通道先后為開坑區和貫穿洞區。Forrestal認為開坑過程發生在4R(R為侵徹體半徑)深度內，侵徹阻力隨侵徹深度增大而線性增大，該區域的力學變化滿足胡克定律[10-11]。當深度大于4R之后，侵徹阻力與深度的關系可采用空腔膨脹理論描述，Forrestal提出整個過程的軸向阻力計算公式為：

(1)

式(1)中，C為侵徹阻力系數，h為侵徹深度，B為動阻力項系數，B=1，fc為材料的單向無側限抗壓強度，N為侵徹體頭部形狀因子，ρ為混凝土密度，v為侵徹體速度，S為侵徹過程中高溫高壓以及高應變率效應對靶體強度的影響。

文獻[12]提出當16.5 MPa≤fc≤176 MPa時，S與fc的轉化關系：

(2)

侵徹阻力系數C和頭部因子N的表達式為:

(3)

N=(8Ψ-1)/(24Ψ2)，

(4)

式(3)、式(4)中，m為侵徹體質量，v0侵徹體初速度，v1為開坑結束時的速度，Ψ為曲徑比。

當侵入深度為4R開坑階段結束，撞擊過程也隨之結束。根據連續性，由式(1)可得，此時的侵徹體速度為：

(5)

聯立式(3)和式(5)，可得侵徹阻力系數C表達式為：

(6)

1.2 摩擦熱力學理論

熱力學是為研究熱現象和力現象兩者之間的關系發展起來的宏觀唯向學理論。摩擦過程熱量場的變化特點是：具有很大的空間和時間梯度，使熱應力以熱沖擊的形式表現出來[1]，將動能轉化為熱能(假設能量無其他形式的損耗)。

摩擦生熱是相互摩擦的物體表面分子相互碰撞生熱的結果：假定一個物體(混凝土靶)靜止，另一物體(侵徹體)相對該物體運動；在此過程中，靜止物體中的分子被撞擊，獲得了運動物體中分子的部分或全部定向動能，使接觸面分子熱運動動能增大，相互碰撞頻率增多；從而導致相互摩擦的物體表面在宏觀上表現為內能增大，溫度升高。一般摩擦產生的熱量用式(7)計算。

Q=Ff·s，

(7)

式(7)中，Ff為兩物體接觸面的滑動摩擦力，s為兩物體的相對位移，計算公式為：

s=v·t，

(8)

式(8)中，v為物體運動的速度，t為物體運動的時間。

2 基于微元法的彈頭摩擦溫度計算

2.1 侵徹體頭部受力分析

由1.1節可知，本文分析彈頭完全撞擊侵入目標靶時的摩擦阻力屬于開坑階段。由式(1)可計算得到開坑階段任何深度侵徹體的阻力大小，再對阻力受力分析，得到彈頭表面所受的摩擦力大小。受力分析如圖1所示。

圖1 受力分析圖Fig.1 Stress analysis diagram

圖1受力分析代表某個時刻摩擦力的大小和方向，圖中F阻為軸向阻力，Fn為法向力，Ff為切向力，即摩擦力。分析可知摩擦力大小取決于該時刻彈頭侵入混凝土靶深度所對應的阻力大小和彈頭的切線方向的變化情況。當彈頭全部侵入混凝土靶板后θ=90°，摩擦力方向不再變化。

由圖1可知，侵徹體撞擊侵入靶板時，彈頭所受摩擦力為：

Ff=F·sinθ。

(9)

聯立式(3)、式(7)和式(9)，可得摩擦力表達式：

(10)

由式(1)和牛頓第二定律式(11)，可得開坑過程中撞擊深度、速度和時間的關系，如式(12)和式(13)所示：

(11)

(12)

(13)

式中h為侵徹體撞擊侵入靶板的深度，vt為某一時刻速度，t為時間，w為無量綱常數，表達式為：

(14)

聯立式(10)和式(12)，可得到摩擦力的表達式如式(15)所示。可知彈頭處的摩擦力是關于時間t和夾角θ的函數。

(15)

2.2 摩擦熱量計算

微元法是一種從部分到整體的思維方法，該法可使一些復雜的物理過程簡化成我們熟悉的規律加以解決。基于微元法的摩擦溫度計算思路是將彈頭表面分解為無數個小長方條，先計算出每個小長方條在整個撞擊過程中的熱量數值，設為Q0，再對整個彈頭表面的Q0積分，即可得到整個彈頭表面的總熱量Q，轉化關系式為：

(16)

式(16)中，r為彈頭侵入靶板的切面圓半徑。

現將總熱量的計算轉換為小長方條的熱量的計算(設摩擦力做功全部轉化為熱量)，式(7)用于計算分析勻速、勻摩擦力條件下的摩擦熱量，而在侵徹過程中速度和摩擦力均處在動態變化的過程中，式(7)結合微元法可得

(17)

式(17)中，Ff為滑動摩擦力，v為速度，t為時間。

由式(12)和式(13)可知，撞擊速度和深度是關于時間t和夾角θ的函數，所以對時間和夾角在對應的區間內積分即可。聯立式(12)、(13)、(14)和式(17)，積分得出Q0的計算公式。分析可知，時間的區間為(0,t)，夾角的區間為(0,π/2)，可得：

(18)

聯立式(16)和式(18)，可得彈頭完全撞擊侵入總熱量Q的表達式：

(19)

2.3 熱量與溫度的轉換

文獻[13]提出高速摩擦條件下瞬態熱量與溫度的轉換方法，對高速移動源熱量被供應到一個固定區域并存在穩態條件時，熱流可以看作流經熱阻。在本文的情況中，彈頭完全侵入混凝土靶時，可將彈靶接觸面的熱量看作供應在固定區域。通過式(19)計算出(0,t)內累積的熱量值Q，對t時刻Q可視為一個穩態值。滿足瞬態熱量與溫度的轉換條件。設彈頭在整個過程的溫度變化量為T，轉換公式為：

(20)

式(20)中，Q為侵徹過程中彈靶接觸面摩擦產生的總熱量，r為撞擊時彈頭完全侵入混凝土靶的接觸半徑，λ為侵徹體導熱系數。

聯立式(19)和式(20)，可得侵徹體侵徹混凝土靶摩擦升溫計算式：

(21)

化簡得

(24)

3 算例分析

3.1 理論計算

設侵徹體半徑和長度分別為1.3 cm和3.9 cm，彈頭形狀為半球型，其他參數如表1所示。混凝土靶參數如表2所示。

表1 侵徹體參數Tab.1 Parameters of penetration body

表2 混凝土靶參數Tab.2 Concrete target parameter

設侵徹體撞擊初速度v0=800 m/s,將表1和表2中的相關參數代入式(22)，可得摩擦溫度變化趨勢圖，如圖2所示。初始溫度為0 ℃彈頭完全侵入混凝土靶時摩擦升溫值對應圖中最高點溫度，此值為1 619 ℃。

圖2 溫度變化趨勢圖Fig.2 Temperature change trend

由圖2可知侵徹體撞擊侵入混凝土靶時,彈頭處摩擦溫度T是關于時間t和夾角θ的時空分布函數。撞擊侵入伊始，彈靶接觸面積小，升溫速率較緩。隨著彈頭不斷侵入，接觸面積變大，升溫速率隨之加快。

分別計算初速度v0為800、600、400和200 m/s下的升溫值可得表3。

表3 不同初速度下的升溫值Tab.3 Heating values at different initial speeds

3.2 仿真計算

有限元軟件ANSYS/LS-DYNA求解非線性動力學等問題功能強大。針對各種情況下的高速碰撞模擬分析、侵徹過程、爆炸過程和坑狀模擬分析等高度非線性瞬態動力學問題優勢巨大。本文使用該軟件完成仿真模擬實驗，仿真類型為熱-固耦合分析，單位制采用cm-g-μs。在建模中彈靶接觸的核心區域劃分為更加密集的網格。在靶板底部施加約束，模擬固定的靶板。由于模型結構和載荷上的對稱性，采用1/4建模的方式。時間步長由式(23)確定[5]，減少無關因素對仿真的影響，有限元模型如圖3所示。

圖3 彈靶有限元模型Fig.3 Finite element model of projectile target

(23)

式(23)中，L為劃分網格單元的特征長度，C1為材料在空氣中傳播時所具有的聲速。

本文針對高速摩擦熱，不考慮塑型變形熱且侵徹體撞擊混凝土靶速度低于800 m/s時，可視為剛體不考慮型變和質量磨損[12]，侵徹體材料定義為Rigid型，在熱固耦合分析中，該型材料不反映升溫值，需通過定義“Contact”和“Eroding_suface_to_suface”等關鍵字，使其具有溫升效應。侵徹體參數如表1所示。混凝土靶定義為HJC型,該模型考慮了應變率、平均正應力以及材料破壞與屈服應力的關系，描述侵徹混凝土具有較高的準確度[14]，混凝土靶參數如表2所示。

本文基于Forrestal半經驗侵徹阻力分析方法，提出侵徹時彈頭處的摩擦升溫計算方法。在侵徹阻力分析式(1)中未包含顯式的摩擦系數項，其原因是該模型中參量S是通過大量的侵徹實驗確定的，摩擦系數隱含其中[15]。而在彈靶有限元仿真模擬中，摩擦系數是一個極為重要的參數，因此需確定一個參考性高的摩擦系數值。

在理論阻力模型中，Forrestal推導的混凝土阻力模型未考慮摩擦系數[16]，結果導致理論值比實驗值偏大，隨后對上述方法進行了修正，提出侵徹體無論侵徹的目標是地質類靶還是金屬靶[17]，摩擦系數取0.1時，理論值與實驗值吻合的更好，基于此結論摩擦系數取0.1。

侵徹過程中存在多種干擾因素，因此對整個仿真過程需做如下假設[3]：1) 摩擦產生的熱量無損失，全部轉換為侵徹體的升溫值；2) 侵徹體材料各物性參數恒定，不隨溫度的變化而變化；3) 選擇彈頭某一網格面，代表彈頭表面溫度變化趨勢。

在ANSYS/LS-DYNA材料參數定義中輸入表1和表2 參數，并完成接觸類型和侵蝕面類型等關鍵字的定義，完成前處理。仿真過程如圖4所示，初始溫度為20 ℃，撞擊初速度為800 m/s，侵入時間5.2 μs。由溫度云圖可知此時侵徹體最高溫度為194.2 ℃。

圖4 仿真過程圖Fig.4 Temperature simulation cloud

本文所用彈靶材料參數與文獻[18]中的參數相同。相同速度下得到的仿真溫度如表4所示。

由表4可知，本文仿真結果與文獻[18]中的結果相差較小，二者存在差距主要原因是彈體幾何形狀的差別，在一定程度上驗證了本文模型的準確性。同時本文仿真得到的侵徹體頭部撞擊侵入靶板的時間與Forrestal侵徹阻力分析公式(式(12))計算得到時間幾乎一致，結果對比如表5所示。

表4 溫度對比Tab.4 Temperature comparison

表5 侵徹時間對比Tab.5 Comparison of penetration time

表5的對比結果，從另一方面也驗證了本文模型的準確性。

3.3 理論與仿真計算結果分析

由式(22)可知侵徹摩擦溫度T是關于時間t夾角θ的函數，而LS-DYNA仿真軟件中只能反映出溫度隨時間的變化關系。為了便于比較，先完成式(22)中變量θ的積分計算，將溫度只表示為時間的函數：

(24)

設初始溫度為20 ℃，侵徹體初速度分別為800、600、400和200 m/s，將表1和表2中的參數帶入式(24)中，理論計算值與LS-DYNA仿真升溫趨勢對比如圖5—圖8所示。

圖5 初速800 m/s撞擊溫度變化趨勢Fig.5 Variation trend of impact temperature with initial velocity of 800 m/s

圖6 初速600 m/s撞擊溫度變化趨勢Fig.6 Variation trend of impact temperature with initial velocity of 600 m/s

圖7 初速400 m/s撞擊溫度變化趨勢Fig.7 Variation trend of impact temperature with initial velocity of 400 m/s

圖8 初速200 m/s撞擊溫度變化趨勢Fig.8 Variation trend of impact temperature with initial velocity of 200 m/s

由圖5—圖8可得出，當侵徹體頭部完全撞擊侵入靶板時，基于微元法的溫度計算結果與LS-DYNA仿真結果吻合較好。圖5—圖8仿真與理論的最高溫度值對比如表6所示。

表6 系列初速下溫度值仿真與理論計算結果對比Tab.6 Comparison of series speed simulation and theoretical calculation results

侵徹體以800、600和200 m/s的初速度撞擊，以仿真結果為參考，相對誤差僅為4.1%、2.2%和1.8%；而以400 m/s的速度撞擊時，二者吻合性較差，相對誤差為8%。侵徹體撞擊初速度為400和200 m/s時，仿真最終的升溫值大于理論計算值。從系列圖中可以看出侵徹體撞擊侵入的初始階段，初速度越低，理論計算結果和仿真結果吻合得更好。

理論和實驗曲線圖存在差異的原因可能是理論公式計算中將彈靶接觸面考慮為光滑面，而建立的彈靶有限元模型是由若干表面不連續的微小網格組成。初始網格數量少，理論計算結果與仿真結果相差較大。隨著彈頭的不斷侵入，網格數量增多，二者結果吻合性也會隨之增高。

侵徹體撞擊初速度為800 m/s時，彈頭處摩擦升溫值超過侵徹體材料的熔點，所以可將800 m/s的撞擊初速度視為侵徹體材料熔化的臨界速度。侵徹體材料的熔化不僅與溫度的高低有關，而且與該溫度持續的時間有關，因此侵徹體材料是否會熔化還有待于進一步論證。

4 結論

本文提出基于微元法的侵徹體彈頭摩擦升溫計算方法。該方法在微元法的基礎上構建了侵徹體軸向阻力-摩擦力-熱量-溫度的函數式，通過設定的侵徹體初速度、靶體密度、導熱系數等物理參數獲得侵徹體頭部的溫升曲線。理論計算結果表明，侵徹體分別以800、600、400和200 m/s的速度侵徹目標時升溫理論計算結果與以仿真結果誤差僅為4.1%、8%、2.2%和1.8%，理論計算結果和仿真結果吻合性較好，更適合用于侵徹體彈頭摩擦升溫理論計算。

本文不足之處：仿真中以選擇的小網格表示整個彈頭面的升溫變化趨勢，存在一定得局限性。由于目前無法通過實測試驗得到侵徹體頭部摩擦升溫值，所以本文只針對靜態條件下的彈靶物理參數，計算獲得理論上的摩擦升溫值。而動態條件下的計算結果，還有待于進一步研究。

仿真結果雖與所引文獻差異較小，但在仿真中仍存在多種影響因素，如彈靶初始網格大小、時間步長等都可能會對仿真結果產生影響，可建立多組模型進一步論證。