Atomic Oxygen Recombination at Surface Defects on Reconstructed

Mar 31, 2015 - (10) obtained a (1×1) pattern on a cut quartz crystal (where the crystal was polished and ...... In a second step, the intermediate D0...
0 downloads 0 Views 7MB Size
Article pubs.acs.org/JPCC

Atomic Oxygen Recombination at Surface Defects on Reconstructed (0001) α‑Quartz Exposed to Atomic and Molecular Oxygen Rubén Meana-Pañeda,† Yuliya Paukku,† Kaining Duanmu,† Paul Norman,‡,§ Thomas E. Schwartzentruber,‡ and Donald G. Truhlar*,† †

Department of Chemistry, Supercomputing Institute, and Chemical Theory Center, University of Minnesota, 207 Pleasant Street SE, Minneapolis, Minnesota 55455-0431, United States ‡ Department of Aerospace Engineering and Mechanics, University of Minnesota, 110 Union Street SE, Minneapolis, Minnesota 55455-0153, United States S Supporting Information *

ABSTRACT: The surface chemistry of silica is strongly affected by the nature of chemically active sites (or defects) occurring on the surface. Here, we employ quantum mechanical electronic structure calculations to study an uncoordinated silicon defect, a non-bridging oxygen defect, and a peroxyl defect on the reconstructed (0001) surface of α-quartz. We characterized the spin states and energies of the defects, and calculated the reaction profiles for atomic oxygen recombination at the defects. We elucidated the diradical character by analyzing the low-lying excited states using multireference wave function methods. We show that the diradical defects consist of weakly coupled doublet radicals, and the atomic oxygen recombination can take place through a barrierless process at defects. We have delineated the recombination mechanism and computed the formation energy of the peroxyl and non-bridging oxygen defects. We found that key recombination reaction paths are barrierless. In addition, we characterize the electronically excited states that may play a role in the chemical and physical processes that occur during recombination on these surface defect sites.

1. INTRODUCTION Silica is a material with a variety of important applications,1−3 many of which are strongly impacted by the properties of its surface; therefore, an accurate characterization of its gas−surface interface is important for the design and modeling of applications. In aerospace applications, silica is used as a major component of both reusable and ablative thermal protection systems for the re-entry of vehicles.4 During high-speed flight through the atmosphere, a shock wave is formed in front of the vehicle, and temperatures between the vehicle and the shock can reach several thousand kelvins. At such extreme temperatures, the gas-phase molecular species dissociate, and under the low-density conditions of high-altitude flight, the resulting atoms can diffuse through the boundary layer and reach the surface where their exothermic recombination increases the surface heating of the vehicle. The oxygen−silica system is relevant to many hypersonic applications because most current high temperature oxidationprotection systems for leading edges on hypersonic vehicles involve silica-based or silica-forming coatings. Furthermore, natural silica oxides form on SiC, Si3N4, and most ultrahigh temperature ceramic (UHTC) composite systems presently under development.5,6 The results of the current study are aimed at increased understanding of reusable ceramic materials and not materials that undergo ablation during flight. The heterogeneous recombination of oxygen is particularly important because © XXXX American Chemical Society

oxygen dissociates at lower temperatures (and therefore over a wider variety of flight conditions) than nitrogen. To improve the design of thermal protection systems, we must understand how heterogeneous reactions are catalyzed on these surfaces. Although there is a large amount of experimental work on the catalysis of oxygen recombination on silica-derived surfaces, as reviewed by Bedra et al.,7 detailed atomic-level information on the interface is not easily accessible to experimental study. The goal of the present work is to use quantum mechanical electronic structure calculations to elucidate the chemically active surface sites on α-quartz silica and study the mechanism of the atomic oxygen recombination that takes place on them. Freshly cleaved silica is highly biotoxic because of the presence of highly reactive surface sites.8,9 Depending on the conditions, the atomic positions of silica and oxygen can form a variety of surface structures that differ in stability. Silica has many forms, the most common of which is quartz, and the current study focuses on α-quartz as a starting structure. The reconstruction of the α-quartz surface by cleavage of the bulk material has been experimentally studied, and, although there are no conclusive experimental results for the atomic positions corresponding to the reconstructed surface, the pattern of the structure is Received: January 6, 2015 Revised: March 29, 2015

A

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

the direct gas-phase recombination of oxygen.18 These simulations were carried out with the ReaxFFGSI SiO potential, which was specifically parametrized to model oxygen−silica gas− surface interactions.19 We identified five types of defects that are potentially catalytic with respect to the direct gas-phase recombination of atomic oxygen; these five types are shown in Figure 1. They do not affect each other, and we already

well-known, and it agrees reasonably well with the data obtained from molecular simulations. Jánossy et al.10 obtained a (1×1) pattern on a cut quartz crystal (where the crystal was polished and then etched in HF to remove damaged layers) using lowenergy electron diffraction (LEED), and their work reveals the reconstruction of the Y-cut, Z-cut, and R planes. Bart et al.11 carried out an exhaustive study of the α-quartz surfaces using LEED, X-ray absorption spectroscopy at the oxygen K edge (XANES), and low-energy electron loss spectroscopy (ELS) to study the long-range, medium-range order, and the electronic structure of the surface band gap, respectively; they concluded that the (1×1) pattern is obtained only after a slight chemical etching in HF, while a (3×1) or (1×3) pattern, which is associated with reconstructed surface regions, is related to a transition phase or a new metastable phase that can be regarded as a (2×1) pattern. More recently, Bart and Gautier12 studied the evolution of the crystallographic structure of the Z-cut surface with respect to the heating temperature, and they found two different patterns, (1×1) and (841/2×841/2). The last pattern emerges at temperatures near 900 K, and it is probably due to the α → β-quartz transition. A more recent study using helium atom scattering (HAS)13 showed that a (2×2) pattern for the same Z-cut can be explained as alternating patches of the (2×1) reconstruction. Rignanese et al.14 investigated the (0001) α-quartz surface using the Car−Parrinello method with the local density approximation (LDA) to Kohn−Sham density functional theory and obtained two stable reconstructions, with (1×1) and (2×1) patterns, respectively. In particular, the pattern (1×1) corresponds to the most stable reconstruction, and it was obtained by heating the cleaved surface to 300 K. This structure, called the dense surface, was more recently identified and studied using density functional calculations and classical molecular dynamics calculations.15,16 A recent study based on Born−Oppenheimer molecular dynamics simulations and static DFT calculations found a reconstruction that led to a 2×2 reoptimized dense surface, which is about 10% lower in energy than the dense surface.17 These patterns agree with their counterparts obtained in the surface reconstruction simulations that we carried out in 19 our previous work18 using the reactive force field ReaxFFGSI SiO , which was specifically parametrized to describe gas surface interactions in silica−oxygen systems. These reconstructed surfaces (shown later in a figure introduced in section 2.1) will be taken as the starting point for the study of the quartz silica surface in this Article. It is well-known that the behavior of silica surfaces is strongly impacted by the nature of defects; therefore, silica surfaces have been characterized by computational as well as experimental means. For instance, the surface and its gas−surface interaction with atomic oxygen have been studied by electronic structure calculations on simple models,19−24 and some experimental data have confirmed the presence of these defects. Silicon dangling bonds, non-bridging oxygen-hole centers, oxygen vacancies, and peroxyl radical defects have been detected by electron spin resonance (ESR), 25,26 electron paramagnetic resonance (EPR),27,28 and reflection energy loss spectroscopy (REELS), and their electronic structures have been computed using semiempirical and density functional methods.29−31 The in situ defects occurring on silica surfaces exposed to atomic oxygen or air plasma have not been experimentally identified. In a previous work, we used molecular dynamics simulations of silica surfaces exposed to atomic oxygen to model which surface structures present a potential catalytic pathway for

Figure 1. Surface defects on the α-quartz surface. Silicon atoms are shown in yellow, and oxygen atoms are in red. Spheres represent the atoms that are part of the defect site, while a wire frame model is used to show the connectivity of the defect with the rest of the crystal.

concluded that they are highly localized.18 The defects can be created by distorting the geometry of a SiO4 center at the silica surface and, for four of them, adding one or two oxygen atoms. They can be divided into two groups (see Figure 1). The first group includes an uncoordinated silicon defect called Si-UC (Figure 1a), a non-bridging oxygen defect denoted NBO I (Figure 1b), and a peroxyl defect (Figure 1c). The second group includes two kinds of non-bridging oxygen defect, one called NBO II (Figure 1d) and another called NBO III (Figure 1e). Within each group, all defects can be obtained from a previous defect by adding an atom of oxygen to the previous defect. In our work reported earlier,18 we found that such defects also form on amorphous silica surfaces. In the present work, we focus on the Si-UC, NBO I, and peroxyl defects, which were found to be the most prevalent defects in our previous study modeling silica surfaces exposed to atomic oxygen.18 In the present work, we use quantum mechanical electronic structure calculations to study these defects. Although the ReaxFFGSI SiO potential is useful for a qualitative understanding of the structures occurring on silica surfaces and the possible reaction pathways for the recombination of atomic oxygen, this interatomic potential is limited to calculations of the potential energy of a system in the ground electronic state. An accurate understanding of the recombination of oxygen on silica surfaces process also involves excited electronic states and transitions between electronic states. For example, the recent work of White et al. found that approximately 10% of the oxygen molecules produced by recombination on fused quartz were in the O2(a1Δ−g ) singlet state.32 This Article will provide the first quantum mechanical calculations relevant to that finding. A common feature of the defects shown in the upper half of Figure 1 is that on the (0001) reconstructed α-quartz surface, B

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C they are diradicals; that is, they not only contain one uncoordinated silicon/oxygen atom on the surface but also involve an uncoordinated oxygen atom in a deeper layer of the crystal, which is shown in Figure 2. Figure 2 further shows that

atomic and molecular oxygen, and their energetic stability in different electronic states. To do this efficiently, we need to understand what size of cluster is adequate for modeling the defects, and we also need to ascertain the preferred electronic spin state of the system. Models that are too small (i.e., clusters that have too few atoms) do not sufficiently reflect the character of the surface and are not reliable for predicting material properties. In a previous work, Lopes et al. pointed out the importance of including both singlet and higher spin states in the calculations and considering the diradical nature of defects.33 Here, we carried out ground- and excited-state calculations for models of various sizes to understand the diradical nature of the defects and to understand to what extent they can be understood by monoradical cluster models. This Article is organized as follows: Section 2 presents the computational details; section 3 contains the computational results and discussion; and section 4 summarizes the main findings.

Figure 2. Close-up picture of the central silicon atom in the nondefective structure and in the Si-UC defect. Silicon atoms are yellow, and oxygen atoms are red. The top of each graphic corresponds to the solid− vacuum interface. The graphic on the left shows the undefected surface, and the graphic on the right shows the silicon atom and oxygen atom radicals in the Si-UC defect as produced by the pyramidal inversion of the central Si atom.

2. COMPUTATIONAL METHODS We use cluster models to study the surfaces. First we describe the way that we created the cluster models of defects; then we describe the electronic structure methods. 2.1. Cluster Models. A reconstructed quartz surface was created by following the recipe described in section 3.1 of our previous article (ref 18). A 2401-atom periodic cluster of the reconstructed quartz surface that contains one central peroxyl defect (see Figure 3) was energy minimized by using the 19 From this structure, we generated a ReaxFFGSI SiO potential. cluster that contains 24 atoms of Si, called here T24 (which denotes that it has 24 silicon atoms, all tetrahedral (T) centers in the nondefective structure); this is the starting point for the generation of smaller clusters explained below. Figure 4a shows the T24 cluster model for the nondefective structure. The T24 cluster for the NBO I defect was generated as follows: three spheres with a radius of 6 Å were created around the three key atoms, the central silicon atom and the two uncoordinated oxygen atoms of the NBO I defect (see Figure 4c). All atoms within the spheres or involved in the 12-membered rings (see Figure 5) that include any of the atoms inside the spheres were

the Si-UC defect can be considered to result from the inversion of the tetrahedron around a Si atom near the surface. In the nondefective structure (left side of Figure 2), the Si atom is bonded to three O atoms above it and one below. The defect can be made by raising this Si atom above the surface plane; this breaks the bond to O below it, so the O atom moves downward, and it creates a diradical (right side of Figure 2) with one radical site on the now tricoordinated Si atom and one on the now singly coordinated O atom below it. The NBO I defect is then made by adding an oxygen atom to the Si-UC defect, and the peroxyl defect is then made by adding an oxygen atom to the NBO I defect. Our goals in this study are to understand and characterize the silica surface defects, how they can be formed, their reaction with

Figure 3. Three different views of a slab of the (1×1) α-quartz reconstructed surface formed by 2401 atoms optimized using ReaxFFGSI SiO . The size of the slab is ∼50 × 40 × 11 Å. (a) View in perspective, (b) top view, and (c) side view. Spheres represent the atoms of the defect surface site. C

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

Figure 4. Nondefective and defective structures of cluster model T24.

retained. The under-coordinated silicon atoms were capped with hydrogen atoms (with a fixed Si−H distance of 1.46 Å34) to restore proper coordination of the model. A partial optimization of the cluster model was then carried out by optimizing all atoms within a radius of 4 Å from any of the three key atoms (see Figure 5). The same atoms and atomic positions were used with the same procedure for generating and energy minimizing the T24 nondefective structure, the Si-UC defect, and the peroxyl defect. Figure 4 shows the T24 structure for the three defects to be studied in this Article: Si-UC (Figure 4b), NBO I (Figure 4c), and peroxyl (Figure 4d) defects. Other smaller clusters that mimic the atomic structure of the defects shown in Figures 6 and 7 have also been considered in this study. We use the notation Tnx where n represents the number of Si atoms in the cluster, and x is a or b. The first two

Figure 5. Cluster model T24 for the NBO I defect where the atomic positions of the atoms represented by spheres are relaxed during the optimization. The bottom shows details of the 12-membered ring.

Figure 6. Various models for the non-bridging oxygen defect for the study of the atomic oxygen recombination on silica surface. D

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

doublet spin ground state. The geometries of D and R2 were partially optimized in their triplet ground spin state. The partial optimizations for all four structures consisted of optimizing the positions of the terminal hydrogen atoms to minimize the energy with the rest of the framework frozen according to the atomic positions of the T24 cluster shown in Figure 5. 2.2. Electronic Structure Methods. Kohn−Sham density functional theory36 was used for geometry optimization of all of the models and to study the mechanism of atomic oxygen recombination on the surface defects. The state-averaged complete active space self-consistent field (CASSCF)37 method followed by the multireference configuration interaction MRCISD+Q method38−46 is used for the calculation of the excited states in the study of the diradical character of the defects and to map out the excited states at some reaction-path geometries. To choose a model chemistry (a model chemistry is a combination of an electronic structure method or exchangecorrelation functional and a basis set) for the study of the clusters described in the previous section, we analyzed the influence of level of theory on the optimized geometry of a test system for which the usage of a larger basis set is affordable; in particular, we studied the silicate Si(OH)3−O−Si(OH)3. The geometry was optimized using three exchange-correlation functionals, M06,47 M06-L,48 and M06-2X,47 that have been successfully used in previous studies of similar systems19 using the basis set MG3S49 (which is the same as 6-311+G(3d2f,2df,2p) for the systems studied here). We note that M06-L is a local functional, and M06 and M06-2X are hybrid functionals; local functionals are less expensive than hybrid functionals for large systems, and a previous test on bond angles, bond lengths, dipole moments, and proton activity of silicaceous materials found M06-L to be the most accurate local functional tested.50 In the present work, the geometrical parameters were compared to their counterparts for α-quartz.51 As shown in Table 1, in all cases the Si−O bond length is in good agreement with the experimental bond lengths in quartz; therefore, we have chosen the least expensive functional, which is M06-L, for optimizing the geometry of the T24 cluster and for structures D, RA, RB, and R2. We used the more broadly applicable47,52 M06 functional for models T1a, T1b, T4a, and T4b. Using the same geometric criteria, several double- and triple-ζ basis sets (see Table 2) containing various numbers of basis functions were evaluated. On the basis of the comparison of Si−O bond length with experiment, for further cluster optimization we have chosen the 6-311G(2d) triple-ζ basis set53 for Si, 6-311G(d) for O, and 6-311G for H; this combination is henceforth denoted as 6-311G(2d,d). Because it is less expensive than MG3S, we used the 6-311G(2d,d) basis set for all our models except for DFT calculations on models T1a, T1b, T4a, and T4b, where the larger MG3S basis set was used. For calculations on the excited states of D, R2, RA, and RB, we performed state-averaged CASSCF37 calculations followed by MRCISD+Q38−42 calculations. The number of states included in the state average is 17 for D and R2 and 6 for RA and RB. In the CASSCF calculations, we placed the 1s, 2s, and 2p orbitals of Si and the 1s and 2s orbitals of O in the inactive space, while the rest of occupied orbitals of the molecules are active; therefore, the CASSCF reference wave function is CASSCF(34,18) for D, CASSCF(31,16) for RA and RB, and CASSCF(30,16) for R2. In the study of excited states for atomic oxygen recombination, we performed restricted active space self-consistent field (RASSCF) calculations54,55 followed by MRCISD+Q calculations.

Figure 7. Models for the study of the nature of the diradical character of the NBO I defect. Top left: Diradical D. Top right: Radical RA. Bottom left: Radical RB. Bottom right: R2, which is a dimer of radicals.

characters correspond to the Tn notation already explained, and the appended a or b distinguishes whether terminal atoms of the cluster are capped by (a) hydrogens or (b) hydroxyl groups. The models in Figure 6, denoted here as T1a, T1b, T4a, and T4b, are monoradicals. The models T1a, T1b, T4a, and T4b have been carved out from the partially optimized geometry of the larger cluster T24 described above (see Figure 4). For the capping groups, the Si−H distance (for T1a and T4a) and the O−H distance (for T1b and T4b) have been fixed to 1.46 Å34 and 0.959 Å,35 respectively. All atomic positions have been optimized, with the exception of the terminal hydrogens for models T1a and T4a and the terminal hydroxyl groups for models T1b and T4b. The models cover a wide size range that extends from the simplest, where the entire crystal had been replaced by an atom of hydrogen (model T0-H, see Supporting Information), to more complex, where we included four tetrahedrons of full SiO4 (model T4b). In the Supporting Information, we include the performance of various model chemistries for the study of binding energy and geometries of model T0-H, and the larger clusters are considered in the main part of this Article. For the study of the excited states, we considered four other structures that are shown in Figure 7. We label them as structure D (a monomer diradical (D) with chemical formula Si3H6O4), RA (a monoradical with chemical formula Si3H7O3, where R denotes radical), RB (a monoradical with the same chemical formula as RA but with a different location of the undercoordinated atom of oxygen), and R2 (a dimer of monoradicals: (SiH3O2)2). The three monomers were cut out from the NBO I defect of the T24 cluster model (see Figure 4), and the dimer was obtained from D by removing the central Si. All of the remaining uncoordinated Si atoms were capped with hydrogens. The geometries of RA and RB were partially optimized in their E

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

Table 1. Comparison of Distances (in Å) and Angles (in deg) of the Optimized Geometry for the Si(OH)3−O−Si(OH)3 Molecule Calculated Using Various Density Functionals with Experimental Measurements

Table 2. Distances (in Å) and Angle (in deg) of the Optimized Geometry of Si(OH)3−O−Si(OH)3 Molecule Calculated Using the Local Density Functional M06-L with Various Basis Sets

a

6-31G(2d) for Si, 6-31G(d) for O, and 6-31G for H. b6-311G(2d) for Si, 6-311G(d) for O, and 6-311G for H. c6-311G(df) for Si, 6-311G(d) for O, and 6-311G for H.

these systems are reported relative to that of the doublet ground state of D0-INT calculated with the same active space (The energy of the doublet ground state of D0-INT is taken as its energy relative to that of SR instead of zero, so we can have a consistent energy scale across all species.) The final active spaces we used are RASSCF(7,7) for SiH3, RASSCF(10,11) for SiH3O, RASSCF(13,15) for SiH3O2, RASSCF(16,19) for doublet states of SiH3O3 (including D0-INT, D0-TS2, Q1-TS4, and Q1-TS4-110), and RASSCF(11,15) for quartet states of SiH3O3. In the calculations carried out to study excited states for atomic oxygen recombination, all RASSCF calculations for molecules containing silicon are state averaged over five states, although MRCISD+Q energies are calculated for all states below 2 eV. The Gaussian 09 software package was used for the density functional calculations,56 and the Molpro package was used for the multireference calculations.57

We needed to consider nine systems in the T1a model: O, O2, SiH3, SiH3O, SiH3O2, D0-INT, D0-TS2, Q1-TS4, and Q1-TS4− 110; the chemical formula of the latter four species is SiH3O3, and a full description of each of these species will be given in sections 3.3 and 3.4. For calculations of doublet spin states of all nine systems and the quartet spin states of the first five, the 1s, 2s, and 2p orbitals of Si and the 1s and 2s orbitals of O are put in the inactive space, while the 3s and 3p orbitals of Si, 1s orbital of H, and 2p orbitals of O are in the active space. To reduce the required memory and computational time while maintaining the accuracy, we employed a restriction that at most two electrons occupy the virtual molecular orbitals. The energies of all species calculated with this active space are the relative energies as compared to SR (SiH3 + O2 + O). For calculations of quartet spin states of the latter four systems with the chemical formula of SiH3O3, smaller active spaces were used to further reduce the memory needed. For the quartet states, the two lowest doubly occupied orbitals and three highest virtual orbitals are excluded from the active space, and the restriction that at most two electrons occupy the virtual molecular orbitals is applied to the rest of the virtual orbitals. The energies of

3. RESULTS AND DISCUSSION 3.1. Defects Analyzed with the T24 Cluster Model. This section presents the results of the electronic structure F

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

nondefective structure in two main ways: (i) the bond between the central silicon atom and one of its corresponding bonded oxygens (located in a deeper layer, see Figures 4b and 2) is broken, and (ii) a pyramidal inversion of the central silicon atom (see Figure 2). Figure 8 shows that the peroxyl defect is calculated to be the most stable structure, and this finding agrees with the fact that the NBO I and peroxyl defects were the most common defects observed in our previous18 simulations; we believe therefore that the peroxyl defect is an important structure for understanding the catalytic recombination of oxygen atoms on the silica surface. To study the dynamics of the atomic oxygen recombination catalyzed on these defects and to propose new models to study it, we need to answer some questions, for example: How strong is the interaction between monoradicals in the diradical structures? What kind of models can be used to study the dynamics of the atomic oxygen recombination? Do we need to consider excited states in addition to the states shown in Figure 8? How does the interaction between radicals influence the energy differences between the excited electronic states? All of these questions will be answered in the following. 3.2. Diradical Character of the NBO I Defect. To gain insight into the diradical character of the NBO I defect, we built the four models shown in Figure 7. Figure 9 shows the structure

calculations for the nondefective structure and the Si-UC, NBO I, and peroxyl defects, as obtained using the T24 cluster. Using the procedures described above, the geometries of the nondefective and defective structures were partially relaxed by M06-L/6-311G(2d,d). All structures of T24 have an even number of electrons; they are closed-shell singlets (nondefective structure) or diradicals (defective structures). For the diradicals, two electronic spin states have been considered, a triplet and a broken-symmetry open-shell singlet. The energetic profile depicted in Figure 8 shows the relative stability of various

Figure 8. Relative energies (in eV) of electronic states of various spin for the nondefective structure and three defects of the silica surface using cluster model T24 at the M06-L/6-311G(2d,d) level of theory. Note that the open-shell singlet (blue lines) and the triplet (red lines) almost coincide for the NBO I and peroxyl defects, so the blue-and-red energy level corresponds to both of the states in those cases.

structures with respect to the nondefective ground electronic state, which corresponds to a closed-shell singlet electronic state. Note that to obtain these relative energies, the energy of one or two atoms of ground-state oxygen (depending on the structure) was added to the energy of the corresponding system, because a meaningful comparison of energies requires that the energies correspond to the same number of atoms with atomic number. The closed-shell singlet ground state for the nondefective structure (plus two ground-state oxygen atoms) has been taken as the zero of energy (see Figure 8). The nondefective structure has a singlet ground state, but we found that all three defective structures have a triplet ground state. In every case, we also tried to optimize the geometry for an open-shell singlet by using the broken-symmetry method.58−60 For the NBO I and peroxyl defects, this led to a stable structure, with a geometry nearly the same as that corresponding to the triplet, and with the singlet energy slightly higher. The spin contamination in the open-shell singlet is large, having the values of 1.0093 and 1.0086 for the expectation value of the square of the total electron spin (⟨S2⟩) for the NBO I and peroxyl defects, respectively; this indicates that these broken-symmetry states are mixtures of the singlet and triplet, but because the energies of the singlet and the triplet are nearly the same, we need not take linear combinations to resolve the singlet component. For the Si-UC defect, the open-shell singlet does not correspond to a stable structure on the potential energy surface. The minimum-energy geometry on the triplet surface of the Si-UC defect is 5.04 eV higher in energy than the singlet energy of the nondefective structure. Its geometry differs from the

Figure 9. Structure and 18 active-space orbitals of D. The orbital number (in order of increasing energy) and orbital label are also shown under each orbital.

and active space orbitals of model D. Orbitals 38−41 are nonbonding pπ orbitals of two radical oxygen atoms, and they are labeled as πx and πy orbitals. Orbitals 34 and 36 are mainly Si−O bonding orbitals, where Si−O is the bond between radical oxygen and silicon, so they are labeled σO orbitals. The triplet states and the singlet states have been labeled as T0, T1, T2, ... and S0, S1, S2, ..., respectively, in order of increasing energy. Table 3 shows the energies and dominant configurations G

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

atoms. T3 is a doubly excited state, which is a combination of the T1 and T2 excitations; T4 and T5 are mainly the excitations from σO bonding orbitals to the local singly occupied πy orbitals; T6 is a doubly excited state, that is, two electrons are excited from 1πx and 2σO to 1πy and 2πy, respectively, which is a combination of T1 and T4, and similarly T7 is a combination of T2 and T5, while T8 is the combination of T4 and T5. Comparison of Table 3 to Figure 8 shows that there are 10 excited singlets and 10 excited triplets of the NBO I structure that lie lower in energy than the Si-UC defective structure with an O atom infinitely far away (5.03 eV above the NBO I ground state; see Figure 4). Thus, we expect that electronically excited states play a role in the recombination of O atoms on defective silica surfaces. Taking into account the dominant configurations and coefficients of the nine lowest-lying singlet states in Table 3, it is clear that each singlet state corresponds to a triplet state. For example, T1 and S1 are both states with 1πx→1πy excitation, and T3 and S3 are both doubly excited states with 1πx→1πy and 2πx→2πy excitations. In all of the systems under consideration, the diradical can therefore be understood as two localized monoradicals that interact only weakly. Combinations of the two local doublets lead to triplet and singlet states. If the two local radicals are weakly coupled, the energy of a triplet state will be nearly degenerate with the energy of a corresponding singlet state. If we define the average energy spitting (ΔE) of triplet and singlet states as

Table 3. MRCISD+Q Energies, Dominant Configurations, and Coefficients for 17 Triplet States and Singlet States of Da

a

state

energy (eV)

dominant configuration

coefficient

T0 T1 T2 T3 T4 T5 T6 T7 T8 T9

0.000 0.233 0.267 0.499 1.993 2.053 2.229 2.326 4.048 4.622

T10

4.876

T11

5.891

T12

6.000

T13

4.914

T14

5.221

T15 T16 S0 S1 S2 S3 S4 S5 S6

6.636 6.359 0.056 0.284 0.332 0.557 2.041 2.119 2.271

S7 S8 S9

2.396 4.106 4.511

S10

4.760

S11 S12

5.723 6.032

S13 S14

6.220 6.580

S15 S16

6.492 6.454

(1πx)2(2πx)2(1πy)1(2πy)1 1πx→1πy 2πx→2πy 1πx→1πy, 2πx→2πy 2σO→2πy 1σO→1πy 2σO→2πy, 1πx→1πy 1σO→1πy, 2πx→2πy 1σO→1πy, 2σO→2πy 9σ→2πy 10σ→2πy 9σ→2πy, 1πx→1πy 10σ→2πy, 1πx→1πy 9σ→1πy, 10σ→1πy 10σ→2πy, 1πx→1πy 10σ→2πy 9σ→2πy 8σ→1πy 10σ→1πy 10σ→1πy, 2σO→2πy 9σ→1πy, 2σO→2πy 9σ→1πy, 1σO→2πy 2π→2πy (1πx)2(2πx)2(1πy)1(2πy)1 1πx→1πy 2πx→2πy 1πx→1πy, 2πx→2πy 2σO→2πy 1σO→1πy 2σO→1πy, 2πx→2πy 2σO→1πy, 2σO→2πy 2πx→1πy, 2πx→2πy 1σO→1πy, 2πx→2πy 1σO→1πy, 2σO→2πy 10σ→2πy 9σ→2πy 10σ→2πy, 1πx→1πy 9σ→2πy, 1πx→1πy 1πy→2πy 1πx→2πy 1σO→2πy, 1πx→1πy 2π→2πy 1σO→1πy, 2σO→2πy 1σO→1πy, 9σ→2πy 8σ→2πy 2π→2πy, 1πx→1πy

0.77 0.80 0.77 0.89 0.72 0.74 0.84 0.82 0.77 0.57 −0.37 0.64 −0.45 0.67 0.36 0.61 0.44 0.40 −0.39 0.46 0.43 0.57 0.62 0.76 0.80 0.80 0.84 0.73 0.73 0.65 −0.40 0.40 0.76 0.69 0.58 0.48 0.60 0.50 0.71 0.57 0.46 0.62 0.55 0.46 0.39 0.70

ΔE =

1 n

n−1

∑ |Ti − Si| i=0

(1)

where Ti and Si represent the energy of a triplet and the energy of a corresponding singlet state, and the index i runs over n states, a small value of ΔE indicates that the two local radicals are weakly coupled. The value of ΔE for the first nine states for model D is only 0.082 eV, which confirms that they interact weakly. To further understand the interaction between the two local monoradicals, we studied the change in the calculated value of ΔE along the rotation of the bottom Si−O bond as shown in Figure 10. This shows that ΔE has a maximum and a minimum

The reference wave function is calculated with CASSCF(34,18).

for the 17 lowest-energy triplet states and the 17 lowest-energy singlets. The ground state (T0) has orbitals 24−39 doubly occupied and orbitals 40 and 41 (both πy) singly occupied. The dominant configurations of excited states, labeled as excitations from the ground-state configuration, and the coefficient of the dominant configuration in the excited state configuration interaction wave function are also shown in Table 3. Our calculations show that the most important configuration for the T1 state is obtained by 1πx→1πy excitation from the ground T0 state, while the most important configuration for the T2 state is 2πx→2πy. They are both local excitations on the radical oxygen

Figure 10. Energy splitting ΔE (in eV) with respect to the torsional angle shown in the picture for the D structure. H

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C of ∼0.1 and ∼0.01 eV, respectively. This difference shows that the weak interaction between the two radical oxygen atoms depends on the distance between the radicals through the dihedral angle along the Si−O bond angle. We also computed the profile of ΔE along the rotation of the same Si−O bond in the R2 dimer (this model was created by removing the center SiH2 group, adding two hydrogen atoms to two oxygen atoms, and optimizing all hydrogen atoms with fixed positions of the remaining heavy atoms). This ΔE profile for the case where the two local radicals are not linked by chemical bonds is shown in Figure 11, where ΔE has a maximum and a minimum value of

along the Si−O bond, but the calculations in Figure 11 show that the interaction of the two radical oxygen atoms takes place through space, and not through chemical bonds connecting them. Because the small energy splittings indicate that the diradical could be understood as two local weakly interacting monoradicals, we further studied two single radical models, that is, RA and RB of Figure 7, where one radical O atom is replaced by an H atom. RA keeps the top oxygen atom, and RB keeps the bottom oxygen atom. The structure and six highest occupied orbitals of RA and RB are shown in Figures 12 and 13, respectively. The energies and dominant configurations of low-lying states for the two models are shown in Table 4.

Figure 11. Energy splitting ΔE (in eV) with respect to the torsional angle for the dimer R2.

Figure 13. Structure and orbitals of RB.

∼0.07 and ∼0.02 eV, respectively. Both values are close to their counterparts for the structure D, which again shows that the two monoradicals may be considered as largely independent. The value of ΔE is still sensitive to the change in the dihedral angle

Table 4. States, Dominant Configurations, and Energies (eV) of Doublet States of RA and RBa state:

D0

D1

D2

D3

D4

D5

dominant configuration:

πy

πx

σO







RA RB

0.000 0.000

0.274 0.277

2.041 1.993

4.120 4.378

5.842 5.511

6.090 6.125

a The dominant configurations are indicated by their singly occupied orbitals.

Now we can compare energies and excitation types of model D with those of models RA and RB. T1 and T2 are obtained by excitations from local πx orbitals to local πy orbitals; similarly, the D1 state is also a state with single excitation from πx orbital to πy orbital, with the energy (0.274 eV for RA and 0.277 eV for RB) close to the energies of T1 and T2 (0.233 and 0.267 eV) and the corresponding S1 and S2 (0.284 and 0.332 eV) states. T3 is a combination of T1 and T2; S3 is a combination of S1 and S2. Therefore, by calculating only D1 states of RA and RB, we can approximately reproduce the excitations and energies of T1, T2, T3, S1, S2, and S3 states of D. Similarly, by calculating D2 states of RA and RB, we can approximately reproduce excitations and energies of excited states of D up to T8 and S8. The result that all nine low-lying states of the diradical molecule D could be simply reproduced by calculating single radical molecules RA and RB further confirms that the two local radicals of D are separable.

Figure 12. Structure and orbitals of RA. I

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

different from those of the models in the previous section; that is, while the ground state of a diradical is a triplet, the ground state of a radical is a doublet. As shown in the previous section, the interaction of the radical sites in the larger cluster is not strong enough to modify the topography of the potential energy surface of the ground electronic state; therefore, we limited our study of the reaction profile for the atomic oxygen recombination to small model clusters depicted in Figure 6. The relative energies calculated by M06/MG3S for varioussize models of the defects are shown in Figure 14. The energies converge nicely as the size of the cluster increases (e.g., the difference in the relative energies between models T4a and T4b is less than 1%), and they are in good agreement with the values obtained for the diradical T24 cluster model (see Figure 8). Comparison of the energies for various structures and model sizes in Figures 14 and 15 (Figure 15 is explained below) shows that T1b, T4a, and T4b are in generally good agreement. This shows that the T1b cluster is already large enough to map out the features of the potential energy surfaces semiquantitatively, and the T1a cluster provides reasonable features of the potential energy surfaces although it is not big enough for this quantitative work. To understand the energetic requirements for the interconversion between the various structures, we carried out nonrelaxed scans (i.e., single-point calculations at a sequence of points along a path) using the M06-L/6-311G(2d,d) model chemistry for the diradical cluster. For the conversion of the Si-UC defect to the NBO I defect, the distance between the central silicon and the impinging oxygen atom was increased over the range from 5.5 to 1.5 Å (see Figure 16a). Scans were also performed for the conversion of NBO I and Si-UC defect to the peroxyl defect, where the O−O distance and the Si−O distance were modified, respectively (see Figure 16b and c). The left panel of Figure 16 involves bringing an O up to the frozen Si-UC structure; the central panel corresponds to bringing O up to a frozen NBO I structure; and the right panel corresponds to

Figure 14. Relative energies (in eV) of defects using various models obtained by the M06/MG3S chemistry model (see models T1a, T1b, T4a, and T4b in Figure 6), where X denotes the silica surface. All atomic coordinates have been optimized, except for the positions of capped hydrogens of models T1a and T4a and capped hydroxyl groups of models T1b and T4b.

This provides a useful strategy for modeling defects in simulations; that is, they can be modeled as localized defects that interact only weakly, and we will use this conclusion in the next section. 3.3. Reaction Profiles for Atomic Oxygen Recombination in the Lowest-Energy States of Each Spin Symmetry. In this section, we present the reaction profiles for atomic oxygen recombination on the silica surface defects. As shown in the previous section, the interaction between radical sites is weak; thus we choose small finite clusters that only contain a single radical site to study the potential energy surface for the reaction (see models in Figure 6). It should be kept in mind that the electronic spin states for the smaller models of this section are

Figure 15. Reaction profile for the atomic oxygen recombination using various models obtained with M06/MG3S model chemistry (see models T1a, T1b, T4a, and T4b in Figure 6). X denotes the silica surface. The energies shown all correspond to the overall system, including separate atoms and diatoms when indicated. The energy levels shown at SR and R are the same for the overall doublet state and the overall triplet state, but after that the paths diverge with the doublet being lower; at P the lowest overall doublet and lowest overall quartet again have the same energy. J

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

Figure 16. Interpolated potential energy curves for the formation of (a) the NBO I from the Si-UC defect, (b) the peroxyl from the NBO I defect, and (c) the peroxyl from the Si-UC defect as calculated with the M06-L/6-311G(2d,d) model chemistry using model T24, where the arrow shows the movement of the atom or atoms to generate the different geometries of the scan.

Table 5. MRCISD+Q Energies (eV) for Low-Lying States of All Molecules in Atomic Oxygen Recombination Reaction with T1a Modela

a

state

SiH3

SiH3O2

Q1-TS4-110

D0-INT

D0-TS2

SiH3O

D0 D1 D2 D3 D4 D5 D6 D7 Q1 Q2 Q3 Q4 Q5

0.00 5.31 5.31 7.75 6.32

−1.67 −1.11 4.13 6.11 7.76

−1.68 −1.68 0.19 0.37 1.20 1.56

0.03 0.07 0.32 0.45 1.22 1.42 1.79

−3.61 −1.34 0.28 0.38 0.88 1.25

−4.03 −3.86 −1.89 0.62 0.93

8.49 8.50 8.62 9.42 9.53

4.20 5.04 6.52 6.65 6.71

−0.12 0.75 1.94 2.37 2.55

−0.79 0.53 0.97 1.34 1.63

0.22 1.12 2.80 4.24 5.12

−3.74 −1.76 −1.54 −0.99 −0.83 −0.68 −0.01 1.00 −1.98 −1.09 0.75 3.05 3.48

state

O2

X3∑−g a1Δ−g b1∑+g

0.00 1.06 1.70

Q1-TS4

state

2.68 4.53 4.64 4.76 4.92 O

3

P 1 D 1 S

0.00 2.15 4.14

The state refers to the state of the subsystem in the column heading, but the energies refer to the overall system as explained in section 3.4.

bringing a frozen O2 up to a frozen Si-UC structure. The interpolated potential energy curves depicted in Figure 16 show that all three processes are barrierless; that is, the NBO I defect can be formed from the uncoordinated silicon defect by adding an atom of oxygen, and the peroxyl defect can be generated from both the NBO I and the Si-UC defects by adding one atom of oxygen and one molecule of oxygen, respectively. In terms of a mechanism for the recombination of atomic oxygen, Figure 16 implies that one peroxyl defect and one atom of oxygen can lead to a NBO I defect and one molecule of oxygen. Figure 15 shows this in more detail by presenting the entire profile of reaction for the models TXa/b (with X = 1 or X = 4) by using the M06/MG3S model chemistry. In a first step, the peroxyl defect forms an intermediate product, labeled here as D0-INT (where D0 denotes that the spin state of the system is a doublet). Our calculations indicate that this reaction is barrierless. In a second step, the intermediate D0-INT produces the molecule of oxygen through the transition state, D0-TS2, which has small barrier height, about 0.06 eV. Adding atomic oxygen to one of the products of reaction (i.e., the NBO I defect), the peroxyl defect can be regenerated. This step closes an autocatalytic cycle that catalyzes the atomic oxygen recombination.

Table 6. MRCISD+Q Energies (eV) for Low-Lying States of SR, R, and P in Atomic Oxygen Recombination Reaction with the T1a Model SR = SiH3 + O2 + O

R = SiH3O2 + O

P = SiH3O + O2

state

energy

state

energy

state

energy

D0 + X3∑−g + 3P D1 + X3∑−g + 3P D2 + X3∑−g + 3P

0.00 1.06 1.70

D0 + 3P D1 + 3P D0 + 1D D1 + 1D

−1.67 −1.11 0.48 1.04

D0 + X3∑−g D1 + X3∑−g D0 + a1Δ−g D1 + a1Δ−g D0 + b1∑+g D1 + b1∑+g D2 + X3∑−g D2 + a1Δ−g D2 + b1∑+g D3 + X3∑−g D4 + X3∑−g D3 + a1Δ−g D4 + a1Δ−g

−4.03 −3.86 −2.97 −2.80 −2.33 −2.16 −1.89 −0.83 −0.19 0.62 0.93 1.68 1.99

We have also included the reaction channel for the first excited state, which is a quartet spin state; this is a concerted reaction passing through the transition state Q1-TS4. K

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

Figure 17. Reaction profile for the atomic oxygen recombination using T1a model obtained with MRCISD+Q/6-311G(2d,d) single point calculations on M06/MG3S geometries. Doublet and quartet states within 2 eV of the ground state of SR are shown in the figure. The lowest doublet and quartet reaction pathways are represented by curved lines and dashed curved lines. A black solid line denotes a doublet spin state, and a red dashed line denotes a quartet spin state; the number in parentheses is the degeneracy of the state, and the asterisk means the state involves O2 in the excited state.

Table 7. Key Geometrical Parameters for All of the Structures in Figure 17a Bond Distances (Å) Si−O1 O1−O2 O2−O3 Bond Angles (deg) Si−O1−O2 O1−O2−O3 Dihedral Angle (deg) Si−O1−O2−O3

SR

R

Q1-TS4-110

NA 1.20 NA

2.45 1.32 NA

NA NA

106.7 NA

107.4 110.0

108.4 110.2

105.9 110.0

NA NA

107.4 168.1

NA

NA

180.0

180.0

174.8

NA

0.0

1.69 1.47 1.57

D0-INT

D0-TS2

1.68 1.52 1.22

1.66 1.73 1.19

P 1.65 NA 1.20

Q1-TS4 1.69 1.47 1.57

a

Note that structure Q1-TS4-110 is obtained from structure Q1-TS4 by changing the O1−O2−O3 bond angle from 168.1° to 110.0°. The structure Q1-TS4 is placed at the far right because it is not located along the low-energy reaction path from SR to R to Q1-TS4-110 to D0-INT to D0-TS2 to P.

3.4. Manifolds of Low-Lying Electronic States. Recent experiments have demonstrated the formation of the molecule of oxygen in its first excited state, O2(a1Δ−g ).32 To understand the mechanism of O2(a1Δ−g ) formation, we next calculated the low-lying excited states using the T1a model. System SR is a combination of three independent subsystems: SiH3, O2, and O; so each state of SR is a combination of three

subsystem states, and the spectrum at SR is formed by merging the three spectra. For example, in the small model, the ground state of SR consists of the D0 state of SiH3, the X3∑−g state of O2, and the 3P state of O, which makes a doublet, a quartet, and a sextet state, and each of them is triply degenerate because the 3P state of O is spatially triply degenerate while the D0 and X3∑−g states have no spatial degeneracy (only spin degeneracy). L

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C Present Address

Similarly, R and P are both combinations of two molecules, so each of their states is also a combination of the states of two subsystems. We list the energies of low-lying states of all molecules in Table 5. In this table, the energies refer to the overall system; for example, the energies listed in the SiH3 column are for SiH3 plus ground-state O2 and ground-state O; the energies listed for SiH3O2 refer to SiH3O2 plus ground-state O; and the energies listed for SiH3O refer to SiH3O plus ground-state O2. The merged energy lists for SR, R, and P are in Table 6. The reaction profile with all low-lying doublet and quartet states (states below 2 eV) for the atomic oxygen recombination reaction is shown in Figure 17. The distance to one of the oxygen atoms is infinite at the peroxyl defect R; as this oxygen atom comes closer to the peroxyl defect, other species such as D0INT, Q1-TS4, and D0-TS2 are formed. The O−O distance in Q1-TS4 is between R and D0-INT; however, the energy of the doublet ground state of Q1-TS4 is not between the energies of the doublet ground states of R and D0-INT because the O−O−O bond angle in Q1-TS4 is quite different from that in D0-INT (see Table 7). To map out the potential energy surface between R and D0-INT, we created a molecule called Q1-TS4-110 by changing the angle to 110°, which is the same as that in D0-INT. The continuous energy changes for both doublet and quartet ground states confirm the smoothness of the energy surfaces. Figure 17 shows a large number of energy levels within 2 eV of the ground state of SR. For the species P, all states that consist of O2 in excited state are marked with asterisks, including states at −2.97, −2.80, −0.83, 1.68, and 1.99 eV involving O2(a1Δ−g ) and states at −2.16, −2.33, and −0.19 eV involving O2(b1∑+g ). These states are all in the doublet spin manifold. With all of the low-lying states in the figure, it is clear that there are many possible doublet potential energy curves leading to the excited states of O2 in the atomic oxygen recombination reaction. Dynamics calculations would be quite interesting, and they could be carried out using a method suitable for a dense manifold of excited states.61

§

Centrus Energy, 350 Centrifuge Way, Oak Ridge, Tennessee 37830, United States. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS We appreciate helpful assistance by Dr. Xuefei Xu and helpful discussions with Dr. Hannah R. Leverentz and Dr. Maxim Makeev. This work was supported in part by the Air Force Office of Scientific Research under grant no. FA9550-12-1-0486. Part of the computation was performed as part of a Computational Grand Challenge project at the Molecular Science Computing Facility (MSCF) in the William R. Wiley Environmental Molecular Sciences Laboratory, a national scientific user facility sponsored by the U.S. Department of Energy’s Office of Biological and Environmental Research and located at the Pacific Northwest National Laboratory, operated for the Department of Energy by Battelle.



(1) Iler, R. K. The Chemistry of Silica; Wiley: New York, 1979. (2) Morey, G. W. The Properties of Glass; Reinhold: New York, 1954. (3) Silica: Physical Behavior, Geochemistry & Materials Applications; Heaney, P. J., Prewitt, C. T., Gibbs, G. V., Eds.; Mineralogical Society of America: Washington, DC, 1996. (4) Cozmuta, I. Molecular Mechanisms of Gas Surface Interactions in Hypersonic Flows, 39th AIAA Thermophysics Conference, Miami, FL, 2007; AIAA-2007−4046. (5) Balat-Pichelin, M.; Badie, J.; Berjoan, R.; Boubert, P. Recombination Coefficient of Atomic Oxygen on Ceramic Materials under Earth Re-entry Conditions by Optical Emission Spectroscopy. Chem. Phys. 2003, 291, 181−194. (6) Alfano, D.; Scatteia, L.; Monteverde, F.; Beche, E.; Balat, M. Microstructural Characterization of ZrB2-SiC based UHTC Tested in the MESOX Plasma Facility. J. Eur. Ceram. Soc. 2010, 30, 2345−2355. (7) Bedra, L.; Balat-Pichelin, M. J. H. Comparative Modeling Study and Experimental Results of Atomic Oxygen Recombination on SilicaBased Surfaces at High Temperature. Aerosp. Sci. Technol. 2005, 9, 318− 328. (8) Davis, G. S. Pathogenesis of Silicosis: Current Concepts and Hypotheses. Lung 1986, 164, 139−154. (9) Murashov, V.; Harper, M.; Demchuk, E. Impact of Silanol Surface Density on the Toxicity of Silica Aerosols Measured by Erythrocyte Haemolysis. J. Occup. Environ. Hyg. 2006, 3, 718−723. (10) Janossy, I.; Menyhard, M. LEED Study of Quartz Crystals. Surf. Sci. 1971, 25, 647−649. (11) Bart, F.; Gautier, M.; Duraud, J. P.; Henriot, M. (0110) α-Quartz Surface: a LEED, XANES and ELS Study. Surf. Sci. 1992, 274, 317−328. (12) Bart, F.; Gautier, M. A LEED Study of the (0001) α-Quartz Surface Reconstruction. Surf. Sci. 1994, 311, L671−L676. (13) Steurer, W.; Apfolter, A.; Koch, M.; Sarlat, T.; Sødergård, E.; Ernst, W. E.; Holst, B. The Structure of the α-Quartz (0001) Surface Investigated using Helium Atom Scattering and Atomic Force Microscopy. Surf. Sci. 2007, 601, 4407−4411. (14) Rignanese, G.-M.; De Vita, A.; Charlier, J.-C.; Gonze, X.; Car, R. First-principles Molecular-Dynamics Study of the (0001) α−Quartz Surface. Phys. Rev. B 2000, 61, 13250−13255. (15) Koudriachova, M. V.; Beckers, J. V. L.; de Leeuw, S. W. Computer Simulation of the Quartz Surface: a Combined ab initio and Empirical Potential Approach. Comput. Mater. Sci. 2001, 20, 381−386. (16) Du, Z. M.; de Leeuw, N. H. A Combined Density Functional Theory and Interatomic Potential-based Simulation Study of the Hydration of Nano-particulate Silicate Surfaces. Surf. Sci. 2004, 554, 193−210.

4. SUMMARY Several defects thought to be important in the recombination of atomic oxygen on silica surfaces have been characterized by using electronic structure calculations. Analysis of the diradical character of these defects leads to the conclusion that they can be modeled as weakly coupled doublet monoradical defect sites for studying atomic oxygen recombination and that several electronically excited states are accessible during this process. We also calculated energies along recombination reaction paths, and these calculations show that the energy profiles are barrierless so atomic oxygen recombination can readily occur through an autocatalytic cycle where the peroxyl defect plays a central role. Furthermore, we find that many potential energy surfaces leading to electronically excited O2 are energetically accessible.



ASSOCIATED CONTENT

S Supporting Information *

Coordinates of all optimized structures. This material is available free of charge via the Internet at http://pubs.acs.org.



REFERENCES

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. M

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

(38) Shavitt, I. In Methods of Electronic Structure Theory; Schaefer, H. F., III, Ed.; Plenum: New York, 1977; pp 189−225. (39) Butscher, W.; Shih, S. K.; Buenker, R. J.; Peyerimhoff, S. D. Configuration Interaction Calculations for the N2 Molecule and its Three Lowest Dissociation Limits. Chem. Phys. Lett. 1977, 52, 457−462. (40) Werner, H.-J.; Knowles, P. J. An Efficient Internally Contracted Multiconfiguration−reference Configuration Interaction Method. J. Chem. Phys. 1988, 89, 5803−5814. (41) Knowles, P. J.; Werner, H.-J. An Efficient Method for the Evaluation of Coupling Coefficients in Configuration Interaction Calculations. Chem. Phys. Lett. 1988, 145, 514−522. (42) Knowles, P. J.; Werner, H.-J. Internally Contracted Multiconfiguration-reference Configuration Interaction Calculations for Excited States. Theor. Chim. Acta 1992, 84, 95−103. (43) Langhoff, S. R.; Davidson, E. R. Configuration Interaction Calculations on the Nitrogen Molecule. Int. J. Quantum Chem. 1974, 8, 61−72. (44) Hirsch, G.; Bruna, P. J.; Peyerimhoff, S. D.; Buenker, R. J. Ab initio CI Study of the Stability and Electronic Spectrum of the HOCl Molecule. Chem. Phys. Lett. 1977, 52, 442−448. (45) Butscher, W.; Shih, S. K.; Buenker, R. J.; Peyerimhoff, S. D. Configuration Interaction Calculations for the N2 Molecule and its Three lowest dissociation limits. Chem. Phys. Lett. 1977, 52, 457−462. (46) Schwenke, D. W.; Steckler, R.; Brown, F. B.; Truhlar, D. G. Estimation of Higher-order Correlation Effects on the Potential Energy Surface for the F+H2 Reaction in the Saddle Point Vicinity. J. Chem. Phys. 1987, 86, 2443−2444. (47) Zhao, Y.; Truhlar, D. G. The M06 Suite of Density Functionals for Main Group Thermochemistry, Thermochemical Kinetics, Noncovalent Interactions, Excited States, and Transition Elements: Two New Functionals and Systematic Testing of Four M06-class Functionals and 12 other Functionals. Theor. Chem. Acc. 2008, 120, 215−241. (48) Zhao, Y.; Truhlar, D. G. A New Local Density Functional for Main-group Thermochemistry, Transition Metal Bonding, Thermochemical Kinetics, and Noncovalent Interactions. J. Chem. Phys. 2006, 125, 194101−18. (49) Lynch, B. J.; Zhao, Y.; Truhlar, D. G. Effectiveness of Diffuse Basis Functions for Calculating Relative Energies by Density Functional Theory. J. Phys. Chem. A 2003, 107, 1384−1388. (50) Zhang, Y.; Li, Z. H.; Truhlar, D. G. Computational Requirements for Simulating the Structures and Proton Activity of Silicaceous Materials. J. Chem. Theory Comput. 2007, 3, 593−604. (51) Le Page, Y. V. O. N.; Donnay, G. Refinement of the Crystal Structure of Low-quartz. Acta Crystallogr. 1976, B32, 2456−2459. (52) Peverati, R.; Truhlar, D. G. Quest for a Universal Density Functional: the Accuracy of Density Functionals Across a Broad Spectrum of Databases in Chemistry and Physics. Philos. Trans. R. Soc. London 2014, 372, 20120476. (53) McLean, A. D.; Chandler, G. S. Contracted Gaussian Basis Sets for Molecular Calculations. I. Second Row Atoms, Z = 11−18. J. Chem. Phys. 1980, 72, 5639−5648. (54) Andersson, K.; Borowski, P.; Fowler, P. W.; Malmqvist, P. Å.; Roos, B. O.; Sadlej, A. J. Electric Properties of the Ozone Molecule. J. Chem. Phys. Lett. 1992, 190, 367−373. (55) Celani, P.; Werner, H.-J. Multireference Perturbation Theory for Large Restricted and Selected Active Space Reference Wave Functions. J. Chem. Phys. 2000, 112, 5546−5557. (56) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; et al. Gaussian 09, revision D.01; Gaussian, Inc.: Wallingford, CT, 2009. (57) Werner, H.-J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Schütz, M. Molpro: A General-purpose Quantum Chemistry Program Package. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2012, 2, 242−253. (58) Nishino, M.; Yamanaka, S.; Yoshioka, Y.; Yamaguchi, K. Theoretical Approaches to Direct Exchange Couplings between Divalent Chromium Ions in Naked Dimers, Tetramers, and Clusters. J. Phys. Chem. A 1997, 101, 705−712.

(17) Malyi, O. I.; Kulish, V. V.; Persson, C. In Search of New Reconstructions of (001) α-Quartz Surface: a First Principles Study. RSC Adv. 2014, 4, 55599−55603. (18) Norman, P.; Schwartzentruber, T. E.; Leverentz, H.; Luo, S.; Meana-Pañeda, R.; Paukku, Y.; Truhlar, D. G. The Structure of Silica Surfaces Exposed to Atomic Oxygen. J. Phys. Chem. C 2013, 117, 9311− 9321. (19) Kulkarni, A. D.; Truhlar, D. G.; Srinivasan, S. G.; van Duin, A. C. T.; Norman, P.; Schwartzentruber, T. E. Oxygen Interactions with Silica Surfaces: Coupled Cluster and Density Functional Investigation and the Development of a New ReaxFF Potential. J. Phys. Chem. C 2013, 117, 258−269. (20) Narayanasamy, J.; Kubicki, J. D. Mechanism of Hydroxyl Radical Generation from a Silica Surface: Molecular Orbital Calculations. J. Phys. Chem. B 2005, 109, 21796−21807. (21) He, Y.; Cao, C.; Wan, Y.-X.; Cheng, H.-P. From Cluster to Bulk: Size Dependent Energetics of Silica and Silica-Water Interaction. J. Chem. Phys. 2006, 124, 024722−5. (22) Rutigliano, M.; Zazza, C.; Sanna, N.; Pieretti, A.; Mancini, G.; Barone, V.; Cacciatore, M. Oxygen Adsorption on β-Cristobalite Polymorph: Ab Initio Modeling and Semiclassical Time-Dependent Dynamics. J. Phys. Chem. A 2009, 113, 15366−15375. (23) Zazza, C.; Rutigliano, M.; Sanna, N.; Barone, V.; Cacciatore, M. Oxygen Adsorption on β-Quartz Model Surfaces: Some Insights from Density Functional Theory Calculations and Semiclassical TimeDependent Dynamics. J. Phys. Chem. A 2012, 116, 1975−1983. (24) Malyi, O. I.; Kulish, V. V.; Persson, C. In Search of New Reconstructions of (001) α-Quartz Surface: A First Principles Study. RSC Adv. 2014, 4, 55599−55603. (25) Hochstrasser, G.; Antonini, J. F. Surface States of Pristine Silica Surfaces: I. ESR Studies of Es′ Dangling Bonds and of CO2− Adsorbed Radicals. Surf. Sci. 1972, 32, 644−664. (26) Griscom, D. L.; Friebele, E. Fundamental Defect Centers in Glass: 29 Si Hyperfine Structure of the Nonbridging Oxygen Hole Center and the Peroxy Radical in a-SiO2. Phys. Rev. B 1981, 24, 4896−4898. (27) Costa, D.; Fubini, B.; Giamello, E.; Volante, M. A Novel Type of Active Site at the Surface of Crystalline SiO2 (α-Quartz) and its Possible Impact on Pathogenicity. Can. J. Chem. 1991, 69, 1427−1434. (28) Warren, W. L.; Kanicki, J.; Rong, F. C.; Poindexter, E. H. Paramagnetic Point Defects in Amorphous Silicon Dioxide and Amorphous Silicon Nitride Thin Films II. J. Electrochem. Soc. 1992, 139, 880−889. (29) O’Reilly, E. P.; Robertson, J. Theory of Defects in Vitreous Silicon Dioxide. Phys. Rev. B 1983, 27, 3780−3795. (30) Baranowski, J. M.; Strzalkowski, I.; Marczewski, M.; Kowalski, M. π-Bonded Model of an Oxygen-vacancy Center in SiO2. J. Appl. Phys. 1987, 61, 2904−2909. (31) Pederson, M. R.; Harrison, J. G.; Klein, B. M. Density Functional Based Studies of Oxygen Vacancies In Crystalline Silicon Dioxide. Mater. Res. Soc. Proc. 1987, 105, 229−234. (32) White, J. D.; Copeland, R. A.; Marschall, J. Singlet Molecular Oxygen Formation by O-Atom Recombination on Fused-Quartz Surfaces. J. Thermophys. Heat Transfer 2015, 29, 24−36. (33) Lopes, P. E. M.; Demchuk, E.; Mackerell, A. D., Jr. Reconstruction of the (011) Surface on α-Quartz: A Semiclassical Ab initio Molecular Dynamics Study. Int. J. Quantum Chem. 2009, 109, 50−64. (34) Ling, S.; El-Sayed, A. M.; Lopez-Gejo, F.; Watkins, M. B.; Afanas’ev, V. V.; Shluger, A. L. A Computational Study of Si−H Bonds as Precursors for Neutral E′ Centres in Amorphous Silica and at the Si/ SiO2 Interface. Microelectron. Eng. 2013, 109, 310−313. (35) Mortier, W. J.; Sauer, J.; Lercher, J. A.; Noller, H. Bridging and Terminal Hydroxyls. A Structural Chemical and Quantum Chemical Discussion. J. Phys. Chem. 1984, 88, 905−912. (36) Kohn, W.; Sham, L. J. Self-consistent Equations Including Exchange and Correlation Effects. Phys. Rev. 1965, 140, A1133−A1138. (37) Knowles, P. J.; Werner, H.-J. An Efficient Second-order MC SCF Method for Long Configuration Expansions. Chem. Phys. Lett. 1985, 115, 259−267. N

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C (59) Bachler, V.; Olbrich, G.; Neese, F.; Wieghardt, K. Theoretical Evidence for the Singlet Diradical Character of Square Planar Nickel Complexes Containing Two o-Semiquinonato Type Ligands. Inorg. Chem. 2002, 41, 4179−4193. (60) Sinnecker, S.; Neese, F.; Noodleman, L.; Lubitz, W. Calculating the Electron Paramagnetic Resonance Parameters of Exchange Coupled Transition Metal Complexes using Broken Symmetry Density Functional Theory: Application to a MnIII/MnIV Model Compound. J. Am. Chem. Soc. 2004, 126, 2613−2622. (61) Valero, R.; Truhlar, D. G. Photochemistry in a Dense Manifold of Electronic States: Photodissociation of CH2ClBr. J. Chem. Phys. 2012, 137, 22A539.

O

DOI: 10.1021/acs.jpcc.5b00120 J. Phys. Chem. C XXXX, XXX, XXX−XXX