- Email: [email protected]
Molecular dynamics analysis of point defects in silicon near solid–liquid interface
Applied Surface Science 159–160 Ž2000. 387–391 www.elsevier.nlrlocaterapsusc
Molecular dynamics analysis of point defects in silicon near solid–liquid interface K. Kakimoto ) , T. Umehara, H. Ozoe Institute of AdÕanced Materials Study, Kyushu UniÕersity, 6-1, Kasuga-Koen, Kasuga, Fukuoka 816-8580, Japan Received 14 October 1999
Abstract Molecular dynamics simulation was carried out to clarify pressure effects on diffusion constants of point defects such as a vacancy and an interstitial atom under constant pressure by using Stillinger–Weber potential. The calculated results indicate that the pressure effect on diffusion of the point defects is small during single crystal growth of silicon, since stress, which was obtained by a global heat and mass transfer model is not enough to modify migration process of the point defects. Activation energy of a vacancy and an interstitial atom was obtained as a function of external pressure. q 2000 Elsevier Science B.V. All rights reserved. PACS: 68.35 y p; 68.60 y p Keywords: Molecular dynamics, point defects; Pressure
1. Introduction Control of distribution of point defects such as a vacancy and an interstitial atom plays an important role to obtain high-quality silicon single crystals. Aggregation of vacancy forms void, which degrades breakdown voltage of LSIs w1x; however, it is hard to control distributions of the point defects. This is mainly attributed to unknown factors such as transport mechanism and thermo-dynamical data of diffusion constants of the defects, although scattered diffusion constants of the defects have been reported, so far w2–16x. Brown et al. w2x reported diffusion constants of defects using molecular dynamics simulation, which was based on constant volume scheme.
Diameter of the crystals is now increasing due to demand of MPU and memories; therefore, temperature gradient in the crystals, especially near a solid– liquid interface of silicon increases. So far, stress effects on diffusion of the point defects have scarcely been studied, although diffusion constants and equilibrium concentrations of the point defects have been discussed. This paper aims to clarify the effects of pressure on diffusion constants of the defects such as a vacancy and an interstitial atom by using molecular dynamics simulation, which is based on constant pressure algorithm. 2. Calculation
) Corresponding author. Tel.: q81-92-583-7836; fax: q81-92583-7838. E-mail address: [email protected] ŽK. Kakimoto..
Global heat and mass transfer calculation was carried out to clarify temperature and stress distribution in silicon crystals during growth furnace. This
0169-4332r00r$ - see front matter q 2000 Elsevier Science B.V. All rights reserved. PII: S 0 1 6 9 - 4 3 3 2 Ž 0 0 . 0 0 1 2 1 - 5
388
K. Kakimoto et al.r Applied Surface Science 159–160 (2000) 387–391
forces. Subsequently, the number was converted into pressure of the system in each time step using Virial theory expressed by Eq. Ž1., 1 1 P s Nk B T q Ýr F V 3V i i i
¦ ; N s žV k T q 31N ¦Ý r F;/ , B
i
Ž 1.
i
i
where V, ri , k B , T and Fi are volume of the system, positions vector of atom i, the Boltzmann constant, temperature and force of ith atom, respectively. The present calculation employed a constant pressure algorithm based on Anderson’s method w19x. Diffusion constants of a vacancy and an interstitial atom were calculated by Eq. Ž2. 1
™` 6 t
D s lim 1
1
¦
N
N
™
™
;
Ý < R i Ž t . y R i Ž 0. < 2 is1
Ž 2.
where N was set to one, since only one atom was migrating during the calculation. Fig. 1. Calculated stress distribution in a silicon crystal by using FEMAG. Unit of the stress indicated in right hand is bar.
calculation was carried out by FEMAG, which was developed by Dupret et al. w17x. The calculation contains heat transfer by convection, conduction and radiation in a growth furnace. Stress distribution in growing crystals can be obtained from the temperature distribution by using elastic model w18x. Molecular dynamics calculation was carried out in the following sequence. Atom positions, volume and velocity were revised in each time step in the calculation using Gear’s fifth-order algorithm w19x. Stillinger–Weber potential w20x was used as a model potential. Total number of silicon atoms used in the calculation was 63 and 215 for a case of a single vacancy, while 65 and 217 for a case of an interstitial atom. Total time step was set to 0.3 = 10y9 s, which corresponds to 1 = 10 6 calculation-steps to obtain diffusion constants, statically. Temperature of the system was controlled by momentum scaling method. Since hydrostatic pressure was applied to the system, isotropic deformation was resulted in the calculation. Virial numbers of the system were calculated in each time-step by summation of two and three body
3. Calculated results Fig. 1 shows stress distribution in a silicon crystal, which was obtained from a calculation of a global heat, and mass transfer analysis by using
Fig. 2. Initial positions of a vacancy and an interstitial atom in silicon crystal. The initial positions of the defects are indicated by reduced units.
K. Kakimoto et al.r Applied Surface Science 159–160 (2000) 387–391
389
Fig. 3. Mean-square displacements of silicon atoms as a function of time for cases of vacancy Ža., interstitial Žb. and perfect crystal Žc..
FEMAG w17x. A gull-shaped solid–liquid interface was resulted in the calculation shown in Fig. 1. The contours show axial component of internal stress with unit of bar, therefore, maximum stress becomes about 20 bar in this case, which is located at a solid–liquid interface. When a diameter of silicon single crystals becomes large: 12 in., maximum pressure might be much higher than the value. Therefore, stress in the present calculation was set up to "50 kbar, which might cover maximum pressure in growing crystals with large diameter. Fig. 2 shows initial positions of vacancy and an interstitial atom in the calculation. The positions are indicated by reduced length, which is based on unit cell of diamond structure of silicon. Initial positions of a vacancy and an interstitial atom were set to Ž1r2, 1r2, 1r2. and Ž9r16, 9r16, 9r16. in reduced length in x, y, z directions, respectively. Fig. 3 shows mean-square displacements of the atoms in the calculating system as a function of duration time, which corresponds to a term of summation in Eq. Ž2. with a vacancy and an interstitial atom, respectively. Mean square displacement of a perfect crystal is also shown in the figure for comparison. Temperature and pressure were imposed at 1600 K and 0 kbar, respectively.
Fig. 4. Arrhenius plots of an interstitial atom. Reported data w2–16x are also indicated.
390
K. Kakimoto et al.r Applied Surface Science 159–160 (2000) 387–391
Fig. 7. Pressure dependence of activation energy for an interstitial atom as a function of external pressure.
Fig. 5. Diffusion constants of vacancy with several data that were reported elsewhere w2–16x.
Mean-square displacements of both cases with a vacancy and an interstitial atom saturated at specific values after 3 = 10y1 0 s, while the value for perfect crystal still decreasing. We can confirm that the mean-square displacement of the defects can be obtained from this calculation statistically; therefore, diffusion constants can be obtained from the data. Figs. 4 and 5 show Arrhenius plots of calculated diffusion-constants of a vacancy and an interstitial atom as a function of temperature. Defect pattern,
which is indicated in the figures means that the data was experimentally obtained from formation process of defect pattern. The external pressure was imposed to 0 kbar. A dashed line and open circles indicate present result, while reported data in the last decades are also listed in the same figure w2–16x. Diffusion constants of a vacancy and an interstitial atom are estimated to 6 = 10y4 cm2rs and 5.8x10y4 cm2rs, respectively, at melting point of silicon in the calculation. Figs. 6 and 7 show activation energies of a vacancy and an interstitial silicon atom as a function of external pressure, respectively. The results indicate that activation energy slightly decreases when pressure increases for a case of a vacancy, while the value gradually increases for a case of an interstitial atom. The above result indicates that diffusion constants are not sensitive to stress within a range from under y1 to 1 kbar, which corresponds to stress range in silicon crystals during growth.
4. Conclusion
Fig. 6. Pressure dependence of activation energy for vacancy as a function of external pressure.
Molecular dynamics simulation was carried out to clarify stress effect on diffusion constants of point defects under constant pressure. The calculation indicates that stress effect is small within a range of stress in silicon crystals during growth, which was obtained by a global heat and mass transfer model. Activation energy of a vacancy slightly decreases with increase in pressure, while that of an interstitial atom increases.
K. Kakimoto et al.r Applied Surface Science 159–160 (2000) 387–391
Acknowledgements This work was conducted as JSPS Research for the Future Program in the Area of Atomic-Scale Surface and Interface Dynamics. The New Energy and Industrial Technology Development Organization ŽNEDO. through the Japan Space Utilization Promotion Center ŽJSUP. supported a part of this work. This work was partially carried out under the support of Grant-in-Aid to the Science Research by the Ministry of Education, Science and Culture.
References w1x M. Itsumi, F. Kiyosumi, Appl. Phys. Lett. 40 Ž1982. 496. w2x R.A. Brown, D. Maroudas, T. Sinno, J. Cryst. Growth 137 Ž1994. 12. w3x P.J. Ungar, T. Halicioglu, W.A. Tiller, Phys. Rev. B 50 Ž1994. 7344. w4x H. Bracht, N.A. Stolwijk, H. Mehrer:, in: Semicond. Silicon, ECS, Pennington, NJ, 1994, p. 593.
391
w5x T.Y. Tan, U. Goesele, Appl. Phys. A 37 Ž1985. 1. w6x H.J. Gossmann, C.R. Rafferty, H.S. Luftmann, Appl. Phys. Lett. 63 Ž1993. 639. w7x G.B. Bronner, J. Cryst. Growth 53 Ž1981. 273. w8x K. Wada, N. Inoue, Defects and properties of semiconductors, in: Defect Engineering, KTK Scientific, 1985, p. 169. w9x R. Habu, A. Tomiura, H. Harada, in: Semicond. Silicon, ECS, Pennington, NJ, 1994, p. 635. w10x B. Leroy, J. Appl. Phys. 50 Ž1979. 7996. w11x C. Boit, J. Cryst. Growth 53 Ž1981. 563. w12x T. Abe, K. Kikuchi, S. Shirai, S. Muraoka, in: H.R. Huff, R.K. Kriegler, Y. Takeishi ŽEds.., Semiconductor Silicon, Electrochem. Soc, Pennnington, 1981, p. 54. w13x H. Zimmermann, H. Ryssel, Appl. Phys. A 55 Ž1992. 121. w14x K. Taniguti, D.A. Antoniadis, Y. Matsushita, Appl. Phys. Lett. 42 Ž1983. 961. w15x M. Yoshida, K. Saito, Jpn. J. Appl. Phys. 6 Ž1967. 573. w16x K. Tempellhoff, Appl. Phys. Lett. 42 Ž1988. 961. w17x F. Dupret, Y. Ryckmans, P. Wouters, M.J. Crochet, J. Cryst. Growth 79 Ž1986. 84. w18x F. Dupret, P. Nicodeme, Y. Ryckmans, J. Cryst. Growth 97 Ž1989. 162. w19x C.W. Gear, Numerical initial value problems in ordinary differential equations’, Prentice Hall, Englewood Cliffs, 1971. w20x F.H. Stillinger, T.A. Weber, Phys. Rev. B 31 Ž1985. 5262.
Molecular dynamics analysis of point defects in silicon near solid–liquid interface K. Kakimoto ) , T. Umehara, H. Ozoe Institute of AdÕanced Materials Study, Kyushu UniÕersity, 6-1, Kasuga-Koen, Kasuga, Fukuoka 816-8580, Japan Received 14 October 1999
Abstract Molecular dynamics simulation was carried out to clarify pressure effects on diffusion constants of point defects such as a vacancy and an interstitial atom under constant pressure by using Stillinger–Weber potential. The calculated results indicate that the pressure effect on diffusion of the point defects is small during single crystal growth of silicon, since stress, which was obtained by a global heat and mass transfer model is not enough to modify migration process of the point defects. Activation energy of a vacancy and an interstitial atom was obtained as a function of external pressure. q 2000 Elsevier Science B.V. All rights reserved. PACS: 68.35 y p; 68.60 y p Keywords: Molecular dynamics, point defects; Pressure
1. Introduction Control of distribution of point defects such as a vacancy and an interstitial atom plays an important role to obtain high-quality silicon single crystals. Aggregation of vacancy forms void, which degrades breakdown voltage of LSIs w1x; however, it is hard to control distributions of the point defects. This is mainly attributed to unknown factors such as transport mechanism and thermo-dynamical data of diffusion constants of the defects, although scattered diffusion constants of the defects have been reported, so far w2–16x. Brown et al. w2x reported diffusion constants of defects using molecular dynamics simulation, which was based on constant volume scheme.
Diameter of the crystals is now increasing due to demand of MPU and memories; therefore, temperature gradient in the crystals, especially near a solid– liquid interface of silicon increases. So far, stress effects on diffusion of the point defects have scarcely been studied, although diffusion constants and equilibrium concentrations of the point defects have been discussed. This paper aims to clarify the effects of pressure on diffusion constants of the defects such as a vacancy and an interstitial atom by using molecular dynamics simulation, which is based on constant pressure algorithm. 2. Calculation
) Corresponding author. Tel.: q81-92-583-7836; fax: q81-92583-7838. E-mail address: [email protected] ŽK. Kakimoto..
Global heat and mass transfer calculation was carried out to clarify temperature and stress distribution in silicon crystals during growth furnace. This
0169-4332r00r$ - see front matter q 2000 Elsevier Science B.V. All rights reserved. PII: S 0 1 6 9 - 4 3 3 2 Ž 0 0 . 0 0 1 2 1 - 5
388
K. Kakimoto et al.r Applied Surface Science 159–160 (2000) 387–391
forces. Subsequently, the number was converted into pressure of the system in each time step using Virial theory expressed by Eq. Ž1., 1 1 P s Nk B T q Ýr F V 3V i i i
¦ ; N s žV k T q 31N ¦Ý r F;/ , B
i
Ž 1.
i
i
where V, ri , k B , T and Fi are volume of the system, positions vector of atom i, the Boltzmann constant, temperature and force of ith atom, respectively. The present calculation employed a constant pressure algorithm based on Anderson’s method w19x. Diffusion constants of a vacancy and an interstitial atom were calculated by Eq. Ž2. 1
™` 6 t
D s lim 1
1
¦
N
N
™
™
;
Ý < R i Ž t . y R i Ž 0. < 2 is1
Ž 2.
where N was set to one, since only one atom was migrating during the calculation. Fig. 1. Calculated stress distribution in a silicon crystal by using FEMAG. Unit of the stress indicated in right hand is bar.
calculation was carried out by FEMAG, which was developed by Dupret et al. w17x. The calculation contains heat transfer by convection, conduction and radiation in a growth furnace. Stress distribution in growing crystals can be obtained from the temperature distribution by using elastic model w18x. Molecular dynamics calculation was carried out in the following sequence. Atom positions, volume and velocity were revised in each time step in the calculation using Gear’s fifth-order algorithm w19x. Stillinger–Weber potential w20x was used as a model potential. Total number of silicon atoms used in the calculation was 63 and 215 for a case of a single vacancy, while 65 and 217 for a case of an interstitial atom. Total time step was set to 0.3 = 10y9 s, which corresponds to 1 = 10 6 calculation-steps to obtain diffusion constants, statically. Temperature of the system was controlled by momentum scaling method. Since hydrostatic pressure was applied to the system, isotropic deformation was resulted in the calculation. Virial numbers of the system were calculated in each time-step by summation of two and three body
3. Calculated results Fig. 1 shows stress distribution in a silicon crystal, which was obtained from a calculation of a global heat, and mass transfer analysis by using
Fig. 2. Initial positions of a vacancy and an interstitial atom in silicon crystal. The initial positions of the defects are indicated by reduced units.
K. Kakimoto et al.r Applied Surface Science 159–160 (2000) 387–391
389
Fig. 3. Mean-square displacements of silicon atoms as a function of time for cases of vacancy Ža., interstitial Žb. and perfect crystal Žc..
FEMAG w17x. A gull-shaped solid–liquid interface was resulted in the calculation shown in Fig. 1. The contours show axial component of internal stress with unit of bar, therefore, maximum stress becomes about 20 bar in this case, which is located at a solid–liquid interface. When a diameter of silicon single crystals becomes large: 12 in., maximum pressure might be much higher than the value. Therefore, stress in the present calculation was set up to "50 kbar, which might cover maximum pressure in growing crystals with large diameter. Fig. 2 shows initial positions of vacancy and an interstitial atom in the calculation. The positions are indicated by reduced length, which is based on unit cell of diamond structure of silicon. Initial positions of a vacancy and an interstitial atom were set to Ž1r2, 1r2, 1r2. and Ž9r16, 9r16, 9r16. in reduced length in x, y, z directions, respectively. Fig. 3 shows mean-square displacements of the atoms in the calculating system as a function of duration time, which corresponds to a term of summation in Eq. Ž2. with a vacancy and an interstitial atom, respectively. Mean square displacement of a perfect crystal is also shown in the figure for comparison. Temperature and pressure were imposed at 1600 K and 0 kbar, respectively.
Fig. 4. Arrhenius plots of an interstitial atom. Reported data w2–16x are also indicated.
390
K. Kakimoto et al.r Applied Surface Science 159–160 (2000) 387–391
Fig. 7. Pressure dependence of activation energy for an interstitial atom as a function of external pressure.
Fig. 5. Diffusion constants of vacancy with several data that were reported elsewhere w2–16x.
Mean-square displacements of both cases with a vacancy and an interstitial atom saturated at specific values after 3 = 10y1 0 s, while the value for perfect crystal still decreasing. We can confirm that the mean-square displacement of the defects can be obtained from this calculation statistically; therefore, diffusion constants can be obtained from the data. Figs. 4 and 5 show Arrhenius plots of calculated diffusion-constants of a vacancy and an interstitial atom as a function of temperature. Defect pattern,
which is indicated in the figures means that the data was experimentally obtained from formation process of defect pattern. The external pressure was imposed to 0 kbar. A dashed line and open circles indicate present result, while reported data in the last decades are also listed in the same figure w2–16x. Diffusion constants of a vacancy and an interstitial atom are estimated to 6 = 10y4 cm2rs and 5.8x10y4 cm2rs, respectively, at melting point of silicon in the calculation. Figs. 6 and 7 show activation energies of a vacancy and an interstitial silicon atom as a function of external pressure, respectively. The results indicate that activation energy slightly decreases when pressure increases for a case of a vacancy, while the value gradually increases for a case of an interstitial atom. The above result indicates that diffusion constants are not sensitive to stress within a range from under y1 to 1 kbar, which corresponds to stress range in silicon crystals during growth.
4. Conclusion
Fig. 6. Pressure dependence of activation energy for vacancy as a function of external pressure.
Molecular dynamics simulation was carried out to clarify stress effect on diffusion constants of point defects under constant pressure. The calculation indicates that stress effect is small within a range of stress in silicon crystals during growth, which was obtained by a global heat and mass transfer model. Activation energy of a vacancy slightly decreases with increase in pressure, while that of an interstitial atom increases.
K. Kakimoto et al.r Applied Surface Science 159–160 (2000) 387–391
Acknowledgements This work was conducted as JSPS Research for the Future Program in the Area of Atomic-Scale Surface and Interface Dynamics. The New Energy and Industrial Technology Development Organization ŽNEDO. through the Japan Space Utilization Promotion Center ŽJSUP. supported a part of this work. This work was partially carried out under the support of Grant-in-Aid to the Science Research by the Ministry of Education, Science and Culture.
References w1x M. Itsumi, F. Kiyosumi, Appl. Phys. Lett. 40 Ž1982. 496. w2x R.A. Brown, D. Maroudas, T. Sinno, J. Cryst. Growth 137 Ž1994. 12. w3x P.J. Ungar, T. Halicioglu, W.A. Tiller, Phys. Rev. B 50 Ž1994. 7344. w4x H. Bracht, N.A. Stolwijk, H. Mehrer:, in: Semicond. Silicon, ECS, Pennington, NJ, 1994, p. 593.
391
w5x T.Y. Tan, U. Goesele, Appl. Phys. A 37 Ž1985. 1. w6x H.J. Gossmann, C.R. Rafferty, H.S. Luftmann, Appl. Phys. Lett. 63 Ž1993. 639. w7x G.B. Bronner, J. Cryst. Growth 53 Ž1981. 273. w8x K. Wada, N. Inoue, Defects and properties of semiconductors, in: Defect Engineering, KTK Scientific, 1985, p. 169. w9x R. Habu, A. Tomiura, H. Harada, in: Semicond. Silicon, ECS, Pennington, NJ, 1994, p. 635. w10x B. Leroy, J. Appl. Phys. 50 Ž1979. 7996. w11x C. Boit, J. Cryst. Growth 53 Ž1981. 563. w12x T. Abe, K. Kikuchi, S. Shirai, S. Muraoka, in: H.R. Huff, R.K. Kriegler, Y. Takeishi ŽEds.., Semiconductor Silicon, Electrochem. Soc, Pennnington, 1981, p. 54. w13x H. Zimmermann, H. Ryssel, Appl. Phys. A 55 Ž1992. 121. w14x K. Taniguti, D.A. Antoniadis, Y. Matsushita, Appl. Phys. Lett. 42 Ž1983. 961. w15x M. Yoshida, K. Saito, Jpn. J. Appl. Phys. 6 Ž1967. 573. w16x K. Tempellhoff, Appl. Phys. Lett. 42 Ž1988. 961. w17x F. Dupret, Y. Ryckmans, P. Wouters, M.J. Crochet, J. Cryst. Growth 79 Ž1986. 84. w18x F. Dupret, P. Nicodeme, Y. Ryckmans, J. Cryst. Growth 97 Ž1989. 162. w19x C.W. Gear, Numerical initial value problems in ordinary differential equations’, Prentice Hall, Englewood Cliffs, 1971. w20x F.H. Stillinger, T.A. Weber, Phys. Rev. B 31 Ž1985. 5262.