H ReaxFF Reactive Potential for Silicon Surfaces Grafted with

Sep 26, 2018 - In this work, we developed Si/C/H ReaxFF force field for the study of the functionalization and decomposition of alkyl monolayers on si...
0 downloads 0 Views 3MB Size
Subscriber access provided by Kaohsiung Medical University

C: Surfaces, Interfaces, Porous Materials, and Catalysis

A Si/C/H ReaxFF Reactive Potential for Silicon Surfaces Grafted with Organic Molecules Federico A. Soria, Weiwei Zhang, Patricia A. Paredes-Olivera, Adri C. T. van Duin, and Eduardo Martin Patrito J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.8b07075 • Publication Date (Web): 26 Sep 2018 Downloaded from http://pubs.acs.org on October 8, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

A Si/C/H ReaxFF Reactive Potential for Silicon Surfaces Grafted with Organic Molecules

Federico A. Soria,a Weiwei Zhang,b Patricia A. Paredes-Olivera,a Adri. C. T. van Duinb and Eduardo M. Patritoc*

a

Departamento de Química Teórica y Computacional and cDepartamento de Fisicoquímica.

Instituto de Investigaciones en Físico Química de Córdoba (INFIQC). Facultad de Ciencias Químicas, Universidad Nacional de Córdoba, X5000HUA Córdoba, Argentina. b

Department of Mechanical and Nuclear Engineering, Pennsylvania State University,

University Park, Pennsylvania 16802, United States.

ACS Paragon Plus Environment

1

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract In this work we developed Si/C/H ReaxFF force field for the study of the functionalization and decomposition of alkyl monolayers on silicon surface. The parameterization was performed based on the main reactions involved in the decomposition of alkyl layers on small silicon clusters. The decomposition mechanisms observed in the molecular dynamics (MD) simulations were validated by the comparison of ReaxFF energy barriers for the elementary steps of the main mechanisms with density functional theory (DFT) calculations. Activation energy barriers obtained from the MD simulations from Arrhenius plots are in excellent agreement with the values calculated from DFT. The trends in the preexponential factor with the alkyl chain length follow the predictions of transition state theory. The results confirm that the main decomposition mechanism of the alkyl chains is the alkene elimination to the gas phase after a β-hydride abstraction by silyl radicals which are formed in a previous step. The ReaxFF force field was used to comparatively investigate the alkyl surface coverage of Si(111), Si(100)−2×1 and “half-flat” Si(100) surfaces as a function of the alkyl chain length, showing good agreement with reported experimental values. Both the DFT and ReaxFF MD calculations predict that decyl monolayers with coverages as high as to 0.8 are thermodynamically stable at moderate temperatures.

ACS Paragon Plus Environment

2

Page 2 of 38

Page 3 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

1. Introduction Silicon is the most important material in modern technology and the functionalization of bulk and nanostructured silicon surfaces is critical for many device applications. The attachment of organic molecules to Si surfaces via Si–C bonds has distinctive advantages due to its high chemical and thermal stability.1 Understanding the structure and reactivity of the grafted silicon surface is required for reliable device manufacturing. High coverage monolayers grafted to silicon form a passivating layer which avoids its oxidation,

2,3,4

an

important step in the development of photovoltaic devices.5 On the nanoscale, when the surface to bulk ratio increases dramatically, the surface functionalization also plays a key role in the chemical stability, electronic structure and optical spectra of silicon nanocrystals.6 The synthesis of water-soluble nanocrystals has opened the route to realistic biomedical applications of functionalized silicon nanoparticles.7 Silicene, the twodimensional allotrope of silicon, is another example for which the functionalization is critical to increase its chemical stability.8 The importance of the surface chemistry of grafted silicon has motivated a series of investigations regarding the mechanisms of decomposition of alkylated silicon surfaces.1 Sung et al.9 studied the thermal stability of octadecyl monolayers in UHV and found that the characteristic vibrational modes of the grafted layer stay unchanged up to about 327 °C. At higher temperature the signals of characteristic modes decrease while the Si−H vibrational modes appear simultaneously. From this observation they concluded that the grafted layer desorbs according to a β-hydride elimination reaction: ≡Si–CH2–CH2–CH2–R → ≡Si–H + H2C=CH–CH2–R(gas) For long alkyl monolayers Yamada et al. observed a loss in the

(1)

νCH

intensity after

annealing above 217°C together with a small shift of the νCH bands which they attributed to an irreversible loss in the organization of the layer.10 Faucheux et al. investigated the thermal stability under oxidizing or reducing atmospheres using IR spectroscopy and found that the layers are thermally stable up to 250 °C.11 In the range of 250–300 °C, the main reaction is alkene desorption accompanied with silicon oxidation. The characteristic

ACS Paragon Plus Environment

3

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 38

desorption temperature is not significantly affected by the chain length from C6 to C18.11 A detailed analysis of the annealing of ≡Si–C2H5 and ≡Si–CH3 surfaces using XPS was carried out by Jaeckel et al. The Si(111)−C2H5 surfaces decompose at temperatures >300°C, whereas Si(111)−CH3 surfaces are stable up to at least 440°C.12 A secondary decomposition mechanism was observed at temperatures above 440°C producing surface methyl groups, for which they proposed the following reaction: ≡Si–CH2–CH2–CH2–R → ≡Si–CH3 + H2C=CH–R(gas)

(2)

The stability of methylated layers was also studied by Hunger et al. who found that the layers are stable up to 447 °C.13 Yang et. al. reported the formation of silicon carbide at temperature higher than 390°C for methylated surfaces.14 The surface coverage of alkylated silicon surfaces depends on the alkyl chain length. Only in the case of methyl groups a complete coverage can be achieved. Nemanick et al. reported that the functionalization of the Si(111) surface with ethyl groups through a two-step chlorination/alkylation procedure produces surfaces with a coverage close to 0.8.15 These results were supported by DFT calculations, which show that the preparation of the ethylterminated Si(111) with an 100% coverage is not limited by thermodynamics but by high reaction barriers.16 For longer alkyl chains the coverage is around 0.5-0.6. Puniredd et al. reported a coverage of around 0.5 with propyl groups.17,18 The surface modification with longer organic chains (C6-C18) also produces coverages around 0.5.11, 15, 17 For hexyl and octyl monolayers DFT calculations show that the most favorable thermodynamic coverage is about 0.516 and the same is found from molecular modeling simulations for octadecyl monolayers which found an optimum coverage of 0.55.19,20 Recently we performed ReaxFF reactive molecular dynamics (MD) simulations to investigate the mechanisms and kinetics of thermal decomposition processes of silicon surfaces grafted with different organic molecules via Si−C bonds at atomistic level.21 We found that the presence of silyl radicals plays a major role in the decomposition of alkyl chains grafted to silicon according to the following reaction: ≡Sia–CH2–CH2–CH2–R + ≡Sib• → ≡Sia• + ≡Sib–H + H2C=CH–CH2–R(gas)

(3)

This reaction involves two silicon atoms (labelled Sia and Sib) and contrasts with Eq. (2) in which desorption of the alkyl chain as an alkene involves a single Si atom.

ACS Paragon Plus Environment

4

Page 5 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Given the importance of the surface chemistry of functionalized silicon, it is essential to clearly understand the mechanisms involved in the thermal desorption of alkyl-grafted surfaces. In this work, we develop a ReaxFF Si/C/H potential which correctly describes the ∆E values of the main reactions involved in the decomposition of alkylated silicon surfaces. The mechanisms of thermal decomposition deduced from the MD simulations were validated by DFT calculations of activation energy barriers of the elementary steps of the main reactions. Next we investigate the functionalization of the bare Si(111), Si(100)−2×1 and “half-flat” Si(100) surfaces with alkyl chains of different lengths in order to evaluate the maximum surface coverages that can be attained on the different crystalline surfaces. We show that even for large alkyl chain lengths, stable monolayers can be produced with coverages higher than 0.5.

2. Theoretical Methods and Surface Modeling 2.1. DFT Calculations Periodic Density Functional Theory (DFT) calculations were performed as implemented in the PWSCF code of the Quantum ESPRESSO (QE) suite.22 Norm-conserving ultrasoft pseudopotentials were used for the atomic species.23 The PBE formulation was employed for the exchange and correlation functional.24 The electron wave functions were expanded in a plane-wave basis set up to a kinetic energy cutoff of 30 Ry (180 Ry for the density). All calculations involving silyl radicals were performed with spin polarization. Reaction pathways and energy barriers were calculated with the climbing-image nudged-elastic-band (CI-NEB) method as implemented in QE.25 A silicon slab with six layers was used to model the (111) face of silicon. The slabs were separated by a vacuum thickness of 10 Å. The dangling bonds of the bottom surface were saturated with hydrogen atoms. We used a0= 5.48 Å for the lattice constant, as determined in a previous work.26 The positions of all of the adsorbate atoms as well as those of the four topmost Si layers were fully optimized. The silicon atoms of the lower bilayer were kept fixed in a bulk configuration. The silicon atoms of the lower bilayer were kept fixed in a bulk configuration. CI-NEB calculations were performed using 2 × 2 and 3 × 3 unit cells and the integration in the first Brillouin zone was performed with a (4×4×1)

ACS Paragon Plus Environment

5

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Monkhorst−Pack mesh.27 Energy profiles along the reaction pathway calculated with the CI-NEB method25 are plotted as a function of a dimensionless reaction coordinate ranging from 0.0 for reactants to 1.0 for products. Selected atomic structures along the reaction coordinate are shown in the panels below each energy curve. Total energy calculations involving different coverages of decyl molecules were performed with 3√3×5 unit cell with 30 top Si atoms and six silicon layers. For this large unit cell, the calculations were performed at the gamma point. Grimme-D2 London dispersion forces were used to describe the van der Waals interactions among the alkyl chains.25 2.2. ReaxFF MD Simulation Setup The (111) surface of silicon was constructed by replicating the unit cell shown in Fig. 1a with dimensions of 30.72 Å × 26.60 Å. It contains 64 top silicon atoms which in Fig. 1a are hydrogenated for clarity. The surface coverage of alkyl chains (R=alkyl chain) is calculated as the ratio of Si−R groups divided by the number of top Si atoms. Fig. 1b shows the unit cell employed to build the Si(100)–2×1 surface, with dimensions of 30.72 Å × 30.72 Å. This reconstruction has the typical rows of Si dimers (one dimer colored in brown) separated by atomically deep channels. There are 32 top Si dimers in Fig. 1b. We also investigated the functionalization of the so called half-flat Si(100) surface, Si(100)-hf. It is experimentally produced when the silicon oxide is removed by successively immersing the Si(100) wafer in NH4F solutions.28,29 This treatment produces a nearly atomically flat surface characterized by missing row structures along [011] directions. Fig. 1c shows that the top Si atoms of this surface are dihydrogenated whereas the Si atoms in the second layer form dimers which are monohydrogenated. The Si(100)-hf surface was constructed by replicating the unit cell in Fig. 1c, with dimensions 30.72 Å × 30.72 Å. It has 32 Si atoms in the first layer (dihydrogenated) and 32 Si atoms in the second layer (forming 16 Si dimers, one dimer is colored in green for clarity). This structure is somewhat different than the experimentally determined structure28,29 in which the elevated rows are terminated by unstrained dihydrides (as in Fig. 1c), whereas the trenches are terminated by strained dihydrides, instead of the monohydride-terminated dimers of Fig. 1c. As the strained dihydrides in the trenches are likely kinetically trapped structures which would convert to dimers at higher temperatures to relieve the strain, we decided for simplicity to use the

ACS Paragon Plus Environment

6

Page 6 of 38

Page 7 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

hydrogenated dimers in our high temperature calculations. This structural difference is not expected to have important effects, as the MD simulations show that only a small number of Si atoms in the second layer are functionalized due to steric effects. The density of surface sites of the Si(111), Si(100)–2×1 and Si(100)–hf surfaces is 7.8, 6.8 and 10.2 sites/nm2, respectively. The three surfaces were grafted with ethyl (C2), propyl (C3), pentyl (C5), and decyl (C10) groups. SiH and Si–alkyl groups were randomly accommodated to yield the desired surface coverage of alkyl groups. The MD simulations of the Si(111) grafted with C10 ( 0.5 to 0.8 surface coverage) were performed at 1600 and 1800 K. MD simulations with the reactive ReaxFF30 force field were performed with the ADF201731 and LAMMPS32 codes. NVT/MD simulations were carried out for 800 ps from 1500 to 2200 K. A Velocity-Verlet algorithm was used with a 0.1 fs time step and a temperature-damping constant of 100 fs. Prior to each simulation, the grafted silicon slab was thermalized at 300 K for 1 ps.

Figure 1. Unit cells used in the ReaxFF simulations: a) Si(111), b) Si(100)–2×1 and c) Si(100)–hf. The surface Si atoms are hydrogenated for clarity. A surface Si dimer of the Si(100)–2×1 and Si(100)–hf surfaces are marked in brown for clarity in b) and c), respectively.

ACS Paragon Plus Environment

7

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 38

3. Results and Discussion 3.1. ReaxFF Si/C/H Force Field Development A Si/C/H/O force field has been reported in a study of the oxidation of silicon carbide by O2 and H2O.33 However, this potential does not correctly describe the dynamics of grafted silicon surfaces as it overestimates the Si−H interactions. This leads to a diffusion of H atoms into the Si bulk thus leaving Si radicals on the surface, which rapidly dehydrogenate the alkyl chains at moderate temperatures as shown in Fig. S1. The force field parameterization in this work was started by combining ReaxFF force fields previously developed for Si/H34 and C/H35 as well as the Si/C/H/O force field.33 The new Si/C/H force field was refined by using a training set containing the main reactions involved in the decomposition of the alkylated Si(111) surface. The training set was constructed using a Si23H30 cluster functionalized with methyl, ethyl and propyl species (Si23H29−CH3, Si23H29−C2H5, Si23H29−C3H7, see Fig. S2 in the Supporting Information). The DFT calculations for the clusters were carried out using the Gaussian 09 suite of programs36 at the PBE / 6-31G (d, p) level of theory. The reactions considered involve the breaking of Si−H, Si−C, C−C and C−H bonds. The training process optimized one parameter at a time in a iterative manner to minimize the sum of the errors.37 The error function for the training set, which minimizes the error values between ReaxFF and QM value, is described by the following equation: Error = ∑ 

, ,

 



(4)

where the Xi,QM represent the energies from quantum mechanics calculations, while Xi,ReaxFF are the energies from ReaxFF, and σi are weighting values assigned to each data point. The ReaxFF force field parameters are listed in the Supporting Information. Table 1 shows the reactions used in the training set together with the ∆E values calculated by DFT and the newly developed Si/C/H force field. Reactions (1) to (3) involve Si–H bond breakage, reactions (4) to (6) involve the breakage of Si–C bonds, whereas reactions (7) to (10) and (11) to (13) involve the breakage of C–H and C–C bonds of the alkyl chains, respectively. Reactions (14) to (18) involve the formation of Si–H bonds on the surface. We chose different weight factors for the data points. For the most important reactions we used the smallest σ

ACS Paragon Plus Environment

8

Page 9 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

values. For example, for reactions (16) to (18) in Table 1 we used σ = 0.05 whereas for reaction (1) we used σ = 2.0.

It can be observed that the force field correctly describes ∆E values for all reactions. This can be better appreciated in Fig. 2a where we plotted the ReaxFF energy as a function of the DFT energy. The points are very close to the line with unity slope. Table 1. Comparison of energies differences for selected reactions from DFT and ReaxFF calculations. ∆E, DFT

∆E, ReaxFF

(kcal/mol)

(kcal/mol)

(1) ≡Si–H → ≡Si• + H•

81.1

78.02

(2) ≡SiH2 → ≡SiH• + H•

83.8

78.4

(3) 2≡Si–H → 2≡Si• + H2

58.06

51.34

(4) ≡Si–CH3 → ≡Si• + CH3•

76.40

76.50

(5) ≡Si–CH2CH3 → ≡Si• + CH2CH3•

70.04

67.41

(6) ≡Si–CH2CH2CH3 → ≡Si• + CH2CH2CH3•

73.64

69.51

(7) ≡Si–CH3 → ≡SiCH2• + H•

102.09

101.78

(8) ≡Si–CH2CH3 → ≡SiCH2CH2• + H•

99.96

94.65

(9) ≡Si–CH2CH2CH3 → ≡Si–CH2CH2CH2• + H•

107.54

96.28

(10) ≡Si–CH2CH2CH3 → ≡Si–CH2CH•CH3 + H•

99.04

92.40

(11) ≡Si–CH2CH3 → ≡SiCH2• + CH3•

89.68

89.81

Reaction

ACS Paragon Plus Environment

9

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 38

(12) ≡Si–CH2CH2CH3 → ≡SiCH2CH2• + CH3•

92.99

85.98

(13) ≡Si–CH2CH2CH3 → ≡SiCH2• + CH2CH3•

88.81

83.48

(14) ≡Si–CH2CH3 →≡SiH + CH2CH2(gas)

30.25

26.82

(15) ≡Si–CH3 + ≡Si• → ≡SiCH2• + ≡SiH

23.87

23.78

(16) ≡Si–CH2CH3 + ≡Si•→ ≡Si–CH2CH2• + ≡SiH

22.25

19.5

(17) ≡Si–CH2CH2CH3 + ≡Si·→ ≡Si–CH2CH•CH3 + ≡SiH

19.23

14.8

(18) ≡Si–CH2CH2CH3 + ≡Si•→ ≡Si– CH2CH2CH2• + ≡SiH

25.6

21.9

Fig. 2b compares the ReaxFF and DFT dissociation curves for Si−H, Si−C, C−C and C−H bonds. The curves were obtained by performing constrained geometry optimizations were the bond length was kept fixed at a specific value and the remaining atoms where optimized. In all cases, DFT geometry optimizations were performed in the singlet and triplet states. The ReaxFF force field correctly predicts the equilibrium distances of the Si−H (1.50 Å), Si−C (1.93 Å), C−C (1.53 Å) and C−H (1.10 Å) bonds. In addition, at long bonding distances, the ReaxFF curve approaches the DFT triplet curve.

ACS Paragon Plus Environment

10

Page 11 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 2. a) ReaxFF energy as a function of the DFT energy. The points are very close to the line with unit slope. b) Comparison of bond dissociation curves calculated with ReaxFF and DFT methods for C−H, Si−H, C−C and Si−C bonds. Note that at long bonding distances, the ReaxFF curves approach the DFT-triplet curves which are the ground state.

ACS Paragon Plus Environment

11

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

3.2. DFT Validation of Elementary Reaction Steps Observed in ReaxFF Simulations 3.2.1. Alkyl Decomposition via Alkene Desorption: One or Two Si Atoms Involved? As mentioned in the introduction, it has been proposed in the literature that desorption of the alkyl chain as an alkene involves only a single Si atom according to Eq. (1). However, this mechanism was never observed in our previous MD simulations21 in which we showed that the decomposition of the surface alkyl group to yield the alkene in the gas phase requires the presence of an adjacent silyl group according to Eq. (3). We therefore calculated the DFT energy profiles along the reaction path for both mechanisms to obtain the corresponding activation energy barriers. Fig. 3 compares the DFT energy profiles in the case of an ethylated surface. It can be observed that the reaction involving a single Si atom has a much higher energy barrier (76.3 kcal/mol, Fig. 3a, black line) than when the decomposition is triggered by an adjacent Si atom (32.3 kcal/mol, Fig. 3a, red line). The panels in Fig. 3b show that in the transition state (panel II) the Si–C bond is broken and an H atom from the terminal methyl group has migrated to the Si atom. Finally, a new Si–H bond is formed, and the ethylene molecule is released to the gas phase (panel III). In the presence of an adjacent silyl radical, the first step involves the H abstraction by the silyl radical (Fig. 3c, panels I’-III’). Then the Si–C bond is broken and the CH2CH2 molecule is released to the gas phase, generating a new silyl radical on the surface with Eact = 11.3 kcal/mol (Fig. 3c, panels III’-V’). This decomposition mechanism is therefore the reverse of the experimentally observed reaction of olefin additions on the hydrogenated Si(111) surface, H−Si(111).38 Cicero et al. clearly demonstrated by STM that the reaction of styrene with H−Si(111) involved the sequential production of silyl radicals in a chain reaction.38 From the relatively low forward and reverse energy barriers in this reaction mechanism it can be inferred that this adsorption-desorption process is reversible. As a can be seen in Fig. 3, the desorption reaction is an endothermic process and therefore the increment in the temperature induces entropic effects that favors the desorbed state.

ACS Paragon Plus Environment

12

Page 12 of 38

Page 13 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 3. a) DFT energy profiles along the reaction coordinate for the decomposition of Si(111)−CH2CH3 involving one (black curve) and two Si atoms (red curve). Panels in b) and c) show the structure of critical points along the reaction path when one and two silicon atoms are involved, respectively. The final state corresponds to desorbed ethylene.

ACS Paragon Plus Environment

13

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 38

Table 2. DFT activation energy barriers for alkyl chain desorption as alkene species through mechanisms that involve one and two Si atoms. Eact / kcal/mol System Two Si atoms

Single Si atom

≡Si–C2H5

32.3

76.3

≡Si–C3H7

31.5

88.5

≡Si–C5H11

39.9

64.7

≡Si–C10H21

40.3

76.6

Table 2 compares the DFT energy barriers for both processes as a function of the alkyl chain length. Much lower barrier energies are observed when two Si atoms are involved. As we shall see in the next section, the breakage of the Si−H bond has a lower barrier than the alkene elimination involving a single Si atom (third column in Table 2). Therefore, before this mechanism can occur under thermal treatment, silyl radicals will be present on the surface to catalyze the alkene elimination with much lower energy barriers (second column in Table 2). The energy barriers for the alkene elimination involving two Si atoms are similar to those found by Takeuchi et al. for the reverse reaction of CH2CH2 with a Si(111)–H surface with a silyl radical on the surface.39 Table 2 shows that the energies barriers for the alkene elimination mediated by a silyl radical are not much affected by the alkyl chain length. This result implies that the thermal stability of organic molecules grafted to Si(111) should not be appreciably modified by the length of the alkyl chain, which is in agreement with the experimental observation that shows that the desorption temperature is not significantly affected by the chain length from C6 to C18.11 In the absence of silyl radicals produced by a thermal treatment, only the single Si atom mechanism is available.

ACS Paragon Plus Environment

14

Page 15 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

The high energy barrier of this mechanism (Table 2) precludes any reactions at room temperature, unless the surface is irradiated with UV light which generates reactive excited states with lower energy barriers for the breakage of Si−C bonds.40 3.2.2. Mechanisms of Silyl Radical Formation: Si−C and Si−H Bond Breakage The MD simulations showed that surface silyl radicals are primarily formed by the breakage of Si–C and or Si−H bonds.21 After the Si–C bond cleavage, the alkyl radical a) may desorb into the gas phase, b) may abstract an H atom from an adjacent SiH group and desorb as an alkane and c) may bind to an adjacent SiH group forming a pentacoodinated Si intermediate.21 The corresponding reactions are: ≡Si–CnH2n+1→ ≡Si• + •CnH2n+1

(5)

≡Si–CnH2n+1 + ≡Si–H → 2≡Si• + CnH2n+2

(6)

≡Si–CnH2n+1 + ≡Si–H → ≡Si• + ≡SiH–CnH2n+1

(7)

Fig. 4a compares the energy profiles for reactions (5)-(7) for Si−CH2CH3 and Si−CH2CH2CH3 surface groups. Figs. 4b-4d show the corresponding structures of reactants, transition states and products for ethane desorption, ethyl desorption and formation of a pentacoordinated Si intermediate, respectively. For these reactions we obtained the following barriers: 61.9 kcal/mol, 68.3 kcal/mol and 65.3 kcal/mol, respectively. These values are lower than the barriers for the alkene desorption involving a single Si atom (Table 2), indicating that during the heating, the formation of silyl radicals will occur in a previous stage. In conclusion, reactions (5)-(6) have similar energy barriers and this validates the ReaxFF MD simulations in which the three mechanisms were observed.21 The MD simulations also showed21 that Si radicals can be produced by the release of H2 according to: ≡Si−H + ≡Si−H → ≡Si· + ≡SiH2

(8)

≡Si−H + ≡Si−H → 2≡Si· + H2(gas)

(9)

However, only some simulations showed the release of only one H2 molecule in the 1000 ps simulation time investigated. The fact that these reactions are less likely to occur is in agreement with the large energy barriers of 66.6 kcal/mol and 55.1 kcal/mol, respectively,

ACS Paragon Plus Environment

15

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

determined from the CI-NEB calculations. These values also indicate that during the heating, the formation of silyl radicals will occur before the alkene desorption involving a single Si atom.

Figure 4. a) DFT Energy profiles along the reaction coordinate for the formation of silyl groups on the Si(111) surface through ethane formation (black line), ethyl radical desorption (red line) and propyl diffusion on the surface. The panels show the structure of reactants, transition states and products for (b) ethane formation, (c) ethyl radical desorption and (d) propyl diffusion to an adjacent site.

ACS Paragon Plus Environment

16

Page 16 of 38

Page 17 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

In fact, the commonly accepted thermal hydrosilylation mechanism requires the heating of the hydrogenated silicon surface at temperatures ranging from 150 °C to 200 °C to homolytically cleave the Si–H bonds to form silyl radicals which subsequently react with alkenes or alkynes.19,41,42,43 As mentioned in the introduction, the thermal stability experiments reported in the literature are performed in a range of temperatures ranging from 200 to 500 °C.9-12 Therefore, at these temperatures the surface will be populated by silyl groups produced by the rupture of Si–H and/or Si–C bonds. These observations further support the fact that the decomposition of the monolayers occurs according to reaction (3) rather than reaction (1) as predicted by the ReaxFF simulations. The high energy barriers in Fig. 4a indicate that formation of silyl radicals is the rate limiting step and this process is also very endothermic. Therefore, high temperatures are required not only to surpass the activation energy barrier but also to shift the reaction to products side. In a previous work we showed that in the temperature range investigated around 10 % of surface Si atoms are silyl radicals and that their concentration remains nearly constant during the simulation.21

3.2.3. Formation of Si–CH3 groups A competing mechanism for the Si−C bond breakage observed in the MD simulations is the C−C bond breakage of the alkyl chains yielding surface methyl groups.21 Experimentally this mechanism was found at temperatures above 440°C for Si(111)–C2H5 surface.12 Therefore, we calculated the corresponding DFT energy profiles (Fig. 5a) during the breaking of the SiCH2−CH3 bond when one and two Si atoms are involved: ≡Si–CH2–CH3 → ≡Si–CH3 + •CH2

(10)

≡Si–CH2–CH3 + ≡Si• → ≡Si–CH3 + ≡Si–CH2

(11)

Figs. 5b and 5c show selected structures along the reaction path for both reactions. For the transition state of reaction (10), Fig. 5b shows that the C–C bond is broken and the H atom from terminal methyl group migrates to SiCH2 surface group. This produces a SiCH3 surface group and the CH2 radical in gas phase. The activation energy barrier for this process is 101.9 kcal/mol. An equivalent calculation in the case of a propyl surface group:

ACS Paragon Plus Environment

17

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

≡Si–CH2–CH2–CH3 → ≡Si–CH3 + CH2CH2 (gas)

Page 18 of 38

(12)

yields a quite close energy barrier of 98.9 kcal/mol. However, when a Si radical is involved according to reaction (11), a lower energy barrier is obtained (83.1 kcal/mol) and the ∆E of the reaction is drastically reduced. The cleavage of the C–C bond occurs concurrently with the formation of the Si−C bond with the adjacent atom (panels II’-III’ in Fig. 5c). For the propylated surface, in the first step the reaction yields the ≡Si–CH3 surface group with a lower barrier (77.5 kcal/mol): ≡Si–CH2–CH2–CH3 + ≡Si• → ≡Si–CH2=CH2 + ≡Si–CH3

(13)

In the next step the ≡Si–CH2=CH2 intermediate yields an ethylene molecule and a silyl radical with a small barrier (11.3 kcal/mol). The barriers for formation of SiCH3 groups induced by adjacent silyl radicals (83.1 kcal/mol and 77.5 kcal/mol for ethylated and propylated surfaces, respectively) must be compared to the barriers when the adjacent silyl radical induces the β-hydride abstraction and alkene elimination of Reaction (3): 32.3 kcal/mol and 31.5 kcal/mol for the ethylated and propylated surfaces, respectively (see Table 2). In conclusion, the formation of surface methyl groups has much higher barriers than for the alkene elimination process mediated by silyl radicals. This result explains why the formation of SiCH3 surface groups is a secondary decomposition mechanism which is observed experimentally only at high temperature12 and validates our ReaxFF MD simulations in which the formation of surface methyl groups according Eq. (13) was much less frequent than alkene desorption according to Eq. (3).

ACS Paragon Plus Environment

18

Page 19 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 5. a) DFT Energy profiles along the reaction coordinate for the formation of Si-CH3 surface groups when one (black line) and two silicon atoms (red line) are involved. Panels in b) and c) show the structure of critical points for the Si−CH3 formation through one atom and two silicon atoms, respectively.

ACS Paragon Plus Environment

19

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 38

3.2.4. Formation of Si−CH2CH2CH2−Si Moieties In the ReaxFF MD simulations on the propylated surface we found that the C3 chains remained on the surface during longer simulation times due to the formation of bicoordinated Si–CH2CH2CH2–Si moieties.21 They are produced by the dehydrogenation of the terminal methyl group according to the overall reaction: ≡Si–CH2CH2CH3 + ≡Si· → ≡Si–CH2CH2CH2–Si–H

(14)

Fig. 6a shows the DFT energy profile as function of the reaction coordinate. The rotation of the alkyl chain from the all trans configuration to a configuration with the terminal methyl group facing the silyl radical (panel I in Fig. 6a) has a small energy barrier of 8.17 kcal/mol (the structure isn’t shown here). The next step is the hydrogen abstraction by the silyl radical: ≡Si–CH2CH2CH3 + ≡Si· → ≡Si–CH2CH2CH2· + ≡Si–H

(15)

which has an energy barrier of 36.6 kcal/mol (panel II, Fig. 6b) and yields the intermediate shown in panel III. In the last step a pentacoordinated Si atom is formed with Eact = 16.8 kcal/mol (Fig. 5b, panels IV-V). The hydrogen abstraction from the terminal methyl group in Eq. (15) has a similar energy barrier than from the abstraction of an H atom from the beta C atom (see Table 2, two Si atoms) which is 31.5 kcal/mol for the propyl group. This explains the competition observed between both processes in the MD simulations.21 There are no experimental results about the thermal stability of ≡Si–C3H7, so we propose that the controlled heating of a surface modified with propyl groups can produce structured surfaces with chains lying down and standing up on the surface.

ACS Paragon Plus Environment

20

Page 21 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 6. a) DFT Energy profile along the reaction coordinate for the formation of Si−CH2CH2CH2−Si moieties. The insert shows the structure of initial structure of the reaction. b) Panels showing the structure of critical points along the reaction coordinate.

ACS Paragon Plus Environment

21

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

3.3. ReaxFF MD simulations 3.3.1. Functionalization of Si(111), Si(100)−2×1 and Si(100)−hf Surfaces In this section we investigate the functionalization process of the bare Si(111), Si(100)−2×1 and Si(100)−hf surfaces to evaluate the maximum attainable surface coverages that can be reached. To speed up the calculations, a gas of alkyl radicals was placed a few Å above the bare surfaces in which the top atoms are silyl radicals. We considered the •CH3, •C2H5, •C3H7, •C5H11 and •C10H21 radicals. Around 200 molecules were added to the simulation box using Packmol44 and the MD calculation was performed at 500 K. As an example, Fig. 7a shows the initial structure in the case of the Si(100) )−2×1 surface immersed in methyl radicals. After 750 ps of simulation time, the surface is completely methylated. The variation of the surface density of alkyl species as a function of time for the three surfaces is shown in Fig. 7b. The horizontal lines in the figure indicate for each case the surface densities corresponding to a coverage of 1.0 (all surface Si atoms alkylated) and 0.5. Fig. 7b shows that in the case of the methyl radical, the maximum surface coverage is readily attained for the Si(111) and Si(100)−2×1 surfaces in the time frame investigated. In the case of the Si(111) surface, it is well documented experimentally that each surface Si atom can be methylated.45,46 The highest density of surface methyl groups is attained on the surface Si(100)−hf as it has the highest density of surface sites (see Fig. 1c). This trend is observed for all the alkyl chains. However on the Si(100)−hf surface the maximum coverage of methyl groups is not observed because the functionalization of the top bicoordinated Si atoms partially precludes the functionalization of the Si dimers which are below (Fig. 1c). For all alkyl chains investigated, around 85 % of functionalized silicon atoms correspond to top bicoordinated Si atoms of the Si(100)–hf surface.

ACS Paragon Plus Environment

22

Page 22 of 38

Page 23 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 7. a) initial (0 ps) and final states (750 ps) for the functionalization process of Si (100)−2×1 surface with methyl radicals. b) Surface density of alkyl species as a function of time for the Si(111), Si (100)−2×1and Si(100)-hf surfaces. The pink horizontal dashed line indicates for each surface the surface densities corresponding to coverages of 1.0 and 0.5. c) Final structure of Si(100)-hf surface functionalized with C5 alkyl chains. Note that almost all functionalized Si atoms correspond to the first layer.

ACS Paragon Plus Environment

23

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 38

For the longer alkyl chains, steric effects limit the number of chains that can be packed on the surface. However, the 0.5 coverage is rapidly reached for C2, C3 and C5. In the timeframe of our simulations, the C2 monolayer reaches a coverage 0.7 on the Si(111) surface and 0.8 on the Si(100)−2×1 surface. These results are consistent with the experimental

observation

that

Si(111)−C2H5

surfaces

prepared

by

the

chlorination/alkylation method are terminated with an estimated Si−C coverage of ∼0.8 of a monolayer on Si top sites of the Si(111) surface.15,45-47 In the case of C3 and C5, the maximum surface density on the Si(100)−hf surface is around 6 and 5 molecules/nm2, respectively. The Si(111) and Si(100)−2×1 have nearly the same surfaces densities: 4.5 molecules/nm2 for C3 and 3.5 molecules/nm2 for C5. These values correspond to surface coverages between 0.5 and 0.6 which are in agreement with experimental observations showing that on both surfaces similar coverages are obtained for organic chains larger than C2.15, 17 In the case of C10, much longer simulation times would be required as many collisions with the surface are non-reactive because the only reactive site of this long alkyl chain molecule is on one end. However, in the time frame of our simulations, we observed that similar surface densities of decyl species are attained on the three surfaces. This is in agreement with recent observations by the Hines's group, where they concluded that the structure, density, and chemical stability of dodecyl (C12H23–) monolayers on Si(111) and Si(100)−hf surfaces are essentially identical.29 This implies that the details of the crystalline surface structure become less important as the chain length increases. 3.3.2. Stability of High Coverage C10 monolayers on Si(111) Fig. 7b shows that it would take very long simulation time to achieve high coverages with the decyl chains. In order to speed up the functionalization process to evaluate the maximum coverage that could be reached for C10 monolayers, we employed the procedure shown in Fig. 8 to build the high coverage structures.

ACS Paragon Plus Environment

24

Page 25 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 8. Procedure for building high coverage structures of C10 layers. a) initial state: slabs of C10 molecules with different densities placed around 5 Å from surface. b) and c) side and top view of the final structure corresponding to a 0.8 of coverage after 100 ps of simulation time.

ACS Paragon Plus Environment

25

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Basically, we placed slabs of C10 radicals with different densities around 5 Å above the surface (Fig. 8a) and we let the system evolve at 300 K during 100 ps. The densest slab of C10 radicals had the same density of solid alkanes where the alkyl chains are separated on average by 4.5 Å. Most decyl radicals rapidly bind to the surface as shown in Figs. 8b and 8c. We reached a maximum coverage of around 0.8 showing that steric effects do not preclude the formation of high coverage monolayers. Next, in order to investigate the thermal stability of these monolayers, unreacted decyl radicals were removed and unreacted silyl groups were hydrogenated. Following this procedure, we built decyl monolayers with a surface coverage ranging from 0.5 to 0.8 in order to evaluate whether higher initial surface coverages could alter the desorption mechanism and eventually have a destabilizing effect on the monolayer structure. The desorption kinetics for the different coverages was characterized by following the amount of surface species (bound decyl molecules, SiH groups and Si–C bonds) as well as 1-decene molecules in the gas phase. Fig. 9 shows the variation in the number of the different species at 1800 K and up to 1000 ps of simulation time.

Figure 9. Number of fragments (Si−C10H21, Si–C bonds, alkene and Si–H bonds) as a function of time for the desorption process. The coverage varies from 0.5 (32 C10) to 0.8 (52 C10). The desorption temperature is 1800 K. From 32 to 52 of the 64 top Si atoms of the simulation box were partially functionalized with decyl chains, yielding surface coverages in the range of 0.5-0.8. Therefore, the initial

ACS Paragon Plus Environment

26

Page 26 of 38

Page 27 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

number of Si−C10H21 species and Si–C bonds coincides and is in the range from 32 to 52. The remaining top Si atoms are hydrogenated. Fig. 9 shows that in all cases the number of Si−C10H21 species decreases with time and after 1000 ps the surface coverage of alkyl chains has dropped to 0.3-0.4. The number of Si−C bonds also decreases with time, but it is always slightly higher than the number of Si−C10H21 species due to the formation of surface species arising by the cleavage of C−C bonds of some alkyl chains.21 The decrease in the number of Si–C10H21 species is correlated to the increase of 1-decene molecules in the gas phase and SiH groups on the surface. This indicates that the decomposition of the monolayers is primarily occurring according to reaction (3). No changes are observed in the reaction mechanism with the variation of the temperature. Fig. S3 in the Supporting Information shows that at 1600 K only the amount of adsorbed and desorbed species changes, but the mechanism is the same. The main features of the decomposition kinetics shown in Fig. 9 also hold under low coverages as shown in Fig. S4 for a 0.25 coverage. In summary, the ReaxFF simulations predict that coverages of up to 0.8 can be achieved for long alkyl chains and that the decomposition of these monolayers follows the same mechanism irrespective of coverage. In order to evaluate the effect of surface coverage on the binding energy and stability of the monolayers, we performed additional DFT calculations to evaluate if steric effects among the alkyl chains for the high coverage structures could eventually destabilize the monolayer. Figs. 10a and 10b show side and top views of the equilibrium structure obtained for a large unit cell in which there are 30 surface sites from which 21 are bound to decyl chains, yielding a surface coverage of 0.7. The stability of the alkyl chains was evaluated from the energy required to remove a decyl chain to the gas phase. London dispersion forces were included in the calculations. Five different molecules were removed from the surface one at a time and the average desorption energy was calculated. For the 0.7 coverage we obtained an average desorption energy of 110.1 kcal/mol whereas an equivalent calculation for a 0.5 coverage yielded 104.8 kcal/mo. The desorption of a single molecule (coverage 1/30) yields a binding energy of 91.9 kcal/mol. In conclusion, the bonding of a decyl chain in not destabilized in the high coverage structure. On the contrary, London dispersion forces have a stabilizing effect.

ACS Paragon Plus Environment

27

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 10. a) side and b) top views of DFT equilibrium structures of Si(111) modified with C10 alkyl chains, in which 21 of 30 surface sites are functionalized, corresponding to a surface coverage of 0.7.

ACS Paragon Plus Environment

28

Page 28 of 38

Page 29 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

3.3.3. Arrhenius Parameters from ReaxFF MD Simulations The MD simulations were used to extract the Arrhenius parameters for the desorption reaction of alkyl chains in which the main products are alkene molecules in the gas phase, Eq. (3). The first order kinetics was followed by computing the Si−alkyl consumption rate. Fig. 11a plots the normalized surface coverage of alkyl chains remaining on the surface as a function of time for C2, C3, C5 and C10 chains in the temperature range from 1500 K to 2200 K. According to a first-order reaction kinetics, the rate constant k at each temperature can be calculated by fitting the profiles in Fig. 11a to [Si–CnHn+1]t/[Si–CnHn+1]0=exp(-kt) where [Si–CnHn+1]0 and [Si–CnH2n+1]t correspond to the concentration of Si−alkyl groups at t = 0 and at an arbitrary time t, respectively. From the Arrhenius plots in Fig. 11b, the Arrhenius parameters for desorption were determined for each alkyl chain. Table 3 collects the activation energy barriers and the pre-exponential factors obtained. Table 3. Activation energy barriers and pre-exponential factors for alkyl chain desorption from Si(111) determined from the ReaxFF MD simulations. Alkyl chain

Eact (kcal/mol)

Pre-exponential Factor (s-1)

C2

32.30

3.74 x 1013

C3

32.40

5.05 x 1013

C5

34.51

7.87 x 1013

C10

37.68

15.2 x 1013

From the comparison of Tables 2 and 3, it can be observed that the activation energy barriers determined from the ReaxFF MD simulations are in excellent agreement with those calculated by DFT using the CI-NEB method (Table 2). For C2 and C3, for example, the ReaxFF and DFT values are around 32 kcal/mol whereas for C5 and C10 both approaches predict barriers which are a few kcal/mol higher. These results imply that the main decomposition mechanism of the alkyl chains grafted to Si involve the desorption of 1alkene molecules mediated by adjacent silyl radicals, according to Eq. (3).

ACS Paragon Plus Environment

29

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 38

Table 3 shows that the pre-exponential factor increases with the alkyl chain length and Fig. 11c shows that such increase is linear. The pre-exponential factor ν can be correlated with the frequency of attempted adsorbed molecules to overcome the desorption energy barrier and become desorbed. The pre-exponential factor includes the change of all rotational, translational and vibrational degrees freedom during desorption. From transition state theory (TST), the pre-exponential factor can be described as: 48 =

   ∗   

(18)

where kb and h are the Boltzman and Planck’s constants, q is partition function of the adsorbed state and q* is the partition function for the desorbed state. For atoms and small molecules both partition functions are similar and therefore  ≈  ⁄ℎ# ≈ 1013 '−1. In the case of large molecules, the desorbed state has many rotational and vibrational degrees of freedom which are not available in the adsorbed state, where only frustrated rotations and vibrations exist. Consequently q* is larger than q which leads to an increase in the pre exponential factor ν.48 Therefore, the ReaxFF MD simulations correctly predict the increase in ν with the alkyl chain length.

ACS Paragon Plus Environment

30

Page 31 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Fig. 11. a) Normalized surface coverage of Si−alkyl groups remaining on the surface as a function of time, in the temperature range from 1500 K to 2200 K, for C2, C3, C5 and C10 chains. b) Arrhenius plot for the thermal desorption of alkyl chains as alkene molecules from Si(111) surface. c) pre-exponential factor obtained from Arrhenius plot as a function of number of C atoms in the alkyl chain.

ACS Paragon Plus Environment

31

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

4. Conclusions In this work we developed a Si/C/H ReaxFF force field for the study of the functionalization of silicon surfaces as well as the decomposition of grafted alkyl monolayers. The main reactions involved in the decomposition of Si(111)−alkyl layers were considered in the parameterization of the force field. The reaction pathways observed in the MD simulations were validated by DFT calculations of activation energy barriers using the CI-NEB method. The main decomposition mechanism of the alkyl chains is the alkene elimination to the gas phase after a β-hydride abstraction by silyl radicals (Eq. (3)). Arrhenius plots obtained from the MD simulations for C2, C3, C5 and C10 chains yielded energy barriers (Table 3) in excellent agreement with the DFT values (Table 2). Even the slight increase in the energy barrier with the chain length is also well reproduced by the ReaxFF force field. Pre-exponential factors increase linearly with the alkyl chain length, in agreement with transition state theory predictions for desorption of large molecules. The Si/C/H potential developed in this work also adequately describes the reactivity of organic compounds having C atoms with sp2 and sp hybridizations, as we showed in a previous work21 in which we considered the stability of Si−CHCHCH3 and Si−CCCH3 monolayers. The ReaxFF force field was employed to investigate the functionalization of Si(111), Si(100)−2×1 and Si(100)−hf surfaces in order to determine the maximum surface coverages that can be reached on each surface. Full coverage methylated surfaces were readily obtained for Si(111) and Si(100)−2×1 whereas a coverage of 0.9 was obtained on the Si(100)−hf, although the surface density of methyl groups is higher on this surface. In the timeframe of our simulations, the C2 monolayer reaches a coverage 0.7 on the Si(111) surface and 0.8 on the Si(100)−2×1 surface whereas for C3 and C5 similar coverages are observed on both surfaces (0.5 and 0.6) in agreement with coverages reported in the literature.11, 13, 15-17 In the case of C10, the details of the crystalline surface structure begin to blur as the same coverage was obtained on all surfaces in agreement with experimental observations.29 Both the DFT calculations and the ReaxFF MD simulations predict that for C10 stable monolayers with coverages of up to 0.8 can be formed, although experimentally coverages around 0.5 have only been reported. Although these high coverages are feasibly

ACS Paragon Plus Environment

32

Page 32 of 38

Page 33 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

from an energetic point of view, they may be difficult to attain in experiments due to steric effects. The structure and dynamics of alkyl chains bonded to surfaces has important implications in reversed-phase liquid chromatography (RPLC) where the separation of organic molecules depends on their interactions with the grafted solid substrate (the stationary phase) and the mobile phase.49 In this context, our ReaxFF potential could be considered a first step to include reactivity in RPLC simulations. In summary, the developed Si/C/H ReaxFF force field correctly describes the details of the functionalization and decomposition mechanisms of alkyl layers on silicon surfaces which paves the way for further studies of functionalized nanostructured silicon systems such as functionalized nanopores, nanotubes, nanowires or two-dimensional silicene.

5. Associated Content The Supporting Information is available free of charge on the ACS Publications website. Molecular clusters used in the training set of the force field, additional MD simulations and the force field parameters. 6. Author Information Corresponding Author E-mail: [email protected] Notes The authors declare no competing financial interest. 7. Acknowledgements EMP thanks Secyt-UNC and Foncyt (PICT-2014-2199) for financial support. ACTvD acknowledges funding from NSF grant MIP/DMR 1539916. This work used computational resources from CCAD – Universidad Nacional de Córdoba (http://ccad.unc.edu.ar/), which is part of SNCAD – MinCyT, República Argentina.

ACS Paragon Plus Environment

33

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

8. References 1. Yang, F. A., P.; Ozanam, F.; Chazalviel, J.-N, Thermal Stability of Organic Monolayers Covalently Grafted on Silicon Surfaces. In Reactions and Mechanisms in Thermal Analysis of Advanced Materials, Tiwari, A.; Raj, B., Eds. John Wiley & Sons: Massachusetts, 2015; pp 3-37. 2. Wong, K. T.; Lewis, N. S., What a Difference a Bond Makes: The Structural, Chemical, and Physical Properties of Methyl-Terminated Si(111) Surfaces. Acc. Chem. Res. 2014, 47 (10), 3037-3044. 3. Scheres, L.; Giesbers, M.; Zuilhof, H., Organic Monolayers onto Oxide-Free Silicon with Improved Surface Coverage: Alkynes versus Alkenes. Langmuir 2010, 26 (7), 47904795. 4. Pujari, S. P.; Filippov, A. D.; Gangarapu, S.; Zuilhof, H., High-Density Modification of H-Terminated Si(111) Surfaces Using Short-Chain Alkynes. Langmuir 2017, 33 (51), 14599-14607. 5. Peng, W.; Rupich, S. M.; Shafiq, N.; Gartstein, Y. N.; Malko, A. V.; Chabal, Y. J., Silicon Surface Modification and Characterization for Emergent Photovoltaic Applications Based on Energy Transfer. Chem. Rev. 2015, 115 (23), 12764-12796. 6. Carroll, G. M.; Limpens, R.; Neale, N. R., Tuning Confinement in Colloidal Silicon Nanocrystals with Saturated Surface Ligands. Nano Lett. 2018, 18 (5), 3118-3124. 7. Cheng, X.; Lowe, S. B.; Reece, P. J.; Gooding, J. J., Colloidal silicon quantum dots: from preparation to the modification of self-assembled monolayers (SAMs) for bioapplications. Chem. Soc. Rev. 2014, 43 (8), 2680-2700. 8. Huey, W. L. B.; Goldberger, J. E., Covalent functionalization of two-dimensional group 14 graphane analogues. Chem. Soc. Rev. 2018, 47 (16), 6201-6223. 9. Sung, M. M.; Kluth, G. J.; Yauw, O. W.; Maboudian, R., Thermal Behavior of Alkyl Monolayers on Silicon Surfaces. Langmuir 1997, 13 (23), 6164-6168. 10. Ryo, Y.; Masato, A.; Hirokazu, T., Temperature Dependence of the Structure of Alkyl Monolayers on Si(111) Surface via Si–C Bond by ATR-FT-IR Spectroscopy. Chem. Lett. 2004, 33 (5), 492-493. 11. Faucheux, A.; Yang, F.; Allongue, P.; Henry de Villeneuve, C.; Ozanam, F.; Chazalviel, J. N., Thermal decomposition of alkyl monolayers covalently grafted on (111) silicon. Appl. Phys. Lett. 2006, 88 (19), 193123. 12. Jaeckel, B.; Hunger, R.; Webb, L. J.; Jaegermann, W.; Lewis, N. S., HighResolution Synchrotron Photoemission Studies of the Electronic Structure and Thermal Stability of CH3- and C2H5-Functionalized Si(111) Surfaces. J. Phys. Chem. C 2007, 111 (49), 18204-18213.

ACS Paragon Plus Environment

34

Page 34 of 38

Page 35 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

13. Hunger, R.; Fritsche, R.; Jaeckel, B.; Jaegermann, W.; Webb, L. J.; Lewis, N. S., Chemical and electronic characterization of methyl-terminated Si(111) surfaces by highresolution synchrotron photoelectron spectroscopy. Phys. Rev. B 2005, 72 (4), 045317. 14. Yang, F.; Roodenko, K.; Hunger, R.; Hinrichs, K.; Rademann, K.; Rappich, J., Near-Ideal Complete Coverage of CD3 onto Si(111) Surfaces Using One-Step Electrochemical Grafting: An IR Ellipsometry, Synchrotron XPS, and Photoluminescence Study. J. Phys. Chem. C 2012, 116 (35), 18684-18690. 15. Nemanick, E. J.; Hurley, P. T.; Brunschwig, B. S.; Lewis, N. S., Chemical and Electrical Passivation of Silicon (111) Surfaces through Functionalization with Sterically Hindered Alkyl Groups. J. Phys. Chem. B 2006, 110 (30), 14800-14808. 16. Nemanick, E. J.; Solares, S. D.; Goddard, W. A.; Lewis, N. S., Quantum Mechanics Calculations of the Thermodynamically Controlled Coverage and Structure of Alkyl Monolayers on Si(111) Surfaces. J. Phys. Chem. B 2006, 110 (30), 14842-14848. 17. Puniredd, S. R.; Assad, O.; Haick, H., Highly Stable Organic Modification of Si(111) Surfaces: Towards Reacting Si with Further Functionalities while Preserving the Desirable Chemical Properties of Full Si−C Atop Site Terminations. J. Am. Chem. Soc. 2008, 130 (29), 9184-9185. 18. Puniredd, S. R.; Assad, O.; Haick, H., Highly Stable Organic Monolayers for Reacting Silicon with Further Functionalities: The Effect of the C−C Bond nearest the Silicon Surface. J. Am. Chem. Soc. 2008, 130 (41), 13727-13734. 19. Sieval, A. B.; Demirel, A. L.; Nissink, J. W. M.; Linford, M. R.; van der Maas, J. H.; de Jeu, W. H.; Zuilhof, H.; Sudhölter, E. J. R., Highly Stable Si−C Linked Functionalized Monolayers on the Silicon (100) Surface. Langmuir 1998, 14 (7), 17591768. 20. Zhang, L.; Wesley, K.; Jiang, S., Molecular Simulation Study of Alkyl Monolayers on Si(111). Langmuir 2001, 17 (20), 6275-6281. 21. Soria, F. A.; Zhang, W.; van Duin, A. C. T.; Patrito, E. M., Thermal Stability of Organic Monolayers Grafted to Si(111): Insights from ReaxFF Reactive Molecular Dynamics Simulations. ACS Appl. Mater. Interfaces 2017, 9 (36), 30969-30981. 22. Giannozzi, P.; Baroni, S.; Bonini, N.; Calandra, M.; Car, R.; Cavazzoni, C.; Ceresoli, D.; Chiarotti, G. L.; Cococcioni, M.; Dabo, I., et al., QUANTUM ESPRESSO: a modular and open-source software project for quantum simulations of materials. J. Phys.: Condens. Matter 2009, 21 (39), 395502. 23. Vanderbilt, D., Soft self-consistent pseudopotentials in a generalized eigenvalue formalism. Phys. Rev. B 1990, 41 (11), 7892-7895. 24. Perdew, J. P.; Burke, K.; Ernzerhof, M., Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77 (18), 3865-3868.

ACS Paragon Plus Environment

35

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

25. Henkelman, G.; Uberuaga, B. P.; Jónsson, H., A climbing image nudged elastic band method for finding saddle points and minimum energy paths. J. Chem. Phys. 2000, 113 (22), 9901-9904. 26. Juarez, M. F.; Soria, F. A.; Patrito, E. M.; Paredes-Olivera, P., Influence of Subsurface Oxidation on the Structure, Stability, and Reactivity of Grafted Si(111) Surfaces. J. Phys. Chem. C 2008, 112 (38), 14867-14877. 27. Monkhorst, H. J.; Pack, J. D., Special points for Brillouin-zone integrations. Phys. Rev. B 1976, 13 (12), 5188-5192. 28. Clark, I. T.; Aldinger, B. S.; Gupta, A.; Hines, M. A., Aqueous Etching Produces Si(100) Surfaces of Near-Atomic Flatness: Strain Minimization Does Not Predict Surface Morphology. J. Phys. Chem. C 2010, 114 (1), 423-428. 29. DeBenedetti, W. J. I.; Li, T. L.; Hines, M. A., Half-flat vs. atomically flat: Alkyl monolayers on morphologically controlled Si(100) and Si(111) have very similar structure, density, and chemical stability. J. Chem. Phys. 2017, 146 (5), 052804. 30. van Duin, A. C. T.; Dasgupta, S.; Lorant, F.; Goddard, W. A., ReaxFF:  A Reactive Force Field for Hydrocarbons. J. Phys. Chem. A 2001, 105 (41), 9396-9409. 31. ADF2017 Molecular Modeling Suite, Scientific Computing and Modeling; Vrije Universiteit, Amsterdam, The Netherlands, http://www.scm.com. 32. Plimpton, S., Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117 (1), 1-19. 33. Newsome, D. A.; Sengupta, D.; Foroutan, H.; Russo, M. F.; van Duin, A. C. T., Oxidation of Silicon Carbide by O2 and H2O: A ReaxFF Reactive Molecular Dynamics Study, Part I. J. Phys. Chem. C 2012, 116 (30), 16111-16121. 34. Psofogiannakis, G.; van Duin, A. C. T., Development of a ReaxFF reactive force field for Si/Ge/H systems and application to atomic hydrogen bombardment of Si, Ge, and SiGe (100) surfaces. Surf. Sci. 2016, 646 (Supplement C), 253-260. 35. Chenoweth, K.; van Duin, A. C. T.; Goddard, W. A., ReaxFF Reactive Force Field for Molecular Dynamics Simulations of Hydrocarbon Oxidation. J. Phys. Chem. A 2008, 112 (5), 1040-1053. 36. 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 B.01. Wallingford CT, 2009. 37. van Duin, A. C. T.; Baas, J. M. A.; van de Graaf, B., Delft molecular mechanics: a new approach to hydrocarbon force fields. Inclusion of a geometry-dependent charge calculation. J. Chem. Soc., Faraday Trans. 1994, 90 (19), 2881-2895.

ACS Paragon Plus Environment

36

Page 36 of 38

Page 37 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

38. Cicero, R. L.; Chidsey, C. E. D.; Lopinski, G. P.; Wayner, D. D. M.; Wolkow, R. A., Olefin Additions on H−Si(111):  Evidence for a Surface Chain Reaction Initiated at Isolated Dangling Bonds. Langmuir 2002, 18 (2), 305-307. 39. Takeuchi, N.; Kanai, Y.; Selloni, A., Surface Reaction of Alkynes and Alkenes with H-Si(111):  A Density Functional Theory Study. J. Am. Chem. Soc. 2004, 126 (48), 1589015896. 40. Reboredo, F. A.; Schwegler, E.; Galli, G., Optically Activated Functionalization Reactions in Si Quantum Dots. J. Am. Chem. Soc. 2003, 125 (49), 15243-15249. 41. Ciampi, S.; Böcking, T.; Kilian, K. A.; James, M.; Harper, J. B.; Gooding, J. J., Functionalization of Acetylene-Terminated Monolayers on Si(100) Surfaces:  A Click Chemistry Approach. Langmuir 2007, 23 (18), 9320-9329. 42. Linford, M. R.; Fenter, P.; Eisenberger, P. M.; Chidsey, C. E. D., Alkyl Monolayers on Silicon Prepared from 1-Alkenes and Hydrogen-Terminated Silicon. J. Am. Chem. Soc. 1995, 117 (11), 3145-3155. 43. Sieval, A. B.; Vleeming, V.; Zuilhof, H.; Sudhölter, E. J. R., An Improved Method for the Preparation of Organic Monolayers of 1-Alkenes on Hydrogen-Terminated Silicon Surfaces. Langmuir 1999, 15 (23), 8288-8291. 44. Martínez, L.; Andrade, R.; Birgin, E. G.; Martínez, J. M., PACKMOL: A package for building initial configurations for molecular dynamics simulations. J. Comput. Chem. 2009, 30 (13), 2157-2164. 45. Webb, L. J.; Rivillon, S.; Michalak, D. J.; Chabal, Y. J.; Lewis, N. S., Transmission Infrared Spectroscopy of Methyl- and Ethyl-Terminated Silicon(111) Surfaces. J. Phys. Chem. B 2006, 110 (14), 7349-7356. 46. Yu, H.; Webb, L. J.; Heath, J. R.; Lewis, N. S., Scanning tunneling spectroscopy of methyl- and ethyl-terminated Si(111) surfaces. Appl. Phys. Lett. 2006, 88 (25), 252111. 47. Yu, H.; Webb, L. J.; Solares, S. D.; Cao, P.; Goddard, W. A.; Heath, J. R.; Lewis, N. S., Scanning Tunneling Microscopy of Ethylated Si(111) Surfaces Prepared by a Chlorination/Alkylation Process. J. Phys. Chem. B 2006, 110 (47), 23898-23903. 48. Winkler, A., Thermal Desorption Of Organic Molecules. In Interface Controlled Organic Thin Films, Al-Shamery, K.; Horowitz, G.; Sitter, H.; Rubahn, H.-G., Eds. Springer Berlin Heidelberg: Berlin, Heidelberg, 2009; pp 29-36. 49. Lindsey, R. K.; Rafferty, J. L.; Eggimann, B. L.; Siepmann, J. I.; Schure, M. R., Molecular simulation studies of reversed-phase liquid chromatography. J. Chromatogr. A 2013, 1287, 60-82.

ACS Paragon Plus Environment

37

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

TOC Graphic

ACS Paragon Plus Environment

38

Page 38 of 38