Thermal processing method for silicon wafer

Abstract

A processing temperature T.sub.S by a rapid thermal processing furnace is 1250° C. or more and 1350° C. or less, and a cooling rate R.sub.d from the processing temperature is in a range of 20° C./s or more and 150° C./s or less, and thermal processing is performed by adjusting the processing temperature T.sub.S and the cooling rate R.sub.d within a range between the upper limit P=0.00207T.sub.S.Math.R.sub.d−2.52R.sub.d+13.3 (Formula (A)) and the lower limit P=0.000548T.sub.S.Math.R.sub.d−0.605R.sub.d−0.511 (Formula (B)) of an oxygen partial pressure P in a thermal processing atmosphere.

Claims

1. A thermal processing method for a silicon wafer, comprising: determining, by a point defect simulation for predicting behavior of vacancies and interstitial silicon atoms in the silicon wafer during thermal processing, combinations each comprising a predetermined value of a processing temperature T.sub.S, a predetermined value of a cooling rate R.sub.d, and a predetermined value of an oxygen partial pressure P, and each determining a value of a width W.sub.DZ of a denuded zone; determining, based on the combinations, a regression equation relating the width W.sub.DZ of the denuded zone to the processing temperature T.sub.S, the cooling rate R.sub.d, and the oxygen partial pressure P; and determining a combination of a specific value of the processing temperature T.sub.S, a specific value of the cooling rate R.sub.d, and a specific value of the oxygen partial pressure P for obtaining a denuded zone with a desired width by using the regression equation, wherein the regression equation is represented by the following formula (C), and the width W.sub.DZ of the denuded zone is made predictable within an error range of ±5 μm
w.sub.DZ=−0.129R.sub.d+15.6P−0.0109T.sub.S.Math.R.sub.d+34.2P/R.sub.d−7.95  (C).

2. The thermal processing method for a silicon wafer according to claim 1, comprising: determining, by a point defect simulation for predicting a behavior of vacancies and interstitial silicon atoms in a silicon wafer during thermal processing, combinations each comprising a predetermined value of a processing temperature T.sub.S, a predetermined value of a cooling rate R.sub.d, and a predetermined value of an oxygen partial pressure P, and each determining a value of a width W.sub.Inter of an intermediate layer is obtained; determining, based on the combinations, a regression equation relating the width W.sub.Inter of the intermediate layer to the processing temperature T.sub.S, the cooling rate R.sub.d, and the oxygen partial pressure P; and determining a combination of a specific value of the processing temperature T.sub.S, a specific value of the cooling rate R.sub.d, and a specific value of the oxygen partial pressure P for obtaining an intermediate layer with a desired width by using the regression equation.

3. The thermal processing method for a silicon wafer according to claim 2, wherein the regression equation is represented by the following Formula (D), and the width W.sub.Inter of the intermediate layer is made predictable within an error range of ±10 μm
W.sub.Inter=−0.141T.sub.S−6.74P+0.00456T.sub.S.Math.R.sub.d+62.7P/R.sub.d+243.1  (D).

4. A thermal processing method for a silicon wafer, comprising: determining, by a point defect simulation for predicting behavior of vacancies and interstitial silicon atoms in the silicon wafer during thermal processing, combinations each comprising a predetermined value of a processing temperature T.sub.S, a predetermined value of a cooling rate R.sub.d, and a predetermined value of an oxygen partial pressure P, and each determining a value of a width W.sub.DZ of a denuded zone; determining, based on the combinations, a regression equation relating the width W.sub.DZ of the denuded zone to the processing temperature T.sub.S, the cooling rate R.sub.d, and the oxygen partial pressure P; determining a combination of a specific value of the processing temperature T.sub.S, a specific value of the cooling rate R.sub.d, and a specific value of the oxygen partial pressure P for obtaining a denuded zone with a desired width by using the regression equation; determining, by a point defect simulation for predicting a behavior of vacancies and interstitial silicon atoms in a silicon wafer during thermal processing, combinations each comprising a predetermined value of a processing temperature T.sub.S, a predetermined value of a cooling rate R.sub.d, and a predetermined value of an oxygen partial pressure P, and each determining a value of a width W.sub.Inter of an intermediate layer is obtained; determining, based on the combinations, a regression equation relating the width W.sub.Inter of the intermediate layer to the processing temperature T.sub.S, the cooling rate R.sub.d, and the oxygen partial pressure P; and determining a combination of a specific value of the processing temperature T.sub.S, a specific value of the cooling rate R.sub.d, and a specific value of the oxygen partial pressure P for obtaining an intermediate layer with a desired width by using the regression equation, wherein the regression equation is represented by the following Formula (D), and the width W.sub.Inter of the intermediate layer is made predictable within an error range of ±10 μm
W.sub.Inter=−0.141T.sub.S−6.74P+0.00456T.sub.S.Math.R.sub.d+62.7P/R.sub.d+243.1  (D).

Description

BRIEF DESCRIPTION OF THE DRAWINGS

(1) FIG. 1A is a conceptual diagram showing a wafer depth direction distribution of the VOx density of a wafer capable of forming a favorable DZ layer and bulk layer.

(2) FIG. 1B is a conceptual diagram showing a wafer depth direction distribution of the VOx density of a wafer in which oxygen precipitates may be formed in a DZ layer.

(3) FIG. 2 is a diagram showing an example of a thermal processing sequence employed in a thermal processing method of a silicon wafer in the present invention.

(4) FIG. 3A is a diagram showing the processing time dependency of the thickness of an oxide film formed when a silicon wafer is thermally processed in oxygen atmospheres at a processing temperature of 1250° C.

(5) FIG. 3B is a diagram showing the processing time dependency of the thickness of an oxide film formed when a silicon wafer is thermally processed in oxygen atmospheres at a processing temperature of 1275° C.

(6) FIG. 3C is a diagram showing the processing time dependency of the thickness of an oxide film formed when a silicon wafer is thermally processed in oxygen atmospheres at a processing temperature of 1300° C.

(7) FIG. 3D is a diagram showing the processing time dependency of the thickness of an oxide film formed when a silicon wafer is thermally processed in oxygen atmospheres at a processing temperature of 1325° C.

(8) FIG. 4 is a diagram showing simulation results of the VOx density.

(9) FIG. 5A is a diagram showing the wafer depth direction distribution of oxygen precipitates at a processing temperature of 1350° C., an oxygen partial pressure of 100%, and a cooling rate of 120° C./s.

(10) FIG. 5B is a diagram showing the wafer depth direction distribution of oxygen precipitates at a processing temperature of 1350° C., an oxygen partial pressure of 100%, and a cooling rate of 50° C./s.

(11) FIG. 5C is a diagram showing the wafer depth direction distribution of oxygen precipitates at a processing temperature of 1350° C., an oxygen partial pressure of 100%, and a cooling rate of 25° C./s.

(12) FIG. 6A is a diagram showing the wafer depth direction distribution of oxygen precipitates at a cooling rate of 25° C./s, an oxygen partial pressure of 0%, and a processing temperature of 1300° C.

(13) FIG. 6B is a diagram showing the wafer depth direction distribution of oxygen precipitates at a cooling rate of 25° C./s, an oxygen partial pressure of 0%, and a processing temperature of 1275° C.

(14) FIG. 6C is a diagram showing the wafer depth direction distribution of oxygen precipitates at a cooling rate of 25° C./s, an oxygen partial pressure of 0%, and a processing temperature of 1250° C.

(15) FIG. 7A is a diagram showing the wafer depth direction distribution of oxygen precipitates at a cooling rate of 120° C./s, an oxygen partial pressure of 0%, and a processing temperature of 1300° C.

(16) FIG. 7B is a diagram showing the wafer depth direction distribution of oxygen precipitates at a cooling rate of 120° C./s, an oxygen partial pressure of 0%, and a processing temperature of 1275° C.

(17) FIG. 7C is a diagram showing the wafer depth direction distribution of oxygen precipitates at a cooling rate of 120° C./s, an oxygen partial pressure of 0%, and a processing temperature of 1250° C.

(18) FIG. 8 is a diagram showing simulation results of the wafer depth direction distribution of the VOx density when the processing temperature was changed under five conditions ranging from 1350 to 11500° C. (at 50° C. pitch), the ramp up rate was set to 75° C./s, the cooling rate was changed under four conditions of 120, 50, 25, and 5° C./s, and the oxygen partial pressure was changed under 10 conditions of 100, 50, 30, 10, 8, 6, 4, 2, 1, and 0%.

(19) FIG. 9 is a diagram showing simulation results of the wafer depth direction distribution of the VOx density variation amount in the wafer depth direction (depth differentiation of the VOx density) in the simulation result shown in FIG. 8.

(20) FIG. 10A is a diagram showing the relationship between the cooling rate and the oxygen partial pressure in a simulation result at a processing temperature of 1250° C. for a silicon wafer having both the integrity of a DZ layer and the gettering ability of a bulk layer.

(21) FIG. 10B is a diagram showing the relationship between the cooling rate and the oxygen partial pressure in a simulation result at a processing temperature of 1300° C. for a silicon wafer having both the integrity of a DZ layer and the gettering ability of a bulk layer.

(22) FIG. 10C is a diagram showing the relationship between the cooling rate and the oxygen partial pressure in a simulation result at a processing temperature of 1350° C. for a silicon wafer having both the integrity of a DZ layer and the gettering ability of a bulk layer.

(23) FIG. 11A is a diagram showing the width of a DZ layer in each simulation result shown in FIGS. 10A to 10C when the processing temperature is 1250° C.

(24) FIG. 11B is a diagram showing the width of a DZ layer in each simulation result shown in FIGS. 10A to 10C when the processing temperature is 1300° C.

(25) FIG. 11C is a diagram showing the width of a DZ layer in each simulation result shown in FIGS. 10A to 10C when the processing temperature is 1350° C.

(26) FIG. 12A is a diagram showing the width of an intermediate layer in each simulation result shown in FIGS. 10A to 10C when the processing temperature is 1250° C.

(27) FIG. 12B is a diagram showing the width of an intermediate layer in each simulation result shown in FIGS. 10A to 10C when the processing temperature is 1300° C.

(28) FIG. 12C is a diagram showing the width of an intermediate layer in each simulation result shown in FIGS. 10A to 10C when the processing temperature is 1350° C.

(29) FIG. 13 is a diagram showing the correlation between the width of a DZ layer by a point defect simulation and the width of a DZ layer predicted by a regression equation (regression equation is prepared from a processing temperature, a cooling rate, and an oxygen partial pressure).

(30) FIG. 14 is a diagram showing the correlation between the width of an intermediate layer by a point defect simulation and the width of an intermediate layer predicted by a regression equation (regression equation is prepared from a processing temperature, a cooling rate, and an oxygen partial pressure).

EMBODIMENTS

(31) A thermal processing method for a silicon wafer (wafer) according to the present invention is a method for performing rapid thermal processing (RTP) by using a ramp annealing furnace capable of rapidly increasing and decreasing temperature in an oxygen atmosphere (oxygen partial pressure is in the range of 1% or more and 100% or less) and for introducing vacancy oxygen complexes (VOx), which are complexes of vacancies and oxygen, into the wafer. As the wafer for thermal processing, a wafer sliced from an ingot grown by the Czochralski method (hereinafter referred to as CZ method) and mirror polished on both sides thereof is used.

(32) As the RTP thermal processing, for example, a sequence shown in FIG. 2 is employed. In this sequence, the temperature is raised at a predetermined ramp up rate R.sub.u from a predetermined temperature T.sub.0 (for example, 700° C.) to a processing temperature T.sub.S, and is held for a predetermined time (for example, 15 seconds) at the processing temperature T.sub.S. Then, the temperature is decreased from the processing temperature T.sub.S to the predetermined temperature T.sub.0 at a predetermined cooling rate R.sub.d.

(33) The VOx density distribution in the wafer depth direction was calculated using a point defect simulation. The point defect simulation is a model of a series of physical phenomena comprising formation of an oxide film on the surface of the wafer by thermal processing in an oxygen atmosphere, generation of interstitial silicon atoms at the interface between the oxide film and the silicon crystal, diffusion of the interstitial silicon atoms into the wafer, and pair annihilation of the interstitial silicon atoms and vacancies. In the following, the point defect simulation will be described in detail, and a method of determining the upper and lower limits of the oxygen partial pressure or predicting the widths of the DZ layer and the intermediate layer by regression analysis from the results of this simulation will be described.

(34) (1) Regarding Point Defect Simulation

(35) (a) Formation of Oxide Film

(36) When RTP is performed in an oxygen atmosphere, an oxide film is formed as described above, and interstitial silicon atoms are injected into the wafer from the interface between the oxide film and the silicon crystal. It is known that the ratio of the concentration of interstitial silicon atoms generated at the interface between the oxide film and the silicon crystal to the thermal equilibrium concentration thereof in the silicon crystal, that is, the degree of supersaturation of interstitial silicon atoms, is a factor which greatly influences the accuracy of point defect simulation, and is known to be determined by the growth rate of the oxide film. Therefore, the growth rate of the oxide film during RTP in an oxygen atmosphere will be examined in detail.

(37) Growth of the oxide film is represented by the Deal-Grove equation (B. Deal, J. Electrochem. Soc. 125 (1978) 576) shown in Formula (1).
[Math. 5]
X.sub.0.sup.2+AX.sub.0=B(t−t.sub.0)  (1)

(38) Here, X.sub.0 is the thickness of the oxide film, t is the time, and t.sub.0 is the time expressing the initial value of the thickness of the oxide film. When the oxide film is sufficiently thick, Formula (1) is approximated as the following Formula (1′), and the thickness X.sub.0 of the oxide film becomes a parabolic function proportional to the root of time.
[Math. 6]
X.sub.2.sup.0=Bt  (1′)

(39) Here, the following Formula (2) is obtained by differentiating Formula (1).

(40) $\begin{array}{cc}\left[\mathrm{Math}.7\right]& \\ {\left(\frac{d{X}_0}{dt}\right)}^{-1}=\frac{1}{B\text{/}A}+\frac{2X_0}{B}=\frac{1}{K_L}+\frac{2X_0}{K_P}& \left(2\right)\end{array}$

(41) Here, B/A is called Linear rate constant (K.sub.L), and B is called Parabolic rate constant (K.sub.P). B/A corresponds to a reaction rate at the oxide film-wafer interface, and B relates to the thickness of the oxide film, and is a term corresponding to a diffusion effect in the oxide film. Although, inherently, Formula (2) is an equation to be applied to the entire oxidation process including the process of increasing the oxide film thickness, Kotani et al. showed that Formula (2) itself is usable also for the initial process of oxidation in which only a thin oxide film is formed such as during RTP by replacing the values of K.sub.L and K.sub.P with values for RTP (Y. Kotani and Y. Omura, Jpn. J. Appl. Phys. 39 (2000) 4549).

(42) Since data during RTP in an oxygen atmosphere obtained by Kotani et al. is up to 1100° C., for the temperature region higher than 1100° C., thickness data of the oxide film when changing the processing temperature, the processing time, and the oxygen partial pressure of RTP was obtained by an experiment. Plots shown in FIGS. 3A to 3D are data obtained by such an experiment. The oxidation rate was calculated from these data and K.sub.L and K.sub.P were determined as a function of temperature and oxygen partial pressure as shown in the following Formulae (3) and (4).
[Math. 8]
K.sub.L=629P.sup.0.6 exp(−0.5437 eV/kT)(A/s)  (3)
[Math. 9]
K.sub.P=1.01×10.sup.11P.sup.0.6 exp(−2.347 eV/kT)(A.sup.2/s)  (4)

(43) Here, P represents the oxygen partial pressure during RTP in an oxygen atmosphere, and T represents the absolute temperature. Calculated values of the thickness of the oxide film obtained from Formulae (2), (3), and (4) are indicated by lines in FIGS. 3A to 3D. These calculated values and experimental data are in good agreement with each other, indicating that the oxidation rate at any processing temperature and oxygen partial pressure can be accurately estimated from the values of K.sub.L and K.sub.P determined by the processing temperature and oxygen partial pressure.

(44) (b) Calculation of VOx Density

(45) The density distribution of VOx formed in the wafer was calculated by the following procedure.

(46) The growth rate of the oxide film formed on the wafer surface can be expressed by the following Formula (2′) from Formula (2).

(47) $\begin{array}{cc}\left[\mathrm{Math}.10\right]& \\ \frac{d{X}_0}{dt}=\frac{1}{1\text{/}{K}_L+2X_0\text{/}{K}_P}& \left(2^{\prime }\right)\end{array}$

(48) The thickness of the oxide film at any time during the RTP process was determined by integrating Formula (2′) over the entire process of RTP while increasing the temperature, while holding the processing temperature, and while decreasing the temperature. In this integration, the film thickness of a natural oxide film was defined as the initial thickness of the oxide film. The values of Formulae (3) and (4) were employed as values of K.sub.L and K.sub.P in Formula (2′). From these, the growth rate of the oxide film at an arbitrary time was obtained using Formula (2′). Using the obtained growth rate of the oxide film, the supersaturation degree of interstitial silicon atoms at the interface between the oxide film and the wafer was obtained by the following Formula (5) proposed by Taniguchi et al. (K. Taniguchi et al., J. Appl. Phys. 65 (1989) 2723).

(49) $\begin{array}{cc}\left[\mathrm{Math}.11\right]& \\ \frac{C_I}{C_I^{\mathrm{eq}}}=\frac{A1\left(d{X}_0\text{/}dt\right)+A2}{{\left(d{X}_0\text{/}dt\right)}^{0.5}+A2}& \left(5\right)\end{array}$

(50) Here, C.sub.I is the concentration of the interstitial silicon atoms, C.sub.I.sup.eq is the thermal equilibrium concentration of the interstitial silicon atoms, and A1 and A2 are parameters related to the oxidation phenomenon. A1 and A2 in a dry oxygen atmosphere were determined to fit to the relationship between the supersaturation degree of interstitial silicon atoms and the growth rate of the oxide film determined experimentally by Dunhum (S. T. Dunhum, J. Appl. Phys. 71 (1992) 685). The density of interstitial silicon atoms generated at the oxidation interface is obtained by multiplying the degree of supersaturation of interstitial silicon atoms obtained from Formula (5) by the thermal equilibrium concentration of interstitial silicon atoms described below.

(51) The vacancy concentration C.sub.V at the interface between the oxide film and the silicon crystal was determined from the following Formula (6) representing a local equilibrium relationship assuming that the relationship of generation/pair annihilation of vacancies V and interstitial silicon atoms I is in a steady state balance.
[Math. 12]
C.sub.VC.sub.I=C.sub.V.sup.eqC.sub.I.sup.eq  (6)

(52) Here, C.sub.V is the concentration of vacancies, and C.sub.V.sup.eq is the thermal equilibrium concentration of vacancies. As can be seen from Formula (6), the C.sub.V at the oxide film-Si crystal interface is determined by the C.sub.I at this interface. The point defect concentration in the wafer is determined by solving the following Formulae (7) and (8) with the C.sub.V and the C.sub.I at the interface as boundary conditions.

(53) $\begin{array}{cc}\left[\mathrm{Math}.13\right]& \\ \frac{\partial C_V}{\partial t}=D_V\frac{{\partial }^2{C}_V}{\partial Z^2}-K_{\mathrm{IV}}\left(C_V{C}_I-C_V^{\mathrm{eq}}{C}_I^{\mathrm{eq}}\right)& \left(7\right)\\ \left[\mathrm{Math}.14\right]& \\ \frac{\partial C_I}{\partial t}=D_I\frac{{\partial }^2{C}_I}{\partial Z^2}-K_{\mathrm{IV}}\left(C_V{C}_I-C_V^{\mathrm{eq}}{C}_I^{\mathrm{eq}}\right)& \left(8\right)\end{array}$

(54) Here, D.sub.V is the diffusion coefficient of vacancies, D.sub.I is the diffusion coefficient of the interstitial silicon atoms, and K.sub.IV is the reaction constant of the pair annihilation reaction of vacancies and interstitial silicon atoms. This reaction constant K can be represented by the following Formula (9).
[Math. 15]
K.sub.IV=4πa.sub.c(D.sub.V+D.sub.I)exp(−ΔG.sub.IV/kT)  (9)

(55) Here, a.sub.c is the critical distance at which the vacancies V and the interstitial silicon atoms I cause pair annihilation, and ΔG.sub.IV is the barrier energy of annihilation. In this calculation, the same numerical values as used by Kissinger et al., or a.sub.c=0.543 nm and ΔG.sub.IV=0 were employed (G. Kissinger et al., J. Electrochem. Soc. 154 (2007) H4547). For the diffusion coefficient and thermal equilibrium concentrations of vacancies and interstitial silicon atoms, the values used by Nakamura et al. shown in Formulae (10)-(13) were employed (K. Nakamura et al., in Semiconductor Silicon, PV2002-2, p. 554, The Electrochemical Society Proceedings Series, Pennington, N.J. (2002)).
[Math. 16]
D.sub.V=3.5×10.sup.−4 exp(−0.3 eV/kT)(cm.sup.2/s)  (10)
[Math. 17]
D.sub.I=2.45×10.sup.−1 exp(−0.9 eV/kT)(cm.sup.2/s)  (11)
[Math. 18]
C.sub.V.sup.eq=3.85×10.sup.26 exp(−3.94 eV/kT)(cm.sup.−3)  (12)
[Math. 19]
C.sub.I.sup.eq=6.12×10.sup.26 exp(−4.05 eV/kT)(cm.sup.−3)  (13)

(56) Next, handling of oxygen atoms and nitrogen atoms in a silicon single crystal will be described. In a silicon single crystal manufactured by the CZ method, a quartz crucible holding a silicon melt at the time of pulling melts into a silicon melt, and oxygen atoms are taken in. Such oxygen atoms are transported to just below the growth interface of the silicon single crystal due to convection of the melt and diffusion in the melt and is mixed into the silicon single crystal from the growth interface. Nitrogen atoms are intentionally added to the silicon melt as a silicon nitride film or the like in order to suppress the size of void defects generated in the silicon single crystal. Normally, these atoms diffuse in the silicon single crystal based on the following Formulae (14) and (15).

(57) $\begin{array}{cc}\left[\mathrm{Math}.20\right]& \\ \frac{\partial C_O}{\partial t}=D_O\frac{{\partial }^2{C}_O}{\partial Z^2}& \left(14\right)\\ \left[\mathrm{Math}.21\right]& \\ \frac{\partial C_N}{\partial t}=D_N\frac{{\partial }^2{C}_N}{\partial Z^2}& \left(15\right)\end{array}$

(58) Here, C.sub.O is the concentration of oxygen atoms and C.sub.N is the concentration of nitrogen atoms. D.sub.O is the diffusion coefficient of oxygen, and D.sub.N is the diffusion coefficient of nitrogen, which are given by the following Formulae (16) and (17), respectively, by Mikkelesen et al. (J. C. Mikkelsen, Jr., in Oxygen, Carbon, Hydrogen and Nitrogen in Crystalline Silicon, ed. J. C. Mikkelsen, Jr. et al. (MRS, Princeton, N.J., 1986) 19), and Y. Itoh et al. (Y. Itoh, T Abe, Appl Phys. Lett.) vol53 (1988) 39).
[Math. 22]
D.sub.O=0.13 exp(−2.53 eV/kT)(cm.sup.2/s)  (16)
[Math. 23]
D.sub.N=2.7×10.sup.3 exp(−2.8 eV/kT)(cm.sup.2/s)  (17)

(59) It is assumed that the concentrations of oxygen and nitrogen atoms at the interface between the oxide film and silicon crystal have thermal equilibrium concentrations represented by the following Formulae (18) and (19), respectively.

(60) $\begin{array}{cc}\left[\mathrm{Math}.24\right]& \\ C_O^{\mathrm{eq}}=2.0\times {10}^{21}\exp \left(-1.03eV\text{/}kT\right)\left({\mathrm{cm}}^{-3}\right)& \left(18\right)\\ \left[\mathrm{Math}.25\right]& \\ C_N^{\mathrm{eq}}=4.5\times {10}^{15}\exp \left(-2.15eV\text{/}k\left(\frac{1}{T}-\frac{1}{1685}\right)\right)\left({\mathrm{cm}}^{-3}\right)& \left(19\right)\end{array}$

(61) Here, C.sub.O.sup.eq is the thermal equilibrium concentration of oxygen and C.sub.N.sup.eq is the thermal equilibrium concentration of nitrogen. As the thermal equilibrium concentration of nitrogen atoms, an energy value of 2.15 eV employed by Voronkov et al. was used (V. V. Voronkov and R. Falster, J. Electrochem. Soc. 149 (2002) G167).

(62) Nitrogen atoms and vacancies in the silicon single crystal easily form complexes (hereinafter referred to as NV complexes) during the cooling process. This density is given by the following Formula (20) according to Voronkov et al.
[Math. 26]
[NV]=C.sub.V/C.sub.V.sup.eq[N.sub.2].sup.1/2 exp(−EneV/kT)(cm.sup.−3)  (20)

(63) Here, the value in [ ] in Formula (20) represents the atomic fraction, and [NV]=C.sub.NV/ρ, and [N.sub.2]=C.sub.N2/ρ. ρ represents the atomic density of the silicon single crystal, and ρ=5×10.sup.22/Cm.sup.3. In this calculation, it is assumed that the NV complexes do not diffuse.

(64) A steady state balance is always established between the concentrations of vacancies V, nitrogen N.sub.2, and NV complexes in the silicon single crystal. In other words, assuming that the total vacancy concentration is C.sub.V.sup.T and its thermal equilibrium concentration is C.sub.V.sup.T.Math.eq, relationships of the following Formulae (21) and (22) are established.
[Math. 27]
C.sub.V.sup.T=C.sub.V+C.sub.NV  (21)
[Math. 28]
C.sub.V.sup.T.Math.eq=C.sub.V.sup.eq+C.sub.NV.sup.eq  (22)

(65) Here, the vacancy concentration C.sub.V and the NV complex concentration C.sub.NV are determined from the total vacancy concentration C.sub.V.sup.T in each element space by using the following Formulae (23) and (24).
[Math. 29]
C.sub.V=(C.sub.V.sup.eq/C.sub.V.sup.T.Math.eq)C.sub.V.sup.T  (23)
[Math. 30]
C.sub.NV=(C.sub.V.sup.eq/C.sub.V.sup.T.Math.eq)C.sub.V.sup.T  (24)

(66) Here, when Formula (23) is modified, the following Formula (25) is obtained.
[Math. 31]
C.sub.V/C.sub.V.sup.eq=C.sub.V.sup.T/C.sub.V.sup.T.Math.eq  (25)

(67) By substituting Formula (25) into Formula (24), the following Formula (26) is obtained.
[Math. 32]
C.sub.NV=C.sub.NV.sup.eq(C.sub.V/C.sub.V.sup.eq)  (26)

(68) Since C.sub.NV.sup.eq is represented by the following Formula (27), Formulae (20) and (26) coincide.
[Math. 33]
C.sub.NV.sup.eq=ρ(C.sub.N2/ρ).sup.1/2 exp(−E.sub.n/kT)  (27)

(69) The Formula (27) represents the thermal equilibrium concentration of the NV complexes when the vacancy concentration is the thermal equilibrium concentration.

(70) In a point defect simulation after RTP in an oxygen atmosphere, first, as the initial state, thermal equilibrium concentrations of vacancies and interstitial silicon atoms corresponding to the initial temperature of RTP are given over the entire depth direction of the wafer. Then, the thermal equilibrium concentrations given by Formulae (5) and (6) are changed as boundary conditions between the oxide film and the silicon crystal, corresponding to each temperature change while the ramp up, while the duration at maximum temperature, and while the cooling in RTP. Along with this, the wafer depth direction distribution of vacancies and interstitial silicon atoms is derived by numerical calculation by using Formulae (7)-(9).

(71) As in the calculation of vacancies and interstitial silicon atoms, firstly, as an initial state, the concentration of oxygen atoms and the concentration of nitrogen atoms corresponding to the initial temperature of RTP are given over the entire depth direction of the wafer. Then, the thermal equilibrium concentrations given in Formulae (18) and (19) are changed as boundary conditions of the concentrations of oxygen atoms and nitrogen atoms between the oxide film and the silicon crystal corresponding to each temperature change while the ramp up, while the duration at maximum temperature, and while the cooling in RTP. Along with this, the wafer depth direction distribution of oxygen and nitrogen is derived by numerical calculation using Formulae (14)-(17).

(72) Regarding the interaction between nitrogen atoms and vacancies, the concentration ([NV]) of the complexes of nitrogen and vacancies is calculated using Formulae (20)-(27). In calculating the vacancy concentration, considering that the nitrogen-vacancy complexes (NV) rapidly diverge (NV.fwdarw.N+V) due to the decrease of independently present single vacancies (V), and becomes a source of single vacancies, the total vacancy concentration was calculated as the sum ([V]+[NV]) of the single vacancy concentration ([V]) and the nitrogen-vacancy complex concentration ([NV]).

(73) During the cooling step of RTP, more specifically at 1000° C. or lower, it was assumed that the total vacancies rapidly react with oxygen atoms present in the wafer at a concentration several orders of magnitude higher than the vacancies, and change to VOx's which are most stable at a low temperature. In other words, in this point defect simulation, it was assumed that the total vacancy concentration and the VOx density correspond one-to-one.

(74) FIG. 4 shows an example of the simulation result of the wafer depth direction distribution of VOx's when the processing temperature, the cooling rate, and the oxygen partial pressure during RTP were varied. The region where the VOx density is less than 1.0×10.sup.12/cm.sup.3 corresponds to the DZ layer, the region where the VOx density is 5.0×10.sup.12/cm.sup.3 or more corresponds to the bulk layer, and the region where the VOx density is 1.0×10.sup.12/cm.sup.3 or more and less than 5.0×10.sup.12/cm.sup.3 corresponds to the intermediate layer.

(75) From the simulation results marked with V in FIG. 4, it is apparent that the higher the cooling rate, the shallower the width of the DZ layer, and the more steeply the VOx density rises. This is because as the cooling rate increases, a sufficient time is not given for vacancies to outwardly diffuse to the front side of the wafer during RTP, and the region having a high VOx density remains near the surface of the wafer.

(76) From the simulation results marked with VI and VII in FIG. 4, it was found that the higher the processing temperature, the higher the VOx density. This is because the higher the temperature, the higher the thermal equilibrium density of vacancies and more VOx's are frozen in the wafer by RTP.

(77) Furthermore, from the comparison of the simulation results marked with V and the simulation results marked with VI and VII in FIG. 4, it was found that the width of the DZ layer was wider when the oxygen atmosphere (here, the oxygen partial pressure is 100%) was used than when the atmosphere containing no oxygen (0% oxygen partial pressure) was used. This is because in the oxygen atmosphere, as described above, interstitial silicon atoms are injected into the wafer surface layer from the interface between the oxide film and the silicon crystal, and pair annihilation of the vacancies and the interstitial silicon atoms decreases the VOx density.

(78) After the RTP was actually performed under the thermal processing conditions described above, thermal processing was performed for revealing oxygen precipitates at 780° C. for three hours followed by 1000° C. for 16 hours. FIGS. 5A to 7C show measurement results of the wafer depth direction distribution of oxygen precipitates using an infrared scattering tomography apparatus (MO441; manufactured by Raytex Corporation). FIGS. 5A to 5C correspond to the simulation result marked with V in FIG. 4, FIGS. 6A to 6C correspond to the simulation result marked with VI in FIG. 4, and FIGS. 7A to 7C correspond to the simulation results marked with VII in FIG. 4. From the comparison of these simulation results and actual measurement results of the wafer depth direction distribution of oxygen precipitates, it was confirmed that the wafer depth direction distribution of oxygen precipitates can be accurately controlled by adjusting the depth direction distribution of the VOx density according to the RTP conditions.

(79) A semiconductor device was formed on each of wafers shown in FIGS. 5A to 5C subjected to RTP under conditions of the cooling rates of 120° C./s (see FIG. 5A), 50° C./s (see FIG. 5B), and 25° C./s (see FIG. 5C), and the electrical characteristics of each semiconductor device were evaluated. As a result, the electrical characteristics of a semiconductor device formed in a region close to the surface (region shallower than about 60 μm) of any of the wafers were considerably favorable. In contrast, when a semiconductor device was formed in a region deeper than about 60 μm of a wafer subjected to RTP at a cooling rate of 120° C./s, when a semiconductor device was formed in a region deeper than about 90 μm of a wafer subjected to RTP at a cooling rate of 50° C./s, and when a semiconductor device was formed in a region (region having a VOx density of 1.0×10′.sup.2/cm.sup.3 or more and denoted as “NG” in FIG. 4) deeper than about 130 μm of a wafer subjected to RTP at a cooling rate of 25° C./s, it was confirmed that the electrical characteristics of the semiconductor device tended to be inferior.

(80) Further, FIG. 8 is a diagram showing simulation results of the wafer depth direction distribution of the VOx density when the processing temperature was changed under five conditions ranging from 1350 to 1150° C. (at 50° C. pitch), the ramp up rate was set to 75° C./s, the cooling rate was changed under four conditions of 120, 50, 25, and 5° C./s, and the oxygen partial pressure was changed under 10 conditions of 100, 50, 30, 10, 8, 6, 4, 2, 1, and 0% (under total of 200 conditions combining all conditions), and FIG. 9 is a diagram showing simulation results of the wafer depth direction distribution of the VOx density variation amount in the wafer depth direction (depth differentiation of the VOx's) in the simulation of FIG. 8.

(81) According to the simulation result of FIG. 8, as described in FIG. 4, under conditions which the cooling rate and the oxygen partial pressure were fixed, although there were some exceptions, the higher the processing temperature, the higher the VOx density. Under conditions which the processing temperature and oxygen partial pressure were fixed, the higher the cooling rate, the higher the VOx density. Further, under conditions which the processing temperature and cooling rate were fixed, the higher the oxygen partial pressure, the lower the VOx density.

(82) As described above, when the region having a VOx density of less than 1.0×10.sup.12/cm.sup.3 is defined as the DZ layer, the region having a VOx density of 1.0×10.sup.12/cm.sup.3 or more and less than 5.0×10.sup.12/cm.sup.3 is defined as the intermediate layer, and the region having a VOx density of 5.0×10.sup.12/cm.sup.3 or more is defined as the bulk layer, whether or not each of these layers is formed by RTP and the thickness of each layer formed by RTP largely change depending on the processing conditions of RTP.

(83) As shown in FIG. 8, for example, when RTP is performed under the conditions of a processing temperature of 1350° C., a ramp up rate of 75° C./s, a cooling rate of 120° C./s, and an oxygen partial pressure of 30%, a DZ layer with a thickness of about 15 μm is formed on the wafer, and an intermediate layer with a thickness of about 50 μm is formed on the wafer. When RTP is performed under conditions of a processing temperature of 1350° C., a ramp up rate of 75° C./s, a cooling rate of 50° C./s, and an oxygen partial pressure of 10%, a DZ layer of about 4 μm is formed on the wafer and an intermediate layer having a thickness of about 59 μm is formed on the wafer.

(84) In contrast, for example, when RTP is performed on a wafer under conditions of a processing temperature of 1350° C., a ramp up rate of 75° C./s, a cooling rate of 120° C./s, and an oxygen partial pressure of 10%, no DZ layer is formed, and when RTP is performed on a wafer under conditions of a processing temperature of 1350° C., a ramp up rate of 75° C./s, a cooling rate of 5° C./s, and an oxygen partial pressure of 30%, no intermediate layer and no bulk layer are formed, in any of the above cases, the wafer is unusable as a substrate for semiconductor devices.

(85) As described above, by performing simulation using the RTP conditions (processing temperature, cooling rate, oxygen partial pressure) as parameters, it is possible to conveniently derive the wafer depth direction distribution of the VOx density which is very difficult to directly evaluate, and it is possible to conveniently determine processing conditions of a wafer having a DZ layer and an intermediate layer having desired thicknesses.

(86) According to the simulation result of FIG. 9, under conditions which the cooling rate and oxygen partial pressure were fixed, although there are some exceptions, the higher the processing temperature, the higher the maximum value of the amount of change (differential value) of the VOx density in the wafer depth direction, and the position of the maximum value moved to the vicinity of the surface of the wafer. Under conditions which the processing temperature and oxygen partial pressure were fixed, the higher the cooling rate, the higher the maximum value of the differential value, and the maximum value moved to the vicinity of the surface of the wafer. Further, under conditions which the processing temperature and the cooling rate were fixed, when the oxygen partial pressure was around 10%, the maximum value of the differential value was the largest, and when the oxygen partial pressure was lower or higher than about 10%, the maximum value was lower. As the oxygen partial pressure increased, the maximum value moved toward the center of the thickness of the wafer.

(87) The rise of the VOx density in the intermediate layer becomes steeper as the maximum value of the differential value becomes larger, so that the DZ layer and the bulk layer can be brought closer to each other via the intermediate layer. In this way, by bringing these layers close to each other, the diffusion length of heavy metals can be shortened, and gettering can be performed efficiently. By setting the position of the maximum value near the surface of the wafer, it is possible to efficiently getter in the same manner as described above even for heavy metals such as molybdenum, tungsten, and cobalt which have a small diffusion coefficient and are difficult to getter.

(88) (2) Regarding Regression Analysis

(89) (a) Determination of Upper Limit and Lower Limit of Oxygen Partial Pressure

(90) From the simulation results shown in FIG. 8, data of a wafer including both a DZ layer with high integrity and a bulk layer with high gettering ability or data of a wafer thermally processed at RTP conditions (processing temperature, cooling rate, and oxygen partial pressure) in which the VOx density of the DZ layer is less than 1.0×10.sup.12/cm.sup.3 and the VOx density of the bulk layer is 5.0×10.sup.12/cm.sup.3 or more was picked up and plotted in FIGS. 10A to 10C. FIG. 10A shows the results when the processing temperature was 1250° C., FIG. 10B shows the results when the processing temperature was 1300° C., and FIG. 10C shows the results when the processing temperature was 1350° C.

(91) Next, a regression equation composed of three parameters of processing temperature T.sub.S, cooling rate R.sub.d, and oxygen partial pressure P was derived in such a manner that the picked-up RTP conditions were included. This regression equation is composed of the equation (A) representing the upper limit of the oxygen partial pressure P and the equation (B) showing the lower limit of the oxygen partial pressure P.
[Math. 34]
P=0.00207T.sub.S.Math.R.sub.d−2.52R.sub.d+13.3  (A)
[Math. 35]
P=0.000548T.sub.S.Math.R.sub.d−0.605R.sub.d−0.511  (B)

(92) By determining the set of three parameters of the processing temperature, the cooling rate, and the oxygen partial pressure so as to satisfy the upper limit equation and the lower limit equation, it is possible to conveniently determine the thermal processing conditions of a wafer including both a DZ layer with high integrity and a bulk layer with high gettering ability.

(93) The simulation results plotted in FIGS. 10A to 10C include both information on the width of the DZ layer and information on the width of the intermediate layer. Therefore, the information on the width of the DZ layer is described in FIGS. 11A to 11C, separately from the information on the width of the intermediate layer, which is described in FIGS. 12A to 12C. The size of the circle of each plot in FIGS. 11A to 12C corresponds to the size of the width of the corresponding DZ layer (FIGS. 11A to 11C) and the size of the width of the intermediate layer (FIGS. 12A to 12C).

(94) The plot of the maximum size in FIGS. 11A to 11C indicates that the width of the DZ layer is about 30 μm or more, and the plot of the minimum size indicates that the width of the DZ layer is about 3 μm or less. From this simulation result, it is apparent that the larger the cooling rate, the lower the oxygen partial pressure, or the higher the processing temperature, the smaller the width of the DZ layer.

(95) The plot of the maximum size in FIGS. 12A to 12C indicates that the width of the intermediate layer is about 91 μm or more, and the plot of the minimum size indicates that the width of the intermediate layer is about 46 μm or less. From this simulation result, it is apparent that the larger the cooling rate, or the higher the processing temperature, the smaller the width of the intermediate layer. It is also found that the oxygen partial pressure does not significantly affect the magnitude of the width of the intermediate layer.

(96) (b) Prediction of Widths of DZ Layer and Intermediate Layer

(97) From the simulation results shown in FIG. 8, data of a wafer including both a DZ layer with high integrity and a bulk layer with high gettering ability or data of a wafer thermally processed at RTP conditions (processing temperature, cooling rate, and oxygen partial pressure) in which the VOx density of the DZ layer is less than 1.0×10.sup.12/cm.sup.3 and the VOx density of the bulk layer is 5.0×10.sup.12/cm.sup.3 or more was picked up, a regression equation was derived from these picked up data, and using the regression equation, the widths of the DZ layer and the intermediate layer predicted from the RTP conditions were derived.

(98) FIG. 13 is a diagram showing the correlation between the width of a DZ layer by a point defect simulation and the width of a DZ layer predicted by a regression equation, and FIG. 14 is a diagram showing the correlation between the width of an intermediate layer by a point defect simulation and the width of an intermediate layer predicted by a regression equation. The regression equation of the width W.sub.DZ of the DZ layer is represented by Formula (C), and the regression equation of the width W.sub.Inter of the intermediate layer is represented by Formula (D).
[Math. 36]
w.sub.DZ=−0.129R.sub.d+15.6P−0.0109T.sub.S.Math.R.sub.d+34.2P/R.sub.d−7.95  (C)
[Math. 37]
W.sub.Inter=−0.141T.sub.S−6.74P+0.00456T.sub.S.Math.R.sub.d+62.7P/R.sub.d+243.1  (D)

(99) By substituting the processing temperature T.sub.S, the cooling rate R.sub.d, and the oxygen partial pressure P into Formulae (C) and (D), it is possible to accurately predict the widths W.sub.DZ and W.sub.Inter of the DZ layer and the intermediate layer (the width of the DZ layer is within the error range of ±5 μm and the width of the intermediate layer is within the error range of ±10 μm), and it is possible to easily determine the RTP conditions for forming desired DZ layer and intermediate layer with the widths W.sub.DZ and W.sub.Inter on the wafer with performing as few as possible experiments involving trial and error.

(100) The above-described embodiment is merely an example, and as long as an object of the present invention, to conveniently determine thermal processing conditions which both the integrity of the denuded zone (DZ layer) and the high gettering ability in the bulk layer can be obtained, is achieved, it is acceptable, for example, to change the parameters of the regression equation.

(101) The present invention can also be applied to crystals to which nitrogen is added at the time of crystal production in order to facilitate the annihilation of cavity defects even in short time thermal processing such as RTP. In such cases, the nitrogen concentration of such a crystal is desirably 1×10.sup.13 or more and 4×10.sup.15 cm.sup.−3 or less, and more desirably 2×10.sup.14 or more and 2×10.sup.15 cm.sup.−3 or less.