Infrared long nanosecond laser pulse ablation of silicon: Integrated two-dimensional modeling and time-resolved experimental study

Applied Surface Science 258 (2012) 7766–7773

Contents lists available at SciVerse ScienceDirect

Applied Surface Science journal homepage: www.elsevier.com/locate/apsusc

Infrared long nanosecond laser pulse ablation of silicon: Integrated two-dimensional modeling and time-resolved experimental study Sha Tao, Yun Zhou, Benxin Wu ∗ , Yibo Gao Department of Mechanical, Materials and Aerospace Engineering, Illinois Institute of Technology, 10 W. 32nd Street, Engineering 1 Building, Chicago, IL 60616, USA

a r t i c l e

i n f o

Article history: Received 6 June 2011 Received in revised form 16 April 2012 Accepted 22 April 2012 Available online 30 April 2012 Keywords: Laser ablation Laser micromachining

a b s t r a c t Nanosecond (ns) laser ablation may provide good solutions to many industrial silicon micromachining applications. However, most of the previous work is on lasers in ultraviolet (UV) or visible spectral ranges, and ns laser ablation of silicon at the infrared (IR) wavelength has not been well understood, particularly for long ns pulses with durations on the order of ∼100 ns. IR ns lasers often have lower costs and less external energy consumption for the same laser energy output than UV or visible lasers, which is desirable for many practical applications. This paper aims to understand the mechanism of IR ns laser ablation of silicon, by combining time-resolved experimental observations with physics-based modeling study. The observation is through a ns-gated intensified charged-coupled devices (ICCD) camera coupled with a microscope tube, while the model is based on two-dimensional (2D) gas dynamic equations for the gaseous phases coupled with the condensed phase heat transfer equation through the Knudsen layer relations. The research shows that the material removal mechanism under the studied laser ablation conditions is surface vaporization in the early stage (yielding a plasma plume above the target), followed by subsequent liquid ejection. The measured plasma front propagation matches reasonably well with the model prediction. The experimentally observed spatial distribution of the plasma radiation intensity is consistent with and has been understood through the model. The study also shows that the observed liquid ejection is induced by the total surface pressure difference between the near-boundary region of the target melt pool and the other remaining region of the pool. The pressure difference is mainly due to the surface vaporization flux drop after laser pulse ends. © 2012 Published by Elsevier B.V.

1. Introduction As a very important semiconductor material, silicon has been extensively used by many applications in photovoltaics, electronics and other industries, where micro-features on silicon are often needed. This may require a micromachining process [1,2], for which laser ablation may provide a good solution in many cases due to the associated advantages, such as high efficiency, good spatial resolution, noncontact (and hence no tool wear), and good flexibility [3–6]. Nanosecond (ns) lasers have a good combination of relatively low cost (as compared to ultrafast lasers) and small heat affected zone, and hence they are often applied. Experimental and theoretical research work has been performed on laser ablation of silicon and reported in literatures (e.g. [3–6]). Laser micromachin-

∗ Corresponding author at: Department of Mechanical, Materials and Aerospace Engineering, Illinois Institute of Technology, 10 W. 32nd Street, Engineering 1 Building, Room 207A, Chicago, IL 60616, USA. Tel.: +1 312 567 3451; fax: +1 312 567 7230. E-mail address: [email protected] (B. Wu). 0169-4332/$ – see front matter © 2012 Published by Elsevier B.V. http://dx.doi.org/10.1016/j.apsusc.2012.04.141

ing (ablation) of other semiconductors (e.g. GaAs [7]) has also been studied. However, most of the previous investigations of ns laser ablation of silicon (such as those in Refs. [5,6]) were mainly carried out in ultraviolet (UV) or visible spectral ranges while those using infrared (IR) lasers are insufficient (particularly for long ns pulses with durations on the order of ∼100 ns), due to the relatively low optical absorption coefficient of silicon at IR wavelength [8]. However, IR ns lasers are often cheaper and consume less external energy than the ns lasers in UV or visible spectral ranges for the same laser energy output. Therefore, IR ns laser ablation of silicon has potential advantages of low cost and energy consumption in practical applications, and the fundamental ablation mechanism deserves lots of further investigations. In the authors’ previous work in Ref. [9], experimental study has been performed for 200-ns pulsed IR (1064-nm) laser ablation of silicon. The research shows that under the studied conditions laserinduced material removal mechanism is surface vaporization in the early stage (which leads to a vapor plasma plume above the target), followed by liquid ejection that takes place at a delay time of a few hundred ns.

S. Tao et al. / Applied Surface Science 258 (2012) 7766–7773

However, despite the authors’ previous work, a fundamental understanding has not been realized for the laser-induced vapor plasma dynamic evolution, and for the underlying mechanism leading to the observed liquid ejection during IR ns laser pulse ablation of silicon. The study in this paper aims to gain this kind of fundamental understanding, by combining time-resolved experimental observations with physics-based modeling studies. The observation is through a ∼ns-gated intensified charged-coupled devices (ICCD) camera that is coupled with a microscope tube. The model simulates laser–silicon interaction, by solving the heat transfer equation in the condensed phase, coupled (through the Knudsen layer relations at the phase interface) with two-dimensional (2D) axisymmetric gas dynamic equations for the gaseous phases (including silicon vapor and the ambient air). This kind of study helps to reveal and clarify the relevant fundamental mechanisms of laser–silicon interactions, which may provide useful guiding information for practical micromachining applications. 2. Experimental setup On the experimental side, as in Ref. [9], the SPI G3.0 laser is applied. The laser operates at 1064 nm with a total pulse duration of 200 ns in this study (the laser pulse has an approximate topflat temporal shape). The IR laser beam is focused by a lens (focal length: 100 mm) in the scan head (ScanLab, HurryScan 14) onto the silicon target surface, where the laser spot has a diameter of about 30 ␮m. The target is a polished and undoped single-crystal silicon wafer with a thickness of ∼680 ␮m, which is placed in air on three-dimensional (3D) linear motion stages (Newport ILS 100PP) that have a spatial resolution of ∼0.5 ␮m. The motion stages are controlled by an XPS-C4 controller from Newport. The process of laser ablation of silicon is observed using an ICCD (Andor iStar DH-734) camera, whose gate width goes down to 10 ns in the experiment in this paper. The ICCD detector area is 13.3 mm × 13.3 mm with 1024 × 1024 pixels. A microscope tube is coupled with the ICCD camera in order to get a high spatial resolution, which consists of a 10× objective and a 2× tube lens. The objective is an infinity-corrected long-working-distance objective from Edmund Optics, which has a focal length of 20 mm and a working distance of 33.5 mm, while the tube lens has a focal length of 400 mm. The ICCD camera and laser pulse are synchronized by an electronic delay generator with high temporal resolutions, so that ICCD images at different delay times after laser pulse starts can be taken for time-resolved observations. Schott glass filters (CG-KG3-2.00-2 from CVI Melles Griot) and neutral density (ND) filters are positioned in front of the microscope tube objective, and the Schott glass filters can block the scattered 1064-nm laser light while the ND filters can reduce the intensity of optical radiation into the ICCD camera to protect the camera from being damaged. 3. The model During the ns laser ablation of silicon, the following are the major physical processes in the early stage: (1) the absorbed laser beam energy raises the silicon target temperature and leads to surface melting, (2) vaporization starts becoming significant at the melted surface, and the velocity of vapor molecules leaving the liquid surface transforms from a non-equilibrium to an equilibrium distribution in a very thin layer (called “Knudsen layer”) immediately above the melted target surface, (3) complicated gas dynamic evolution occurs in silicon vapor and the ambient air above the target (the vapor may have high temperatures and hence in a partially ionized state, i.e. the plasma state).

7767

Fig. 1. Schematic of the model setup (KL, Knudsen layer).

All of the above important physical processes have been considered in the model (see Fig. 1 for the schematic). In the target condensed phase (solid and liquid, at z ≤ 0), the thermal transport process is simulated by solving the 2D axisymmetric heat transfer equation, while the 2D axisymmetric gas dynamics equations are solved for the gaseous phase (silicon vapor and the ambient air, at z ≥ 0). The equations in the condensed and gaseous phases are coupled through Knudsen layer (KL) relations (at z = 0) (KL relations have been studied in literatures, e.g., in Refs. [10–15]). The model will be compared with the time-resolved experimental observations to verify the model and understand the experiments. The goal is to fundamentally understand the laser-induced vapor plasma dynamic evolution, and reveal the underlying mechanism for the previously observed liquid ejection [9] from the silicon target at a few hundred ns after the laser pulse starts. As the silicon vapor plume expands, it will compress the ambient air and generate a shock wave (see Fig. 1). The model calculation shows that the thickness of the compressed air layer between the shock wave front and the silicon vapor front is relatively thin during the early stage, although the thickness does increase with time. 3.1. 2D axisymmetric heat transfer equation During laser ablation of silicon, the laser energy absorption may lead to surface melting of the silicon target, where vaporization may become significant. Because of material removal due to surface vaporization, the target surface may recede in −z direction (see Fig. 1) at a velocity of us (called vaporization velocity). In a coordinate system where the original point is fixed at the target surface, the target material will move in +z direction at us if vaporization occurs, and the thermal transport process in the silicon target is governed by the following heat transfer equation [9–11,16]: 

∂h 1 ∂ ∂h + us = r ∂r ∂t ∂z

∂I = ˛Si I ∂Z



rk

∂T ∂r



+

∂ ∂z



k

∂T ∂z



+

∂I ∂z

(1a)

(1b)

where t is time, z and r are spatial coordinates,  is density, h is enthalpy, which is a function of temperature, k is thermal conductivity, T is temperature, us is the vaporization velocity (see Section 3.3 Knudsen layer relations for the procedure of calculating us during surface vaporization), I is the laser beam intensity (which is assumed to have a Gaussian spatial profile in the r direction), and ˛Si is the optical absorption coefficient of silicon, which depends on temperature [8]. Solving Eq. (1) is an approximate approach to describe heat transfer process in the condensed phase, and the approximation should be reasonable when the vaporization velocity (us ) and the vaporization depth are relatively small as in this study. The required material properties are taken from literatures [6,8,17–21]. In particular, it should be noted that the silicon optical surface reflectivity will increase significantly upon melting to more than ∼0.7 [6,20].

7768

S. Tao et al. / Applied Surface Science 258 (2012) 7766–7773

3.2. 2D axisymmetric gas dynamic equations

˛ei =

The gas dynamic evolutions of the ambient air and the silicon vapor ejected from the target surface are governed by the 2D axisymmetric gas dynamic equations [12,22–25]: 1 ∂(rv ur ) ∂(v uz ) ∂v + + =0 r ∂t ∂r ∂z

(2a)

1 ∂(ra ur ) ∂(a uz ) ∂a + + =0 r ∂t ∂r ∂z

(2b)

∂P 1 ∂(ru2r ) ∂(ur uz ) ∂ur + + =− r ∂t ∂r ∂z ∂r

(2c) (2d)

1 ∂[r(Et + P)ur ] ∂[(Et + P)uz ] ∂Et + + = ˛I r ∂t ∂r ∂z

(2e)

where v and a are the density of silicon vapor and the ambient air, respectively,  = a + v , ur and uz are the velocity components in r and z directions, respectively, P is the pressure, ˛ is the optical absorption coefficient of the plasma, and Et is the total energy per unit volume including both the kinetic and internal energy of the gaseous phase. The silicon vapor–air inter-diffusion is not considered in Eq. (2), because diffusion is typically a slow process, while the studied time-scale in this paper is only around a few hundred ns or less. The gas dynamic equations need to be supplemented by suitable equations of state (EOS). When the gaseous material is not ionized, the ideal gas EOS is applied. If ionization occurs, the electron number density can be calculated by solving the Saha equation [26], based on which the gas pressure and internal energy can be obtained for a given temperature and density. This will provide all the EOS information needed by the gas dynamic equations. The simulation results show that under the studied laser irradiation conditions, the ionization degree is typically less than 1.0% for most of the spatial points in the plasma during most of the simulated time stages. Hence, the effect of ionization on the plume gas dynamic process can be neglected. The major possible optical absorption mechanisms of laser beam by the plasma plume are photon-ionization and inverse Bremsstrahlung, and hence, the plasma optical absorption coefficient ˛ is given by [24,26–30]:

 ∞

Nn vn

(4)

n=n0

where Nn is the number density at the energy level n, and  vn is the cross-section for the capture of a photon by an electron at the energy level n. n0 is the lower limit in the summation in Eq. (4), and is determined by the condition that the photon energy should be larger than the binding energy of the electron in the atom. It should be noted that because the 1064-nm laser photon energy is much smaller than the ionization potential of silicon [18], the photon-ionization absorption coefficient is very small. The inverse Bremsstrahlung absorption coefficient is given by [24,26–30]: ˛IB = ˛ei + ˛ea

− e2 ne = na Q ea mc2





ne ni

(5b)

8kB T m

(5c)

where ˛IB is the total electron inverse Bremsstrahlung absorption coefficient, ˛ei the electron–ion inverse Bremsstrahlung absorption coefficient, ˛ea is the electron–atom inverse Bremsstrahlung absorption coefficient, v is the laser frequency, T is the plasma temperature, h is the Planck constant, kB is the Boltzmann constant, ne is the electron number density, ni is the ion number density, na is the atom number density, m is the electron mass, c is the light speed, and Q ea is electron–atom collision cross section. Based on Eqs. (3–5) and the calculated plasma temperature, density and electron density, etc., in the simulation, it has been found that the plasma optical absorption coefficient is on the order of 1.0 m−1 or smaller under the studied conditions in this paper, while the plasma plume size is only around a few hundred ␮m or less. Therefore, when laser beam propagates through the plasma plume, its absorption by the plume is negligible. The simulations also show that for the studied cases the compressed air layer between the silicon vapor front and the shock wave front has relatively low temperatures and is not ionized (i.e. not in the plasma state). 3.3. Knudsen layer relations The velocities of vapor molecules leaving the target surface transform from a non-equilibrium to an equilibrium distribution within a thin layer, which is the so called Knudsen layer (KL) [10–15]. Usually the thickness of KL is only several mean free paths, and can be treated as a discontinuity in a thermal model, across which the material state parameters (e.g. temperature) may have a jump. The Knudsen layer (KL) conditions during surface vaporization have been studied and developed in literatures [10–15,31–33], which relate the material state parameters on the two sides of KL (the condensed phase side right below the layer and the vapor side immediately above the layer) [10–15,31]: T1 = Ts

(3)

where ˛PI is the absorption coefficient due to photon-ionization and ˛IB is the absorption coefficient due to the inverse Bremsstrahlung process. The photon-ionization absorption coefficient is given by [24,26–30]: ˛PI =



−

∂(u2z ) 1 ∂(rur uz ) ∂P ∂uz + + =− r ∂t ∂r ∂z ∂z

˛ = ˛PI + ˛IB

˛ea



3.7 × 108 h 1 − exp − √ kB T 3 T

(5a)

1 = s

1+



+

Ts T1



  − 1 s 2  +1 2

s2 +

1 2



√  −1 s −   +1 2

2

s exp(s2 )erfc(s) − √ 

√ 1 Ts [1 − s exp(s2 )erfc(s)] 2 T1

(6a)



(6b)



M, M is the Mach number where  is the specific heat ratio, s = 2 of the vapor at the edge of the KL, and erfc(s) is the complementary error function. In Eq. (6), T1 and 1 are the temperature and density on the vapor side of KL, while Ts is the temperature on the condensed phase side of KL (i.e. the condensed target surface temperature), and s is the corresponding saturation vapor density at Ts . The saturation vapor pressure Ps at the temperature Ts can be calculated through the Clausius–Clapeyron equation [31,34]: Ps (Ts ) = Pb exp

L  ev Rs Tb

1−

Tb Ts



(7)

where Pb is the saturation vapor pressure at Tb (in this work Pb = 1.01 × 105 Pa and Tb = 3514 K [6] are used), Lev is the latent heat of vaporization [6], and Rs is the specific gas constant [26]. From Eq. (7) Ps can be calculated, based on which s can be obtained.

S. Tao et al. / Applied Surface Science 258 (2012) 7766–7773

The pressure ratio P/Ps determines if surface vaporization or vapor re-condensation is dominant [10] (where P is the vapor pressure just above the KL and is determined by solving Eq. (2)). If the ratio is less than 1.0, then the surface vaporization process dominates re-condensation. If the ratio is larger than 1.0, then vapor re-condensation dominates. The pressure ratio can be related to the Mach number by the following relations [10]: M=−

M=

1 − (Ps /P) 1 − (I0 )

1 − (P/Ps ) 1 − Ab

␣1

␣1

,

(P/Ps > 1)

(8a)

b

,

(otherwise)

(8b)

where I0 and ˛1 are parameters that depend on ln(T1 /Ts ) (where T1 is the temperature on the vapor side of KL), and the dependence relations are taken from Ref. [10], and A and b are constants [10]. When vapor particles come back to the silicon surface during recondensation, some of them can be reflected [31]. In this paper, it has been found that the amount of silicon vapor re-condensed to the target surface is very small during the studied time range, and hence the reflected amount is even smaller and can be neglected during the simulation without causing any significant errors. The gas dynamic equations and the condensed phase heat transfer equation are coupled through the Knudsen layer relations, and solved numerically. The numerical time step is on the order of 0.1 ns. In each time step, the heat transfer equation is solved based on a non-uniform mesh with finer meshes near the target surface, and the gas dynamic equations are solved using the finite difference method based on the Essentially Non-Oscillatory (ENO) scheme from Ref. [35], which is fast, robust and relatively easy to apply. During surface vaporization, the vaporization velocity us is determined through the following procedure at each time step: (1) The surface temperature of silicon condensed phase, Ts , can be obtained by solving Eq. (1), and then the corresponding saturation vapor pressure at this temperature, Ps , is calculated through Eq. (7). (2) The pressure above the KL, P, is obtained by solving the gas dynamic equations (Eq. (2)). (3) Then the pressure ratio P/Ps is calculated. (4) Based on the pressure ratio and Eqs. (6) and (8), the Mach number M, and the vapor temperature T1 and density 1 above the KL can be obtained. Then the vapor velocity (in z direction) above the KL, u1 , can be determined from the relation [11,12]: M = (u1 /cs ) = (u1 /( kB T1 /m)), where cs is the local sound speed,  is the ratio of specific heats, kB is the Boltzmann constant, and m is the molecular mass. (5) Then the vaporization velocity us can be determined based on the mass conservation relation: us = 1 u1 , where  is the silicon condensed phase density. In the developed model, the interface between the silicon vapor and the ambient air is not tracked explicitly, which is extremely difficult for two-dimensional (2D) compressible flow. Instead, two continuity equations are solved: one for the silicon vapor and the other for the ambient air (see Eq. (2)). Based on the solutions, the spatial distributions of the silicon vapor density v and the air density a can be obtained, from which the air–vapor interface can be determined. Eq. (2) theoretically implies a sharp air–vapor interface. However, when it is solved numerically, there will always be some numerically induced air–vapor inter-diffusion. It has been found that the diffusion layer is very thin and hence does not cause any big problem. In the very thin air–vapor numerical diffusion layer, v and a are both non-zero, and the effective gas heat capacity ratio is obtained by calculating the weighted average of the heat

7769

capacity ratios of silicon vapor and air based on their molecular number density fractions.

4. Results and discussions 4.1. Dynamic evolution of vapor plasma induced by laser ablation of silicon Fig. 2 shows the ICCD images for laser ablation of silicon taken at different delay times after the laser pulse starts, and also the model-predicted gas phase temperature and silicon vapor density at the same delay times. The 200-ns laser pulse begins at t = 0 and completes at about t = 200 ns, and it has a fluence of 30 J/cm2 . The ICCD images show that during the early stage (at delay times of less than ∼200 to 300 ns), laser-induced material removal from the silicon target is mainly realized through surface vaporization, which yields a vapor plasma plume above the target. Because laser radiation at 1064 nm has been blocked by the Schott glass filters, the relative brightness in the ICCD images reflects the relative plasma plume radiation intensity at different spatial locations. One interesting phenomenon has been found: the radiation intensity of the plume is not uniform, and instead two high-radiation-intensity (HRI) regions exist. One is the small spot immediately above the silicon target surface, and the other is the region right behind the plume expanding front. It is expected that the HRI regions have relatively higher temperatures and/or densities. This point has been confirmed by the physics-based modeling results that are also shown in Fig. 2. Numerical modeling studies have been reported in literatures on laser-induced plume evolution (e.g. Refs. [36–39]). Although these previous studies are not about the plume induced by infrared nanosecond laser ablation of silicon, they do show that numerical models provide a powerful tool for understanding laserinduced plume evolution. In Fig. 2, the simulation results predict two regions with relatively high temperature and vapor density. One is the small spot immediately above the target surface, and the other is the region near and behind the vapor expanding front. Between these two regions, the temperature and vapor density are relatively low. The model-predicted temperature and vapor density distribution structures are very similar to the plume brightness profiles in the ICCD images. The agreement is good considering that the model has no free adjustable variables while the involved processes are complicated. To better understand the formation mechanisms for the above temperature and density distribution structures, Fig. 3 shows the model-predicted pressure contour plot together with the velocity vector plot for the gaseous phase at t = 150 ns. The model calculation suggests that due to the large vaporization flux from the target surface, the small region immediately above the target surface has relatively high density and pressure. As the vapor goes away from the target surface to the plume middle core region, the pressure drops, while due to the pressure gradient the velocity increases as shown in Fig. 3. However, due to the confinement effect of the ambient air, the vapor velocity significantly decreases when approaching the plume front, where vapor accumulation occurs and relatively high densities and pressures are developed. With the vapor compression by the high pressure and the reduction of its velocity (and hence kinetic energy), the vapor internal energy and hence temperature will increase. Therefore, the region behind the vapor plasma front also has relatively high temperatures in addition to high densities. The above also shows that even without any significant direct absorption of laser beam energy by the vapor plume, the vapor temperature can still be significantly increased due to the complicated gas dynamic process.

7770

S. Tao et al. / Applied Surface Science 258 (2012) 7766–7773

Fig. 2. Model-predicted gas phase temperature and silicon vapor density distribution compared with ICCD images from the experiment (target: silicon; laser fluence: 30 J/cm2 ; the 200-ns and 1064-nm laser pulse begins at t = 0 and completes at about t = 200 ns; for each simulation plot and ICCD image, the domain size is 665 ␮m × 665 ␮m; please see the paper’s electronic version if some ICCD images are not clear in the hardcopy).

Fig. 4 compares the model-predicted vapor plasma front locations with measurements using the ICCD camera. The good agreement between simulations and experiments, as shown in Figs. 2 and 4, has verified the developed model. On the other hand, the model has provided a fundamental understanding of the experimentally observed vapor plasma dynamic evolution. In the future work, we plan to perform time-resolved measurements of the plasma temperature through emission spectroscopy [18], which will be used to further verify and improve the developed model. 4.2. Liquid ejection from silicon target and the underlying mechanism Fig. 5 shows the ICCD images for laser ablation of silicon after the early stage (i.e. after t = ∼300 ns). The 200-ns laser pulse ends at t = 200 ns, after which the plasma plume can not obtain energy (mainly indirectly) from the laser beam, and hence the plume

radiation intensity drops significantly (and the plume may even become unobservable) as can be seen from the ICCD images. On the other hand, after a delay time of roughly around ∼300 ns, the ICCD images show that a cloud of liquid is ejected from the target surface, whose morphology evolution is totally different from the vapor plasma plume in the early stage, and the liquid cloud front propagates at a much slower speed (on the order of ∼50 to 100 m/s). This suggests that target material removal through liquid ejection has occurred. Our previous work [9] has shown that the liquid ejection is very unlikely due to phase explosion (PE), because the target near-surface temperature does not satisfy the required conditions for PE to occur [6]. It has also been found through experiments that the laser ablation rate of silicon increases continuously and smoothly from the ablation threshold to the fluence of 30 J/cm2 without any jumps in material removal. This gives another strong evidence supporting that PE does not occur under the studied conditions, because it is known that if PE occurs, the laser ablation rate will have a jump at the PE threshold [6,40].

S. Tao et al. / Applied Surface Science 258 (2012) 7766–7773

7771

Fig. 6. The model-predicted total pressure exerted on the surface of the target melt pool (laser irradiation conditions are the same as Fig. 2; the radial position of the target solid–liquid interface is at r = ∼13.5 ␮m at t = 200 ns, and r = ∼11.5 ␮m at t = 500 ns; in this figure, the curves are all truncated at r = ∼11.5 ␮m). Fig. 3. The model-predicted pressure and velocity for the gaseous phase at t = 150 ns (laser irradiation conditions are the same as Fig. 2; the arrow length and direction represent the velocity magnitude and direction, respectively; the domain size is 665 ␮m × 665 ␮m).

Fig. 4. The comparison of the model-predicted vapor plasma front locations above the target surface with those measured from experiments (laser irradiation conditions are the same as Fig. 2).

Fig. 6 shows the model-predicted total pressure exerted on the target surface (i.e. the surface of the target melt pool). The total pressure felt by the liquid at the melt pool surface includes both the surface vaporization-induced recoil pressure (which refers to v u2z in this paper) and the vapor pressure immediately above the target (P). The simulation results show that at the end of the laser pulse (t = 200 ns), the total surface pressure is on the order of ∼107 Pa at the melt pool center, and decreases very quickly along the r direction. After t = ∼200 to 300 ns, Fig. 6 shows that the total pressure on the melt pool drops very quickly, except the near-boundary region of the pool where the pressure has even slightly increased. This leads to a surface pressure difference between the melt pool nearboundary region and the other remaining region of the pool roughly after t = ∼300 ns. Fig. 6 shows that the pressure difference is very large (the maximum over minimum pressure ratio is up to over ∼10 to 100). The simulations have helped to reveal the mechanism of the observed liquid ejection. It is expected that this large pressure difference will drive significant motion of the material in the target melt pool, and some of the liquid may be ejected out from the surface, which has explained the liquid ejection observed in the experiments. This explanation is further supported by the fact that the time when this pressure difference occurs (after a delay time of roughly around ∼300 ns) is consistent (at least in the semi-quantitative or qualitative sense)

Fig. 5. ICCD images of laser–silicon interaction after t ≥ 300 ns (laser fluences are given in the plot and the other laser irradiation conditions are the same as Fig. 2; the ejected liquid cloud is indicated by the white arrows; please see the paper’s electronic version if some ICCD images are not clear in the hardcopy).

7772

S. Tao et al. / Applied Surface Science 258 (2012) 7766–7773

certain pool regions as shown in Fig. 7b), and the previously vaporized material quickly moves away from the melt pool surface in z and r direction, leaving behind a low-pressure zone. However, for the melt pool near-boundary region, the vaporization flux is very low both before and after laser pulse ends, and hence the vaporization-induced recoil pressure (v u2z ) does not change significantly after laser pulse ends. On the other hand, the simulation shows that the vapor motion in r direction (away from r = 0) has increased the vapor pressure (P) immediately above the melt pool near-boundary region after t = ∼200 to 300 ns. This has explained why the total surface pressure (P + v u2z ) has even slightly increased for the melt pool near-boundary region, while the pressure has significantly dropped for the other remaining region of the pool after t = ∼200 to 300 ns. It is expected that this pressure difference will induce the observed liquid ejection from the target. 5. Conclusions A comprehensive study, combining physics-based modeling and time-resolved in situ experimental observations, has been performed on IR ns laser ablation of silicon. The study shows:

Fig. 7. Model-predicted vapor flux (v uz ) distribution along r direction at the target surface: (a) before the laser pulse ends, and (b) after the laser pulse ends (laser irradiation conditions are the same as Fig. 2).

with the time when the experimentally observed liquid ejection starts. Another evidence, which supports that the liquid ejection is due to the higher surface pressure at the melt pool near-boundary region instead of higher pressure at the pool center, is that: in the latter case, much less amount of ejected liquid would be observed in the region vertically above the pool center (that is, the region near r = 0), which would contradict the ICCD images (Fig. 5). Certainly, more studies are still needed to future test the above explanation, in particular to perform fast shadowgraphy [6] observations, in order to eliminate the possible effect of plasma radiation for a more accurate determination of the time when liquid ejection starts. Then the next question is: why the total pressure exerted on the surface of target melt pool drops significantly after t = ∼200 to 300 ns except for the pool near-boundary region? Fig. 7 shows the vapor flux (v uz ) at the target surface before and after the laser pulse ends. It can be seen that within the laser pulse duration (Fig. 7a), the vapor flux (v uz ) is very large at the melt pool center, which yields both a high vaporization-induced recoil pressure (v u2z ) and a high vapor pressure (P) immediately above the target surface. Therefore, the total surface pressure (P + v u2z ) at the melt pool center is very high. Because surface vaporization is very sensitive to the target surface temperature that decreases along the r direction, the vapor flux decreases along the r direction, and so does the total surface pressure of the melt pool. After the laser pulse ends at t = 200 ns, because the target surface temperature drops, the vaporization flux decreases significantly, and so does the total surface pressure of the melt pool (except for the pool near-boundary region that will be discussed later). The pressure may drop below the ambient air pressure, because vaporization flux drops significantly (slight condensation even occurs in

(1) Under the studied laser irradiation conditions, the material removal mechanism for IR ns laser ablation of silicon is surface vaporization in the early stage, followed by liquid ejection that occurs at a delay time of roughly around ∼300 ns. Surface vaporization leads to a plasma plume above the silicon target, which has two high-radiation-intensity (HRI) zones during the laser pulse: one is immediately above the silicon target, and the other is right behind the plume expansion front. (2) The model-predicted plasma front propagation agrees reasonably well with experiments. The model-predicted hightemperature and density regions are consistent with (and hence provide a fundamental understanding of) the two HRI zones in the plasma observed in experiments. The HRI zone immediately above the target surface is due to the large vaporization flux from the target, while the other HRI zone is due to the vapor accumulation and compression behind the vapor expansion front induced by the confinement effect of the ambient air. (3) The study shows that the observed liquid ejection is induced by the total surface pressure difference between the nearboundary and central regions of the target melt pool. The pressure difference is due to the surface vaporization flux drop after laser pulse ends, and the vapor motion in radial direction away from r = 0. The study has revealed a new general ablation mechanism for ns laser ablation: the liquid ejection due to the total surface pressure difference. It has also provided useful information for the practical applications of silicon ablation using infrared ns lasers, which typically have much lower cost and energy consumption as compared to UV or visible ns lasers. In the future we will perform fast shadowgraphy observations to more accurately determine the time when liquid ejection starts, and also extend the developed model to simulate the liquid ejection process from the silicon condensed phase. References [1] K.F. Ehmann, A synopsis of U.S. micro-manufacturing research and development activities and trends, in: 4M 2007 – Proceedings of the 3rd International Conference on Multi-Material Micro Manufacture, 3–5 October 2007, Borovets, Bulgaria, 2007. [2] K.F. Ehmann, D. Bourell, M.L. Culpepper, T.J. Hodgson, T.R. Kurfess, M. Madou, K. Rajurkar, R. DeVor, Micromanufacturing: International Research and Development, Springer, AA Dordrecht, The Netherlands, 2007.

S. Tao et al. / Applied Surface Science 258 (2012) 7766–7773 [3] P.R. Miller, R. Aggarwal, A. Doraiswamy, Y.J. Lin, Y. Lee, R.J. Narayan, Laser micromachining for biomedical applications, Journal of the Minerals Metals and Materials Society 61 (2009) 35–40. [4] S. Zhu, M.H. Hong, M.L. Koh, Y.F. Lu, Laser ablation of Si in water and ambient air, Proceedings of SPIE 4426 (2002) 39–42. [5] D.M. Karnakis, High power single-shot laser ablation of silicon with nanosecond 355 nm, Applied Surface Science 252 (2006) 7823–7825. [6] J.H. Yoo, S.H. Jeong, R. Greif, R.E. Russo, Explosive change in crater properties during high power nanosecond laser ablation of silicon, Journal of Applied Physics 88 (2000) 1638–1649. [7] X. Zhang, J. Wen, C. Sun, Thermal and carrier transport originating from photon recycling and non-radiative recombination in laser micromachining of GaAs thin films, Applied Physics A: Materials Science and Processing 76 (2002) 261–267. [8] G.E. Jellison, D.H. Lowndes, Optical absorption coefficient of silicon at 1.152 ␮ at elevated temperatures, Applied Physics Letters 41 (1982) 594–596. [9] Y. Zhou, B. Wu, S. Tao, A. Forsman, Y. Gao, Physical mechanism of silicon ablation with long nanosecond laser pulses at 1064 nm through time-resolved observation, Applied Surface Science 257 (2011) 2886–2890. [10] A.V. Gusarov, I. Smurov, Thermal model of nanosecond pulsed laser ablation: analysis of energy mass transfer, Journal of Applied Physics 97 (2005) 014307. [11] S.H. Jeong, R. Greif, R.E. Russo, Laser heating of metal targets including the influence of the background pressure, Proceedings of the ASME Heat Transfer Division 351 (1997) 63–73. [12] A.V. Gusarov, A.G. Gnedovets, I. Smurov, Gas dynamics of laser ablation: influence of ambient atmosphere, Journal of Applied Physics 88 (2000) 4352–4364. [13] C.J. Knight, Transient vaporization from a surface into vacuum, American Institute of Aeronautics and Astronautics Journal 20 (1982) 950–954. [14] C.J. Knight, Theoretical modeling of rapid surface vaporization with back pressure, American Institute of Aeronautics and Astronautics Journal 17 (1979) 519–523. [15] S.I. Anisimov, Vaporization of metal absorbing laser radiation, Soviet Physics JETP 27 (1968) 182–183. [16] F.P. Incropera, D.P. Dewitt, T. Bergman, A.S. Lavine, Fundamentals of Heat and Mass Transfer, Wiley, Hoboken, NJ, 2007. [17] A.C. Diebold, Handbook of Silicon Semiconductor Metrology, Marcel Dekker, New York, 2001. [18] H.C. Liu, X.L. Mao, J.H. Yoo, R.E. Russo, Early phase laser induced plasma diagnostics and mass removal during single-pulse laser ablation of silicon, Spectrochimica Acta Part B 54 (1999) 1607–1624. [19] D.V. Makhov, L.J. Lewis, Isotherms for the liquid–gas phase transition in silicon from NPT Monte Carlo simulations, Physical Review B 67 (2003) 153202. [20] H. Kawamura, H. Fukuyama, M. Watanabe, T. Hibiya, Normal spectral emissivity of undercooled liquid silicon, Measurement Science and Technology 16 (2005) 386–393. [21] H.M. van Driel, Kinetics of high-density plasmas generated in Si by 1.06- and 0.53-␮m picosecond laser pulses, Physical Review B 35 (1987) 8166–8176.

7773

[22] S. Gurlui, M. Agop, P. Nica, M. Ziskind, C. Focsa, Experimental and theoretical investigations of a laser-produced aluminum plasma, Physical Review E 78 (2008) 026405. [23] J.C. Tannehill, D.A. Anderson, R.H. Pletcher, Computational Fluid Mechanics and Heat Transfer, 2nd ed., Taylor and Francis, Philadelphia, PA, 1997. [24] B. Wu, Y.C. Shin, Modeling of nanosecond laser ablation with vapor plasma formation, Journal of Applied Physics 99 (2006) 084310. [25] Y. Zhou, S. Tao, B. Wu, Backward growth of plasma induced by long nanosecond laser pulse ablation, Applied Physics Letters 99 (2011) 051106. [26] Y.B. Zel’dovich, Y.P. Raizer, Physics of Shock Waves and High-temperature Hydrodynamic Phenomena, Academic Press, New York/London, 1966/1967. [27] B. Wu, Y.C. Shin, A self-closed thermal model for laser shock peening under the water confinement regime configuration and comparisons to experiments, Journal of Applied Physics 97 (2005) 113517. [28] L. Spitzer, Physics of Fully Ionized Gases, Interscience Publishers, Inc., New York, 1956. [29] R.J. Harrach, Theory for Laser-induced Breakdown over Vaporizing Target Surface, UCRL-52389, Lawrence Livermore National Laboratory, Livermore, CA, 1977. [30] M. Mitchner, C.H. Kruger, Partially Ionized Gases, Wiley, New York, 1973. [31] A. Peterlongo, A. Miotello, R. Kelly, Laser-pulse sputtering of aluminum: vaporization, boiling, superheating, and gas-dynamic effects, Physical Review E 50 (1994) 4716–4727. [32] R. Kelly, On the dual role of the Knudsen layer and unsteady, adiabatic expansion in pulse sputtering phenomena, Journal of Chemical Physics 92 (1990) 5047–5056. [33] R. Kelly, Gas dynamic of the pulsed emission of a perfect gas with applications to laser sputtering and to nozzle expansion, Physical Review A 46 (1992) 860–874. [34] Q. Lu, S.S. Mao, X. Mao, R.E. Russo, Delayed phase explosion during highpower nanosecond laser ablation of silicon, Applied Physics Letters 80 (2002) 3072–3074. [35] X.D. Liu, S. Osher, Convex ENO high order multi-dimensional schemes without field by field decomposition or staggered grids, Journal of Computational Physics 142 (1998) 304–330. [36] J.R. Ho, C.P. Grigoropoulos, J.A.C. Humphrey, Computational study of heat transfer and gas dynamics in the pulsed laser evaporation of metals, Journal of Applied Physics 78 (1995) 4696–4709. [37] G. Callies, H. Schittenhelm, P. Berger, H. Hügel, Modeling of the expansion of laser-evaporated matter in argon, helium and nitrogen and the condensation of clusters, Applied Surface Science 127–129 (1998) 134–141. [38] A.V. Bulgakov, N.M. Bulgakova, Gas-dynamic effects of the interaction between a pulsed laser-ablation plume and the ambient gas: analogy with an underexpanded jet, Journal of Physics D: Applied Physics 31 (1998) 693–703. [39] D.C. Lobao, A. Povitsky, Interaction of plume with shock waves in laser ablation, American Institute of Aeronautics and Astronautics Journal 43 (2005) 595–607. [40] C. Porneala, D.A. Willis, Observation of nanosecond laser-induced phase explosion in aluminum, Applied Physics Letters 89 (2006) 211121.