激光輻照下衛星筒體部分多物理建模及毀傷效應分析

孫銘遠，張昊春，曲博巖，金 亮

（哈爾濱工業大學 能源科學與工程學院，黑龍江 哈爾濱 150001）

引言

激光武器是一種利用定向發射的激光束直接毀傷目標或使之失效的定向能武器，具有準確度高、延遲低、不受電磁干擾等優點[1]。近年來激光武器作為發展最為迅速且已進入實戰校驗階段的新概念武器被人們所重視，因此，針對高能激光毀傷效應的研究具有十分重要的意義[2]。

Cai等[3]建立二維軸對稱氣體動力學模型，模擬了包括逆韌致輻射、熱輻射、熱傳導和對流過程在內的激光能量傳遞過程。Du等[4]基于有限元法，通過仿真研究了激光功率和掃描速度對粉末床熱行為的影響，建立了SLM過程中AlSi10mg熔池溫度場的預測模型，并通過實驗驗證了模型的正確性。Nguyen等[5]使用ANSYS軟件，基于有限元分析的方式，模擬了激光輻射下0.5 mm×22 mm密封邊緣表面0.2 mm厚膠層的溫度場和熱流密度，研究介紹了最佳溫度、激光束產生的熱流密度以及膠粘劑的保持性能。張昊春等[6-7]建立了飛機紅外成像仿真模型，求解其蒙皮溫度與紅外輻射亮度；同時采用建模仿真的方法，基于尾流區海平面溫度分布特征，結合三維坐標變換和投影映射方法實現熱尾流目標的紅外成像仿真過程，得到在相同探測條件下熱尾流的輻射能量隨波段的分布。Luo等[8]采用有限元法分析了峰值功率密度為50 kw/cm2時，激光光斑速度分別為0 m/s、4 mm/s和10 mm/s時探測器的熱效應，基于熱傳導理論，建立了10.6 mm連續波激光輻照碲化鎘探測器的理論模型。李欣榮等[9]使用2種激光器對鋁、鋼、銅和鈦靶材進行燒蝕試驗，得到不同材料在不同能量密度下的燒蝕率。郭芳等[10]采用電子束熱蒸發方法制備了SiO2薄膜，并在相同實驗條件下制備的薄膜加以不同能量的激光輻照，研究在激光輻照前后樣片的透射率、折射率、消光系數、膜厚、表面形貌及激光損傷閾值（LIDT）的變化。王剛等[11]利用熱彈性理論對CO2激光器輻照K9玻璃材料進行研究，建立激光輻照材料溫升及熱應力分布二維平面模型，通過解析計算得到由激光輻照半導體材料引起的溫度場和應力場的瞬態分布。宋乃秋[12]等建立了典型工況下瞬態一維溫度場模型，通過比較分析解來驗證數值解的精確性，在此基礎上，考慮了相變引起的熱導率變化和熔化潛熱對溫度場分布的影響，求解了相應工況下的溫度場分布。此外，Yilbas等[13]針對某些復雜結構件的斜向激光沖擊情況，建立了不同沖擊角和不同激光偏振狀態下的沖擊波壓力模型。Qiao等[14]采用有限元方法，分析了LSP后表面和深度方向的殘余應力分布，建立了不同沖擊角和不同激光偏振狀態下的沖擊波壓力模型。Bibhu等[15]通過測量熱應變對模型進行了分析驗證，研究了掃描速度和激光功率對熱殘余應力發展的影響。Chen等[16]利用格林函數法和拉普拉斯變換，得到了兩種加熱方式的溫度場和應力場解。Boustie等[17]計算了強激光輻照下靶材產生的反沖沖量以及材料在高速沖擊波沖擊下產生的破壞效應，進一步確定材料的破壞閾值。

高能激光的實際毀傷過程十分復雜，影響其毀傷效應的因素包括天氣狀況、跟瞄誤差、大氣湍流、等離子效應等等[18-19]。考慮各種因素，為獲得不同條件和工況下的多物理場，需要大量的實驗測量數據，工作量比較大，成本高，且很難在某一外部條件下進行重復實驗。通過數值模擬的方法對衛星筒體部分的毀傷效應進行分析，可以極大地降低研究成本，對實際實驗過程有一定的指導作用。本文基于Fourier的非穩態導熱方程，對典型工況激光輻照下衛星筒體的毀傷過程進行了數值模擬，獲得了不同功率下靶材的溫度場和應力場分布，研究了靶材的熔化曲線及應力耦合效應。

1 激光輻照下衛星筒體的數學模型

高功率激光在發射后經過大氣耗散、湍流效應、等離子效應等多種作用的衰減后，輻射到人造衛星上，同時由于靶材內部的熱傳導及自身高溫特性的影響，形成了不均勻的溫度場和應力場。靶材的大致結構和毀傷的具體情況可參考圖1所示。

圖1 激光輻射下筒體靶材毀傷過程示意圖Fig.1 Schematic diagram of damage process of cylinder target under laser irradiation

當激光輻射到靶材表面時，一部分激光被反射，另一部分激光則被靶材吸收，隨著激光輻射過程的進行，目標表面的溫度不斷升高，導致目標熔化，甚至汽化，最終達到毀傷靶材的作用。激光毀傷人造衛星的過程涉及了激光與靶材之間復雜的物理和化學反應，為了便于建模分析，做如下假設：

1）假設在激光輻射過程中，靶材的熱物性及光學特性不隨溫度升高而改變；

2）忽略靶材與周圍環境的熱交換作用；

3）假設激光的功率密度在其輻射區域分布均勻：

4）假設激光光斑直徑d0遠大于板厚L，遠小于熱擴散長度。

激光在輻射靶材的過程中，靶材表面不斷吸收熱量，其吸收的熱量又會以熱傳導的形式傳遞到靶材內部，進而形成溫度場。為了獲得靶材溫度場的數學模型表達式，必須根據能量守恒方程和傅里葉定律來建立導熱微分方程：

式中：ρ表示密度，單位kg/m3；c表示比熱容，單位kJ/（kg·K）；λ表示導熱系數，單位kW/（m·K）；Q表示其他源項。人造衛星在遙遠的太空中，其初始溫度隨其位置不斷變化，但激光毀傷所需的時間非常短，因此可以認為靶材在激光輻射過程中其溫度和熱邊界條件是不變的。其初始條件和邊界條件為

為了對靶材各位置的溫度進行求解，需要對網格進行離散化劃分。現將網格劃分為均勻的n個子區域，其中各子區域長度為Δx，各節點之間距離為δx，這里取Δx=δx=h（空間步長）。那么其內部節點、左右邊界與相鄰節點之間的關系可表示為

式中：P為所研究的控制節點；W為左節點；E為右節點；r為網格比，反應了網格的質量，Δt和Δx分別為每次計算的時間和空間步長；ρ和c分別是靶材的密度和比熱容；A為靶材表面的吸收系數；qinc表示輻照目標表面的激光功率密度。

（3）式為推導的邊界迭代方程，其截斷誤差均為0（Δx），對于內部節點的迭代方程，其截斷誤差均為0（Δx2），這對提高運算精度是十分不利的。為保證整個網格的截斷誤差一致，在左右邊界各設一個虛點，即在左邊界處設i=0的節點，在右邊界處設i=n+1的節點。這時左右邊界迭代方程的截斷誤差為0（Δx2），與中間節點相同。對于整個網格來說，其中各節點截斷誤差都為0（Δx2），保證了運算的精度和穩定性。

那么內部節點、左右邊界與相鄰節點之間的關系可改進為

根據（4）式中各節點的溫度關系，在給定具體邊界初始條件后，便可以獲得最基本的溫度場分布。但激光在實際毀傷過程中，靶材溫度往往會在極短的時間內達到其熔點，甚至汽化溫度。在前面溫度場分析的過程中，并沒有考慮熔化和汽化潛熱，以及熔化前后導熱系數變化的影響，這些復雜情況將會在使用編程時分別考慮計算，具體流程圖如圖2所示。

圖2 靶材的溫度場計算流程圖Fig.2 Flow chart of target temperature field calculation

薄殼圓柱體作為人造衛星主體的一種常見結構類型，在激光輻射的過程中，圓柱表面會對靶材的吸收功率產生一定衰減作用，建立的激光輻射圓柱模型如圖3和圖4所示。

圖3和圖4分別為兩組坐標系，其中圖3是以柱心為原點的中心坐標系，θ為激光光軸與圓筒軸之間的夾角。圖4是以柱面上激光輻射中心為原點的柱面坐標系，X2、Y2（環向距離）坐標軸在圓筒外表面，設β為XY平面離開X2軸角度，則β=D/R。

圖3 中心坐標系Fig.3 Central coordinate system

圖4 柱面坐標系Fig.4 Cylindrical coordinate system

于是圓筒殼表面的入射激光功率可表示為

式中：θ為激光光軸與圓筒軸之間的夾角（rad）；β為XY平面離開X2軸角度（rad）；Q0為入射的激光功率（W/m2）；Qy為筒體結構靶材表面實際吸收的激光功率（W/m2）。

靶材經過激光輻射后形成溫度場，又由于溫度場溫度分布不均勻形成應力場。假設圓筒殼體兩端自由且絕熱，筒內無內熱源，對于圓筒殼結構，在無外部載荷時可以不考慮切向應力。利用邊界條件、溫度場求得積分常數，那么由于溫差引起的熱應力在壁厚方向的分布為[20]

總應力為

式中：σr為徑向應力（MPa）；σθ為環向應力（MPa）；σz為軸向應力（MPa）；α為熱膨脹系數（℃-1）；E為彈性模量（GPa）；a為圓筒內徑（m）；b為圓筒外徑（m）；r為圓筒所求位置處半徑（m）；D1，D2為積分常數；υ為泊松比；ΔT為內外壁溫差（℃）。

2 典型工況下靶材的多物理場分布及毀傷效應分析

2.1 靶材溫度場的建模與計算

Al及其合金廣泛應用于人造衛星的制造中。現取Al作為靶材來進行溫度場的計算，其熱物性參數如表1[21]所示。由于金屬熔化前后導熱系數發生了變化，所以網比r也發生了變化，但變化前后都應滿足0≤r≤0.5。考慮到r可能取最大值，設空間步長Δx=0.1 cm，Δt=0.005 s，數值法計算工況分別取材料厚度d=2 cm，激光輻射時間tf=1.5 s，迭代次數n=300次。金屬材料的吸收系數與激光的波長等因素有關，對于表面光滑保存良好的鋁材，若激光波長為1.315 μm時，靶材表面的吸收系數在0.1～0.2之間，鋁的表面吸收系數A=0.2，輻射到靶材的激光功率密度為10 kW/cm2，則可計算得到一維溫度場隨時間和空間變化的數值解，人造衛星外殼的溫度場分布如圖5所示。

表1 鋁的熱物性參數Table 1 Thermophysical property parameters of aluminum

從圖5可看出，在初始階段靶材的溫度場各位置溫度平穩上升，在達到熔化溫度時由于相變潛熱的影響，溫度將在一定的時間內維持在600 ℃左右，之后由于熔化前后靶材導熱系數發生了變化，靶材溫度的上升速率有了明顯的提升。同時還可以看出，在毀傷時間內，靶材底部溫度并未超出其熔化溫度，因此靶材并未因溫度過高而發生損毀。

圖5 靶材的溫度場分布Fig.5 Temperature field distribution of target

為了研究激光功率密度對靶材溫度場的影響，在保證其他條件一致的情況下，分別取激光的功率密度為8 kW/cm2、10 kW/cm2、12 kW/cm2、15 kW/cm2，討論靶材表面溫度在不同激光功率密度下溫度隨時間變化的差異，具體情況如圖6所示。

圖6 不同功率密度下靶材表面溫度分布Fig.6 Surface temperature distribution of target at different power densities

隨著激光功率密度的升高，金屬達到相變點溫度后，激光在相變點處停留的時間不斷減小，同時隨著激光功率密度的升高，曲線的斜率（溫升速率）也不斷升高，靶材達到相變點所需時間不斷減少。這是由于相同金屬從一定溫度到完全融化過程中，所需熱量是一定的，激光功率密度不斷升高，靶材表面單位時間吸收的熱量也隨之升高，故金屬毀傷過程所需的時間更短。

為了更直觀地反映靶材毀傷速率與激光功率密度的關系，本文取激光輻照下達到或超過靶材熔點的材料最大深度作為靶材毀傷深度。在所得溫度場的基礎上獲得靶材的毀傷深度曲線，如圖7所示。從圖7可以看出，相較于靶材的溫度分布，該曲線對激光功率的敏感性并不明顯，因此，這里取激光功率密度分別為10 kW/cm2、20 kW/cm2、40 kW/cm2、80 kW/cm2進行對比分析。

圖7 不同功率密度激光的擊穿曲線Fig.7 Breakdown curves of lasers with different powerdensities

由圖7可看出，在相同時間內，隨著激光功率密度的不斷增大，靶材的毀傷深度也隨之增加。曲線的斜率為靶材的毀傷速率，可以看出高功率的激光在初始時刻毀傷速率高，隨著毀傷效應的不斷進行，各功率的毀傷速率趨近相同。當激光功率密度一定時，隨著時間的增加靶材毀傷速率在某一時間內會有所降低，隨后基本保持不變，約為10 cm/s。這是由于激光在毀傷過程中表層金屬發生汽化現象所影響。

2.2 靶材應力場的建模與計算

為了獲得靶材的應力場分布，在（6）式的基礎上，需要知道靶材的力學性能參數和靶材的內外壁溫差，這里仍以鋁作為基體材料，其參數可通過查表獲得，如表2所示。靶材的內外壁溫差可以根據前面獲得的溫度場數據得出，如表3所示。

表2 鋁的力學性能參數[21]Table 2 Mechanical property parameters of aluminum

表3給出了激光功率密度為20 kW/cm2時靶材的內外表面溫差。其分布規律與外壁面溫度相似，這是由于內壁面相比于外壁面溫度變化并不明顯，從而導致了內外表面溫差隨時間不斷變大。根據（6）式，表2和表3的壁面溫差對圓筒壁面的徑向應力、環向應力和軸向應力進行計算，以鋁的熱物性參數和力學性能參數作為基本參數，整理可得3個方向分應力及總應力隨時間和空間的變化趨勢，如圖8～圖11所示。

表3 靶材的部分內外壁溫差數據Table 3 Partial inner and outer surface temperature difference data of targets

圖8 靶材徑向應力分布Fig.8 Target radial stress distribution

由圖8可知，在相同時間內，靶材沿厚度方向徑向應力逐漸減小，呈近似線性關系，在靶材底層處徑向應力達到最小值。同時，徑向應力的大小明顯受溫差曲線的影響，整體應力在1 000 MPa～2 000 MPa之間，應力在一段時間內保持不變，有明顯的停頓。

圖9為靶材的環向應力沿壁厚方向隨時間的變化情況。其分布規律與徑向應力近似相同，靶材沿厚度方向徑向應力逐漸減小，在靶材底層處達到最小值，呈近似線性關系。

圖9 靶材環向應力分布Fig.9 Target hoop stress distribution

圖10為靶材的軸向應力沿壁厚方向隨時間的變化情況。與徑向和環向應力不同，靶材的軸向應力沿壁厚方向逐漸增大，靶材表面軸向應力近似為零，近視成錐型，而應力在靶材底層處達到最大值，且靶材的軸向應力大小明顯要比徑向和環向應力小得多，僅為幾百MPa。

圖10 靶材軸向應力分布Fig.10 Target axial stress distribution

對3個應力進行矢量求和，可以得到靶材沿壁厚方向隨時間的變化情況，如圖11所示。靶材的總應力分布，明顯受徑向與環向應力分布的影響，其分布規律也近似相同，在靶材內部溫度達到熔點處，其應力大小也出現短暫停頓現象，隨后應力迅速升高，直至超出靶材許用應力大小，進而發生損毀。

圖11 靶材總應力分布Fig.11 Target total stress distribution

靶材總應力相同條件下，在靶材表面處最大，也是最容易超出許用應力的位置。同時，靶材由于應力耦合效應，其內部許用應力也隨時間不斷下降，如表4所示。為了更加直觀地從應力角度對激光的毀傷效應進行分析，取靶材表面處總應力與許用應力進行分析，靶材表面總應力和許用應力隨時間的變化曲線如圖12所示。

表4 Al的許用應力[20]Table 4 Allowable stress of Al

從圖12可以看出，靶材在強激光的輻照下表面溫度不斷升高，受到的總應力也在不斷變大，但其許用應力卻在不斷降低。在0.003 s左右靶材表面的總應力已經大于許用應力，此時靶材的安全性已經無法保證，隨時可能發生力學破壞。

圖12 靶材表面總應力和許用應力分布Fig.12 Target surface total stress and allowable stress distribution

3 結論

本文建立了激光輻照下衛星圓型筒體的計算模型，實現了靶材的溫度場、應力場建模仿真，研究了激光不同功率密度對毀傷效果的影響，得出如下結論：

1）建立了激光輻射下各節點之間的迭代方程，結合初始條件和邊界條件，建立了衛星的溫度場模型。在此基礎上考慮了靶材的熔化和汽化現象，對溫度場進一步進行完善。最后基于所得溫度場，獲得了靶材的熔化速率及位移曲線，可以從溫度方面對激光毀傷效應進行評估。隨著激光功率密度的不斷增加，初始時刻靶材的毀傷速率也在不斷增加，然而隨著毀傷過程的不斷進行，到達相變處其熔化速率趨于穩定，約為10 cm/s。

2）建立了激光輻射下衛星的應力場模型。在應力場建模計算中結合衛星結構，從圓筒結構出發進行應力場的建模計算，并考慮應力耦合效應，分析靶材的毀傷效果，從應力方面對激光的毀傷效應進行評估。在典型工況下，靶材內部應力在0.003 s時便超出了許用應力，其安全性能無法保證。