Molecular Dynamics Simulation of Silicon Dioxide Etching by Hydrogen Fluoride Using the Reactive Force Field.

In this study, we develop a reactive force field (ReaxFF) for a Si/O/H/F system to perform etching simulations of SiO2 with an HF etchant. Quantum mechanical (QM) training sets from density functional theory calculations, which contain structures of reactant/product and energies with bond dissociation, valence angle distortions, and reactions between SiO2 clusters and SiO2 slab with HF gases, are used to optimize the ReaxFF parameters. Structures and energies calculated using the ReaxFF match well with the QM training sets. Using the optimized ReaxFF, we conduct molecular dynamics simulations of the etching process of SiO2 substrates with active HF molecules. The etching yield and number of reaction products with different incident energies of the HF etchant are investigated. These simulations show that the developed ReaxFF offers insights into the atomistic surface reaction of the SiO2 etching process.

Introduction

Owing to the continuous improvement in the performance of electronic devices, it has become more important to understand plasma–surface interactions at the atomic scale in semiconductor processing. Dry etching processes have been developed and widely used to achieve both high aspect ratios and selectivities in the etching process. In dry etching techniques using remote plasma sources, defects and structure quality degradation can be avoided.[1] Fluorine-based etching gases such as SF6, CF4, and SiF4 are typically used in plasma processing for etching because F atoms are the most reactive among all halogens on Si and produce volatile products such as SiF4 after reacting with Si atoms.[2−6]

For plasma–surface interactions, computational studies have been widely performed to analyze the mechanisms of various plasma processes, such as etching and deposition. Hoshino and Nishioka[7] theoretically suggested the etching reaction mechanism of SiO2 caused by HF molecules using ab initio quantum chemical calculations. Kang and Musgrave[8] conducted a detailed quantum chemical investigation on the chemical etching mechanism of SiO2 with HF and H2O etchants. Several studies based on using molecular dynamics (MD) simulations have been performed to analyze etching processes involving a system comprising Si and SiO2 with halogen/fluorine etchant gases. Abrams and Graves[9] presented a new empirical interatomic potential for Si/C/F systems to simulate their etching reactions with Si and CF3+. Barone and Graves[10] reported MD simulation results of energetic F+ and Cl+ impacting Si surfaces at normal incidence. Ohta and Hamaguchi[11] performed MD simulations for Si and SiO2 etching using energetic halogens (F and Cl) with sets of interatomic potentials. Iwakawa et al.[12] conducted Si etching using Cl-based plasmas, including both high-energy Cl+ ions and low-energy-neutral Cl radicals, using MD simulations. Brichon et al.[13] performed MD simulations of low-energy Cl+ and Cl2+ bombardment on Si(100) surfaces to investigate the effect of plasma dissociation on Si etch applications. Miyake et al.[14] investigated the mechanisms of SiN and SiO2 etching by fluorocarbon or hydrofluorocarbon plasmas using MD simulations. Nakazaki et al.[15] performed MD simulations for Cl+ and Br+ ions incident on Si(100) surfaces with Cl and Br neutrals, respectively, to obtain a better understanding of surface reaction kinetics during Si etching. Numazawa et al.[16] investigated the adsorption mechanisms of F radicals on Si, SiO2, and Si3N4 surfaces during the reactive ion etching using MD simulations and an extended Langmuir model.

Although significant studies based on MD simulations have been reported, it is necessary to study the dynamics of the formation and breaking of bonds to further understand the dynamics of surface etching reactions during dry etching. Hence, the reactive force field (ReaxFF), which allows the formation and breaking of bonds in chemical reactions, should be used as a force field in MD simulations.[17] The ReaxFF method combines a bond order/bond distance relationship with a polarizable charge description using the electronegativity equalization method (EEM)[18] and bond-order-dependent three- and four-body interactions. Owing to the combination of bond/nonbond-order terms, the ReaxFF is applicable to a wide range of materials, including covalent,[19,20] metallic,[21,22] and multicomponent metal hydride/oxide/carbide systems.[23−25] When optimizing the ReaxFF parameters, substantial quantum mechanical (QM) training sets that contain energies of the corresponding structures and reactions are used. As parameters are trained against QM data, MD simulations using the ReaxFF guarantee the accuracy of density functional theory (DFT) calculations with a larger simulation system size and lower computational cost compared with DFT calculations.

In this study, we conducted MD simulations of SiO2 substrates using active HF molecules. First, we optimized the ReaxFF parameters of the Si/O/H/F system against QM training sets using DFT calculations. Next, we conducted MD simulations of the etching process of an α-cristobalite SiO2 substrate using HF molecules and our optimized ReaxFF parameters. Finally, simulation results with different incident energies of the HF molecules were analyzed for a systematic study of the system.

Computational Methods

ReaxFF Parameter Optimization

In the ReaxFF simulation, the total energy of the system is expressed aswhere Ebond is the bond energy calculated from the bond order, Eover is the correction term for overcoordination, Eunder is the correction term for undercoordination, Elp is the lone-pair energy, Eval is the valence angle energy, EvdW is the nonbonded van der Waals interaction energy, and ECoulomb is the nonbonded Coulomb potential associated with the charge of the system. The van der Waals interaction energy EvdW contains short-range Pauli repulsion and long-range dispersion terms. The Coulombic potential ECoulomb of the system was calculated using a geometry-dependent charge distribution from the EEM.

The system of interest includes etching of the SiO2 layer using high-energy HF gas. Therefore, the ReaxFF should capture chemical reactions between Si/O/H/F atoms for the etching simulation, where F ions participate in reactions as anions in the etching process. Previously, the ReaxFF, which included Si/O/H/F atoms, was developed to simulate Kapton polyimide, polyhedral oligomeric silsesquioxane (POSS), silica, and Teflon to investigate the surface chemistry of these materials when exposed to oxygen.[26] However, to utilize the previously developed ReaxFF, the parameters for Si/O/F must be modified to capture the etching of SiO2 and the formation of SiF. Based on the previously developed ReaxFF parameters for Si/O/H/F atoms, we reoptimized parameters for Si/O/F and bond parameters for Si–F and H–F, and the valence angle parameters for F–Si–F and F–Si–O against QM resulted in the training sets, which included the structures and energies of reactions using the covariance matrix adaptation evolution strategy (CMA-ES) technique.[27,28] Each CMA-ES step was iteratively conducted to improve a multivariate normal distribution in the parameter space to identify a distribution that minimizes the objective function or cost function.

Figure shows a comparison of the energies of the SiF4 and SiO2 clusters (Si4O10H4 containing two or three Si–F bonds) with respect to the bond length and valence angle using DFT and our optimized ReaxFF. DFT calculations were conducted using the Amsterdam density functional program[29,30] with a generalized gradient approximation in the Perdew–Burke–Ernzerhof Grimme DFT-D3 dispersion exchange–correction functional form and triple zeta with two polarization function basis sets for all atoms. In Figure a,b, the energy differences from the equilibrium using DFT and the ReaxFF between Si and F in SiF4 are shown as functions of Si–F distance and valence angle, respectively. The bond distance between Si and F changed from 1.3 to 2.3 Å, and the valence angle of F–Si–F changed from 89.5 to 129.5°. Energy differences using DFT and the ReaxFF near the equilibrated bond lengths of SiO2 clusters containing SiF2 and SiF3 are shown in Figure c,d, respectively. In general, the energy difference plots using the ReaxFF showed good agreement with the results from DFT calculations.

Figure 1

Minimum energy difference calculated using DFT and ReaxFF of (a) bond dissociation between Si and F for SiF4, (b) valence angle distortion between (F–Si–F) for SiF4, (c) bond dissociation between Si and F for SiO2 cluster containing two Si–F bonds, and (d) bond dissociation between Si and F for SiO2 cluster containing three Si–F bonds.

In optimizing the ReaxFF parameters, reaction paths including transition states and their associated energies (enthalpy changes and activation energies of reactions using DFT) that produce Si–F bonds or SiF gases were included in the training sets. The main reaction between SiO2 and HF molecules that produce SiF4 gas is as follows

To simulate the etching mechanism, we assumed that four Si–F bonds were subsequently formed to produce SiF4 from the previously suggested subsequent reaction mechanism of SiO2 etching.[7,8] Using the previously explained SiO2 cluster model (Si4O10H4 containing a number of Si–F bonds), we calculated the enthalpy changes and activation energies of four subsequent reactions forming SiF4, and the molecular illustration of each state (initial, transition, and final state in each reaction) is presented in Figure . The transition state location and activation energy for each reaction were calculated using both the nudged elastic band (NEB) method and dimer method.[31,32] Using the NEB method, one can easily determine the intermediate states of the reaction because the NEB method conducts energy calculations in many molecular images. The transition state search process is as follows: first, 30 iterations of geometry relaxation were conducted using the generated images from linear interpolation. Next, the saddle point was searched using the climbing image NEB method.[33] When the energy diverged, we converged the force in the image of the latest step using the dimer method. To implement the NEB and dimer methods, the Vienna Ab-initio Software Package was used[34−37] with the projector augmented wave pseudo-potential and the exchange–correlation function described by the Perdew–Burke–Ernzerhof functional. The plane-wave energy cutoff was set to 400 eV. The structures were fully relaxed until the energy difference in the self-consistent field cycle was lower than 10–7 eV and maximum force acting on each atom was lower than 0.02 eV/Å in ionic steps.

Figure 2

Four subsequent reactions forming SiF4 using the SiO2 cluster model: HF reacting with (a) Si4O10H4, (b) F–Si4O9H3, (c) F2–Si4O9H4, and (d) F3–Si4O9H5.

For further optimization, we calculated the enthalpy change and activation energy of the additional reactions using the SiO2 slab model shown in Figure S1. Additional reactions occurred are in the relatively early stage of the etching process, where HF reacted with a slab containing either no or one Si–F bond. For the reaction of HF with a SiO2 slab containing one Si–F bond, two different paths are possible depending on the location of the additional Si–F bond from the reaction. For the reaction where an additional Si–F bond is formed at the same Si atom containing the Si–F bond, additional O–H bonds are formed at neighboring Si atoms, as shown in Figure S1b. Meanwhile, when an additional Si–F bond is formed at the neighboring Si atom, H2O is produced, as shown in Figure S1c. Despite the mismatches of enthalpy change and activation energy in some reactions, both the ReaxFF and DFT results demonstrated consistency in general.

MD Simulation Using Optimized ReaxFF

The box for the MD simulation was set to 2.1 nm × 1.1 nm × 6.5 nm, which included an α-cristobalite SiO2 substrate measuring 2.1 nm × 1.1 nm × 2.4 nm. The initial geometry of the SiO2 substrate was optimized using a limited Broyden–Fletcher–Goldfarb–Shanno minimization method and a convergence criterion of 1.0 (kcal/mol)/Å, followed by a relaxation process with an NVT ensemble for 100 ps with a time step of 0.25 fs. After the SiO2 substrate was relaxed, the SiO2 etching process was performed for 300 ps, which is enough time to investigate various etching properties in our simulation. For the etching process, incident HF molecules were added to the system with incident energies of 20, 30, 40, and 80 eV. Each HF molecule was added from the top surface of the simulation box with randomly selected x- and y-coordinates every 250 fs, and 200 molecules were inserted over 50 ps. In each insertion step, the initial velocity of the inserted HF molecule contained only the z-direction component with the associated incident energy to simplify the zero-incident-angle system. After 200 HF molecules were inserted, additional simulations were conducted for 50 ps without further insertion of the HF molecule. Subsequently, gas molecules produced such as H2, H2O, O2, and SiF were removed to simulate the purge process. This entire simulation procedure was repeated three times, resulting in a total simulation time of 300 ps. All MD simulations with an optimized ReaxFF were conducted at T = 308.15 K using a Berendsen thermostat with a damping constant of 100 fs.

Results and Discussion

Figure shows images of the etching process simulation with an incident energy of 40 eV at t = 0, 100, and 300 ps. At t = 0 in Figure a, the stable structure of the SiO2 film is clearly depicted. However, incident HF molecules collided with the SiO2 surface, and the H–F bond of the HF molecule dissociated and created a new Si–F bond at the SiO2 surface, as shown in the inset of Figure b at t = 1.075 ps. The O–H bond was created with a left H atom near the previously formed Si–F bond (inset of Figure b at t = 6.075 ps), as suggested in the NEB calculation results. When a sufficient number of Si–F bonds were created in a single Si atom, SiF4 was formed, as shown in the inset of Figure b at t = 39.425 ps, after the dangling Si–O bond was dissociated when an additional O–H bond was formed or the incident HF molecule collided with the dangling Si–O bond. At t = 100 ps, the produced gas molecules were purged and additional insertion of active HF was repeated, as previously explained. At the final stage of the simulation at t = 300 ps, a single layer of SiO2 was removed, and it is expected that more SiO2 layers will be etched as the simulation time progresses.

Figure 3

Images of the etching process simulation with an incident energy of 40 eV at (a) t = 0, (b) t = 100 ps, and (c) t = 300 ps. Insets of (b) show images of Si–F bond formation at t = 1.075 ps, O–H bond formation at t = 6.075 ps, and SiF4 formation at t = 39.425 ps.

For the systematic study of incident energy dependence on etching, we conducted additional simulations with incident energies of 20, 30, and 80 eV. Images of the final structure of all the simulations are presented in Figure . As shown, the amount of etched SiO2 increased with the incident energy. To analyze the effect of the incident energy of each HF molecule on the local chemical reaction at the SiO2 surface, we calculated the dissociation fraction as a function of the incident energy, as shown in Figure a. The dissociation fraction is the probability of incident HF dissociating at the SiO2 surface upon collision, and it is calculated by counting undissociated HF molecules during the purge process. As shown in Figure a, the dissociation fraction is an incremental function of the incident energy of HF. For a better understanding, we calculated the etching initiation time as a function of the incident energy, as shown in Figure b. The etching initiation time was defined as the moment when the first S–F bond was formed at the SiO2 surface. Based on Figure b, we can conclude that the earlier initiation of etching with higher incident energy resulted in more chemical reactions at the SiO2 surface during the etching process.

Figure 4

Final images of each etching simulation with different incident energies: (a) 20, (b) 30, (c) 40, and (d) 80 eV.

Figure 5

(a) Dissociation fraction and (b) etching initiation time with respect to the incident energy of the HF molecule.

To provide a comprehensive explanation of the etching process for each simulation with different incident energies of HF molecules, we counted the number of SiF (x = 4, 5, and 6) and H2O molecules produced in the systems with different incident energies with respect to time, as shown in Figure a,b, respectively. Although SiF4 was the primary volatile gas among the SiF gases that contributed to the etching process, SiF5– and SiF62– were also produced from unreacted active HF reacting with SiF4. It is noteworthy that the numbers of SiF and H2O molecules from t = 50 to 100 ps, t = 150 to 200 ps, and t = 250 to 300 ps did not show significant changes because no incident HF molecules were added to the system in the corresponding time steps for all cases with different incident energies. As shown in Figure a,b, more SiF and H2O molecules were produced through the etching reaction with higher incident energy of HF molecules. Additionally, undesired byproduct molecules such as HF3OSi, HF5Si, and HF4Si, other than SiF and H2O, were present prior to the purge process. However, the amounts of these byproducts were relatively low, and some of these byproducts dissociated naturally prior to the purge process, as shown in Figure c; hence, they did not significantly affect the entire etching process.

Figure 6

(a) Numbers of generated SiF, (b) H2O, and (c) byproduct molecules with respect to time for four incident energies.

Finally, we plotted the etching yield, which refers to the number of removed SiO2 molecules from the initial SiO2 substrate per incident HF molecule, as shown in Figure . We discovered that the etching yield of the high-energy HF on the SiO2 substrate increased with the incident energy, which is consistent with previous explanations regarding the effect of incident energy on the dissociation fraction and etching initiation time. During etching, physical sputtering of the substrate occurred when the incident energy exceeded 100 eV. Therefore, because the incident energy of the etchant was less than that of the physical sputtering regime, we believe that our optimized ReaxFF can well describe the incident energy behavior of SiO2 dry etching with incident HF molecules. Under the incident energy level where physical sputtering can be disregarded, the simulation results suggest that more etching reactions occurred when the incident energy of the etchant molecules increased. We believe that our simulation study will enable the systematic study of the dry etching process and that the model used can be expanded to other chemical etching processes for various applications in the semiconductor industry.

Figure 7

Etching yield of SiO2 by an HF molecule as a function of ion incident energy.

Conclusions

In summary, we conducted an intermediate-size etching MD simulation using newly optimized ReaxFF parameters for a Si/O/H/F system. The ReaxFF was developed by training against DFT data for energies of geometries and associated reactions in SiO2 clusters and SiO2 slabs (surface) for the etching process using HF molecules as an etchant. To validate the newly optimized ReaxFF parameters, we compared the calculated results using the ReaxFF and DFT, which indicated good agreement. Using this optimized ReaxFF, MD simulations of the etching process of the SiO2 surface using HF molecules were conducted. Because physical sputtering was accompanied by etching reactions with incident energies exceeding 100 eV, simulations with etchant incident energies of less than 100 eV were performed to investigate the chemical reactions between the SiO2 substrate and HF molecules. In the SiO2 etching process with an HF incident energy of 40 eV and a simulation runtime of 300 ps, we discovered that a single layer of SiO2 substrate was removed and SiF gases were produced from etching reactions. By performing a systematic study based on different incident energies of 20, 30, and 80 eV, we observed a faster initiation of the etching reaction and higher dissociation fractions with higher incident energies. The SiO2 etching simulation with an HF etchant using the ReaxFF potential developed in this study facilitates a more comprehensive understanding of the computational chemical modeling of the etching process using Si-related surfaces with halogen/fluorine etchants.