Ind. Eng. Chem. Res. 2008, 47, 4791–4796
4791
Nonlinear Model Identification for Temperature Control in Single Wafer Rapid Thermal Processing Wonhui Cho and Thomas F. Edgar Department of Chemical Engineering, UniVersity of Texas, Austin, Texas 78712
Jietae Lee* Department of Chemical Engineering, Kyungpook National UniVersity, Taegu 702-701, Korea
Single wafer rapid thermal processing (RTP) has become one of the key technologies in semiconductor manufacturing due to its faster wafer processing with reduced thermal budget and precise control of processing conditions. As the standard size of the silicon wafer grows and integration of integrated circuits increases, better control to improve the processing time, uniformity, and repeatability of processing is required. In RTP, identification and control are complicated because of high nonlinearity, drift, and time varying nature of the wafer dynamics. Physical models for the wafer dynamics have been available, but they are not utilized fully for identification and control. Here, based on a physical model of wafer, a practical identification method and an analytic linearizing control method are proposed. 1. Introduction Single wafer rapid thermal processing is used for various processes in the fabrication of semiconductor devices, such as rapid thermal annealing, oxidation, nitration, and chemical vapor deposition. It has three steps (1) rapid heating to the desired operating temperature, (2) processing for a prescribed time at constant operating temperature, and (3) rapid cooling to an ambient condition. A control system is required to follow the desired temperature trajectory while maintaining temperature uniformity on the wafer. This control becomes more stringent as the size of integrated circuits become smaller and the wafer diameters grow. In control of rapid thermal processing (RTP), there are several difficulties such as high nonlinearity of the heating process, drift, time-varying nature of wafer processing, and noncontact measurements of temperatures. A review of RTP control has been provided by Edgar et al.1 Various control methods have been applied to RTP. They include decoupling control,2 iterative learning control,3,4 adaptive control,5 internal model control,6 model-based control,2 and nonlinear model predictive control.7 Huang et al.8 used a proportional-double integral-derivative (PI2D) controller to eliminate offset error during the heating step. To design and analyze control systems for RTP, models based on first principles are available.7,9 However, these physical nonlinear models have not been utilized in identifying RTP models and designing control systems. Because temperature uniformity across the wafer is maintained, conductive heat transfer can be ignored and the model becomes a dynamical system with nonlinear diagonal elements and full static input gain matrix, allowing model complexity to be reduced considerably. Consequently, simpler identification and control for RTP can be obtained. Practical methods that use this simple physical nonlinear model in designing the control system and identification to tune control systems are investigated here. To show that the RTP dynamics can be described by a simple nonlinear model, a nonlinear model is identified experimentally. Identification of the nonlinear model of the wafer dynamics is * To whom correspondence should be addressed. Tel.: +82-53-9505620. Fax: +82-53-950-6615. E-mail:
[email protected].
very difficult because there is drift due to unmodeled dynamic components. To overcome this difficulty, a set of local linear models is identified and then a nonlinear model that matches local linear models is constructed. Bauman and Rugh10 proposed this concept in designing control systems for nonlinear processes and Lee et al.11 used this approach in obtaining nonlinear first order plus time delay models for SISO (single-input singleoutput) processes. In obtaining local linear models, a closedloop identification method with a controller having integral action is used to reduce effects of drift and keep the operating condition constant.12 Local linear models for various operating temperatures and a corresponding nonlinear model have been obtained. On the basis of this nonlinear model, a linearizing control system is proposed. It linearizes the system first, and then, a linear control law is applied. In this way, a simple analytic control law whose responses are predictable can be constructed. 2. Rapid Thermal Processing System This work was carried out with experimental RTP equipment developed by Texas Instruments to fabricate 6 in. wafers. The main chamber is shown in Figure 1.12 The edge, middle, and center lamp zones consist of 24, 12, and 1 tungsten-halogen lamps, respectively. They are used for power sources P1, P2, and P3, to control the wafer temperature. The reflector is covered with gold-plated material to enhance the light reflection. The lamp housing and chamber wall are cooled by water continuously during the experiment. A quartz plate (10 in. in diameter and 1/2 in. thick) is placed to divide the lamp zone from the 6 in. wafer. The 6 in. wafer is equipped with chromel-alumel thermocouples (type K, temperature range 0-1100 °C) to measure wafer temperatures. Thermocouples are bonded with high temperature ceramic cement (composed of 60% SiO2 and 40% Al2O3) onto the wafer surface. Locations of the thermocouples are shown in Figure 1. Physical models for the wafer dynamics of different complexity have been published.7,8,13 Here, a lumped model as in the work of Huang et al.8 is used.
10.1021/ie071378+ CCC: $40.75 2008 American Chemical Society Published on Web 06/13/2008
4792 Ind. Eng. Chem. Res., Vol. 47, No. 14, 2008
FVicP
dTi conv ) qrad + qcond + qwall + qlamps , i + qi i i i dt
) -Aihi(Ti - Tgas) qconv i
i)
1, 2, ... , N (1) where i is the index of the radial N cells on the wafer with volume Vi, mass density F, specific heat capacity cP, and temperature Ti [K]. The specific heat capacity of the wafer is dependent on temperature,7 but it is assumed to be constant because its error is within (10% for temperatures between 500 and 1000 °C. The radiation term can be described as10 N
qrad i ) Ai
∑a
eff eff 4 eff 4 i Γijej σTj - 2Aiei σTi
(2)
where hi is the heat transfer coefficient. The convective heat transfer coefficient depends on the gas properties and flow pattern.14 The conduction term depends on the gradient of wafer temperature across the radial direction and is known to be negligible compared to other radiation and convective terms. Because the wafer temperature is controlled so that the gradient is small, ignoring this conduction term works well. The heat exchange term qiwall is heat addition due to the heated wall. Its dynamics is about ten times slower than that of the wafer. Usually the wall temperature is not measured and it is considered as a slow drift. Heat addition from the lamp power is given as
j)1
where Ai is the surface area of the ith cell, aieff is the effective absorptivity, and eieff is the effective emissivity. The first term in eq 2 Γ describes the fraction of radiation emitted from the surface element j. These radiative exchange factors depend on the geometry of the RTP chamber and the optical properties of the wafer surface and can require complicated computation such as Monte Carlo simulations.9 Schaper et al.6 showed that the off-diagonal terms in Γ are small and can be ignored for the purpose of control system design. The convective heat transfer term is given as
(3)
N
) Aiaeff qlamp i i
∑L P
(4)
ij j
j)1
where Lij consists of the direct and reflective view factors from the lamp j to the annular zone i on the wafer. In summary, we have FVicP
N dTi 4 eff ) -Ai(hiTi + 2eeff σT ) + A a LijPj + di(t), i i i i dt j)1
Figure 1. Experimental single wafer rapid thermal processing equipment developed by Texas Instruments.12
∑
1, 2, ... , N
i) (5)
Ind. Eng. Chem. Res., Vol. 47, No. 14, 2008 4793
j , the For a given operating condition, T1s ) T2s ) T3s ) T system of eq 5 with N ) 3 can be linearized as
[ ] [ ] [
T1(s) P1(s) T2(s) ) G(s) P2(s) ) T3(s) P3(s) 1/(τ1s + 1) k11 k12 k13 P1(s) 0 0 s + 1) k 1/(τ 0 0 2 21 k22 k23 P2(s) s + 1) k 1/(τ 0 0 3 31 k32 k33 P3(s)
][
][ ]
(6)
Because it is hard to maintain the operating condition due to drift, a closed-loop identification method is used. To obtain transfer functions in the first column of G(s), a proportional-integral (PI) controller is designed for the first loop and the loop is closed (Figure 2). The closed-loop transfer functions are
[ ]
[ ]
T1(s) g11(s) c1(s) T2(s) ) g21(s) R (s) 1 + g11(s)c1(s) g (s) 1 T3(s) 31
Figure 2. System for closed loop identification.
(7)
A step change of size ∆T is introduced in the set point for the first loop. The step set point response will be T1(t) ) L-1
(
)
g11(s)c1(s) ∆Ts + d˜(t) 1 + g11(s)c1(s) s
(8)
where L-1(.) means the inverse Laplace transformation. For a given c1(s), assuming d˜(t) to be a linear function of time, the model parameters of g11(s), k11, and τ1 are estimated by the nonlinear least-squares method minimizing the following objective function: J1 )
∫
(
( (
T1(t) - L-1
)
))
g11(s)c1(s) ∆Ts + d˜0 + d˜1t 1 + g11(s)c1(s) s
2
dt (9)
In a similar manner, model parameters of g21(s) and g31(s) can be estimated from the responses of T2(t) and T3(t). For the inverse Laplace transformation and nonlinear optimization, MATLAB routines are used. The SISO (single-input singleoutput) controller c1(s) for this identification is designed for g11(s). Because g11(s) is a first order transfer function without time delay, c1(s) is not difficult to design and robust. A typical design for c1(s) is 1 (10) 2.5s The same procedure can be applied to the second and third loops to obtain transfer functions in the second and third columns of G(s). Identification experiments were performed for various temperature levels. Experimental temperature responses at levels of 400, 500, 600 and 700 °C are shown in Figure 3. PI controllers do not produce oscillatory responses as in Figure 3. These oscillatory responses will be due to dynamic lags of the lamps and the sampling time of 1 s. Time constants identified contain these lags. The experimental procedure for each temperature level is as follows. First, wafer temperatures are moved near the given identification temperature level with roughly tuned multiloop PI controllers. Then, the second and third loops are opened and the set point of the first loop is increased by a specified amount (∆Ts, mostly 5 °C). After the step set point response of the first loop reaches near the steady state, wafer temperatures are returned to the given identification temperature level by closing all three multiloop controllers. Second, responses for the step set point change of the second
(
c1(s) ) 2.5 1 +
Figure 3. Experimental responses for closed loop identification.12
where di(t) includes drift and effects of all components ignored. 3. Identification of Local Linear Models Because parameters such as hi, eieff and Lij in eq 5 are very hard to estimate theoretically, they are usually identified from experimental dynamic data. A large perturbation that will cover the whole operating region is required to obtain the nonlinear model accurately. However, this approach to identify the nonlinear model may not be applied for this process because of the nonlinear time-varying drift. Here, we use a different approach in which local linear models at various wafer temperature levels are first identified and then a nonlinear model is constructed from the local linear models.
)
4794 Ind. Eng. Chem. Res., Vol. 47, No. 14, 2008 Table 1. Steady State Gains and Time Constants for Linearized Models K
τ [s]
0.5519 0.3226 0.0591 0.4366 0.3622 0.0842 0.4371 0.3304 0.1260
13.1728 12.7944 9.4286
0.4449 0.3116 0.0566 0.3564 0.3902 0.0853 0.3463 0.3488 0.1311
8.8518 10.0599 8.246
0.3200 0.2192 0.0512 0.2519 0.2938 0.0631 0.2497 0.2552 0.0988
6.0098 6.7375 4.9344
0.2234 0.1761 0.0269 0.1790 0.2199 0.0437 0.1739 0.1990 0.0783
4.5327 4.9591 3.8760
temp [°C] 400
500
600
700
Figure 4. Fitting results for the NM1 model.
data. From the first column in eqs 12 and 13, we estimate a11, a14, b11, b12, and b13 that minimize
( (
a11 + 4a14T1s 1 wj k b1j 1j j)1
Nonlinear Model 1 (NM1): From eq 5, with assuming that cP, hi, eieff, and aieff are all constant, we have ˙ T 1 ˙ ) -f(T) + BP + d(t) ) T 2 ˙ T
[
3
a11T1 + a14T14
- a21T2 + a24T24 a31T3 + a34T34
] [ ][ ]
b11 b12 b13 P1 + b21 b22 b23 P2 + d(t) b31 b32 b33 P3
(11)
K) b11/(a11 + 4a14T1s3) b12/(a11 + 4a14T1s3) b13/(a11 + 4a14T1s3) b21/(a21 + 4a24T2s3) b22/(a21 + 4a24T2s3) b23/(a21 + 4a24T2s3)
]
b31/(a31 + 4a34T3s3) b32/(a31 + 4a34T3s3) b33/(a31 + 4a34T3s3) (12)
and time constants are
[
T1s
+ w4
1 - (a11 + τ1
))
4a14T1s3)
2
(14)
The weights are set to {wj}){1 1 0.01 100} so that each term has similar magnitude. In a similar manner, we estimate other elements in eq 11. Fitting results are shown in Figure 4, and the minimized cost J2 is shown in Table 2. The model obtained is ˙ ) -f (T) + B P + d(t) T NM1 NM1 T ) (T1,T2,T3)T
[
P ) (P1,P2,P3)T
3.6170 × 10-2T1 + 4.7203 × 10-11T14
fNM1(T) ) 3.4455 × 10-2T2 + 4.2441 × 10-11T24
Steady state gains at the steady state temperatures T1s ) T2s ) T3s are
[
∑∑
(15)
where
4. Nonlinear Model Construction
[]
J2 )
) (
3 2
3
loop are obtained by opening the first and third controllers. The step set point change for the third loop follows. Closed-loop identification makes longer experiments possible in the presence of temperature drifts. Table 1 shows the model parameters identified. Time constants obtained coincide closely to those in the literature.2,8 Drifts make an accurate estimation of linearized models difficult. These drifts are mainly due to the thermal capacitance of the window, dividing the lamp, and wafer zones, which shows slow dynamics. To reduce these drifts, Kim and Cho15 measure the window temperature and include its dynamics in their model.
1/(a11 + 4a14T1s3)
τ ) 1/(a21 + 4a24T2s3) 1/(a31 + 4a34T3s ) 3
]
(13)
In the previous section, we obtained the steady state gains and time constants at various temperature levels T1s ) T2s ) T3s. A nonlinear model (eq 11) can be obtained by fitting these
[
4.5144 × 10-2T3 + 5.5126 × 10-11T34
]
]
4.9823 × 10-2 3.6034 × 10-2 0.6401 × 10-2 BNM1 ) 3.6167 × 10-2 4.1290 × 10-2 0.8654 × 10-2 4.6108 × 10-2 4.8056 × 10-2 1.8020 × 10-2 d(t) ) initial and drift term We call this model NM1. We can see in Figure 4 that fits of the steady state gains are excellent but those of the time constants are somewhat poor. This is due to poor experimental estimation of time constants in the previous section. Nonlinear Model 2 (NM2): In eq 5, coefficients for the term Ti4 are Aieeff 2eeff i σ i σ ) (16) FVicP FzcP where z is the thickness of the wafer. Assuming these coefficients to be constant and independent of the zone (i.e., κi ) κ), we can have a nonlinear model expressed as κi ≡
[] [ ][
][ ]
˙ T a11T1 + κT41 b11 b12 b13 P1 1 ˙ ) - a T + κT4 + b21 b22 b23 P2 + d(t) T 21 2 2 2 4 b31 b32 b33 P3 ˙ a T + κT T 31 3 3 3 (17)
Ind. Eng. Chem. Res., Vol. 47, No. 14, 2008 4795 a
Table 2. Estimated Models
a
model
sum of minimized cost J2
estimation of B
NM1
2.90
x
x
NM2
3.10
|BNM2BNM1-1 -I|i2 ) 0.129
NM3
3.29
|BNM3BNM1-1 - I|i2 ) 0.169
1.00 0 0 RGA(BNM1BNM2-1) ) 0 1.00 0 0 0 1.00 1.00 0 0 RGA(BNM1BNM3-1) ) 0 1.00 0 0 0 1.00
|A|i2 and RGA(A) are the induced 2-norm and the relative gain array of a matrix A, respectively.
Figure 5. Fitting results for the NM2 model.
The coefficient κ is dependent on physical properties of the wafer available in the literature. Here, instead of identifying it from experimental responses, its theoretical estimate is used. Lee et al.14 used F ) 2330 [kg/m3], z ) 6.75 × 10-4 [m], cP ) 748.38 + 0.1678T [W s/(kg K)], eieff ) 0.6, and σ ) 5.67 × 10-8 [W/m2 K4]. In this case, the value κ is between 4.75 × 10-11 and 5.03 × 10-11 for T ∈ [400, 700 °C]. When data in Dassu et al.8 are used, it is between 3.91 × 10-11 and 5.91 × 10-11. The identified coefficients in fNM1(T) of eq 15 are well within these analytic estimations. Here, we used a constant κ ) 5.00 × 10-11 [1/s K3], instead of letting unknown parameters for a simpler identification. Coefficients ai1 and bij are estimated from the set of linear models. We solve the minimization problem of the cost function of eq 14 with a14 ) a24 ) a34 ) κ ) 5.00 × 10-11. The model obtained is ˙ ) -f (T) + B P + d(t) T NM2 NM2 where
[
3.4136 × 10-2T1 + 5.00 × 10-11T1
fNM2(T) ) 2.8803 × 10-2T2 + 5.00 × 10-11T2
[
(18)
]
5.0210 × 10-2T3 + 5.00 × 10-11T3 5.1461 × 10-2 3.7243 × 10-2 0.6612 × 10-2 BNM2 ) 3.9662 × 10-2 4.5393 × 10-2 0.9502 × 10-2 4.4117 × 10-2 4.5899 × 10-2 1.7780 × 10-2
]
Figure 5 shows the fitting results, and Table 2 shows the sum of minimized cost J2. The nonlinear part of this model is predetermined. Hence, it is very useful in identification and controller design. Nonlinear Model 3 (NM3): We further simplify the nonlinear model by assuming a11 ) a21 ) a31 ) a. Because ai1
Figure 6. Fitting results for the NM3 model.
contains the convective heat transfer coefficient hi that depends on the gas properties and flow pattern, it is dependent on wafer zones. However, time constants in the works of Huang et al.8 and Cho12 show that their effects are negligible. Thus the assumption of a11 ) a21 ) a31 ) a will be made. Coefficients a and bij are estimated from the set of linear models. We solve the minimization problem of the cost function of eq 14 and the model obtained is ˙ ) -f (T) + B P + d(t) T NM3 NM3 where
[
3.7130 × 10-2T1 + 5.00 × 10-11T14
fNM3(T) ) 3.7130 × 10-2T2 + 5.00 × 10-11T24
[
(19)
]
3.7130 × 10-2T3 + 5.00 × 10-11T34 5.2402 × 10-2 3.7907 × 10-2 0.6733 × 10-2 BNM3 ) 4.1739 × 10-2 4.7685 × 10-2 0.9991 × 10-2 4.0880 × 10-2 4.2641 × 10-2 1.6520 × 10-2
]
Figure 6 shows the fitting results, and Table 2 shows the sum of minimized cost J2. We can see in Figure 6 that there are considerable mismatches in time constant estimations. Estimations of steady state gains are still excellent. 5. External Linearization Control System and Discussion When a model is given, model-based control systems can be designed. A simple analytic control system designed by the external linearization method is tested: P(t) ) P + B-1(f(T) - RT + V(t))
(20)
4796 Ind. Eng. Chem. Res., Vol. 47, No. 14, 2008
the closed-loop identification method with loosely tuned PI controllers is used to maintain the operating point and suppress drift. From the set of local linear models, three nonlinear models are obtained by imposing different assumptions on the physical model. Model parameters closely match those in the literature and estimates based on physical property data. The second and third models use analytic nonlinear terms whose coefficients are not identified but estimated analytically through physical property data. They can provide a simpler identification method where local linear models at different temperature levels are not needed. A simple linearizing controller is designed and yields excellent closed-loop performances. The proposed model can also be used to test control systems for RTP.
Figure 7. Linearizing control system.
Acknowledgment This work was supported by the Korea Research Foundation Grant funded by the Korean Government (MOEHRD, Basic Research Promotion Fund; KRF-2006-D00099-100554). Literature Cited
Figure 8. Closed-loop responses for the process of NM1 with drifts (controllers are designed from NM2 and NM3).
j is the initial value of P(t). Then, the system becomes where P linear, and the transfer functions between T(s) and V(s) are all 1/(s + R). A PI controller can be used for this linearized system: R I (21) s 3 Figure 7 shows this linearizing control system. For R ) 0.2 and kc ) 3.33, control performances are shown in Figure 8. In the simulation, the nonlinear model NM1 is used for the process and nonlinear models NM2 and NM3 are used to design controllers. Artificial drifts are introduced as shown in Figure 8. Control performances are excellent as expected. In simulations, saturation in input power is not considered. Thus, control performances in real applications can be worse than those in Figure 8 due to physical limitations in input power (i.e., negative inputs are not allowed because heat removal via lamps are not possible). The control law of eq 20 uses the inversion of B. Hence, accurate estimation of B is important. Relative errors of B among models NM1, NM2, and NM3 are shown in Table 2. There are over 10% relative errors in the induced 2-norm. However, as shown in the relative gain array of BNM1BNM2-1 and BNM1BNM3-1, the control law of eq 20 can decouple the process well. Because the decoupled system is first order without time delay, the control system can be robust.
(
PI(s) ) kc 1 +
)
6. Conclusion A method to obtain nonlinear models which use a physical model of the single wafer rapid thermal processing and a set of local linear models is proposed. In obtaining local linear models,
(1) Edgar, T. F.; Butler, S. W.; Campbell, W. J.; Pfeiffer, C.; Bode, C.; Hwang, S. B. Automatic Control in Microelectronics Manufacturing: Practices, Challenges, and Possibilities. Automatica 2000, 36, 1567. (2) Balakrishnan, K. S.; Edgar, T. F. Model-based Control in Rapid Thermal Processing. Thin Solid Films 2000, 365, 322–333. (3) Choi, J. Y.; Do, H. M. A. Learning Approach of Wafer Temperature Control in a Rapid Thermal Processing System. IEEE Trans. Semicond. Manuf. 2001, 14, 1. (4) Cho, M.; Lee, Y.; Joo, S.; Lee, K. S. Semi-empirical Model-based Multivariable Iterative Learning Control of an RTP System. IEEE Trans. Semicond. Manuf. 2005, 18, 430. (5) Choi, J. Y.; Do, H. M.; Choi, H. S. Adaptive Control Approach of Rapid Thermal Processing. IEEE Trans. Semicond. Manuf. 2003, 16, 621. (6) Schaper, C. D.; Kailath, T.; Lee, Y. J. Decentralized Control of Wafer Temperature for Multizone Rapid Thermal Processing Systems. IEEE Trans. Semicond. Manuf. 1999, 12, 193. (7) Dasssu, E.; Grosman, B.; Lewin, D. R. Modeling and Temperature Control of Rapid Thermal Processing. Comput. Chem. Eng. 2006, 30, 686. (8) Huang, C. J.; Yu, C. C.; Shen, S. H. Selection of Measurement Locations for the Control of Rapid Thermal Processor. Automatica 2000, 36, 705. (9) Kersch, A.; Schafbauer, T. Thermal Modeling of RTP and RTCVD Processes. Thin Solid Films 2000, 365, 307. (10) Bauman, W. T.; Rugh, W. J. Feedback Control of Nonlinear Systems by Extended Linearization. IEEE Trans. Automatic Control 1986, 13, 40. (11) Lee, J.; Cho, W.; Edgar, T. F. Control System Design Based on a Nonlinear First Order Plus Time Delay Model. J. Process Control 1997, 7, 65–73. (12) Cho, W. Temperature Control and Modeling of the Rapid Thermal Processing Chamber. PhD Dissertation, University of Texas at Austin, 2005. (13) Lee, K. S.; Lee, J.; Chin, I.; Choi, J.; Lee, J. H. Control of Wafer Temperature Uniformity in Rapid Thermal Processing Using an Optimal Iterative Learning Control Technique. Ind. Eng. Chem. Res. 2001, 40, 1661. (14) Campbell, S. A.; Ahn, K. H.; Knutson, K. L.; Liu, B. Y. H.; Leighton, J. D. Steady-state thermal uniformity and gas flow patterns in a rapid thermal processing chamber. IEEE Trans. Semicond. Manuf. 1991, 4, 14. (15) Kim, S. J.; Cho, Y. M. Optimal Design of a Rapid Thermal Processor via Physics-based Modeling and Convex Optimization. Control Eng. Practice 2003, 10, 1199–1210.
ReceiVed for reView October 15, 2007 ReVised manuscript receiVed April 4, 2008 Accepted April 30, 2008 IE071378+