Toward a Process-Based Molecular Model of SiC ... - ACS Publications

Article pubs.acs.org/JPCC

Toward a Process-Based Molecular Model of SiC Membranes. 1. Development of a Reactive Force Field Saber Naserifar,† Lianchi Liu,‡ William A. Goddard, III,§ Theodore T. Tsotsis,† and Muhammad Sahimi*,† †

Mork Family Department of Chemical Engineering & Materials Science, University of Southern California, Los Angeles, California 90089-1211, United States ‡ Computational Chemistry Group, School of Chemistry and Chemical Engineering, Shanghai Jiao Tong University, Min Hang, Shanghai 200240, China § Materials and Process Simulation Center, California Institute of Technology, Pasadena, California 91125, United States ABSTRACT: A broad class of important materials, such as carbon molecular sieves, silicon carbide (SiC), and silicon nitride, are fabricated by temperature-controlled pyrolysis of preceramic polymers. In particular, the fabrication of SiC membranes by pyrolysis of a polymer precursor that contains Si is quite attractive for separation of hydrogen from other gases. It has been quite difficult to extract atomistic-scale information about such SiC membranes since they are amorphous. In principle, ab initio quantum mechanics (QM) can provide information about the structure of the amorphous systems. However, to determine the structure of the SiC membrane layer one should capture in the simulations the various reactive processes involved in forming the layer. This requires QM simulations on systems with about 3000 atoms per cell at temperature of 1200 K for microseconds, which are far beyond the current QM capabilities. Instead, this paper extends the ReaxFF reactive force field, validated for high temperature reactions of other materials, to describe the processes involved in the thermal decomposition of hydridopolycarbosilane (HPCS) to form SiC nanoporous membranes. First, we carry out QM calculations on models meant to capture important reaction steps and structures. Then, we develop a model of the HPCS polymer and utilize ReaxFF to describe the thermal degradation and decomposition of the polymer as the system is heated in the molecular dynamics (MD) simulations. Analysis of the pyrolysis studies and their results leads to various quantities that can be compared with experimental data. Good agreement is found between the data and the results of the MD simulations.

1. INTRODUCTION Many important materials are fabricated by temperaturecontrolled pyrolysis of preceramic polymers. Examples include carbon molecular sieves, silicon carbide (SiC), silicon nitride (Si3N4), and aluminum nitride (AlN).1−6 Such materials have important applications in the design and synthesis of fuels and chemicals as separation and reactive media, such as adsorbents, membranes, and catalysts as well as various types of sensors. The fabrication of such materials involves a series of tailored reactions via which the original polymer precursor is converted into the final ceramic. Consider, for example, fabrication of nanoporous membranes. Growing interest in the hydrogen economy has motivated research on inorganic, hydrogen-permselective membranes for use in the processes related to H2 production that take place at high temperatures and pressures. Due to its many unique properties, such as high thermal conductivity,7 thermal shock resistance,8 biocompatibility,9 resistance to acidic and alkali environments,10 chemical inertness, and high mechanical strength,11,12 a promising material for fabricating such inorganic membranes is SiC, for which two distinct approaches have been developed. One involves the use of chemical-vapor deposition/ chemical-vapor infiltration techniques, while the second © 2013 American Chemical Society

approach is based on a dip-coating technique combined with the pyrolysis of a polymer precursor. Both approaches have been pursued,1−5,13,14 but the first one suffers from certain difficulties13,14 that make it less attractive. The second method utilizes the pyrolysis of a polymer precursor, hydridopolycarbosilane (HPCS), or a partially allyl-substituted hydridopolycarbosilane (AHPCS). The role of the allyl group is to slightly increase the carbon content of the final ceramic and enable the polymer to cross-link at lower temperatures. The AHPCS is curable in the presence of argon, rather than oxygen that introduces Si−O−C bonds in the polycarbosilane (PCS)-derived ceramics that are thermally and hydrothermally unstable. The primary pyrolysis product of both the HPCS and AHPCS is a SiC ceramic with Si−C bonds as well as some Si−Si and C−C ones (see below). To fabricate SiC nanoporous membrane using the dip-coating technique and the pyrolysis of the AHPCS,1−5 porous SiC support tubes are prepared using uniaxial cold-pressing of β-SiC powder, together with the appropriate sintering aids. The pores of the support are, however, too large to be effective for Received: August 6, 2012 Revised: January 23, 2013 Published: January 28, 2013 3308

dx.doi.org/10.1021/jp3078002 | J. Phys. Chem. C 2013, 117, 3308−3319

The Journal of Physical Chemistry C

Article

chemistry of the AHPCS and a hyperbranched polycarbosilane. The results indicated that H2 was the main product of the pyrolysis whose loss caused the polymers to cross-link above 525 K. Amorphous SiC was produced after the polymer was heated up to 1275 K, which was then transformed into a nanocrystalline SiC at about 1875 K. A quite similar process was reported by Interrante et al.20 for the pyrolysis of the HPCS. Its related and compositionally equivalent linear polymer analog, polysilaethylene (PSE),[SiH2CH2]n, was studied by Liu et al.21 Both polymers produced H2 above 575 K by a (1,1)-elimination process to form a cross-linked network. The process, also called unimolecular elimination, is a model for explaining a particular type of chemical elimination reaction that consists of ionization and deprotonation. More H2 and an amorphous SiC were produced between 725 and 1275 K. Mass spectrometric analysis of the gas production during the pyrolysis of the PCS, such as the HPCS, indicated that although most of the hydrogen production occurs below 1075 K, a significant amount is also produced21 above this temperature. The loss of hydrogen results in further condensation of the amorphous material. A reaction pathway was postulated for cross-linking and pyrolysis of the HPCS in which both (1,1)-H2 elimination and intramolecular H-transfer reactions lead to highly reactive silylene intermediates, producing eventually Si−Si bonds that rapidly rearrange to Si−C at such temperatures and form Si−C interchain cross-links. To obtain a fundamental understanding of such pyrolysis processes we carried out detailed atomistic-scale simulations. Due to the large size of the polymers, modeling of their pyrolysis using ab initio methods is computationally prohibitive. Instead, we develop the ReaxFF reactive force field to match our quantum mechanical (QM) calculations on small systems and then use ReaxFF to carry out reactive molecular dynamics (MD) simulation of the pyrolysis processes. There are currently several MD methods that can simulate reactive systems. These are reviewed by Farah et al.22 Three of the better known of such methods are as follows. ReaxFF reactive force field was developed by van Duin, Goddard and co-worker23 and applied to a large number of reactions in metallic, ceramic, semiconductor, and polymer systems that we study here. Another method is the reactive empirical bond-order potentials developed by Brenner and co-workers24 and applied to carbon materials. Nyden and Noid25 used a conventional nonreactive FF but dynamically varied it to describe chemical reactions. Chemical reactions can also be simulated by the tight binding and density functional methods, but the focus of this paper is on MD simulation of chemical reactions. ReaxFF bridges the gap between the computational methods based on quantum chemical (QC) and empirical force fields. Although the QC methods are, in principle, applicable to all materials, their computational expense makes them inapplicable for large systems (with more than a few hundred atoms). ReaxFF, on the other hand, enables one to carry out MD simulations of large-scale reactive chemical systems with accuracy near QC at a fraction of the cost. It also allows for accurate description of bond breaking and bond formation, and, therefore, it may be used in the MD simulation of pyrolysis of a polymer. The present paper develops the ReaxFF aimed at simulating pyrolysis of the HPCS. ReaxFF has been previously developed for hydrocarbons and hydrocarbon oxidation,23,26 metaloxides27 and metal hydrides,28 Si/Si oxides,29 aluminum/ aluminum oxides,30 nitramines,31 and a variety of complex chemical systems.32−35 There has also been a number of studies for simulating pyrolysis and thermal degradation of polymers

separating gasesous (and even liquid) mixtures. To endow the support with high selectivity for the separation of fluid mixtures, one deposits on the support (via dip-coating) a thin layer of the preceramic polymer precursor, which leaves behind a thin nanoporous SiC film after pyrolysis.1−5 During pyrolysis bonds in the polymer break, releasing various atoms and radicals that, in turn, react with each other and produce such gases as H2, CH4, and SiH4. Understanding equilibrium and nonequilibrium phenomena, such as adsorption, reaction, and transport of fluid mixtures in SiC membranes entails having an accurate model of the materials’ pore space, which can be utilized for designing the optimal pore structure for the separation process. One of the most important ingredients for developing such a model is a fundamental understanding of the pyrolysis of the polymer precursors, such as the HPCS and AHPCS, and the evolution of the structure of the thin film deposited on the support. It is, therefore, instructive to briefly review the experimental results for pyrolysis of the PCS to SiC ceramics, as they will be compared with the results of the molecular simulations presented in this paper and its sequel. Hasegawa and Okamura15 studied pyrolysis of several types of the PCS using various techniques and suggested six stages for their conversion to SiC ceramic. In the first stage, up to about 625 K, gas generation is very small and degradation of low-molecular weight molecules occurs. The average molecular weight of the PCS begins to change during the second stage from 625 to 825 K, but the gas production is still small. Dehydrogenation of the Si−H bonds occurs to some extent, and, due to cross-linking, the Si−Si bonds also appear and some H2 is produced. Depending on the type of the PCS, some alkanes may form that are the result of dehydrocarbonation condensation. In the third stage from 825 K to about 1075 K the Si−H and C−H bonds break readily, resulting in the formation of many free hydrogen radicals that bond together and generate H2. For the PCS that have methyl group as side chains, CH4 may also be produced as a result of the decomposition of the methyl group and the subsequent reaction with free hydrogen radicals. With increasing temperature the connectivity of the Si−C−Si bonds increases as a result of dehydrogenation of the Si−H and C−H bonds and demethanation of SiCH3, producing eventually an inorganic amorphous residue. At about 975 K the residue includes CC bonds that are produced by dehydrogenation of the −CH−CH bonds16,17 as well as hydrogen atoms. The fourth and fifth stages happen between 1075 and 1475 K, with the former known as a transition stage in which the structure is still amorphous. The remaining hydrogen atoms decompose with increasing temperature to 1275 K, above which nucleation of Si−C begins and the amorphous SiC is transformed into a crystalline structure. The carbon content of the system in this stage is more than its Si content and results in a nonstoichiometric crystalline structure with crystal sizes of 2−3 nm. Finally, above 1475 K crystal growth occurs in which a stoichiometric ratio between Si and C is obtained by burning the extra carbon. Tang et al.18 pyrolyzed the PCS in H2 atmosphere and reported the initiation of the process by cleavage of the Si−H bonds below 825 K and of the Si−CH3 bonds that break significantly above 1175 K. The properties of SiC membranes and fibers at high temperatures are affected negatively by the presence of extra Si. Therefore, Tang et al. used a hydrogen atmosphere to achieve a Si/C ratio close to one in the final ceramic fiber via reaction of the free radicals with H2 and generation of CH4. Interrante and Moraes19 studied the pyrolysis 3309



Article

using ReaxFF. Chenoweth et al.35 used ReaxFF to investigate the thermal decomposition of Poly(dimethylsiloxane). Desai et al.36 studied the pyrolysis of phenolic polymer in the presence of carbon nanotubes, while Saha and Schatz studied37 thermal decomposition of polyacrylonitrile in the carbonization process. The rest of this paper is organized as follows. In the next section we describe the essentials of ReaxFF. Section 3 explains the QC computations for estimating the parameters of ReaxFF for the pyrolysis process of HPCS. In Section 4 we describe the generation of a molecular model of the HPCS, while Section 5 provides the details of minimization of the energy of HPCS’ molecular model and its pre-equilibration, before it is heated up. Section 6 describes the results of MD simulation of thermal decomposition of the polymer. The paper is summarized in the last section.

energy, however, contributes weakly to ReaxFF, and, thus, only 9 terms out of 26 are used in the ReaxFF that is developed here.

3. QUANTUM CHEMICAL CALCULATIONS The parameters of ReaxFF were estimated such that they correctly reproduce the first principles quantum mechanical interactions in the HPCS and provide a way of extending the first-principles’ accuracy to large length- and time-scale simulations. This was achieved by adding the SiC crystal chemistry and the QC data for systems relevant to the HPCS to the ReaxFF training set for Si materials.29 The QC data for nonperiodic systems were obtained from the density-functional theory (DFT) calculations, carried out using Jaguar (version 5.5),38 B3LYP functional,38,39 and Pople’s 6-311G** basis set.41,42 The Mulliken charges were also computed using the 6-31G** basis set developed by Pople and co-workers.43,44 To obtain the equation of state for the β-SiC crystal structure, the DFT calculations were carried out with Quantum-Espresso.45 The generalized gradient approximation (GGA), proposed by Perdew et al.,46 was used for the exchange-correlation energy, and ultrasoft pseudopotentials were utilized to replace the core electrons. We used the Perdew’s implementation of the GGA47 with a kinetic energy cutoff of 320 eV. The Monkhorst-Pack scheme48 was used to generate the k-space grid with a 0.1 Å−1 spacing. In general, the accuracy of the results are similar to or better than that provided by the parametrized model number 3 (PM3),49 a semiempirical method for quantum mechanical computation of the electronic structure of molecules. ReaxFF is about 100 times faster than PM3, which in turn is about 100 times faster than the QC calculations.

2. THE REACTIVE FORCE FIELD REAXFF The original ReaxFF was developed23,29 to make practical MD simulation of large-scale reactive compounds with thousands of atoms. The conventional empirical force fields used in the standard MD simulations assume rigid connectivity of the atoms and cannot describe the complex chemistry of real materials, but ReaxFF may be used to examine such chemical processes as polymer decomposition. It uses a general relationship between bond distance and bond order, on the one hand, and between bond order and bond energy, on the other hand, which lead to proper dissociation of the bonds and atoms. Its development is based on the fundamental assumption that the bond order between a pair of atoms can be obtained directly from the interatomic distance rij between atoms i and j. Other valence terms present in ReaxFF are also defined in terms of the same bond order, so that all the terms vanish smoothly as the bonds break. In addition, ReaxFF contains Coulomb and van der Waals potentials for describing nonbond interactions between all the atoms. They are shielded at short-range, so that they become constant as the distance between two atoms vanishes. Similar to nonreactive empirical force fields, ReaxFF computes the total energy E of a system as the sum of partial energies, each representing a particular contribution to E: E = E bond + Eover + Eunder + Eval + Econ + Etor + Econj + EvdWaals + ECoulumb

(1)

Here, Ebond represents the bond energy, and Eover is an overcoordination penalty term to account for a degree of overcoordination number that may remain in the molecule since we know, for example, that the bond order of C should not exceed 4. Similarly, Eunder represents the energy contribution due to the resonance of the π-electron between attached undercoordinated atomic centers. The contribution of the valence angle energy is given by Eval, Econ, and Econj represent, respectively, the three- and four-body conjugation energies, and Etors is the torsional rotation barriers energy. Each partial energy term is parametrized by a number of parameters. Their parametrized functional forms are long and given elsewhere29 and, thus, are not repeated here. The overall number of parameters for the HPCS is 438, of which 39 are general, 50 are computed specifically in this paper for the HPCS, while the remaining 349 parameters are those computed previously23 for ReaxFF for poly(dimethylsiloxane). Since we have three types of atoms in the HPCS, namely, Si, C, and H, then there are 26 combinations of torsion angles. The torsion

Figure 1. Comparison of the energies computed by ReaxFF and QC for (a) Si−C and (b) double-bond dissociations. 3310



Article

Figure 2. Comparison of distortion energies computed by ReaxFF and QC for the angles in (a) Si−C−Si, (b) Si−C−H, (c) C−Si−H, (d) C−C−Si, (e) C−Si−C, and (f) C−Si−Si triplets in the HPCS.

interest fixed, in order to compute the angle-bending energies relative to the optimal geometry. The results obtained by the DFT computations are compared in Figure 2 with those that ReaxFF yielded. The charge distributions in ReaxFF were calculated using the electronegativity equalization method (EEM).49 The Mulliken charge distributions obtained from the DFT calculations were used to optimize the EEM charge parameters. The partial charges for SiH3CH2SiH3, CH3SiH2CH3, and C−(SiH2−CH2)3 (six-membered ring), computed by the DFT and ReaxFF, are compared in Figure 3. They all indicate good agreement between the two sets of results, indicating the accuracy of ReaxFF. To sample the possible reactions in Si−C systems relating to the decomposition of the HPCS, the DFT- and ReaxFFcomputed reaction energies were computed and compared. The reactions include not only the energetically favorable clusters but also the strained and energetically less favorable clusters that test the ability of ReaxFF to represent both cases. Cleavage of the Si−H and C−H bonds are important in the degradation process of the HPCS, when exposed to heat. The energies computed by ReaxFF and DFT for these systems are listed in Table 1. As pointed out earlier, our goal is to develop ReaxFF in order to carry out MD simulation of the structure that results from

The DFT calculations were utilized to compute the geometries and energies of the species and complexes in the HPCS as well as the products of the pyrolysis. The ReaxFF parameters were determined by combining the Si training set with the results from the DFT calculations, followed by optimization of the parameters, in order to minimize the differences between the energies calculated by the DFT and ReaxFF. 3.1. Optimization of ReaxFF. To optimize the Si−C bond energy, we used the dissociation pathway for the bond in H3Si− CH3, shown in Figure 1(a), and the SiC bond in H2SiCH2, depicted in Figure 1(b). To optimize the ReaxFF parameters, the DFT results for the singlet case were used around the equilibrium bond distance. The valence angle parameters were optimized to represent the angle-bending energies (between three atoms), obtained from the DFT calculations on various clusters of atoms. There are six valence angles in the HPCS that need to be accounted for, namely, Si−C−Si, Si−C−H, C−Si−H, C−C−Si, C−Si−C, and C−Si−Si. For the first three the triplets H3Si− CH2−SiH3 and for the rest H3C−CH2−SiH3, H3C-SiH2−CH3, and H3Si-SiH2−CH3 were used, respectively, as the representative molecules for valence angles. For each angle, the geometry of the representative molecule was optimized with the angles of 3311



Article

Figure 3. Comparison of the Mulliken charge distributions computed by the DFT (red) and ReaxFF (blue).

Figure 4. The equations of state (compression and expansion) for the three crystal structures of SiC. (a) β-SiC, (b) 4H-SiC, and (c) α-SiC, as computed by the DFT and ReaxFF.

Table 1. Comparison of the Reaction Energies E (in kcal/ mol), Computed by the DFT and ReaxFF reaction SiH2(CH3)2 + CH2(SiH3)2→ c-(SiH2−CH2)3 + 2H2 SiH2(CH3)2 + CH2(SiH3)2 → (CH3)2SiH−CH(SiH3)2 + H2 CH3(SiH2CH2)2SiH3 → c-(SiH2−CH2−SiH2−CH2−CH)− SiH3 + H2 CH3(SiH2CH2)2SiH3 → c-(SiH2−CH2−SiH2−CH)− SiH2CH3 + H2 CH3(SiH2CH2)2SiH3 → c-(SiH2−CH2−SiH2−CH2−SiH)− CH3 + H2

E E (DFT) (ReaxFF) 17.15 7.86 27.84

12.75 5.15 22.08

24.43

25.81

13.00

17.88

two dimensions but differ in the third one. Thus, the ability of ReaxFF for representing the equation of state of the SiC crystal was tested against the aforementioned crystalline structures of SiC. The Crystal Builder module52 in Ceruis2 and the unit cell parameters reported by others52,53 were used to construct the SiC crystals. For each crystalline form of SiC, the quantum energies were computed for a broad range of volumes (densities), describing both expansion and compression. Figure 4 compares the ReaxFF results against those obtained by the DFT calculations, indicating that ReaxFF correctly describes the equation of state for the SiC structures. There is, however, a rather large difference between the energies computed by ReaxFF and DFT for 4H-SiC at the expansion part of α-SiC. This does not, however, cause a problem since the most relevant and physically interesting parts are the equilibrium volumes from which the relative phase stabilities are deduced. Optimization of β-SiC crystal with

pyrolysis of the HPCS. The structure is a condensed phase of SiC that may be either in crystalline or amorphous form, particularly at high pressures and temperatures. There are, in fact, more than 250 crystalline structures of SiC and three major SiC polytypes known51 as (β) 3C-SiC, (α) 6H-SiC, and52 4H-SiC. Polytypes refer to a family of similar crystalline structures and represent variations of the same chemical compound that are identical in 3312



Article

Table 2. Atom Parameters of ReaxFF, Computed for the HPCSa

a

atom

χ (eV)

η (eV)

γ (Å)

pov/uv

pval3

pboc3

pboc4

pboc5

C H Si

5.725 3.819 2.912

6.923 9.883 7.427

0.871 0.762 0.629

−2.898 −15.768 −5.849

4.782 3.350 2.285

33.563 3.246 8.0583

5.713 3.865 4.511

11.996 1.000 0.838

The first three from the right are the EEM parameters, while the rest are the bond-order correction coefficients.

Table 3. Bond Energies De (in kcal/mol) and Bond-Order Parameters of ReaxFF, Computed for the HPCSa

a

pair

De,σ

De,π

pbe1

pbo1

pbo2

C−H Si−H C−Si

159.8520 238.3492 100.3210

0.000 0.000 25.714

−0.4646 −0.7467 0.1000

−0.0098 −0.0438 −0.1100

8.5954 6.8851 6.0000

pbo3

−0.2900

pbo4

rσ

rπ

8.800

1.0386 1.3266 1.640

1.570

The parameter rσ is in Å.

Table 4. van der Waal Parameters of ReaxFF, Computed for the HPCS

Table 6. Torsion Parameters of ReaxFF for Quadrupletes, Computed for the HPCSa

pair

σ (Å)

ε (kcal/mol)

α

γvdW

quadruplet

V1

V2

V3

ptor1

pcot1

C−H Si−H C−Si

1.7207 1.4723 1.9000

0.0431 0.0640 0.0600

10.3632 13.8391 11.0000

3.2757 2.6229 1.7008

X−Si−C−X C−C−C−C C−C−C−H H−C−C−H H−Si−Si−H H−Si−Si−Si

0.0000 0.0000 0.0000 0.0000 0.0000 0.0000

0.6675 38.9174 49.1001 34.0265 0.0000 0.0000

0.0000 0.3649 0.2713 0.3804 0.0640 0.1587

−8.2352 −8.2931 −8.5284 −6.3917 −2.4426 −2.4426

0.0000 −2.0127 −1.5309 −0.9965 0.0000 0.0000

ReaxFF yielded a density of 3.16 g/cm3, only 1.6% smaller than the experimental density 54 of 3.21 g/cm 3, measured at 300 K. Values of the ReaxFF parameters that we computed for the atoms, bond energies and bond orders, van der Waals, angle, and torsion contributions to the HPCS and SiC crystals are given, respectively, in Tables 2−6. They may be used directly in any MD simulator, and in particular those that are freely available such as the LAMMPS package.55 One goal of the development of ReaxFF is to estimate parameters that are transferable to other Si-based materials and polymers. Since the complete training set has been developed for Si-containing polymers, the new parameters of ReaxFF that we have computed, in addition to be relevant to the HPCS, make the FF general and capable of describing other polymers that contain other combinations of Si with C and H. We point out that the important point in using ReaxFF is to enable chemical reactions to proceed with appropriate energy barriers. Thus, ReaxFF uses a bond energy-to-bond order-tobond distance relation that can capture double and single bonds as well as lower-order bonds during reactions. As a result, the compressibility curves computed by ReaxFF are generally not as accurate as one would like them to be. For covalent systems calculations of accurate second-order properties are best carried out by using ReaxFF to obtain the structure and then using a nonreactive FF to analyze the distorions near an equilibrium structure.

a

Vi is in kcal/mol.

4. ATOMISTIC MODEL OF THE POLYMER To test the accuracy of ReaxFF we study thermal decomposition of the HPCS. The goal is to see whether the reaction products that MD simulation using reaxFF determines agree with what is known experimentally, described earlier. To do so, we first generate an atomistic model of the HPCS. A modified version of the self-avoiding walk (SAW) method of Theodorou and Suter56−58 was used to generate the initial structure of the HPCS. A cubic simulation cell was used in which three connected atoms of the HPCS’ backbone were placed in random orientations. The polymer was then constructed by adding one atom at a time in a stepwise fashion, using the modified SAW algorithm. The allowed rotational states of the successive bonds between adjacent atoms were determined from the probability distribution functions that were governed by energy considerations (see below). To suppress the surface effects, periodic boundary conditions were used. The polymer-consistent force field (PCFF) was utilized in order to generate the atomistic model of the HPCS in which the total potential energy E of the polymer is given by E = Es + Eθ + Eϕ + Enonbond (2)

Table 5. Parameters for the Triplets in ReaxFF, Computed for the HPCSa

a

triplet

θ0,0

pval1

pval2

pval4

pval7

pcoa1

ppen1

C−C−Si C−Si−C Si−Si−C Si−C−Si Si−C−H C−Si−H

72.52 69.33 69.33 69.33 72.59 72.59

22.3583 18.9270 19.6964 18.9270 13.8347 14.8347

2.0393 2.1333 2.0703 2.0703 2.4952 2.4952

1.040 1.040 1.040 1.040 1.040 1.040

1.0031 1.0031 1.0031 1.0031 1.0000 1.0000

0.000 0.000 0.000 0.000 0.000 0.000

0.000 0.000 0.000 0.000 0.000 0.000

The angles are in degrees. 3313



Article

(Other methods of simulating Coulombic interactions are also available.62) It would perhaps be somewhat more accurate to also include corrections due to the introduction of the cutoffs, but for the main purpose of the present work, namely, creating the initial atomistic model of the HPCS to be used with ReaxFF for modeling its pyrolysis, the cutoffs without any corrections seemed to suffice (see below). The parameters associated with Enonbond are listed in Table 10. We emphasize that the PCFF was used only to generate the initial structure of the HPCS.

with 4

Es =

4

∑ K i(l − l0)i ,

Eθ =

i=1

∑ Hi(θ − θi)i , i=2

3

Eϕ =

∑ Di(1 − cos iϕ)

(3)

i=1

and Enonbond =

∑ i,j

qiqj rij

+

⎡ ⎛ ⎞9 ⎛ ⎞6 ⎤ ⎢2⎜ σij ⎟ − 3⎜ σij ⎟ ⎥ ε ∑ ij⎢ ⎜ ⎟ ⎝ rji ⎠ ⎥⎦ ij ⎣ ⎝ rij ⎠

Table 10. Values of the Size and Energy Parameters of the Nonbond Part of the PCFF (4)

Here, Es is the energy associated with bond stretching, with l0 being a bond’s equilibrium length and l its actual length at any time during the simulations, and Ki is a stretching-force constant. Values of Ki and the equilibrium lengths l0 are listed in Table 7 for

l0 (Å)

K2

K3

K4

C−C C−H C−Si Si−H

1.417 1.101 1.899 1.478

470.836 345.000 189.653 202.779

−627.618 −691.890 −279.421 −305.360

1327.635 844.60 307.513 280.268

all the atoms. Eθ represents the energy associated with the changes in the bonds’ angles, where θ0 is the equilibrium angle of a bond, θ is its angle at any given time during the simulations, and Hi is the corresponding force constant. Numerical values of Hi and θ0 are listed in Table 8. The contribution of the torsional Table 8. Values of the Force Constants Hi and the Equilibrium Angles θ0 for Various Triplets, Used with the PCFF triplet

θ0 (deg)

H2

H3

H4

C−Si−C C−Si−H H−C−H H−Si−H S−C−H

113.185 112.098 107.660 108.605 112.035

36.206 36.483 39.641 32.541 28.772

−20.394 −12.809 −12.921 −8.316 −13.952

20.017 0.000 −2.432 0.000 0.000

σi (Å)

εi (kcal/mol)

C H Si

4.010 2.995 4.450

0.054 0.020 0.190

To obtain the equilibrated atomistic model of the polymer, energy minimization and MD simulations were utilized. The density of the polymer at the beginning of the energy minimization was set to be 0.15 g/cm3, so as to avoid ring spearings and catenations. The total energy E of the polymer was then minimized using the conjugate-gradient method.63 The resulting structure was then compressed using MD simulation in the (NPT) ensemble at pressures between 0.5 and 0.6 GPa, in order to increase the polymer density and bring it to as close to the experimental value of 0.98 g/cm3 as possible. Energy minimization and short-duration MD simulations in the (NVT) ensemble at 1000 K were then followed, in order to further relax the material. Then, MD simulations were carried out in the (NPT) ensemble for several ns at 1 atm and 300 K. The time step was 1 fs. The pressure was controlled by the Andersen method,64 while the temperature was held fixed using the Nosé-Hoover thermostat.65,66 We used W = 20 amu for masslike parameter of the Andersen barostat. For the Nosé-Hoover thermostat, the parameter that was used to scale the fictitious mass was 1. Figure 5 presents the change with the time of the potential energy of the HPCS at 300 K during MD simulations in the

Table 7. Values of the Force Constants Ki and the Equilibrium Lengths l0 for Various Pairs, Used with the PCFF pair

atom

forces is represented by Eϕ, with Di being a force constant, and ϕi being the dihedral angle. Values of the parameters are given in Table 9 for all the quadruplets. The parameters σi and εi are, Table 9. Values of the Force Constants Di for the Quadruplets, Used with the PCFF quadruplet

D1

D2

D3

C−Si−C−H H−C−Si−H

0.00 0.00

0.00 0.00

−0.0657 −0.0657

respectively, the Lennard-Jones (LJ) size and energy parameters for atom i, and qi is the partial charge of i. The nonbond interactions were cut off at an interatomic distance of 12.5 Å, while a spline switching function was used between 9.5 and 12.5 Å. The cutoffs were selected carefully, based on extensive preliminary simulations, as well as our previous extensive experience with modeling of Coulombic interactions using the multipole expansion59 and the Ewald summation60,61 techniques.

Figure 5. Time-dependence of the potential energy of the HPCS during (NPT)-MD simulation at 300 K. The force field used was the PCFF.

(NPT) ensemble. We note that while the constant value of the total energy of the polymer may be a good indication of whether it has attained its true equilibrium configuration, it might not be sufficient under certain conditions, and, therefore, we also computed its radial distribution function g(r). Figure 6 presents 3314



Article

Figure 6. The computed radial distribution function of the HPCS.

the results at 300 K. The computed g(r) and, in particular, the locations of the various peaks that correspond to the three types of the covalent bonds in the polymer are in very good agreement with the experimental data,67 hence confirming the accuracy of the equilibrium configuration of the HPCS. Figure 7 presents the

Figure 7. The amorphous structure of the HPCS, as generated by the PCFF.

structure of the HPCS at the end of MD simulation, which is composed of 1928 atoms within a cell of dimensions of 28.9 × 28.9 × 28.9 (Å)3. The cell contained 4 separate chains, each composed of 482 atoms. The computed density of the polymer is 0.973 gr/cm3, which should be compared with the experimental value of 0.980 gr/cm3. To further check whether the generated polymer structure is accurate, we also computed the bond distances as well as the triple and dihedral angles that exist in the polymer. The results for Si−C bonds, Si−C−Si angles, and H−C−Si−C dihedral are presented in Figure 8. The most probable Si−C distance was computed to be 1.9 Å, very close to the actual distance of 1.89 Å. The most probable angle in the Si−C−Si triplets is 109°, which is in the expected range for such triplets. As for the dihedral angles, the most probable angles for the cis and trans configurations were computed to be 58° and 180°, respectively, which are again very close to the theoretical values of 60° and 180°. To our knowledge, there are no experimental data for the angles in our polymer. Similar results were obtained for other bond distances, as well as triplet and dihedral angles and, thus, are not shown.

Figure 8. Distribution of (a) Si−C distances, (b) Si−C−Si angles, and (c) H−C−Si−C dihedral angles in the HPCS polymer, as generated by the PCFF.

The atomistic model of the HPCS was then utilized with ReaxFF to carry out MD simulation of its thermal decomposition. Before heating the polymer up, energy minimization and pre-equilibration simulations using Reaxff were utilized to further relax the atomistic structure of the polymer. All the simulations in this section were implemented using the MD simulator GRASP - General Reactive Atomistic Simulation Program - which can be utilized in parallel computations. The total energy E of the polymer was minimized with relaxation in all the directions. The maximum displacement allowed in any direction was 0.01 Å, while the time step was 0.1 fs. Figure 9 presents the change with the time of the total energy and pressure of the system during minimization. It is clear that the energy of the polymer has reached its minimum after 0.025 ps. To further check whether the true equilibrium was reached, the atoms were numbered, and their energies were monitored in 3315



Article

temperature was kept close to its target value, while the pressure fluctuated very little after a few ps.

5. COOK-OFF SIMULATION After determining the parameters of ReaxFF and generating a molecular model of the HPCS, cook-off simulation of the polymer was carried out to obtain an overview of its thermal degradation as a result of heating it up. Molecular dynamics simulations were carried out in the (NVT) ensemble using a heating rate of 0.1 K/fs, during which the temperature was increased from 50 to 5000 K. The variations of the pressure during the simulation are shown in Figure 11. The parameters of the simulations were the same as before. Figure 9. Time-dependence of the energy of the HPCS and the system’s pressure during minimization using the ReaxFF.

the last time frame of the minimization simulation. The results confirmed that all the energies are constant for each atomic type. To further equilibrate the polymer structure after energy minimization, MD simulation was carried out in the (NVT) ensemble for 5 ps at 50 K. The minimum absolute temperature deviation for which the rescaling of the velocities was performed was 10 K, and after each 10 time steps the temperature was compared with its target value. The rescaling coefficient was 0.7, so selected carefully to rescale the velocities to attain the target temperature. The dependence on time of the temperature and pressure of the system is shown in Figure 10, which indicate that

Figure 11. Time-dependence of the temperature and pressure during the cook-off simulation of the HPCS, using ReaxFF.

For temperatures less than but close to 1000 K, the Si−H and C−H bonds were very mobile and continuously formed and broke before the initiation of the polymer decomposition. The simulations indicated that the degradation of the HPCS is initiated at about 1000 K by releasing a few H atoms from the Si− H and C−H bonds; see Figure 12. The Si−H bond is about 20 kcal/mol weaker than the C−H bond. Therefore, its cleavage initiated the decomposition process and happened more frequently, consistent with the experimental results.67 At higher temperatures, extensive hemolytic Si−H and C−H bond breaking occurred, which led to the formation of a large number of H radicals, shown in Figure 13. A secondary reaction involving bond formation between the H radicals led to the formation of hydrogen gas, the main pyrolysis product. The loss of H from SiHn (n = 1,2) led to the formation of reactive silylene (Si:) intermediates, resulting in the formation of the Si−Si bonds due to the initial chain cross-linking. But, above 1000 K it quickly broke because there is tendency for Si atoms to form more stable bonds with C, hence explaining the large fluctuations for the number of the Si−Si bonds shown in Figure 12. Cleavage of the Si−C bonds initiated at about 2000 K and accelerated at 2400 K and higher, which led to the formation of CH2, SiH2, and CH4Si radicals, shown in Figure 13. [In Part 2 (DOI 10.1021/jp307799p) we carry out the pyrolysis simulation at same temperature of 2000 K and explain in more detail why this particular temperature was used to generate the amorphous SiC.] They reacted with the H radicals to form CH3, SiH3, and CH5Si radicals and eventually methane, silane, and methylsilane gases. In addition, various C2 hydrocarbons, including C2H6 and C2H4, and such radicals as C2H5 and C2H3 as well as C2H2 were

Figure 10. Time-dependence of the temperature of the HPCS and the system’s pressure during pre-equilibration at 50 K, using ReaxFF. 3316



Article

Figure 12. The evolution of the number of bonds during the cook-off simulation of the HPCS, using ReaxFF.

Figure 13. The evolution of the number of the hydrogen radicals during the cook-off simulation of the HPCS, using ReaxFF.

Figure 14. The evolution of the C2 hydrocarbons during the cook-off simulation of the HPCS, using ReaxFF.

first formed, as shown in Figure 14. The time evolution of the gaseous compounds is shown in Figure 15. Table 11 lists all the products that are obtained after heating up the HPCS to 3600 K. Analysis of the gaseous products by mass spectroscopy and gas chromatography has indicated that the

simulation results are consistent, both quantitatively and qualitatively, with the experimental data.67 Reference 67 reported that most of the produced gas was hydrogen, followed in much smaller quantities by methane, methyl silanes, and C2 hydrocarbons. These are also what we obtained as well as other gases 3317



Article

6. SUMMARY The reactive force field ReaxFF was extended in order to carry out molecular dynamics simulation of pyrolysis of hydripolycarbosilane (HPCS), which is a crucial step in the fabrication of nanoporous SiC membranes. The parameters of ReaxFF that are specifically associated with the HPCS were estimated using quantum chemical computations based on density-functional theory. An atomistic model of the HPCS was then generated and utilized with ReaxFF in MD simulation of its thermal decomposition. According to the MD simulations the degradation of the HPCS is initiated at about 1000 K. Cleavage of the Si−H bonds initiates the decomposition. At higher temperatures a large number of H radicals are formed, along with a secondary reaction that involves bond formation between the radicals that produces hydrogen gas, the main pyrolysis product. Reactive silylene intermediates are also produced, which lead to the formation of the Si−Si bonds that quickly break above 1000 K, as Si tends to form more stable bonds with C. Cleavage of the Si−C bonds is initiated at about 2000 K and accelerated at 2400 K and higher, which led to the formation of some intermediate gaseous reaction products that eventually form methane, silane, and methylsilane gases as well as some C2 hydrocarbons. The results are consistent, both quantitatively and qualitatively, with the experimental data.64 The next step is to carry out MD simulation of pyrolysis of the polymer under the conditions that closely resemble those in the fabrication of SiC nanoporous membranes, in order to generate the SiC material and to compare its properties with the experimental data. This will be studied in the sequel to this paper.

Figure 15. The evolution of various compounds during cook-off simulation of the HPCS, using ReaxFF.

Table 11. Pyrolysis Products after Heating up the Polymer to 3600 Ka N

product

MW

N

product

MW

122 26 3 11 2 3 7 5 10 1 4 14 16 10 5 1 1 1 3 3 1 1 2 4 1 3 3 4 1 2 1 1 2

H H2 CH3 CH4 C2H4 HSi H2Si H3Si H4Si CSi CH2Si CH3Si CH4Si CH5Si CH6Si CH7Si C2H2Si C2H3Si C2H5Si C2H6Si C2H7Si CH4Si2 CH5Si2 CH6Si2 CH7Si2 C2H5Si2 C2H6Si2 C2H7Si2 C2H8Si2 C2H9Si2 C2H10Si2 C3H8Si2 C3H9Si2

1.00 2.00 15.03 16.04 28.05 29.09 30.10 31.11 32.12 40.09 42.11 43.12 44.13 45.13 46.14 47.15 54.12 55.13 57.14 58.15 59.16 72.21 73.22 74.23 75.24 85.23 86.24 87.25 88.25 89.26 90.27 100.27 101.27

1 1 3 1 1 3 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 1 1 1 1

C2H5Si3 C2H7Si3 C2H10Si3 C2H11Si3 C5H9Si2 C3H8Si3 C3H9Si3 C3H10Si3 C4H12Si3 C4H13Si3 C2H12Si4 C3H10Si4 C3H11Si4 C4H12Si4 C4H13Si4 C5H14Si4 C3H15Si5 C4H13Si5 C4H14Si5 C4H15Si5 C5H15Si5 C5H19Si5 C6H16Si5 C7H16Si5 C5H15Si6 C6H19Si6 C6H19Si7 C8H25Si7 C9H20Si8 C14H37Si12 C14H39Si13 C17H42Si13

113.32 115.33 118.36 119.36 125.29 128.35 129.36 130.37 144.39 145.40 148.46 158.45 159.46 172.48 173.49 186.50 191.58 201.57 202.58 203.59 215.60 219.63 228.62 240.63 243.69 259.73 287.81 317.88 352.94 542.47 572.57 611.63

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS The support of the Office of Energy Sciences, Chemical Sciences, Geosciences & Biosciences Division of the Department of Energy is gratefully acknowledged. T.T.T. also acknowledges the support of the National Science Foundation. W.A.G. received support from the Department of Transportation.

■

REFERENCES

(1) Elyassi, B.; Sahimi, M.; Tsotsis, T. T. J. Membr. Sci. 2007, 288, 290. (2) Elyassi, B.; Sahimi, M.; Tsotsis, T. T. J. Membr. Sci. 2008, 316, 73. (3) Elyassi, B.; Kim, T. W.; Sahimi, M.; Tsotsis, T. T. Mater. Chem. Phys. 2009, 118, 259. (4) Elyassi, B.; Sahimi, M.; Tsotsis, T. T. In Encyclopedia of Chemical Processing; Lee, K. B., Ed.; Taylor & Francis: London, 2009; p 1. (5) Elyassi, B.; Sahimi, M.; Tsotsis, T. T. Submitted for publication. (6) Kotani, M.; Nishiyabu, K.; Matsuzaki, S.; Tanaka, S. J. Ceram. Soc. Jpn. 2011, 119, 563. (7) Takeda, Y.; Nakamura, K.; Maeda, K.; Matsushita, Y. J. Am. Ceram. Soc. 1987, 70, 266. (8) Schulz, K.; Durst, M. Filtr. Sep. 1994, 31, 25. (9) Rosenbloom, A. J.; Sipe, D. M.; Shishkin, Y.; Ke, Y.; Devaty, R. P.; Choyke, W. J. Biomed. Microdevices 2004, 6, 261. (10) Li, Z.; Kusakabe, K.; Morooka, S. J. Membr. Sci. 1996, 118, 159. (11) Zorman, C. A.; Fleischman, A. J.; Dewa, A. S.; Mehregany, M.; Jacob, C.; Nishino, S.; Pirouz, P. J. Appl. Phys. 1995, 78, 5136. (12) Kenawy, S. H.; Nour, W. M. N. J. Mater. Sci. 2005, 40, 3789.

a

N is the number of the produced product, and MW is its molecular weight.

that may not be possible to measure experimentally with any precision (due to their trace amounts). In fact, this is an advantage that ReaxFF offers. 3318



Article

(13) Chen, F.; Mourhatch, R.; Tsotsis, T. T.; Sahimi, M. Chem. Eng. Sci. 2008, 63, 1460. (14) Mourhatch, R.; Tsotsis, T. T.; Sahimi, M. J. Membr. Sci. 2010, 356, 138. (15) Hasegawa, Y.; Okamura, K. J. Mater. Sci. 1983, 18, 3633. (16) Soraru, G. D.; Babonneau, F.; Mackenzie, J. D. J. Non-Cryst. Solids. 1988, 106, 256. (17) Soraru, G. D.; Babonneau, F.; Mackenzie, J. D. J. Mater. Sci. 1990, 25, 3886. (18) Tang, X.; Zhang, L.; Tu, H.; Gu, H.; Chen, L. J. Mater. Sci. 2010, 45, 5749. (19) Interrante, L. V.; Moraes, K. Ceram. Trans. 2002, 144, 125. (20) Interrante, L. V.; Whitmarsh, C. W.; Sherwood, W.; Wu, H.-J.; Lewis, R.; Maciel, G. In NATO ASI Series E; Harrod, J. F., Laine, R. M., Eds.; Kluwer Academic: Amsterdam, 1995; p 173. (21) Liu, Q.; Wu, H.-J.; Lewis, R.; Maciel, G. E.; Interrante, L. V. Chem. Mater. 1999, 11, 2038. (22) Farah, K.; Müller-Plathe, F.; Böhm, M. C. ChemPhysChem 2012, 13, 1127. (23) van Duin, A. C. T.; Dasgupta, S.; Lorant, F.; Goddard, W. A. J. Phys. Chem. A 2001, 105, 9396. (24) Brenner, D. W. Phys. Rev. B 1990, 42, 9458. (25) Nyden, M. R.; Noid, D. W. J. Phys. Chem. 1991, 95, 940. (26) Chenoweth, K.; van Duin, A. C. T.; Goddard, W. A. J. Phys. Chem. A 2008, 112, 1040. (27) Chenoweth, K.; van Duin, A. C. T.; Persson, P.; Cheng, M. J.; Oxgaard, J.; Goddard, W. A. J. Phys. Chem. C 2008, 112, 14645. (28) Cheung, S.; Deng, W. Q.; van Duin, A. C. T.; Goddard, W. A. J. Phys. Chem. A 2005, 109, 851. (29) van Duin, A. C. T.; Strachan, A.; Stewman, S.; Zhang, Q.; Xu, X.; Goddard, W. A. J. Phys. Chem. A 2003, 107, 3803. (30) Zhang, Q.; Cagin, T.; van Duin, A.; Goddard, W. A.; Qi, Y.; Hector, L. G. Phys. Rev. B 2004, 69, 045423. (31) Strachan, A.; van Duin, A. C. T.; Chakraborty, D.; Dasgupta, S.; Goddard, W. A. Phys. Rev. Lett. 2003, 91, 098301. (32) van Duin, A. C. T.; Zeiri, Y.; Dubnikova, F.; Kosloff, R.; Goddard, W. A. J. Am. Chem. Soc. 2005, 127, 11053. (33) Nielson, K. D.; van Duin, A. C. T.; Oxgaard, J.; Deng, W. Q.; Goddard, W. A. J. Phys. Chem. A 2005, 109, 493. (34) Liu, L.; Jaramillo-Botero, A.; Goddard, W. A.; Sun, H. J. Phys. Chem. A 2012, 116, 3918. (35) Chenoweth, K.; Cheung, S.; van Duin, A. C. T.; Goddard, W. A.; Kober, E. M. J. Am. Chem. Soc. 2005, 127, 7192. (36) Desai, T. G.; Lawson, J. W.; Keblinski, P. Polymer 2011, 52, 577. (37) Saha, B.; Schatz, G. C. J. Phys. Chem. B 2012, 116, 4684. (38) Jaguar 5.5; Schrödinger, LLC: Portland, Oregon, 1991−2003. (39) Becke, A. D. J. Chem. Phys. 1993, 98, 1372. (40) Lee, C.; Yang, W.; Parr, R. G.; Pople, J. A. Phys. Rev. B 1988, 37, 785. (41) Krishnan, R.; Binkley, J. S.; Seeger, R.; Pople, J. A. Chem. Phys. 1980, 72, 650. (42) Clark, T.; Chandrasekhar, J.; Spitznagel, G. W.; Schleyer, P. V. R. J. Comput. Chem. 1983, 4, 294. (43) Haharan, P. C.; Pople, J. A. Theor. Chim. Acta 1973, 28, 213. (44) Francle, M. M.; Pietro, W. J.; Hehre, W. J.; Binkley, J. S.; Gordon, M. S.; DeFrees, D. J.; Pople, J. A. J. Chem. Phys. 1982, 77, 3654. (45) P. Giannozzi, P.; et al. J. Phys.: Condens. Matter. 2009, 21, 395502. (46) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. Rev. Lett. 1996, 77, 3865. (47) Perdew, J. P.; Chevary, A.; Vosko, S. A.; Jackson, K. A.; Pederson, M. R.; Singh, D. J.; Fiolhais, C. Phys. Rev. B 1992, 46, 6671. (48) Monkhorst, H. J.; Pack, J. D. Phys. Rev. B 1976, 13, 5188. (49) Stewart, J. J. P. J. Comput. Chem. 1989, 10, 209. (50) Mortier, J.; Ghosh, S. K.; Shankar, S. J. Am. Chem. Soc. 1986, 108, 4315. (51) Taylor, A.; Jones, R. M. In Silicon Carbide: A High Temperature Semiconductor; O’Connor, J. R., Smiltens, J., Eds.; Pergamon Press: Oxford, 1960; p 147.

(52) Goldberg, Yu.; Levinshtein, M. E.; Rumyantsev, S. L. In Properties of Advanced Semiconductor Materials GaN, AlN, SiC, BN, SiC, SiGe; Levinshtein, M. E., Rumyantsev, S. L., Shur, M. S., Eds.; Wiley: New York, 2001; p 93. (53) Cerius2, version 4.0; Accelyrs: San Diego, CA, 1999. (54) Patnaik, P. Handbook of Inorganic Chemicals; McGraw-Hill: New York, 2002. (55) Plimpton, S. J. Comput. Phys. 1995, 117, 1. (56) Theodorou, D. N.; Suter, U. W. Macromolecules 1985, 18, 1467. (57) Tocci, E.; Hofmann, D.; Paul, D.; Russo, N.; Drioli, E. Polymer 2001, 42, 521. (58) Ostwal, M. M.; Tsotsis, T. T.; Sahimi, M. Phys. Rev. E 2009, 79, 061801. (59) Mehrabi, A. R.; Sahimi, M. Phys. Rev. Lett. 1999, 82, 735. (60) Kim, N.; Harale, A.; Tsotsis, T. T.; Sahimi, M. J. Chem. Phys. 2008, 127, 224701. (61) Mehrabi, A. R.; Sahimi, M. J. Chem. Phys. 2008, 128, 234503. (62) Wolf, D.; Keblinski, P.; Phillpot, S. R.; Eggebrecht, J. J. Chem. Phys. 1999, 110, 8254. (63) Sahimi, M. Heterogeneous Materials II, Nonlinear and Breakdown Properties and Atomistic Modeling; Springer: New York, 2003; Chapter 8. (64) Andersen, H. C. J. Chem. Phys. 1980, 72, 2384. (65) Nosé, S. J. Chem. Phys. 1984, 81, 511. (66) Hoover, W. Phys. Rev. A 1985, 31, 1695. (67) Interrante, L. V.; Whitmarsh, C. W.; Sherwood, W. Mater. Res. Soc. Symp. Proc. 1994, 346, 593.

3319


Toward a Process-Based Molecular Model of SiC ... - ACS Publications

Recommend Documents