DETAILED DESCRIPTION

[0040] A typical vertical thermal reactor 100 is shown in FIG. 1. The vertical thermal reactor 100 includes a long quartz or silicon carbide process tube 110 delimiting a process region. A batch of wafers 152, accommodated in a wafer boat 150, placed on a pedestal 151 for support and thermal isolation, are inserted into the process tube 110. The process tube 110 includes an inlet 111 and an outlet 112 for process gas. The process tube is surrounded by a heating element 120 having multiple zone electric heating coils 121 to 125. Each zone has one or more temperature sensors. In FIG. 1, each zone has a spike ThermoCouple (TC) 130 and a “profile” or paddle ThermoCouple (TC) 140. The spike TC produces a spike TC signal corresponding to a spike temperature. The paddle TC produces a paddle TC signal corresponding to a paddle temperature. The spike TCs 130 are located outside the process tube 110 relatively near the heating element and the paddle TCs are located inside the tube 110 relatively near the wafers. The vertical reactor system 100, using the resistive heating element 120 to control temperature, is an inherently non-linear system because an electric heating element can only generate, not absorb, heat. Further, due to the large physical mass of the heating element 120, process tube 110, and wafer batch 152, and a correspondingly high thermal mass or heat capacity, the vertical thermal reactor 100 exhibits long time constants or delay times. This means that after increasing the power input of one or more of the heating coils 121-125 , it takes a relatively long time before a new steady-state at a higher temperature is achieved. When the reactor temperature is higher than the control setpoint, cooling occurs at a rate that depends on the reactor design and its temperature, not on the controller. In the case of overshoot, it takes a relatively long time before the system is cooled down to the control setpoint again, especially at lower temperatures where the natural cooling rate is lower.

[0041] Historically, the temperature in thermal reactors was controlled by a spike TC control loop, using a PID control algorithm. By profiling the furnace in a static mode, using a paddle TC, the relation between paddle TC and spike TC under static conditions was established and stored in a profile table. Such a profiling procedure was performed at regular intervals or after maintenance. Because the paddle TC gives a more relevant reading for the actual wafer temperature, there has been a desire to use a paddle TC control loop that would make the time-consuming profiling procedure unnecessary. A control configuration 200 employing such a paddle TC control loop is shown in FIG. 2. In the configuration 200, an adder 210 computes an error signal Es_Pd from the paddle control setpoint Pdset and the actual paddle temperature Pd. Based on the error signal Es_Pd, the PID controller 220 generates a power output signal Pw that is provided, via a thyristor unit, not shown, to the heating elements 230. The tube and wafers are indicated by 240. A feedback loop 250 that includes a digital filter 251, provides the actual paddle TC signal to the adder 210. However, such a paddle TC control loop, as shown in FIG. 2, has such a strong non-linear behavior and long time constants that it is difficult or impossible to achieve a stable control loop with acceptable performance under dynamic conditions.

[0042] A so-called cascade controller 300, as shown in FIG. 3, uses an inner spike temperature control loop PID controller, with the spike temperatures as inputs, to control the heating elements. An additional outer loop, with the paddle TC temperatures as inputs, is used to generate time-dependent setpoints given to the spike temperature control loop. In the cascade controller 300, the spike TC temperatures are controlled by a PID controller 340 in an inner loop 370. An adder 330 provides an error signal Es_Sp to the PID controller 340. The error signal Es_Sp is based on a spike TC control setpoint Sp_set and the actual spike temperature Sp. The PID controller 340 provides a power output signal Pw that is provided via a thyristor unit, not shown, to the heating elements 350 to heat a tube and wafers 360. A second PID controller 320 in an additional, outer, control loop 380 generates a spike control setpoint Sp_set. The second controller 320 receives a paddle error signal Es_Pd from the adder 310, calculated from the paddle control setpoint Pdset and the actual paddle temperature Pd. Feedback loops 370 and 380 include digital filters 371 and 381 respectively to remove spurious data.

[0043] A typical thermal process starts at a standby temperature at which the wafers are loaded into the thermal reactor. After loading, the thermal reactor heats up to the desired process temperature for oxidation, annealing, drive, or CVD. After performing the process, the thermal reactor cools to the stand-by temperature again and unloads the wafers. If the standby temperature, ramp up/down rate, and process temperatures are set in reasonable ranges, acceptable temperature control performance can be achieved during the process by using the cascade PID controller 300. However, optimizing the performance of the cascade PID controller 300 often requires significant off-line time for tuning of the controller parameters such as the PID parameters. Tuning of a cascade controller is often more of an art than a science and usually very time consuming. The best choice of tuning parameters depends on a variety of factors including the dynamic behavior of the controlled process, the controller's objectives, and the operator's understanding of the tuning procedures. For a cascade PID controller, the inner and outer loop tuning is strongly coupled, which adds to the tuning complexity. Besides the tuning of PID parameters, for dealing with long time delay, the outer or “profile” PID control loop of a cascade PID controller still needs a “profiling table” to provide constraints. Generating the profiling table involves a procedure that requires many hours of off-line equipment time. Off-line time cannot be used for useful wafer processing and is thus very expensive.

[0044] With the advances in modern control technology and system identification, more advanced control systems, such as, for example Model-Based Predictive Controllers (MBPC), have been developed, but these more advanced control methods are often computationally complex, typically requiring matrix inversion during online processing. FIG. 4 shows a single loop control system 400 that includes a MBPC 410. The MBPC 410 receives as input a paddle control setpoint Pdset. The controller 410 also receives as input the actual spike temperatures Sp and actual paddle temperatures Pd via feedback loops 440 and 450, respectively. The controller 410 generates a power output signal Pw that is provided via a thyristor unit, not shown, to heating elements 420 that heat a tube and wafers 430. Feedback loops 440 and 450 include digital filters 441 and 451, respectively. The use of spike temperatures is optional. The MBPC 410 uses a complex dynamic model of the controlled process to compute the predictive control signals by minimizing an objective function to provide on-line optimization control. In the MBPC controller 410, tuning is relatively easy as compared to a cascade PID controller, and the MBPC 410 can compensate for delay times. In the MBPC 410, treatment of constraints in the system to be controlled is conceptually simple and multivariable control is conceptually straightforward, though computationally complex. The MBPC 410 typically relies on dynamic models, which are not required for the traditional PID controllers. The MBPC 410 typically requires a relatively large amount of computational resources, particularly when constraints are considered. This can cause problems in practice when the control processor cannot provide the necessary computational resources.

[0045] FIG. 5A shows a hybrid cascade controller having an outer loop controller 501 in an outer control loop and an inner loop controller 502 in an inner control loop. The outer loop controller 501 provides, to a non-inverting input of an adder 530, a control setpoint for the inner loop controller 502. An output of the adder 530 is an error signal that is provided to a control input of the inner loop controller 502. A control output of the inner loop controller 502 is provided to a plant 509. A sensor output 505 from one or more inner loop sensors of the plant 509 is provided to an input of a filter 507. An output from the filter 507 is provided to an inverting input of the adder 530. A sensor output 506 from one or more outer loop sensors of the plant 509 is provided to an input of a second filter 508. An output of the second filter 508 is provided to an input of the outer loop controller 501. The filters 507 and 508 can be omitted. In one embodiment, the inner loop sensors (corresponding to the sensor output 505) tend to measure one or more operating parameters of a first portion 503 of the plant 509. In one embodiment, outer loop sensors (corresponding to the sensor output 506) tend to measure one or more operating parameters of a second portion 504 of the plant 509. In one embodiment, the first portion 503 responds to changes in the control output of the inner loop controller relatively more rapidly than the second portion 504.

[0046] In one embodiment, the sensor output 505 corresponds to sensors that tend to respond relatively more quickly but relatively less accurately to certain desired parameters than the sensors corresponding to the sensor output 506. In one embodiment, the sensor output 506 corresponds to sensors that tend to respond relatively less quickly but relatively more accurately to certain desired parameters than the sensors corresponding to the sensor output 505. Thus, in one embodiment, the inner loop controller 502, using the sensor output 505 is able to respond relatively quickly to certain changes in the plant 509 but relatively less accurately. The outer loop controller 501, using the sensor output 506 is able to respond relatively less quickly to certain changes in the plant 509 but relatively more accurately to certain desired parameters. In one embodiment, the outer loop controller 501 is configured to produce a setpoint for the inner loop controller 502 to improve the controlled properties of the plant. In one embodiment the inner loop controller 502 includes a conventional controller, such as a PID controller. In one embodiment, the outer loop controller includes a predictive-type controller. In one embodiment, the outer loop controller includes a MBPC controller.

[0047] FIG. 5B shows a hybrid cascade thermal process plant controller 500 that uses a MBPC 520 in an outer control loop and a conventional controller 540 in an inner control loop. The cascade controller 500 is one embodiment of the hybrid cascade control system shown in FIG. 5A. The MBPC 520 receives as input a paddle control setpoint temperature Pdset and the actual paddle temperature Pd received via feedback loop 580. In one embodiment, the feedback loop 580 includes a digital filter 581. The MBPC 520 receives as input for the model calculations the actual spike temperatures via an input 572. The MBPC 520 calculates as output a spike TC control setpoint Spset. An adder 530 calculates a spike error signal Es_Sp by subtracting the actual spike TC temperature Sp, received via feed-back loop 570, provided by a digital filter 571, from the spike TC control setpoint. The spike error signal Es_Sp is provided to the conventional controller 540 which generates a power output signal Pw that is provided via a thyristor unit, not shown, to heating elements 550 to heat a tube and wafers 560. The sampling time ts1 in the inner control loop is preferably shorter than the sampling time ts2 in the outer MBPC control loop. In one embodiment, ts1 is in the range of 1 second and ts2 is in the range of 4 seconds. The conventional controller 540 can be implemented using a PID controller, an H controller, etc.

[0048] As the name implies, the MBPC 520 is based on a model of the controlled plant. In a thermal process reactor, the MBPC 520 typically relies on several models corresponding to the different thermal zones of the vertical thermal reactor. The simplest model is a static model, describing the relation between spike TC temperature and paddle TC temperatures under steady-state conditions according to Equation (3). In one embodiment, the static model is based on 4^th order polynomial models representing the relation between spike TC and paddle TC temperatures over a specified temperature range. Dynamic models describe the dynamic response of the system. A dynamic paddle model gives the paddle TC temperature as a function of spike TC temperature according to Equation (1) and a dynamic spike model gives the spike TC temperature as a function of power output and paddle TC temperature, according to Equation (2). By dividing the temperature range in a plurality of temperature sub-ranges, a set of linear dynamic models can be obtained. This simplifies the required calculations. The various models are acquired experimentally from a measurement procedure as described below.

[0049] FIG. 6 shows a hybrid control system 501. The control system 501 is an embodiment of the control system 500. In the control system 501, the conventional controller 540 is based on a PID controller 542.

[0050] FIG. 7 shows a vertical thermal reactor system 700 where the vertical thermal reactor 100 is controlled by the hybrid cascade control system 501. In the system 700, the process tube 110 is surrounded by the heating element 120, comprising multiple zone electric heating coils. Each zone has the spike TC 130 and the “profile” or paddle TC 140. The spike TC is located outside the process tube relatively near the heating element and the paddle TC is located inside the tube relatively near the wafers. A paddle control setpoint Pd_set and the actual paddle temperatures Pd are provided to an MBPC controller 720 (corresponding to an embodiment of the MBPC controller 520), which generates a spike control setpoint Sp_set. An adder 730 computes a spike TC error signal using the spike control setpoint Sp_set and the actual spike temperatures, provided to the adder 730 via an inverter 732. A PID controller 740 uses the spike error signal to generate a power output signal that is provided to a power actuator 750 to provide power to control the heating element 120.

1. Experiment Design and Data Acquisition for Model Identification

[0051] The identification test design plays an important role in a successful model identification and MBPC design. Current practices of identification methods are to use single variable step or finite impulse tests for MBPC model identification. The tests are carried out manually. The advantages of these methods are that the system dynamic responses are described in an intuitive manner. One drawback with the step or finite impulse tests is that the data from these tests may not contain enough information about the multivariable characteristics of the process because the step or pulse signals may not induce enough dynamic behavior of the process. A second drawback of these step or impulse tests is that they can be very time consuming. To avoid the above-mentioned drawbacks of common identification methods, during the identification procedure for the MBPC 720, the inner-loop PID control loop 740 is actively used. The PID constants used during the model identification are based on previous control experiences of the vertical thermal reactor. Then, using this inner-loop PID controller 740, the system identification can be carried out automatically by using a model identification and data acquisition recipe. The signals inducing dynamic behavior of the system are real control signals, and the conditions are similar to real process conditions. In this case, the PID controller 740 helps to keep the spike TC temperatures within their limits. The models based on these inner closed-loop data enhance the performance and stability margins of the system 700.

[0052] In one embodiment, the modeling data acquisition is achieved by using the control scheme shown in FIG. 8. In FIG. 8, a temperature control setpoint for model identification is provided to a trajectory planning unit 800, which creates a time-dependent spike TC control setpoint temperature Sp_set such as, for example, a controlled ramp-up rate. An adder 810 calculates a spike error signal Es, using the spike control setpoint Sp_set and the actual spike temperatures, received via a feedback loop 850, provided with a digital filter 851. A PID controller 820 uses the spike error signal Es to calculate a power output Pw signal that is provided via a thyristor unit, not shown, into the heating elements 830 to heat a tube and wafers 840.

[0053] An example identification control process sequence (i.e., a “recipe”) for a vertical thermal reactor starts at room temperature and ramps up at a ramp rate (that is, a time rate of change) of 10° C./min, stabilizing the spike temperature for approximately 45 minutes at 200, 400, 600, 800 and finally 1000° C., using a PID controller in the configuration of FIG.8. In order to prevent slip, the ramp rate during the last step is 5° C./min. During execution of the model identification and data acquisition recipe, the resulting actual spike and paddle temperature signals (shown in FIGS. 9-11, respectively) and the power output signals are recorded and stored for modeling.

[0054] The data collected from the recipe is divided into five data subsets, corresponding to the five stabilization temperatures of the model identification and data acquisition recipe. Each subset starts at the beginning of a ramp up period and ends just before the beginning of the next ramp up period.

[0055] In contrast to prior art model identification methods where data is acquired in open loop control, in the present invention the data is acquired in closed loop control as shown in FIGS. 7 and 8. Therefore, the collected data is typically free of spurious data, so little or no data pre-processing is required to remove spurious data points.

2. Model Structure

2.1 Static Models

[0056] For each zone, a static model is derived. For each stabilization temperature in the identification recipe, at least one value is extracted both for the spike TC temperature and the paddle TC temperature. These pairs of values, 5 in the example shown in FIGS. 9-11, are used to estimate the parameters for the static model using a parameter estimation technique, such as, for example, a polynomial fit, least-squares fit, etc. In one embodiment, the model equation for the static model is: 1

[0057] where n is the zone number, s_nq are the static models parameters to be determined, and q is the order of the static model. The model according to equation (4) gives an adequate description of the relationship between the spike temperature Sp and the paddle temperature Pd over the desired temperature range.

2.2 Dynamic Models

[0058] For each zone, a dynamic model is derived from each data subset for the spike TC and the paddle TC. This means that in the case of five zones and five temperature sub-ranges, according to the present example, 25 dynamic paddle linear models and 25 dynamic spike linear models are used. It will be clear to one of ordinary skill in the art that any number of thermal zones and temperature sub-ranges can be selected, depending on the circumstances. In one embodiment, a linear least-squares algorithm is used. The model equations for the dynamic linear models are: 2

[0059] where n is the zone number, r is the temperature range number, t is the discrete time index, a_nrl and b_nrm are paddle model parameters, l and m are model orders of the paddle models, and c_nrx, d_nry and p_nrz are spike model parameters, Further, x, y, and z are model orders of the spike models, and e_nrp and e_nrs are model errors or disturbances. Typically, a first or second order approximation (m=1 or 2) is adequate for the model of equation (5) whereas the order for the model of equation (6) is typically 20 or more (e.g., l=28) for an adequate description.

[0060] Model validation can be provided by visual comparison of the measured and calculated model output, simulation, residual analysis, cross-correlation analysis, etc.

3. Mathematics of Model Extraction

[0061] The methods used for parameter estimation for both the static models and the dynamic models involve solving a linear least-squares problem (LLS). This can be done via the so-called normal equations, or via a QR-decomposition. The method via QR-decomposition typically requires more calculations, but tend to be numerically more stable. For most cases, solving the normal equations gives good results, but for higher orders it is safer to use the QR-decomposition.

[0062] The method of solving the LLS problem is described below. Given a system of equations defined as A·x=b, where A is a matrix and x and b are vectors, the linear least-squares problem is to find a vector x that minimizes

φ(x)=∥A·x−b∥ (7)

[0063] If the matrix A is non-singular (i.e., if the inverse of A exists), this problem has a unique solution.

3.1 Solving LLS via the Normal Equations

[0064] The direct method of solving a LLS problem is via the normal equations. If a function has a minimum, its derivative at that minimum is zero. Thus φ(x) has a minimum where d φ(x)/dx=0 or:

A^TA ·x−A^T·b=0 (8)

[0065] where A^T is the inverse matrix of A. This results in the so-called Gaussian normal equations:

A^TA ·x=A^T·b (9)

[0066] From these equations, x can be calculated as

x=(A< sup>TA)⁻¹A^T< highlight>b (10)

3.2 Solving LLS via OR-decomposition

[0067] Another method of solving a LLS problem is via QR-decomposition. With this algorithm the matrix A is expressed as the product of an orthogonal matrix Q and an upper triangular matrix R.

A=QR (11)

[0068] The QR-decomposition of the matrix A can be computed by calculating the Householder reflection H for each column of A. The Householder reflection H_k of the k^th column can be calculated as:

s_k= {square root}{square root over (a_k,k²+a_k+1,k²+ . . . +a_n,k ²)} (12)

[0069] 3 H_k=I−2w_k w< highlight>_k^T (14)

[0070] where a_i,j is the element of matrix A at the i^th row and j^th column, and n is the number of columns and I is the unity matrix.

[0071] The matrices Q and R can be calculated as

Q=H₁H₂ . . . H_n−1 (15)

R=H_n−1H_n−2 . . . H₁A (16)

[0072] For a matrix A of m rows and n columns, where m is greater than n, the last m-n rows of the upper-triangular matrix R are completely zero.

[0073] Once the matrix A is expressed as A=QR, the LLS problem can be written as: 4$\begin{array}{cc}\{\begin{array}{c}\mathrm{Qy}=b\\ \mathrm{Rx}=y\end{array}& \left(\mathrm{1\; 7}\right)\end{array}$

[0074] Q is orthogonal, so

y=Q^Tb (18)

[0075] With y and R known, x can be calculated. Since R is upper-triangular, this is simply done via backward substitution. By way of example, the solution is shown for a 3 by 3 matrix in the equation below. 5

[0076] This system can be solved via backward substitution starting at the bottom row.

r₃₁ ·x₁< highlight>=y₃x₁=y₃/r₃₁

r_{21< /sub>}·x ₁+r₂₂·x₂< /highlight>=y₂x₂=y₂/r₂₂−r₂₁·x1=y₂/r ₂₂−r_{21< /sub>}·y ₃/r₃₁

r_{11< /sub>}·x ₁+r₁₂·x₂< /highlight>+r₁₃·x ₃=y₁

x₃=y₁/r₁₃−r₁₁·x< sub>1−r₁₂·x₂=y₁ /r ₁₃−r< /italic>₁₁·y₃/r₃₁−r₁₂·(y₂< /highlight>/r₂₂−r₂₁ ·y₃/r ₃₁) (20)

[0077] When using the QR-decomposition for solving a LLS problem, the orthogonal matrix Q is usually not explicitly calculated. Instead, R and y are calculated in a recursive algorithm, with initial conditions

R=A

y=b (21)

[0078] Next for each column of R, the vector w is calculated. With this vector, R and y are updated according to the following formulas.

R_k= R_k− 2(w_k ^TR_k)w_k

y=y−2(< italic>w_k^Ty< /italic>)w_k (22)

[0079] where R_k is the k^th column of R. This is repeated for all columns of R. The solution of the LLS problem can then be calculated via backward substitution as described before.

4. Extraction of the Dynamic Models

[0080] The dynamic models (shown in Equations (5) and (6)) for the MBPC can be represented by the following equations:

{circumflex over (P)}d_nr(t)=φ< sub>nrp^T (t)·θnrp (23)

Âp_nr(t)=φnrs^T( t)·θnrs (24)

[0081] where {circumflex over (P)}d_nr and Ŝp_nr are model predictive outputs, and 6

[0082] With the model structure defined above, given the model orders and a set of input and output data, the parameters of the model are found by minimizing a so-called loss function. An often-used loss function is the summed squared error: 7

[0083] where N is the number of input and output samples and ε is the prediction error vector, defined as the difference between measured and predicted output:

ε_pd(t|θ)= Pd_nr(t)−{circumflex over (P)}d_nr(t|θ)= Pd_nr(t)−φ_nrp^T(t)θ_nrp (26)

ε _sp(t|θ)=Sp _nr(t)−Ŝp_nr(t |θ)=Sp< highlight>_nr(t)−φ_nrs ^T(t)θ_nrs (27)

[0084] For simplification, the matrix Φ and vector Y are used, defined as: 8