- Email: [email protected]
Modeling of long pulse laser thermal damage of Silicon using analytical solutions
G Model IJLEO-58424; No. of Pages 8
ARTICLE IN PRESS Optik xxx (2016) xxx–xxx
Contents lists available at ScienceDirect
Optik journal homepage: www.elsevier.de/ijleo
Original research article
Modeling of long pulse laser thermal damage of Silicon using analytical solutions Guibo Chen ∗ School of Science, Changchun University of Science and Technology, Changchun, 130022, PR China
a r t i c l e
i n f o
Article history: Received 8 September 2016 Accepted 6 November 2016 Keywords: Long pulse laser Thermal damage Analytical solutions Modeling
a b s t r a c t In this paper, we present the modeling of long pulse laser thermal damage of Silicon using analytical solutions. Firstly, the two-dimensional axisymmetric physical model of the temperature rise problem for long pulse laser heating is established. By using the integral transformation method, analytical solutions of the heat conduction equation are obtained. Then, temperatures in the surface and inside of the Silicon are modeled, and the change characteristics and rules of the radial and axial distributions of temperature are studied. Finally, an analytical expression for the melting damage threshold of Silicon is given, and the threshold of different radial distances is analyzed. © 2016 Elsevier GmbH. All rights reserved.
1. Introduction The long pulse laser has been widely used in many fields such as industry, military and so on, because of its high energy density and good heating effects. The modeling of long pulse laser-matter interaction is an important method to obtain the heating mechanism of long pulse laser, it is able to reduce the experimental cost and minimizes the experimentation time [1,2]. Moreover, modeling studies can provide results of laser heating for sufficient conditions even if in the environment that traditional experiment cannot achieve. The modeling method of laser heating can be divided into numerical modeling and analytical modeling. Analytical modeling establishes a direct functional relationship between the parameters and the laser heating process, which can provide very useful information for revealing mechanism of the irradiation effects and parameters optimization of the laser [3–10]. Considerable research studies were carried out to solve the laser heating using analytical modeling. The absorption mechanism of the material on the laser energy is studied and 1-D analytical solution to the temperature rise problem of the laser heating is given by Ready [11]. El-Adawi studied on the 1-D analytical solutions of heat conduction in homogeneous material and two layered material respectively [12,13]. A closed-form solution for temperature rise inside solid substrate due to time exponentially varying laser pulse was obtained using Laplace transform by Yilbas [14]. Gospavic studied the 2-D axisymmetric analytical solution of the thermal effect problem of laser irradiation [8]. Yilbas introduced analytical simulation methods and the related results of laser heating in the application field systematically [1]. In this paper, the Silicon material is selected as the example bulk. It is organized as follows. In Section 2 we introduce the physical model and analytical solutions of heat transfer problems induced by a long pulse laser heating Silicon bulk with body absorption. Results of temperature distributions for different cases are presented in Section 3. Our main conclusions are summarized in Section 4.
∗ Correspondence to: No.7089 of WeiXing Road, ChaoYang District, ChangChun City, JiLin Province, PR China. E-mail address: [email protected] http://dx.doi.org/10.1016/j.ijleo.2016.11.009 0030-4026/© 2016 Elsevier GmbH. All rights reserved.
Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
G Model
ARTICLE IN PRESS
IJLEO-58424; No. of Pages 8
G. Chen / Optik xxx (2016) xxx–xxx
2
Laser beam
0
H
z
r
R Axis of symmetry
Fig. 1. Schematic diagram of laser irradiation material surface.
2. Physical model The classical Fourier heat transfer equation for a long pulse laser heating with a 2-D axisymmetric form can be written as [15]: 2
2
∂ T (r, z, t) 1 ∂T (r, z, t) ∂ T (r, z, t) Q (r, z, t) 1 ∂T (r, z, t) = + + + r ˛ k ∂r 2 ∂r ∂z 2 ∂t
(1)
where, k is the thermal conductivity, ˛ = k/c the thermal diffusivity, the mass density and c the heat capacity of the material. The temperature T is defined here as a function of (r, z, t), and variable ranges of the positional arguments r, z are 0 < r ≤ R,0 < z ≤ H respectively (as shown in Fig. 1). In Eq. (1), Q (r, z, t) represents the source function of laser heating. If we assume that the laser intensity is Gaussian distribution, and the energy gain mechanism of the material to the laser is the body absorption, then Q (r, z, t) can be expressed as:
Q (r, z, t) = I0 1 − rf ı exp −ız exp −r 2 /r02 g (t)
(2)
where, I0 is the laser power density,rf is the reflection coefficient, ı is absorption coefficient of the material, r0 is the waist radius of laser. g(t) is the temporal distribution function of the laser intensity, for the single millisecond pulse laser, g(t) yields:
g (t) =
1,
0 ≤ t ≤ tp
0,
t > tp
(3)
where, tp is the pulse width of the incident laser beam. It is assumed that the boundary conditions of Eq. (1) are adiabatic as: −k
∂T (r, z, t) =0 | ∂r r=R
(4)
−k
∂T (r, z, t) ∂T (r, z, t) = −k =0 | | ∂z ∂z z=0 z=H
(5)
and the initial condition yields: T (r, z, t) |t=0 = T0
(6)
In order to solve the Eqs. (1)–(6), forward and inverse transform pairs for T (r, z, t) about r are introduced based on integral transformation method [15] and written as: T (r, z, t) =
∞ J0 (n r) · T¯ (n , z, t) n=1
(7)
C1 (n )
R
r J0 n r T r , z, t dr
T¯ (n , z, t) =
(8)
0
where, J0 (n r) is the zero order Bessel function of the first kind, n (n ≥ 0, n = 1, 2, 3, · · ·) is the n-th root of equation J 1 (n R) = 0, C1 (n ) is the normalization integral as:
R
2
r J0 n r
C1 (n ) = 0
dr =
R2 /2n = 0 R2 J02 (n R) /2n = / 0
(9)
Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
G Model
ARTICLE IN PRESS
IJLEO-58424; No. of Pages 8
G. Chen / Optik xxx (2016) xxx–xxx
3
The Eq. (8) is used to transform every term in Eq. (1) and the initial condition (6), one can obtain that: 2
∂ T¯ (n , z, t) Q¯ (n , z, t) 1 ∂T¯ (n , z, t) + = ˛ k ∂z 2 ∂t
−2n T¯ (n , z, t) + and,
R
(10)
r J0 n r T 0 dr
T¯ (n , z, t) |t=0 = T¯ 0 =
(11)
0
where,
R
r J0 n r Q r , z, t dr .
Q¯ (n , z, t) =
(12)
0
Next, the forward and inverse transform pairs for T¯ (n , z, t) about z are written as: ∞ cos (m z)
T¯ (n , z, t) =
C2 (m )
m=1
H
T˜¯ (n , m , t) =
· T˜¯ (n , m , t)
(13)
cos m z T¯ n , z , t dz
(14)
0
where, m (m ≥ 0, m = 1, 2, 3, · · ·) is the m-th root of equation sin (m H) = 0, C2 (m ) is the normalization integral as:
H
C2 (m ) =
cos
2
m z
dz =
Hm = 0
(15)
/ 0 H/2m =
0
The Eq. (14) is used to transform each term of the Eqs. (10) and (11), and considering the boundary condition (5), we can obtain that:
dT˜¯ (n , m , t) ˛ + ˛ 2n + 2m T˜¯ (n , m , t) = Q˜¯ (n , m , t) dt k and,
T˜¯ (n , m , t) |t=0 = T˜¯ 0 (n , m ) =
H
R
(16)
r J0 n r cos m z T0 dr dz 0
0
.
(17)
T0 R sin (m H) J1 (n r) = m n The Eq. (16) is an ordinary differential equation about T˜¯ (n , m , t), which is easy to be solved. We substitute the solution of Eq. (16) into Eq. (13), solutions of T¯ (n , z, t) can be obtained. Furthermore, using T¯ (n , z, t) and Eq. (7), solutions of T (r, z, t) can be written as: T (r, z, t) =
˛ k
∞ ∞ J0 (n r) · cos (m z)
C1 (n ) · C2 (m )
n=1 m=1
t
exp ˛
2n
+ 2m
· exp −˛ 2n + 2m t ×
˜ t · Q¯ n , m , t dt + T˜¯ 0 (n , m )
0
(18)
Using Eq. (12), double integral Q˜¯ n , j , t in Eq. (18) can be written as:
Q˜¯ n , j , t =
H
R
r J0 n r cos m z Q r , z , t dr dz 0
0
(19)
= I0 1 − rf ıg (t ) ˝1 ˝2 where,
0
=
H
cos m z exp −ız dz
˝1 =
ı + e−ıH m sin ıH − ı cos (m H)
(20)
2m + ı2
Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
G Model
ARTICLE IN PRESS
IJLEO-58424; No. of Pages 8
G. Chen / Optik xxx (2016) xxx–xxx
4 Table 1 Laser parameters and the Silicon material parameters. Parameters
Values
laser power density/(W/m2 ) waist radius/mm laser pulse width/ms density/(kg/m3 ) specific heat/(J/kg·K) coefficient of heat conductivity/(W/m·K) reflection coefficient absorption coefficient/m melting point/K initial temperature/K thickness of material/mm radius of material/mm
109 3 1 2330 712 156 0.33 5000 1683 300 5 10
1200
N=5 N=10 N=20 N=30 FEM
1000
T/K
800 600 400 200 0.0
0.2
0.4
0.6
0.8
1.0
t/ms Fig. 2. Temperature distributions change with laser irradiating time (r = 0 mm, z = 0 mm).
R
r J0 n r I p r D dr
˝2 =
(21)
0
Then, the closed form of (18) can be written as: T (r, z, t) =
∞ ∞ J0 (n r) · cos (m z) n=1 m=1
I0 1 − rf ı˝1 ˝2
2
C1 (n ) · C2 (m )
k˛ 2n + m
e˛(
2 2 n +m
· exp −˛ 2n + 2m t ×
(22)
)t − 1 + T0 R sin (m H) J1 (n r) m n
It notes that ˝2 needs to be calculated by numerical integration, the Gauss–Kronrod quadrature method is used to handle it [16]. 3. Results and discussions In this part, we will carry out the modeling of long pulse laser heating using the above analytical solutions. The laser parameters and the Silicon material parameters used in the modeling are shown in Table 1, here it is assumed that the physical property parameters of the Silicon do not change with temperature. 3.1. Method validation We note that in Eq. (22), only the finite N and M terms in the series of n and m need to be evaluated, while the rest of the terms can be omitted. In order to validate the presented method, temperature distributions for three cases are given in Figs. 2–4 using the analytical solutions and the finite element method (FEM) of the published literature [2]. In order to facilitate the comparison and validation, the results obtained from different N and M in Eq. (22) are also given. It can be seen that the presented analytical solutions converge quickly with the gradual increase of N and M, and the calculated results of the analytical solutions and the FEM are very close when N and M are equal to 30. Therefore, we can conclude that the presented method is correct and effective. Taking into account the calculation accuracy and computational efficiency, the following examples in this paper are calculated based on N = M = 30.
Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
G Model
ARTICLE IN PRESS
IJLEO-58424; No. of Pages 8
G. Chen / Optik xxx (2016) xxx–xxx
5
1200
N=5 N=10 N=20 N=30 FEM
1000
T/K
800 600 400 200
0
2
4
6
8
10
r/mm Fig. 3. Temperature distributions change with radial distances (t = 1 ms, z = 0 mm).
1200
N=5 N=10 N=20 N=30 FEM
1000
T/K
800 600 400 200 0
1
2
3
4
5
z/mm Fig. 4. Temperature distributions change with axial depths (t = 1 ms, r = 0 mm).
1200
z=0.0mm z=0.1mm z=0.2mm z=0.5mm z=1.0mm
1000
T/K
800 600 400 200
0
2
4
6
8
10
r/mm Fig. 5. Temperature distributions change with radial locations for different axial depths (t = 1 ms).
3.2. Temperature results Next, temperature distributions in different cases are calculated using the above analytical solutions. Fig. 5 gives the temperature changes with radial locations for different axial depths. It can be seen that with the increase of the radial distance, the temperatures gradually decrease. This is due to the maximum intensity of the Gaussian beam at the radial center, and the laser intensity decreases with the increase of the radial distance. When the radial location is fixed, the temperature of the material surface (z = 0) is higher than the temperature below surface (z > 0). With the increase of the axial depth, the temperature gradually decreased. This indicates that in the thin layer of Silicon surface vicinity, temperature rises due to absorption of irradiated energy dominates over the conduction energy transport from surface to vicinity. As the depth increases, the absorption decreases, whereas the conduction enhances, but the temperature rise induced by conduction is much smaller than direct absorption of laser energy. Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
G Model
ARTICLE IN PRESS
IJLEO-58424; No. of Pages 8
G. Chen / Optik xxx (2016) xxx–xxx
6
r=0.0mm r=1.0mm r=2.0mm r=3.0mm r=4.0mm
1200 1000
T/K
800 600 400 200 0.0
0.2
0.4
0.6
0.8
1.0
t/ms Fig. 6. Temperature distributions change with time for different radial distances (z = 0 mm).
1200
z=0.0mm z=0.1mm z=0.2mm z=0.3mm z=0.4mm
1000
T/K
800 600 400 200 0.0
0.2
0.4
0.6
0.8
1.0
t/ms Fig. 7. Temperature distributions change with time for different axial depths (r = 0 mm).
1200
r=0.0mm r=1.0mm r=2.0mm r=3.0mm r=4.0mm
1000
T/K
800 600 400 200
0
1
2
3
4
5
z/mm Fig. 8. Temperature distributions change with axial locations for different radial distances (t = 1 ms).
Fig. 6 gives the surface temperatures change with laser irradiating time for different radial distances. It can be seen that the temperatures rise gradually with the increase of laser irradiating time. At the same irradiating time, the temperatures decrease with the increase of the radial distance. This is due to the reason that the intensity of the Gaussian beam becomes smaller with the increase of the radial distance. Fig. 7 gives the radial center temperature changes with laser irradiating time for different axial depths. It also can be seen that the temperatures rise gradually with the increase of laser irradiating time. And at the same irradiating time, the temperatures decrease with the increase of the axial depth. This is due to the different energy gain mechanisms for different axial depths. Fig. 8 gives the temperature changes with axial locations for different radial distances when the laser irradiating time is 1 ms. It can be seen that in the thin layer of Silicon surface vicinity, temperatures fall rapidly. And as the depth increases, temperature changes gradually become slight. This is because that the temperature rise is mainly dependent on the direct absorption of laser energy in the thin layer of surface vicinity, and the laser energy is e exponential decay in the z direction, Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
G Model
ARTICLE IN PRESS
IJLEO-58424; No. of Pages 8
G. Chen / Optik xxx (2016) xxx–xxx
r=0.0mm r=1.0mm r=2.0mm r=3.0mm r=4.0mm
8
2
10
Ep/ J/m
7
7
10
6
10
0.0
0.2
0.4
0.6
0.8
1.0
t/ms Fig. 9. Thermal damage thresholds change with time for different radial distance.
as shown in Eq. (2). When the axial depth is fixed, the temperatures decrease with the increase of the radial distance, due to the Gaussian distribution of laser intensity. Next, we analyze the thermal damage threshold of Silicon. When the surface temperature of the Silicon material reaches its melting point, the melting damage will happen. We define the laser energy density at this time as the thermal damage threshold. Using the expression (22), the thermal damage threshold can be written as:
∞ ∞ J ( r )·cos(m z) T R sin H J r ( ) ( ) n 1 m n 0 Tm − · exp −˛ 2n + 2m t · 0 m n C ( )·C ( )
Eth = tp
n=1 m=1
1
n
2
m
∞ ∞
(1−rf )ı˝1 ˝2 ˛(2 +2 )t J (n r )·cos(m z) 0 n m · exp −˛ 2n + 2m t · e −1 2 2 C
n=1 m=1
1
(n )·C2 (m )
(23)
k˛(n +m )
where tp is the pulse width of laser. Fig. 9 gives the thermal damage thresholds change with laser irradiating time for different radial distance. It can be seen that with the increase of laser irradiating time, the thermal damage threshold of Silicon is gradually reduced. When the laser irradiating time is fixed, the greater the radial distance, the higher the thermal damage threshold is. 4. Conclusions (1) Temperature profiles obtained from the analytical solutions agree well with the existing finite element method. It can provide 2-D axisymmetric modeling of long pulse laser heating materials with body absorption. (2) The temperature rises with the increase of laser irradiating time. At the same irradiating time, the temperature gradually decreases with the increase of the radial distance or the axial depth. (3) The thermal damage threshold of Silicon is gradually reduced with the increase of laser irradiating time. When the laser irradiating time is fixed, the greater the radial distance, the higher the thermal damage threshold is. Acknowledgement This work was supported by the Changchun University of Science and Technology Young Scientists Project “Numerical modeling and inverse problems of heat transfer in laser induced materials” under Grant No. XQNJJ-2014-03. References [1] B.S. Yilbas, Laser Heating Applications: Analytical Modeling, Elsevier, Amsterdam, 2012. [2] J. Bi, X.H. Zhang, X.W. Ni, G.Y. Jin, C.L. Li, Numerical simulation of thermal damage effects on Gallium Arsenide (Silicon) induced by a 0.53 m wavelength long pulsed laser, Lasers Eng. 22 (2011) 37–46. [3] J. Bi, G.Y. Jin, X.W. Ni, X.H. Zhang, Z.J. Yao, Analysis of 532 nm long pulse laserinduced thermal decomposition damage to Silicon by semi-analytical method, Acta Phys. Sin. 61 (2012) 244209. [4] J. Bi, G.B. Chen, G.Y. Jin, X.H. Zhang, Analytical modelling of annular millisecond laser-induced damage in Silicon, Lasers Eng. 29 (2014) 175–187. [5] J. Bi, G.B. Chen, Characteristic analysis of 532 nm millisecond pulse laser-induced surface damage in the shape of thermal decomposition to Silicon by means of a Semi-analytical method, Lasers Eng. 32 (2015) 99–117. [6] Y. Qin, Y.B. Chen, X.W. Ni, Z.H. Shen, J. Bi, X.H. Zhang, Axisymmetric numerical simulation of plastic damage in aluminum alloy induced by long pulsed laser, Opt. Laser Eng. 48 (2010) 361–367. [7] P.K. Srivastava, A.P. Singh, A. Kapoor, Theoretical analysis of pit formation in silicon surfaces in picosecond and femtosecond laser ablation regimes, Opt. Laser Technol. 38 (2006) 649–653. [8] R. Gospavic, M. Sreckovic, V. Popov, Modelling of laser-material interaction using semi-analytical approach, Math. Comput. Simulat. 65 (2004) 211–219.
Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
G Model IJLEO-58424; No. of Pages 8
8
ARTICLE IN PRESS G. Chen / Optik xxx (2016) xxx–xxx
[9] Z.W. Li, H.C. Zhang, Z.H. Shen, X.W. Ni, Time-resolved temperature measurement and numerical simulation of millisecond laser irradiated silicon, J. Appl. Phys. 114 (2013) 033104. [10] M. Kalyon, B.S. Yilbas, Repetitive laser pulse heating analysis: pulse parameter variation effects on closed form solution, Appl. Surf. Sci. 252 (2006) 2242–2250. [11] J.F. Ready, Effects due to absorption of laser radiation, J. Appl. Phys. 36 (1963) 462–470. [12] M.K. El-Adawi, E.F. Elshehawey, Heating a slab induced by a time-dependent laser irradiance-an exact solution, J. Appl. Phys. 60 (1986) 2250–2255. [13] M.K. El-Adawi, M.A. Abdel-Naby, S.A. Shalaby, Laser heating of a two-layer system with constant surface absorption: an exact solution, Int. J. Heat Mass Transfer 38 (1995) 947–952. [14] B.S. Yilbas, A closed form solution for temperature rise inside solid substrate due to time exponentially varying pulse, Int. J. Heat Mass Transfer 45 (2002) 1993–2000. [15] M.N. Ozisik, Heat Conduction, John Wiley and Sons, New York, 1993. [16] D.P. Laurie, Calculation of Gauss-Kronrod quadrature rules, Math. Comput. 66 (1997) 1133–1145.
Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
ARTICLE IN PRESS Optik xxx (2016) xxx–xxx
Contents lists available at ScienceDirect
Optik journal homepage: www.elsevier.de/ijleo
Original research article
Modeling of long pulse laser thermal damage of Silicon using analytical solutions Guibo Chen ∗ School of Science, Changchun University of Science and Technology, Changchun, 130022, PR China
a r t i c l e
i n f o
Article history: Received 8 September 2016 Accepted 6 November 2016 Keywords: Long pulse laser Thermal damage Analytical solutions Modeling
a b s t r a c t In this paper, we present the modeling of long pulse laser thermal damage of Silicon using analytical solutions. Firstly, the two-dimensional axisymmetric physical model of the temperature rise problem for long pulse laser heating is established. By using the integral transformation method, analytical solutions of the heat conduction equation are obtained. Then, temperatures in the surface and inside of the Silicon are modeled, and the change characteristics and rules of the radial and axial distributions of temperature are studied. Finally, an analytical expression for the melting damage threshold of Silicon is given, and the threshold of different radial distances is analyzed. © 2016 Elsevier GmbH. All rights reserved.
1. Introduction The long pulse laser has been widely used in many fields such as industry, military and so on, because of its high energy density and good heating effects. The modeling of long pulse laser-matter interaction is an important method to obtain the heating mechanism of long pulse laser, it is able to reduce the experimental cost and minimizes the experimentation time [1,2]. Moreover, modeling studies can provide results of laser heating for sufficient conditions even if in the environment that traditional experiment cannot achieve. The modeling method of laser heating can be divided into numerical modeling and analytical modeling. Analytical modeling establishes a direct functional relationship between the parameters and the laser heating process, which can provide very useful information for revealing mechanism of the irradiation effects and parameters optimization of the laser [3–10]. Considerable research studies were carried out to solve the laser heating using analytical modeling. The absorption mechanism of the material on the laser energy is studied and 1-D analytical solution to the temperature rise problem of the laser heating is given by Ready [11]. El-Adawi studied on the 1-D analytical solutions of heat conduction in homogeneous material and two layered material respectively [12,13]. A closed-form solution for temperature rise inside solid substrate due to time exponentially varying laser pulse was obtained using Laplace transform by Yilbas [14]. Gospavic studied the 2-D axisymmetric analytical solution of the thermal effect problem of laser irradiation [8]. Yilbas introduced analytical simulation methods and the related results of laser heating in the application field systematically [1]. In this paper, the Silicon material is selected as the example bulk. It is organized as follows. In Section 2 we introduce the physical model and analytical solutions of heat transfer problems induced by a long pulse laser heating Silicon bulk with body absorption. Results of temperature distributions for different cases are presented in Section 3. Our main conclusions are summarized in Section 4.
∗ Correspondence to: No.7089 of WeiXing Road, ChaoYang District, ChangChun City, JiLin Province, PR China. E-mail address: [email protected] http://dx.doi.org/10.1016/j.ijleo.2016.11.009 0030-4026/© 2016 Elsevier GmbH. All rights reserved.
Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
G Model
ARTICLE IN PRESS
IJLEO-58424; No. of Pages 8
G. Chen / Optik xxx (2016) xxx–xxx
2
Laser beam
0
H
z
r
R Axis of symmetry
Fig. 1. Schematic diagram of laser irradiation material surface.
2. Physical model The classical Fourier heat transfer equation for a long pulse laser heating with a 2-D axisymmetric form can be written as [15]: 2
2
∂ T (r, z, t) 1 ∂T (r, z, t) ∂ T (r, z, t) Q (r, z, t) 1 ∂T (r, z, t) = + + + r ˛ k ∂r 2 ∂r ∂z 2 ∂t
(1)
where, k is the thermal conductivity, ˛ = k/c the thermal diffusivity, the mass density and c the heat capacity of the material. The temperature T is defined here as a function of (r, z, t), and variable ranges of the positional arguments r, z are 0 < r ≤ R,0 < z ≤ H respectively (as shown in Fig. 1). In Eq. (1), Q (r, z, t) represents the source function of laser heating. If we assume that the laser intensity is Gaussian distribution, and the energy gain mechanism of the material to the laser is the body absorption, then Q (r, z, t) can be expressed as:
Q (r, z, t) = I0 1 − rf ı exp −ız exp −r 2 /r02 g (t)
(2)
where, I0 is the laser power density,rf is the reflection coefficient, ı is absorption coefficient of the material, r0 is the waist radius of laser. g(t) is the temporal distribution function of the laser intensity, for the single millisecond pulse laser, g(t) yields:
g (t) =
1,
0 ≤ t ≤ tp
0,
t > tp
(3)
where, tp is the pulse width of the incident laser beam. It is assumed that the boundary conditions of Eq. (1) are adiabatic as: −k
∂T (r, z, t) =0 | ∂r r=R
(4)
−k
∂T (r, z, t) ∂T (r, z, t) = −k =0 | | ∂z ∂z z=0 z=H
(5)
and the initial condition yields: T (r, z, t) |t=0 = T0
(6)
In order to solve the Eqs. (1)–(6), forward and inverse transform pairs for T (r, z, t) about r are introduced based on integral transformation method [15] and written as: T (r, z, t) =
∞ J0 (n r) · T¯ (n , z, t) n=1
(7)
C1 (n )
R
r J0 n r T r , z, t dr
T¯ (n , z, t) =
(8)
0
where, J0 (n r) is the zero order Bessel function of the first kind, n (n ≥ 0, n = 1, 2, 3, · · ·) is the n-th root of equation J 1 (n R) = 0, C1 (n ) is the normalization integral as:
R
2
r J0 n r
C1 (n ) = 0
dr =
R2 /2n = 0 R2 J02 (n R) /2n = / 0
(9)
Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
G Model
ARTICLE IN PRESS
IJLEO-58424; No. of Pages 8
G. Chen / Optik xxx (2016) xxx–xxx
3
The Eq. (8) is used to transform every term in Eq. (1) and the initial condition (6), one can obtain that: 2
∂ T¯ (n , z, t) Q¯ (n , z, t) 1 ∂T¯ (n , z, t) + = ˛ k ∂z 2 ∂t
−2n T¯ (n , z, t) + and,
R
(10)
r J0 n r T 0 dr
T¯ (n , z, t) |t=0 = T¯ 0 =
(11)
0
where,
R
r J0 n r Q r , z, t dr .
Q¯ (n , z, t) =
(12)
0
Next, the forward and inverse transform pairs for T¯ (n , z, t) about z are written as: ∞ cos (m z)
T¯ (n , z, t) =
C2 (m )
m=1
H
T˜¯ (n , m , t) =
· T˜¯ (n , m , t)
(13)
cos m z T¯ n , z , t dz
(14)
0
where, m (m ≥ 0, m = 1, 2, 3, · · ·) is the m-th root of equation sin (m H) = 0, C2 (m ) is the normalization integral as:
H
C2 (m ) =
cos
2
m z
dz =
Hm = 0
(15)
/ 0 H/2m =
0
The Eq. (14) is used to transform each term of the Eqs. (10) and (11), and considering the boundary condition (5), we can obtain that:
dT˜¯ (n , m , t) ˛ + ˛ 2n + 2m T˜¯ (n , m , t) = Q˜¯ (n , m , t) dt k and,
T˜¯ (n , m , t) |t=0 = T˜¯ 0 (n , m ) =
H
R
(16)
r J0 n r cos m z T0 dr dz 0
0
.
(17)
T0 R sin (m H) J1 (n r) = m n The Eq. (16) is an ordinary differential equation about T˜¯ (n , m , t), which is easy to be solved. We substitute the solution of Eq. (16) into Eq. (13), solutions of T¯ (n , z, t) can be obtained. Furthermore, using T¯ (n , z, t) and Eq. (7), solutions of T (r, z, t) can be written as: T (r, z, t) =
˛ k
∞ ∞ J0 (n r) · cos (m z)
C1 (n ) · C2 (m )
n=1 m=1
t
exp ˛
2n
+ 2m
· exp −˛ 2n + 2m t ×
˜ t · Q¯ n , m , t dt + T˜¯ 0 (n , m )
0
(18)
Using Eq. (12), double integral Q˜¯ n , j , t in Eq. (18) can be written as:
Q˜¯ n , j , t =
H
R
r J0 n r cos m z Q r , z , t dr dz 0
0
(19)
= I0 1 − rf ıg (t ) ˝1 ˝2 where,
0
=
H
cos m z exp −ız dz
˝1 =
ı + e−ıH m sin ıH − ı cos (m H)
(20)
2m + ı2
Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
G Model
ARTICLE IN PRESS
IJLEO-58424; No. of Pages 8
G. Chen / Optik xxx (2016) xxx–xxx
4 Table 1 Laser parameters and the Silicon material parameters. Parameters
Values
laser power density/(W/m2 ) waist radius/mm laser pulse width/ms density/(kg/m3 ) specific heat/(J/kg·K) coefficient of heat conductivity/(W/m·K) reflection coefficient absorption coefficient/m melting point/K initial temperature/K thickness of material/mm radius of material/mm
109 3 1 2330 712 156 0.33 5000 1683 300 5 10
1200
N=5 N=10 N=20 N=30 FEM
1000
T/K
800 600 400 200 0.0
0.2
0.4
0.6
0.8
1.0
t/ms Fig. 2. Temperature distributions change with laser irradiating time (r = 0 mm, z = 0 mm).
R
r J0 n r I p r D dr
˝2 =
(21)
0
Then, the closed form of (18) can be written as: T (r, z, t) =
∞ ∞ J0 (n r) · cos (m z) n=1 m=1
I0 1 − rf ı˝1 ˝2
2
C1 (n ) · C2 (m )
k˛ 2n + m
e˛(
2 2 n +m
· exp −˛ 2n + 2m t ×
(22)
)t − 1 + T0 R sin (m H) J1 (n r) m n
It notes that ˝2 needs to be calculated by numerical integration, the Gauss–Kronrod quadrature method is used to handle it [16]. 3. Results and discussions In this part, we will carry out the modeling of long pulse laser heating using the above analytical solutions. The laser parameters and the Silicon material parameters used in the modeling are shown in Table 1, here it is assumed that the physical property parameters of the Silicon do not change with temperature. 3.1. Method validation We note that in Eq. (22), only the finite N and M terms in the series of n and m need to be evaluated, while the rest of the terms can be omitted. In order to validate the presented method, temperature distributions for three cases are given in Figs. 2–4 using the analytical solutions and the finite element method (FEM) of the published literature [2]. In order to facilitate the comparison and validation, the results obtained from different N and M in Eq. (22) are also given. It can be seen that the presented analytical solutions converge quickly with the gradual increase of N and M, and the calculated results of the analytical solutions and the FEM are very close when N and M are equal to 30. Therefore, we can conclude that the presented method is correct and effective. Taking into account the calculation accuracy and computational efficiency, the following examples in this paper are calculated based on N = M = 30.
Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
G Model
ARTICLE IN PRESS
IJLEO-58424; No. of Pages 8
G. Chen / Optik xxx (2016) xxx–xxx
5
1200
N=5 N=10 N=20 N=30 FEM
1000
T/K
800 600 400 200
0
2
4
6
8
10
r/mm Fig. 3. Temperature distributions change with radial distances (t = 1 ms, z = 0 mm).
1200
N=5 N=10 N=20 N=30 FEM
1000
T/K
800 600 400 200 0
1
2
3
4
5
z/mm Fig. 4. Temperature distributions change with axial depths (t = 1 ms, r = 0 mm).
1200
z=0.0mm z=0.1mm z=0.2mm z=0.5mm z=1.0mm
1000
T/K
800 600 400 200
0
2
4
6
8
10
r/mm Fig. 5. Temperature distributions change with radial locations for different axial depths (t = 1 ms).
3.2. Temperature results Next, temperature distributions in different cases are calculated using the above analytical solutions. Fig. 5 gives the temperature changes with radial locations for different axial depths. It can be seen that with the increase of the radial distance, the temperatures gradually decrease. This is due to the maximum intensity of the Gaussian beam at the radial center, and the laser intensity decreases with the increase of the radial distance. When the radial location is fixed, the temperature of the material surface (z = 0) is higher than the temperature below surface (z > 0). With the increase of the axial depth, the temperature gradually decreased. This indicates that in the thin layer of Silicon surface vicinity, temperature rises due to absorption of irradiated energy dominates over the conduction energy transport from surface to vicinity. As the depth increases, the absorption decreases, whereas the conduction enhances, but the temperature rise induced by conduction is much smaller than direct absorption of laser energy. Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
G Model
ARTICLE IN PRESS
IJLEO-58424; No. of Pages 8
G. Chen / Optik xxx (2016) xxx–xxx
6
r=0.0mm r=1.0mm r=2.0mm r=3.0mm r=4.0mm
1200 1000
T/K
800 600 400 200 0.0
0.2
0.4
0.6
0.8
1.0
t/ms Fig. 6. Temperature distributions change with time for different radial distances (z = 0 mm).
1200
z=0.0mm z=0.1mm z=0.2mm z=0.3mm z=0.4mm
1000
T/K
800 600 400 200 0.0
0.2
0.4
0.6
0.8
1.0
t/ms Fig. 7. Temperature distributions change with time for different axial depths (r = 0 mm).
1200
r=0.0mm r=1.0mm r=2.0mm r=3.0mm r=4.0mm
1000
T/K
800 600 400 200
0
1
2
3
4
5
z/mm Fig. 8. Temperature distributions change with axial locations for different radial distances (t = 1 ms).
Fig. 6 gives the surface temperatures change with laser irradiating time for different radial distances. It can be seen that the temperatures rise gradually with the increase of laser irradiating time. At the same irradiating time, the temperatures decrease with the increase of the radial distance. This is due to the reason that the intensity of the Gaussian beam becomes smaller with the increase of the radial distance. Fig. 7 gives the radial center temperature changes with laser irradiating time for different axial depths. It also can be seen that the temperatures rise gradually with the increase of laser irradiating time. And at the same irradiating time, the temperatures decrease with the increase of the axial depth. This is due to the different energy gain mechanisms for different axial depths. Fig. 8 gives the temperature changes with axial locations for different radial distances when the laser irradiating time is 1 ms. It can be seen that in the thin layer of Silicon surface vicinity, temperatures fall rapidly. And as the depth increases, temperature changes gradually become slight. This is because that the temperature rise is mainly dependent on the direct absorption of laser energy in the thin layer of surface vicinity, and the laser energy is e exponential decay in the z direction, Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
G Model
ARTICLE IN PRESS
IJLEO-58424; No. of Pages 8
G. Chen / Optik xxx (2016) xxx–xxx
r=0.0mm r=1.0mm r=2.0mm r=3.0mm r=4.0mm
8
2
10
Ep/ J/m
7
7
10
6
10
0.0
0.2
0.4
0.6
0.8
1.0
t/ms Fig. 9. Thermal damage thresholds change with time for different radial distance.
as shown in Eq. (2). When the axial depth is fixed, the temperatures decrease with the increase of the radial distance, due to the Gaussian distribution of laser intensity. Next, we analyze the thermal damage threshold of Silicon. When the surface temperature of the Silicon material reaches its melting point, the melting damage will happen. We define the laser energy density at this time as the thermal damage threshold. Using the expression (22), the thermal damage threshold can be written as:
∞ ∞ J ( r )·cos(m z) T R sin H J r ( ) ( ) n 1 m n 0 Tm − · exp −˛ 2n + 2m t · 0 m n C ( )·C ( )
Eth = tp
n=1 m=1
1
n
2
m
∞ ∞
(1−rf )ı˝1 ˝2 ˛(2 +2 )t J (n r )·cos(m z) 0 n m · exp −˛ 2n + 2m t · e −1 2 2 C
n=1 m=1
1
(n )·C2 (m )
(23)
k˛(n +m )
where tp is the pulse width of laser. Fig. 9 gives the thermal damage thresholds change with laser irradiating time for different radial distance. It can be seen that with the increase of laser irradiating time, the thermal damage threshold of Silicon is gradually reduced. When the laser irradiating time is fixed, the greater the radial distance, the higher the thermal damage threshold is. 4. Conclusions (1) Temperature profiles obtained from the analytical solutions agree well with the existing finite element method. It can provide 2-D axisymmetric modeling of long pulse laser heating materials with body absorption. (2) The temperature rises with the increase of laser irradiating time. At the same irradiating time, the temperature gradually decreases with the increase of the radial distance or the axial depth. (3) The thermal damage threshold of Silicon is gradually reduced with the increase of laser irradiating time. When the laser irradiating time is fixed, the greater the radial distance, the higher the thermal damage threshold is. Acknowledgement This work was supported by the Changchun University of Science and Technology Young Scientists Project “Numerical modeling and inverse problems of heat transfer in laser induced materials” under Grant No. XQNJJ-2014-03. References [1] B.S. Yilbas, Laser Heating Applications: Analytical Modeling, Elsevier, Amsterdam, 2012. [2] J. Bi, X.H. Zhang, X.W. Ni, G.Y. Jin, C.L. Li, Numerical simulation of thermal damage effects on Gallium Arsenide (Silicon) induced by a 0.53 m wavelength long pulsed laser, Lasers Eng. 22 (2011) 37–46. [3] J. Bi, G.Y. Jin, X.W. Ni, X.H. Zhang, Z.J. Yao, Analysis of 532 nm long pulse laserinduced thermal decomposition damage to Silicon by semi-analytical method, Acta Phys. Sin. 61 (2012) 244209. [4] J. Bi, G.B. Chen, G.Y. Jin, X.H. Zhang, Analytical modelling of annular millisecond laser-induced damage in Silicon, Lasers Eng. 29 (2014) 175–187. [5] J. Bi, G.B. Chen, Characteristic analysis of 532 nm millisecond pulse laser-induced surface damage in the shape of thermal decomposition to Silicon by means of a Semi-analytical method, Lasers Eng. 32 (2015) 99–117. [6] Y. Qin, Y.B. Chen, X.W. Ni, Z.H. Shen, J. Bi, X.H. Zhang, Axisymmetric numerical simulation of plastic damage in aluminum alloy induced by long pulsed laser, Opt. Laser Eng. 48 (2010) 361–367. [7] P.K. Srivastava, A.P. Singh, A. Kapoor, Theoretical analysis of pit formation in silicon surfaces in picosecond and femtosecond laser ablation regimes, Opt. Laser Technol. 38 (2006) 649–653. [8] R. Gospavic, M. Sreckovic, V. Popov, Modelling of laser-material interaction using semi-analytical approach, Math. Comput. Simulat. 65 (2004) 211–219.
Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009
G Model IJLEO-58424; No. of Pages 8
8
ARTICLE IN PRESS G. Chen / Optik xxx (2016) xxx–xxx
[9] Z.W. Li, H.C. Zhang, Z.H. Shen, X.W. Ni, Time-resolved temperature measurement and numerical simulation of millisecond laser irradiated silicon, J. Appl. Phys. 114 (2013) 033104. [10] M. Kalyon, B.S. Yilbas, Repetitive laser pulse heating analysis: pulse parameter variation effects on closed form solution, Appl. Surf. Sci. 252 (2006) 2242–2250. [11] J.F. Ready, Effects due to absorption of laser radiation, J. Appl. Phys. 36 (1963) 462–470. [12] M.K. El-Adawi, E.F. Elshehawey, Heating a slab induced by a time-dependent laser irradiance-an exact solution, J. Appl. Phys. 60 (1986) 2250–2255. [13] M.K. El-Adawi, M.A. Abdel-Naby, S.A. Shalaby, Laser heating of a two-layer system with constant surface absorption: an exact solution, Int. J. Heat Mass Transfer 38 (1995) 947–952. [14] B.S. Yilbas, A closed form solution for temperature rise inside solid substrate due to time exponentially varying pulse, Int. J. Heat Mass Transfer 45 (2002) 1993–2000. [15] M.N. Ozisik, Heat Conduction, John Wiley and Sons, New York, 1993. [16] D.P. Laurie, Calculation of Gauss-Kronrod quadrature rules, Math. Comput. 66 (1997) 1133–1145.
Please cite this article in press as: G. Chen, Modeling of long pulse laser thermal damage of Silicon using analytical solutions, Optik - Int. J. Light Electron Opt. (2016), http://dx.doi.org/10.1016/j.ijleo.2016.11.009