- Email: [email protected]
Temperature distribution in multilayers covered by liquid layer and processed by focused laser beam
~ii~i~i~i~ii~iii~i~i~!~iiiii!ii~i!~i~ii!ii!i!~ii!i~iii~ii!~i!i applied surface science ELSEVIER
Applied Surface Science 106 (1996) 422-428
Temperature distribution in multilayers covered by liquid layer and processed by focused laser beam Zs. Geretovszky a,*, L. Kelemen b, K. Piglmayer c Research Group on Laser Physics of the Hungarian Academy of Sciences, Ddm t~r 9, H-6720 S=eged, Hungary b Department of Optics and Quantum Electronics, JA TE Unit'el:s'io:, D6m tgr 9, H-6720 Szeged, Hungary c Ange~,andte Physik, Johannes Kepler Uni~'ersit~it Lin=, A-4040 Lin:, Austria Received 17 September 1995; accepted 31 December 1995
Abstract
As the experimental determination of surface temperatures on micrometer scale can only be carried out in very special cases, the estimation of spatial temperature distributions requires numerical model calculations. In the present contribution experimental results obtained by scanning a focused beam of an Ar-ion laser over a tungsten covered glass slide immersed in water, are compared with the results of temperature calculations. The steady-state temperature distribution is calculated by solving numerically the three-dimensional heat diffusion equation using the second-order finite difference method. During the discretization procedure the upwind scheme is applied to ensure the stability of the equation system at scanning speeds as high as few tens of meter per second. The temperature dependence of the thermophysical properties of the materials is explicitly included. The effects of processing parameters, in particular laser power, focus diameter and scanning speed, are described and found to be consistent with the experimental results.
1. Introduction
In spite of continuous efforts to reduce the temperature during both processing and operating conditions [1-3], thermal damage of the previously formed (heat sensitive) patterns through intermixing, surface degradation or recrystallization can never be completely avoided. Thus the exact knowledge of temperature distribution is essential in process optimization. Although under strictly optimized conditions sub0.5 /zm resolution can be achieved by industrial
:: Corresponding author. Tel.: +36-62-454274; fax: +36-62322529; e-mail: [email protected].
scaled equipment, the experimental monitoring of the temperature distribution has only been realized for at least two to three orders of magnitude higher lateral scales [4,5] and experimental studies regarding its time evolution have not been published to date to the best knowledge of the authors. Thus direct experimental results are usually replaced by various model calculations [6-13]. Since the analytical solution of the heat diffusion equation [6,7,9-11] is restricted usually to the simplest idealized cases, numerical approaches should be preferred [12,13] to obtain more realistic results. The temperature dependence of the materials' thermophysical properties, such as heat conduction coefficient, mass density, specific heat, optical absorption coefficient, and the effect of the incident intensity distribution, especially in the
0169-4332/96/$15.00 Copyright © 1996 Elsevier Science B.V. All rights reserved. PII S 0 1 6 9 - 4 3 3 2 ( 9 6 ) 0 0 4 2 1 - 7
Zs. Geretocszkv et al. / Applied Surface Science 106 (1996) 422-428
case of multilayered substrates can only be evaluated using a numerical approach. Though most of the related papers devoted their attention to the investigation of gas phase systems [14], liquids or solutions as precursor materials has already been demonstrated to be a promising, cheap alternative for surface micro modification [15,16], especially when the contamination of the processed surface is not critical (e.g. in the field of sensors or micro machining). The aim of the present contribution is to show that results of a thermal analysis are worth to be considered as a starting point for a deeper understanding of the complex nature of the underlying heat transfer phenomena and prediction of the temperature distribution in a thin liquid film covered bilayered substrate.
2. Model system and caiculational procedure The model system, used for temperature calculations, consists of three layers, as shown in Fig. 1. While a thin tungsten layer covered glass slide is taken as a heated bilayer substrate, the top layer is assumed to have the thermophysical properties of water. This structure, lying in the x-y plane of the orthogonal coordinate system, is perpendicularly illuminated with a focused circular Gaussian laser beam and scanned towards the positive x direction. Only the midlayer tungsten is assumed to absorb the incident laser radiation at the processing wavelength.
423
As the liquid motion - - and its effect on the energy transport - - will be incorporated into an improved future work, the results presented here refer to temperature conditions in which the fluid layer is considered to be 'rigid' and participating in the heat conduction only through its thermal properties. As the time evolution of the temperature distribution was not of present interest, the most general non-linear form of the heat equation was simplified to
V( K(T)VT) - p ( T ) c ( T ) v V T = - Q ( x, y, z, T),
(1) where Q(x, y, z, T) is the source term (the supplied heat per unit time and unit volume at the point x, y, z of temperature T), c' is the scanning speed of the laser spot relative to the (static) sample, •(T), p(T) and c(T) are the temperature dependent heat conduction coefficient, mass density and specific heat of the material, respectively. As the source term, we have considered the in-depth absorption of the incident laser light in the tungsten layer alone. Thus the absorbed intensity as a function of temperature can be written as P(1 - R ( T ) )
l ( x , v, z , T ) =
W27/.
X exp
(
x2+y2 w2
foaC~d z)
(2) Laser beam
~
~
~ ~ .... water
./;t111 // .
100 gm
. z
<
.
_ 20 ,W
l
// / ,/
V
t
Vscanning
Fig. 1. Schematic g e o m e t r y of the scanning laser illuminated system, used for model calculations. The lateral ( x y ) dimensions are given in focal diameter units, and the thickness of the layers in micrometers.
where P is the total incident power of the beam impinging on the water-tungsten interface, w is the ( l / e ) radius of the focal spot, R(T) is the temperature dependent reflectivity of the HzO-W interface and a is the optical absorption coefficient of tungsten, which is approximately temperature independent in the visible spectral domain. Other energy generation terms, e.g. latent heat of phase changes, were not considered here. A second order finite difference method was used to solve the heat equation describing the heat transport in the structure. The discretization was performed in a Cartesian coordinate system for which a 3D mesh of non equidistant nodes was set up within the structure, with boundary nodes lying at the bor-
424
Z.s. Geretot'szky et al. / Applied Surface Science 106 (1996) 4 2 2 - 4 2 8
der of the structure. The nodes were defined more densely where a greater heat gradient was expected (e.g. within the spot diameter and in the tungsten layer) and the node distance was gradually increased outside these critical areas. The boundary conditions were the following: (i) At the interfaces of the layers the heat flux and the temperature is continuous, ensuring the energy conservation
the structure, except in the x - z plane at y = 0 (where T ( x i, y 1, zk) is replaced by T ( x i, Y+l, zk)), and setting the temperature outside the structure to the ambient temperature, TO (e.g. T(x,,+L, Yj, Zk ) = T ( X _ l , y~, Z k ) = T 0 ). As an example, the boundary condition applied at the far side of the structure in the x direction is TO- T ( x,,, y~, zk) ~3"n.n-
I
Xn+ 1 -- Xn
T( xi' Y.i' zk+,) - T ( x i , Yi' zk) -a,~-~
v(x,,, y j,
Zk+ I -- Zk
at zk
=
Ik
gk+ l
-k
-
-
I .~
at zk = Zw gl....
~
Yi, Zk -- Zk
]
(3)
where Kw , K_T and KP denote the heat conduction coefficient of'the water, tungsten and glass, respectively, estimated between the ith and jth nodes, as an arithmetic average of heat resistance (i.e. l / K ) [12,17]. (ii) Escape of energy into the surroundings is prohibited by setting zero heat flux at the border of
Although the central difference discretization procedure and the above described boundary conditions are widely used [12,13], it could lead to mathematically unstable and physically incorrect result, i.e. when the convective term is becoming dominant over the conductive one in the heat equation [17]. This is obviously the case when relatively high scanning speeds are applied. Since in liquid environment, the higher the scanning speed (typically exceeding 1 m m / s ) the better the morphology of the deposited stripes [15,16], it is important to make model calculations at high scanning speeds. At high speed a weak point in the central difference formulation is the assumption that the temperature at an interface perpendicular to the direction of the motion is the average of the temperatures in the nodes lying on its sides. The rationale behind the
Table 1 The applied values of the thermophysical properties for the water, tungsten and glass layer Property
Assumed value or function
KH~o ( W i n -I K - I ) Kw-(Wm I K-I) Kglass ( W m I K - t) CH,o(Jkg ~ K-I) C w ( J k g -I K - I ) cgl~ (J k g - I K - I )
- 0.58180 + 6.357044 X 10- 3T - 7.9662523 X 10- 6 T 2 42.645 + 1.8985 X 104T - i _ 1.898 x 106T 2 0.9094 + 1.422 x 10-3T - 2 . 4 2 1 3 9 X 10 S T 3 + 2 . 6 8 5 3 6 X 10-ST 2 - 9 . 6 8 1 3 7 X 1 0 - 3 T + 2 . 1 3 9 7 4 2.5767 X 1 0 - 2 T + 1.24735 × 102 0.275T + 587.5 9.982 x l02 19.3 x 103 2.65 x 103 0.31 6.6164 x 107
PH,O (kg m - 3 ) pw- (kg m -3 ) P~l,~ (kg m -3) R, H 2 0 - W interface a w @514.5 nm ( m - i)
0.
(4)
Zk
yj,
K.G
L, Yj,
Zk -- Z k - 1
ZH20-W
-k,~.l
-r(x,,
I,.
Xn -- X n - 1
T ( x i ' Yi' z k + , ) - T ( x i , Yi' zk)
=
/4%
r ( x i , Yi' z k ) - T ( x i , vj, z k - , )
= KT -~
--
Zs. Geretot'szlo' et al./ Applied Surface Science 106 (1996) 422-428
upwind scheme [17] is that at high speeds the temperature at an interface is much closer to the value of temperature at the node on the upwind side of the face than that on the other side. It means that one can get a better estimate if one substitutes all the T(xi+ 1, Yj, zk) formulae with T(xi, Y,i' Zk) in the x terms of each discretized equation and calculates the first order difference in the convective term of each equation in two consecutive nodes rather than on three. It is also necessary to make the far end side of the structure (Fig. 1) thermally 'open', i.e. applying non-fixed temperatures for the nodes lying outside the structure and thus allowing a non zero heat flux at that boundary. Consequently using upwind formulation a new boundary condition will substitute for Eq. (4) in the x direction:
3. Result and discussion
3.1. Scanning speed Since in the previous experimental studies [15,16] it was shown that the scanning speed plays a fundamental role in pyrolytic liquid phase deposition we start with the description of scanning speed related effects.
380
(a)
-
K):n I. n
Xn
-
-
=Fluxou,(xn, Yj, zk).
Wa~r:
3 Fm above
360
the •~
,o\/
340
m
v( x,,, yj, zk) -V(x,,_,, yj, zk) -
425
~. 320
,\\~/../100
E
Xn- 1
I-
(5)
Although the temperature variation of the thermal properties causes nonlinearity in the heat diffusion equation, their neglect can result in unrealistic temperature distributions. Thus the temperature dependence of the material parameters has accurately been taken into account. The thermal conductivity and heat capacity of water, tungsten and glass, as a function of temperature, used in this work, are listed in Table 1. Although the temperature dependence of reflectivity and absorption coefficient of W could also be incorporated into the calculation procedure, as they are practically constant in our case, temperature independent values were used (cf. Table 1). Finally, the thermal dilatation of each material have been neglected, as the variation of their density with temperature is less than a few percent, even for the most critical water layer. The heat diffusion problem was mathematically reduced to a set of non-linear algebraic equations, which was solved by using the Gauss-Seidel iterative method [18]. This procedure can be further improved by SOR (successive overrelaxation) [18]. In the next section the effects of practically important process parameters, namely scanning speed, laser power and focal diameter, on the calculated temperature distribution will be discussed in detail.
300 280 -15
-10
-5
0
5
10
15
-15 38O
-10
-5
0
5
10
15
380 360 ,,-1 340 Q. 320
E I,-
3OO
280
(c)
~,,ss:
360
2 i~m below
i 340
10.
320
/1~ the tungsten-glasl tl/l~~ interface
,o "W // I '\, O.i~/7 /
300
-15
-10
-5
/" J,,,,"' 0
':>¢~-~-)~/' ~°°0 5
10
5
Normalized x position
Fig. 2. Typical temperature profiles in the scanning (x) direction calculated for five different scanning speeds (given in c m / s units) in the water (a), tungsten (b) and glass (c) layer, respectively. The laser power is 10 roW. The x coordinates are given in focal diameter units (2w = ~ - x 9/zm).
Zs. GeretoL,szky et al. / Applied Surface Science 106 (1996) 422-428
426
Fig. 2a, b and c show representative temperature distributions as a function of normalized x coordinate (defined as x . . . . . lized = X/W) at five different scanning speeds from 1 m m / s up to 10 m / s in the water (3 /xm above the H 2 0 - W interface), tungsten (in the middle of the 90 nm thin film) and glass (2 /xm below the W-glass interface) layer, respectively. Although the temperature profiles deviate from a symmetrical distribution even at 1 m m / s or 1 c m / s scanning speed, they become strongly distorted only above I c m / s , as seen in each graph in Fig. 2. (i) In the direction of scanning the ascending edge of the profiles becomes steeper, while the slope of the descending edge decreases, as the speed is increasing. For the same scanning speed the temperature profile is the narrowest in the tungsten layer, which is the best heat conductor within the system. (ii) The maximum temperature (referred to as peak-temperature) decreases with increasing scanning speed. (iii) A slight shift of the peak-temperature in the x direction (i.e. in the direction of scanning) is also observed. This shift is more pronounced in the less
O Glass 1 m/s
•
conducting water and glass and can result in more than ten micrometer lagging of the peak temperature behind the center of the illuminated area. It should be noted, that most curves, plotted in Fig. 2, could not be calculated without the upwind scheme when scanning speeds higher than a few c m / s are applied. This can be shown to be a consequence of the fact, that above this speed the convective term within the glass layer becomes dominant over the conductive one. As the absorbed energy of the incoming laser light is spreading through the system, heat also diffuses laterally, thus the motion of the laser spot relative to the sample will also modify the symmetry of the temperature distribution in the y direction. In general, the farther we are (in the y direction) from the illumination center, the smaller is the actual-peak temperature (causing an increase in the width of the temperature profile) and the bigger the displacement of the actual-peak temperature, as it is emphasized by the highlighted trace on the 3D graph of Fig. 3. On the inset in Fig. 3 the positions of the actual-peak temperatures are projected to the x - y plane for the
Water 1 m/s
A
Glass 10 cm/s
•
Water 10 cm/s ]
i j
315
6 8 ed y position
20
15
]0
5
0
-5
-]0
-]5
j
10
-20
Normalized x position
Fig. 3. The temperature distribution in the x-y plane, 3 /~m above the water-tungsten interface. The highlighted curve connects the actual-peak temperature points. The inset shows the position of actual-peak temperature points in the x-y plane for water and glass, at l0 cm/s and l m/s scanning speeds, respectively.
Z.s. Geretot,szky et al./ Applied Surface Science 106 (1996) 422-428
10 c m / s and 1 m / s speeds, to show a parabolic behavior of the peak shift. For comparison, similar curves calculated for the glass layer are also presented with empty symbols. Nevertheless the actualpeak temperature as a function of normalized y position resembles the shape of the temperature distribution for ~' = 0.
'
For comparing the calculated results with experimental data the following experiment was designed and carried out. For an At-ion laser (A = 514.5 rim) based liquid phase deposition system, detailed elsewhere [16], the variation of the focal diameter was performed by the lateral displacement of the focusing microscope objective, perpendicularly to the horizontally aligned sample surface. Scanning the beam over a bilayer system consisting of a 90 nm thick vacuum evaporated tungsten covered glass slide immersed into water, the threshold power for the onset of bubble formation due to boiling has been detected at three different focal diameters. In Fig. 5 the obtained experimental results are compared with the results of numerical temperature
I
I
I
E
I
I
I
'
I
I
I
'
I
I
I
..~.=
360 34°
- 320 P
300 280
i
2
3.3. Focal diameter
'
380
3.2. Laser p o w e r
The power of the laser beam is also an important process parameter. The variation of the total incident power, as an input parameter of our calculational procedure, shows that a few mW laser power are sufficient to increase the peak temperature in the tungsten layer and consequently in the water film, being in contact with it, close to the boiling point of the water. Within the 93°C wide temperature domain (ranging from 280 K to 373 K) the peak temperature in all layers was found to be directly proportional to the total incident laser power (cf. Fig. 4). The slope of the peak temperature versus power function can be modified by changing the focal diameter. In general, the wider the illuminated area (i.e., greater diameter) the smaller the temperature increment, belonging to a unit power increase. The peak temperature at the water-tungsten interface as a function of incident power are given in Fig. 4 for three different spot radii.
I
427
4
6
8
10
I
12
Incident power [mW] Fig. 4. The peak temperature of the water-tungsten interface as a function of incident power for three different 1 / e spot radii (t, = 1 m m / s ) . The calculated slopes are 12.1 K / m W , 9.6 K / r o W and 7.4 K / r o W for 9 / z m (©), ~ - × 9 / z m ( O ) and 18 /~m ( D ) spot diameters, respectively.
calculations according to the method described in this work. Although, both the calculated and the experimentally determined threshold power are directly proportional to the l / e diameter of the laser spot, the calculated values obtained were always far below the measured ones, indicating that some important energy consuming process has not been included into the model. This is not surprising, since the present model does not take into account phase changes and effects relevant to the dynamics of the bubble formation. The fact that the ratio of the experimentally and theoretically determined threshold power increases with increasing focal diameter (from 1.6 at 6.4 /J,m up to 2.1 at 12.7 /zm) implies
25 • .E. 20
Calculated i
/ /
[] Measured I
t" / / [ 3 ~ /
eJ 15 "0
,.I
/
/-
/"
oe-
/
/" f
el"-
6
8
10
12
14
l i e f o c a l d i a m e t e r [gm]
Fig. 5. Measured ([23) and calculated ( 0 ) threshold power of bubble formation as a function of focal diameter using a scanning speed of 1 m m / s .
428
Zs. Geretocszky et al. / Applied Surface Science 106 (1996) 422-428
the existence of a process, which plays an increasing role in the removal of the heat from the vicinity of the hottest region as a greater liquid volume is heated. The vortical liquid motion may be such a process, as it is reported to have an important role in the heat transport of convective fluid motion. Thus the incorporation of fluid dynamics into our model may improve the compatibility of calculated results with the experiments.
4. Conclusion Steady-state temperature distributions calculated for a system, consisting of a water film covered bilayer sample are obtained by solving numerically the three-dimensional non-linear heat diffusion equation. Since liquid phase processes have shown experimentally the necessity of the application of high scanning speeds, the second-order central finite difference method was improved by applying the upwind scheme during the discretization procedure. The effects of practically important process parameters, in particular scanning speed, laser power and focal diameter, on the temperature profiles were given. The comparison of the calculated threshold power versus spot diameter function with the experimental results, as obtained from scanning a focused beam of an Ar-ion laser over a tungsten covered glass slide immersed in water, reveals that although our thermal analysis describes qualitatively the experiments, the incorporation of fluid dynamics and its related thermal effects may be needed to obtain a better quantitative model, being able to account for finer details of the complex heat conduction problem.
Acknowledgements The partial support of the B.6 Hungarian-Austrian Intergovernmental Scientific and Technological Program and the OTKA W 015239 is greatly acknowledged.
References [l] I.W. Boyd and R.B. Jackman, Photochemical Processing of Electronic Materials, I st ed. (Academic Press, London, 1992). [2] DJ. Elliott, Microlithography: Process Technology for IC Fabrication, 1st ed. (McGraw-Hill, New York, 1986). [3] I.P. Herman, in: Laser Chemical Processing for Microelectronics, Eds. K.G. Ibbs and R.M. Osgood (Cambridge University Press, 1989) p. 61. [4] T.T. Kodas, T.H. Baum and P.B Comita, J. Appl. Phys. 61 (1987) 2749. [5] MJ.J. Theunissen, R.P.M.L.C. van de Nieuwenhoff and H.K. Kuiken, J. Appl. Phys. 68 (1990) 806. [6] M. von Allmen, Laser-Beam Interactions with Materials, 1st ed. (Springer, Berlin, 1987). [7] DJ. Sanders, Appl. Opt. 23 (1984) 30. [8] J. Ascough, Opt. Lasers Eng. 6 (1985) 137 [9] A.H.M. Holtslag, J. App. Phys. 66 (1989) 1530. [10] J.E. Moody and R.H. Hendel, J. Appl. Phys. 53 (1982) 4364. [11] Y,I. Nissim, A. Lietoila, R.B. Gold and J.F. Gibbons, J. Appl. Phys. 51 (1980) 274. [t2] T.T. Rantala and J. Levoska, J. Appl. Phys. 65 (1989) 4475. [13] C. Garrido, B. Le6n and M. P6rez-Amor, J. Appl. Phys. 69 (1991) 1133. [14] G, Auvert, Appl. Surf. Sci. 86 (1995) 466. [15] S. Pfliiger, M. Wehner, F. Jansen, Th. Kruck and F. Lupp, Appl. Surf. Sci. 86 (1995) 504. [16] Zs. Geretovszky, L. Kelemen, K. Bali and T. Sz~Sr6nyi, Appl. Surf. Sci. 86 (1995) 495. [17] S.V. Patankar, Numerical Heat Transfer and Fluid Flow, 1st ed. (Hemisphere, Washington, DC, 1980). [18] G,D. Smith, Numerical Solution of Partial Differential Equations: Finite Difference Methods, 2nd ed. (Oxford University Press, New York, 1978).
Applied Surface Science 106 (1996) 422-428
Temperature distribution in multilayers covered by liquid layer and processed by focused laser beam Zs. Geretovszky a,*, L. Kelemen b, K. Piglmayer c Research Group on Laser Physics of the Hungarian Academy of Sciences, Ddm t~r 9, H-6720 S=eged, Hungary b Department of Optics and Quantum Electronics, JA TE Unit'el:s'io:, D6m tgr 9, H-6720 Szeged, Hungary c Ange~,andte Physik, Johannes Kepler Uni~'ersit~it Lin=, A-4040 Lin:, Austria Received 17 September 1995; accepted 31 December 1995
Abstract
As the experimental determination of surface temperatures on micrometer scale can only be carried out in very special cases, the estimation of spatial temperature distributions requires numerical model calculations. In the present contribution experimental results obtained by scanning a focused beam of an Ar-ion laser over a tungsten covered glass slide immersed in water, are compared with the results of temperature calculations. The steady-state temperature distribution is calculated by solving numerically the three-dimensional heat diffusion equation using the second-order finite difference method. During the discretization procedure the upwind scheme is applied to ensure the stability of the equation system at scanning speeds as high as few tens of meter per second. The temperature dependence of the thermophysical properties of the materials is explicitly included. The effects of processing parameters, in particular laser power, focus diameter and scanning speed, are described and found to be consistent with the experimental results.
1. Introduction
In spite of continuous efforts to reduce the temperature during both processing and operating conditions [1-3], thermal damage of the previously formed (heat sensitive) patterns through intermixing, surface degradation or recrystallization can never be completely avoided. Thus the exact knowledge of temperature distribution is essential in process optimization. Although under strictly optimized conditions sub0.5 /zm resolution can be achieved by industrial
:: Corresponding author. Tel.: +36-62-454274; fax: +36-62322529; e-mail: [email protected].
scaled equipment, the experimental monitoring of the temperature distribution has only been realized for at least two to three orders of magnitude higher lateral scales [4,5] and experimental studies regarding its time evolution have not been published to date to the best knowledge of the authors. Thus direct experimental results are usually replaced by various model calculations [6-13]. Since the analytical solution of the heat diffusion equation [6,7,9-11] is restricted usually to the simplest idealized cases, numerical approaches should be preferred [12,13] to obtain more realistic results. The temperature dependence of the materials' thermophysical properties, such as heat conduction coefficient, mass density, specific heat, optical absorption coefficient, and the effect of the incident intensity distribution, especially in the
0169-4332/96/$15.00 Copyright © 1996 Elsevier Science B.V. All rights reserved. PII S 0 1 6 9 - 4 3 3 2 ( 9 6 ) 0 0 4 2 1 - 7
Zs. Geretocszkv et al. / Applied Surface Science 106 (1996) 422-428
case of multilayered substrates can only be evaluated using a numerical approach. Though most of the related papers devoted their attention to the investigation of gas phase systems [14], liquids or solutions as precursor materials has already been demonstrated to be a promising, cheap alternative for surface micro modification [15,16], especially when the contamination of the processed surface is not critical (e.g. in the field of sensors or micro machining). The aim of the present contribution is to show that results of a thermal analysis are worth to be considered as a starting point for a deeper understanding of the complex nature of the underlying heat transfer phenomena and prediction of the temperature distribution in a thin liquid film covered bilayered substrate.
2. Model system and caiculational procedure The model system, used for temperature calculations, consists of three layers, as shown in Fig. 1. While a thin tungsten layer covered glass slide is taken as a heated bilayer substrate, the top layer is assumed to have the thermophysical properties of water. This structure, lying in the x-y plane of the orthogonal coordinate system, is perpendicularly illuminated with a focused circular Gaussian laser beam and scanned towards the positive x direction. Only the midlayer tungsten is assumed to absorb the incident laser radiation at the processing wavelength.
423
As the liquid motion - - and its effect on the energy transport - - will be incorporated into an improved future work, the results presented here refer to temperature conditions in which the fluid layer is considered to be 'rigid' and participating in the heat conduction only through its thermal properties. As the time evolution of the temperature distribution was not of present interest, the most general non-linear form of the heat equation was simplified to
V( K(T)VT) - p ( T ) c ( T ) v V T = - Q ( x, y, z, T),
(1) where Q(x, y, z, T) is the source term (the supplied heat per unit time and unit volume at the point x, y, z of temperature T), c' is the scanning speed of the laser spot relative to the (static) sample, •(T), p(T) and c(T) are the temperature dependent heat conduction coefficient, mass density and specific heat of the material, respectively. As the source term, we have considered the in-depth absorption of the incident laser light in the tungsten layer alone. Thus the absorbed intensity as a function of temperature can be written as P(1 - R ( T ) )
l ( x , v, z , T ) =
W27/.
X exp
(
x2+y2 w2
foaC~d z)
(2) Laser beam
~
~
~ ~ .... water
./;t111 // .
100 gm
. z
<
.
_ 20 ,W
l
// / ,/
V
t
Vscanning
Fig. 1. Schematic g e o m e t r y of the scanning laser illuminated system, used for model calculations. The lateral ( x y ) dimensions are given in focal diameter units, and the thickness of the layers in micrometers.
where P is the total incident power of the beam impinging on the water-tungsten interface, w is the ( l / e ) radius of the focal spot, R(T) is the temperature dependent reflectivity of the HzO-W interface and a is the optical absorption coefficient of tungsten, which is approximately temperature independent in the visible spectral domain. Other energy generation terms, e.g. latent heat of phase changes, were not considered here. A second order finite difference method was used to solve the heat equation describing the heat transport in the structure. The discretization was performed in a Cartesian coordinate system for which a 3D mesh of non equidistant nodes was set up within the structure, with boundary nodes lying at the bor-
424
Z.s. Geretot'szky et al. / Applied Surface Science 106 (1996) 4 2 2 - 4 2 8
der of the structure. The nodes were defined more densely where a greater heat gradient was expected (e.g. within the spot diameter and in the tungsten layer) and the node distance was gradually increased outside these critical areas. The boundary conditions were the following: (i) At the interfaces of the layers the heat flux and the temperature is continuous, ensuring the energy conservation
the structure, except in the x - z plane at y = 0 (where T ( x i, y 1, zk) is replaced by T ( x i, Y+l, zk)), and setting the temperature outside the structure to the ambient temperature, TO (e.g. T(x,,+L, Yj, Zk ) = T ( X _ l , y~, Z k ) = T 0 ). As an example, the boundary condition applied at the far side of the structure in the x direction is TO- T ( x,,, y~, zk) ~3"n.n-
I
Xn+ 1 -- Xn
T( xi' Y.i' zk+,) - T ( x i , Yi' zk) -a,~-~
v(x,,, y j,
Zk+ I -- Zk
at zk
=
Ik
gk+ l
-k
-
-
I .~
at zk = Zw gl....
~
Yi, Zk -- Zk
]
(3)
where Kw , K_T and KP denote the heat conduction coefficient of'the water, tungsten and glass, respectively, estimated between the ith and jth nodes, as an arithmetic average of heat resistance (i.e. l / K ) [12,17]. (ii) Escape of energy into the surroundings is prohibited by setting zero heat flux at the border of
Although the central difference discretization procedure and the above described boundary conditions are widely used [12,13], it could lead to mathematically unstable and physically incorrect result, i.e. when the convective term is becoming dominant over the conductive one in the heat equation [17]. This is obviously the case when relatively high scanning speeds are applied. Since in liquid environment, the higher the scanning speed (typically exceeding 1 m m / s ) the better the morphology of the deposited stripes [15,16], it is important to make model calculations at high scanning speeds. At high speed a weak point in the central difference formulation is the assumption that the temperature at an interface perpendicular to the direction of the motion is the average of the temperatures in the nodes lying on its sides. The rationale behind the
Table 1 The applied values of the thermophysical properties for the water, tungsten and glass layer Property
Assumed value or function
KH~o ( W i n -I K - I ) Kw-(Wm I K-I) Kglass ( W m I K - t) CH,o(Jkg ~ K-I) C w ( J k g -I K - I ) cgl~ (J k g - I K - I )
- 0.58180 + 6.357044 X 10- 3T - 7.9662523 X 10- 6 T 2 42.645 + 1.8985 X 104T - i _ 1.898 x 106T 2 0.9094 + 1.422 x 10-3T - 2 . 4 2 1 3 9 X 10 S T 3 + 2 . 6 8 5 3 6 X 10-ST 2 - 9 . 6 8 1 3 7 X 1 0 - 3 T + 2 . 1 3 9 7 4 2.5767 X 1 0 - 2 T + 1.24735 × 102 0.275T + 587.5 9.982 x l02 19.3 x 103 2.65 x 103 0.31 6.6164 x 107
PH,O (kg m - 3 ) pw- (kg m -3 ) P~l,~ (kg m -3) R, H 2 0 - W interface a w @514.5 nm ( m - i)
0.
(4)
Zk
yj,
K.G
L, Yj,
Zk -- Z k - 1
ZH20-W
-k,~.l
-r(x,,
I,.
Xn -- X n - 1
T ( x i ' Yi' z k + , ) - T ( x i , Yi' zk)
=
/4%
r ( x i , Yi' z k ) - T ( x i , vj, z k - , )
= KT -~
--
Zs. Geretot'szlo' et al./ Applied Surface Science 106 (1996) 422-428
upwind scheme [17] is that at high speeds the temperature at an interface is much closer to the value of temperature at the node on the upwind side of the face than that on the other side. It means that one can get a better estimate if one substitutes all the T(xi+ 1, Yj, zk) formulae with T(xi, Y,i' Zk) in the x terms of each discretized equation and calculates the first order difference in the convective term of each equation in two consecutive nodes rather than on three. It is also necessary to make the far end side of the structure (Fig. 1) thermally 'open', i.e. applying non-fixed temperatures for the nodes lying outside the structure and thus allowing a non zero heat flux at that boundary. Consequently using upwind formulation a new boundary condition will substitute for Eq. (4) in the x direction:
3. Result and discussion
3.1. Scanning speed Since in the previous experimental studies [15,16] it was shown that the scanning speed plays a fundamental role in pyrolytic liquid phase deposition we start with the description of scanning speed related effects.
380
(a)
-
K):n I. n
Xn
-
-
=Fluxou,(xn, Yj, zk).
Wa~r:
3 Fm above
360
the •~
,o\/
340
m
v( x,,, yj, zk) -V(x,,_,, yj, zk) -
425
~. 320
,\\~/../100
E
Xn- 1
I-
(5)
Although the temperature variation of the thermal properties causes nonlinearity in the heat diffusion equation, their neglect can result in unrealistic temperature distributions. Thus the temperature dependence of the material parameters has accurately been taken into account. The thermal conductivity and heat capacity of water, tungsten and glass, as a function of temperature, used in this work, are listed in Table 1. Although the temperature dependence of reflectivity and absorption coefficient of W could also be incorporated into the calculation procedure, as they are practically constant in our case, temperature independent values were used (cf. Table 1). Finally, the thermal dilatation of each material have been neglected, as the variation of their density with temperature is less than a few percent, even for the most critical water layer. The heat diffusion problem was mathematically reduced to a set of non-linear algebraic equations, which was solved by using the Gauss-Seidel iterative method [18]. This procedure can be further improved by SOR (successive overrelaxation) [18]. In the next section the effects of practically important process parameters, namely scanning speed, laser power and focal diameter, on the calculated temperature distribution will be discussed in detail.
300 280 -15
-10
-5
0
5
10
15
-15 38O
-10
-5
0
5
10
15
380 360 ,,-1 340 Q. 320
E I,-
3OO
280
(c)
~,,ss:
360
2 i~m below
i 340
10.
320
/1~ the tungsten-glasl tl/l~~ interface
,o "W // I '\, O.i~/7 /
300
-15
-10
-5
/" J,,,,"' 0
':>¢~-~-)~/' ~°°0 5
10
5
Normalized x position
Fig. 2. Typical temperature profiles in the scanning (x) direction calculated for five different scanning speeds (given in c m / s units) in the water (a), tungsten (b) and glass (c) layer, respectively. The laser power is 10 roW. The x coordinates are given in focal diameter units (2w = ~ - x 9/zm).
Zs. GeretoL,szky et al. / Applied Surface Science 106 (1996) 422-428
426
Fig. 2a, b and c show representative temperature distributions as a function of normalized x coordinate (defined as x . . . . . lized = X/W) at five different scanning speeds from 1 m m / s up to 10 m / s in the water (3 /xm above the H 2 0 - W interface), tungsten (in the middle of the 90 nm thin film) and glass (2 /xm below the W-glass interface) layer, respectively. Although the temperature profiles deviate from a symmetrical distribution even at 1 m m / s or 1 c m / s scanning speed, they become strongly distorted only above I c m / s , as seen in each graph in Fig. 2. (i) In the direction of scanning the ascending edge of the profiles becomes steeper, while the slope of the descending edge decreases, as the speed is increasing. For the same scanning speed the temperature profile is the narrowest in the tungsten layer, which is the best heat conductor within the system. (ii) The maximum temperature (referred to as peak-temperature) decreases with increasing scanning speed. (iii) A slight shift of the peak-temperature in the x direction (i.e. in the direction of scanning) is also observed. This shift is more pronounced in the less
O Glass 1 m/s
•
conducting water and glass and can result in more than ten micrometer lagging of the peak temperature behind the center of the illuminated area. It should be noted, that most curves, plotted in Fig. 2, could not be calculated without the upwind scheme when scanning speeds higher than a few c m / s are applied. This can be shown to be a consequence of the fact, that above this speed the convective term within the glass layer becomes dominant over the conductive one. As the absorbed energy of the incoming laser light is spreading through the system, heat also diffuses laterally, thus the motion of the laser spot relative to the sample will also modify the symmetry of the temperature distribution in the y direction. In general, the farther we are (in the y direction) from the illumination center, the smaller is the actual-peak temperature (causing an increase in the width of the temperature profile) and the bigger the displacement of the actual-peak temperature, as it is emphasized by the highlighted trace on the 3D graph of Fig. 3. On the inset in Fig. 3 the positions of the actual-peak temperatures are projected to the x - y plane for the
Water 1 m/s
A
Glass 10 cm/s
•
Water 10 cm/s ]
i j
315
6 8 ed y position
20
15
]0
5
0
-5
-]0
-]5
j
10
-20
Normalized x position
Fig. 3. The temperature distribution in the x-y plane, 3 /~m above the water-tungsten interface. The highlighted curve connects the actual-peak temperature points. The inset shows the position of actual-peak temperature points in the x-y plane for water and glass, at l0 cm/s and l m/s scanning speeds, respectively.
Z.s. Geretot,szky et al./ Applied Surface Science 106 (1996) 422-428
10 c m / s and 1 m / s speeds, to show a parabolic behavior of the peak shift. For comparison, similar curves calculated for the glass layer are also presented with empty symbols. Nevertheless the actualpeak temperature as a function of normalized y position resembles the shape of the temperature distribution for ~' = 0.
'
For comparing the calculated results with experimental data the following experiment was designed and carried out. For an At-ion laser (A = 514.5 rim) based liquid phase deposition system, detailed elsewhere [16], the variation of the focal diameter was performed by the lateral displacement of the focusing microscope objective, perpendicularly to the horizontally aligned sample surface. Scanning the beam over a bilayer system consisting of a 90 nm thick vacuum evaporated tungsten covered glass slide immersed into water, the threshold power for the onset of bubble formation due to boiling has been detected at three different focal diameters. In Fig. 5 the obtained experimental results are compared with the results of numerical temperature
I
I
I
E
I
I
I
'
I
I
I
'
I
I
I
..~.=
360 34°
- 320 P
300 280
i
2
3.3. Focal diameter
'
380
3.2. Laser p o w e r
The power of the laser beam is also an important process parameter. The variation of the total incident power, as an input parameter of our calculational procedure, shows that a few mW laser power are sufficient to increase the peak temperature in the tungsten layer and consequently in the water film, being in contact with it, close to the boiling point of the water. Within the 93°C wide temperature domain (ranging from 280 K to 373 K) the peak temperature in all layers was found to be directly proportional to the total incident laser power (cf. Fig. 4). The slope of the peak temperature versus power function can be modified by changing the focal diameter. In general, the wider the illuminated area (i.e., greater diameter) the smaller the temperature increment, belonging to a unit power increase. The peak temperature at the water-tungsten interface as a function of incident power are given in Fig. 4 for three different spot radii.
I
427
4
6
8
10
I
12
Incident power [mW] Fig. 4. The peak temperature of the water-tungsten interface as a function of incident power for three different 1 / e spot radii (t, = 1 m m / s ) . The calculated slopes are 12.1 K / m W , 9.6 K / r o W and 7.4 K / r o W for 9 / z m (©), ~ - × 9 / z m ( O ) and 18 /~m ( D ) spot diameters, respectively.
calculations according to the method described in this work. Although, both the calculated and the experimentally determined threshold power are directly proportional to the l / e diameter of the laser spot, the calculated values obtained were always far below the measured ones, indicating that some important energy consuming process has not been included into the model. This is not surprising, since the present model does not take into account phase changes and effects relevant to the dynamics of the bubble formation. The fact that the ratio of the experimentally and theoretically determined threshold power increases with increasing focal diameter (from 1.6 at 6.4 /J,m up to 2.1 at 12.7 /zm) implies
25 • .E. 20
Calculated i
/ /
[] Measured I
t" / / [ 3 ~ /
eJ 15 "0
,.I
/
/-
/"
oe-
/
/" f
el"-
6
8
10
12
14
l i e f o c a l d i a m e t e r [gm]
Fig. 5. Measured ([23) and calculated ( 0 ) threshold power of bubble formation as a function of focal diameter using a scanning speed of 1 m m / s .
428
Zs. Geretocszky et al. / Applied Surface Science 106 (1996) 422-428
the existence of a process, which plays an increasing role in the removal of the heat from the vicinity of the hottest region as a greater liquid volume is heated. The vortical liquid motion may be such a process, as it is reported to have an important role in the heat transport of convective fluid motion. Thus the incorporation of fluid dynamics into our model may improve the compatibility of calculated results with the experiments.
4. Conclusion Steady-state temperature distributions calculated for a system, consisting of a water film covered bilayer sample are obtained by solving numerically the three-dimensional non-linear heat diffusion equation. Since liquid phase processes have shown experimentally the necessity of the application of high scanning speeds, the second-order central finite difference method was improved by applying the upwind scheme during the discretization procedure. The effects of practically important process parameters, in particular scanning speed, laser power and focal diameter, on the temperature profiles were given. The comparison of the calculated threshold power versus spot diameter function with the experimental results, as obtained from scanning a focused beam of an Ar-ion laser over a tungsten covered glass slide immersed in water, reveals that although our thermal analysis describes qualitatively the experiments, the incorporation of fluid dynamics and its related thermal effects may be needed to obtain a better quantitative model, being able to account for finer details of the complex heat conduction problem.
Acknowledgements The partial support of the B.6 Hungarian-Austrian Intergovernmental Scientific and Technological Program and the OTKA W 015239 is greatly acknowledged.
References [l] I.W. Boyd and R.B. Jackman, Photochemical Processing of Electronic Materials, I st ed. (Academic Press, London, 1992). [2] DJ. Elliott, Microlithography: Process Technology for IC Fabrication, 1st ed. (McGraw-Hill, New York, 1986). [3] I.P. Herman, in: Laser Chemical Processing for Microelectronics, Eds. K.G. Ibbs and R.M. Osgood (Cambridge University Press, 1989) p. 61. [4] T.T. Kodas, T.H. Baum and P.B Comita, J. Appl. Phys. 61 (1987) 2749. [5] MJ.J. Theunissen, R.P.M.L.C. van de Nieuwenhoff and H.K. Kuiken, J. Appl. Phys. 68 (1990) 806. [6] M. von Allmen, Laser-Beam Interactions with Materials, 1st ed. (Springer, Berlin, 1987). [7] DJ. Sanders, Appl. Opt. 23 (1984) 30. [8] J. Ascough, Opt. Lasers Eng. 6 (1985) 137 [9] A.H.M. Holtslag, J. App. Phys. 66 (1989) 1530. [10] J.E. Moody and R.H. Hendel, J. Appl. Phys. 53 (1982) 4364. [11] Y,I. Nissim, A. Lietoila, R.B. Gold and J.F. Gibbons, J. Appl. Phys. 51 (1980) 274. [t2] T.T. Rantala and J. Levoska, J. Appl. Phys. 65 (1989) 4475. [13] C. Garrido, B. Le6n and M. P6rez-Amor, J. Appl. Phys. 69 (1991) 1133. [14] G, Auvert, Appl. Surf. Sci. 86 (1995) 466. [15] S. Pfliiger, M. Wehner, F. Jansen, Th. Kruck and F. Lupp, Appl. Surf. Sci. 86 (1995) 504. [16] Zs. Geretovszky, L. Kelemen, K. Bali and T. Sz~Sr6nyi, Appl. Surf. Sci. 86 (1995) 495. [17] S.V. Patankar, Numerical Heat Transfer and Fluid Flow, 1st ed. (Hemisphere, Washington, DC, 1980). [18] G,D. Smith, Numerical Solution of Partial Differential Equations: Finite Difference Methods, 2nd ed. (Oxford University Press, New York, 1978).