Nano-Scale Corrugations in Graphene: A Density Functional Theory

Subscriber access provided by SUNY DOWNSTATE

Article

Nano-Scale Corrugations in Graphene: a Density Functional Theory Study of Structure, Electronic Properties and Hydrogenation Antonio Rossi, Simone Piccinin, Vittorio Pellegrini, Stefano de Gironcoli, and Valentina Tozzini J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/jp511409b • Publication Date (Web): 18 Mar 2015 Downloaded from http://pubs.acs.org on March 25, 2015

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 free 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 accessible to all readers and 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.

The Journal of Physical Chemistry C 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 25

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

Nano-Scale Corrugations in Graphene: A Density Functional Theory Study of Structure, Electronic Properties and Hydrogenation Antonio Rossia,b , Simone Piccininc,d, Vittorio Pellegrinib,e, Stefano de Gironcolid,c, Valentina Tozzinie,f∗ a

Dipartimento di Fisica ‘E. Fermi’, Università di Pisa Largo B. Pontecorvo 3-56127 Pisa, Italy, b Graphene Labs, iit Istituto Italiano di Tecnologia, Via Morego, 30 16163 Genova, Italy, c CNR-IOM DEMOCRITOS c/o SISSA, Via Bonomea 265, 34136 Trieste, Italy, d SISSA, Scuola Internazionale Superiore di Studi Avanzati, via Bonomea 265, 34136 Trieste, Italy, e Istituto Nanoscienze, Cnr, Piazza San Silvestro 12, 56127 Pisa, Italy, f NEST- Scuola Normale Superiore, Piazza San Silvestro 12, 56127 Pisa, Italy

Abstract Graphene rippled at the nano-scale level is gathering attention for advanced applications, especially in the field of nano-electronics and hydrogen storage. Convexity enhanced reactivity towards H was demonstrated on naturally corrugated graphene grown by Si evaporation on SiC, which makes this system a platform for fundamental studies on the effects of rippling. In this work, we report a Density Functional Theory study on a model system specifically designed to mimic graphene on SiC. We first study the supercell geometry and configuration that better reproduce the corrugated monolayer. The relatively low computational cost of this model system allows a systematic study of the dependence of stability, structure and electronic properties of graphene subject to different levels of stretching and corrugation. The most representative structure is then progressively hydrogenated, imitating the exposure to atomic hydrogen, and stability, structural and electronic properties are evaluated as a function of hydrogenation. Our results quantitatively reproduce the measured evolution of electronic properties as a function of hydrogenation, offering the possibility of evaluating the coverage by means of STS measurements. The dependence of hydrogen binding energy on coverage extends our previous results on reactivity of corrugated graphene, including the effect of H clustering. This work reports quantitative results directly comparable with experimental measurements performed on epitaxial graphene on SiC and reveals the quantitative interplay between local structure, electronic properties and reactivity to hydrogen, which could be used to design devices for flexible nano-electronics and for H storage. Keywords: nano-structured 2D materials, ab initio molecular dynamics, structure optimization, reactivity

∗ Corresponding author: Valentina Tozzini, Istituto Nanoscienze – Cnr, NEST-SNS, Piazza San Silvestro 12, 56127 Pisa Italy. Tel: +39 050 509433. Fax: +39 050 509 417. E-mail: [email protected]

ACS Paragon Plus Environment


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

Introduction The interaction between hydrogen and graphene has recently received much attention, due to its potential interest for different technological applications. While graphene is a high mobility conductor1, its alkane counterpart, graphane2 obtained covalently bonding one hydrogen atom to each carbon site, is a wide band gap insulator 3 . Partially hydrogenated graphene displays intermediate properties4: stripped graphane/graphene hybrids are semiconductors with tunable band gaps 5,6,7,8,9; different hydrogen decorations, such as graphane islands10, have potentially interesting transport properties and possible applications in nano-electronics. Clearly, exploiting these properties requires a control of the hydrogenation at the nano-scale. Graphene hydrogen interaction is clearly also interesting for H storage applications. Due to its low molecular weight, graphene is potentially a storage mean with high gravimetric capacity11. However, molecular hydrogen is very weakly physisorbed by means of van der Waals interactions12, and although these are stronger in nanostructured13,14 or multilayered graphene15, low temperatures are necessary for stable storage 16 , 17 . Chemisorption was alternatively considered: the covalent bond of hydrogen to graphene is robust, leading to stable storage up to high temperature. But chemi(de)sorption of H2 are high barrier processes (~1–1.5 eV/atom18) implying slow loading/release kinetics. Chemical (e.g. with Pd19) or “alternative” catalysis (e.g. electric fields20 or N-substitutional doping21) were also suggested. We recently proposed that the two-dimensionality and structural/mechanical properties of graphene could be also exploited to control chemisorption: the local curvature enhances carbon reactivity 22 , specifically to hydrogen 23 , due to the “pyramidalization” of the carbon sites on convexities, and consequent extrusion of an unpaired sp3-like orbitals. This suggests the possibility of (de)hydrogenating graphene by controlling the local curvature. The proof of this concept was given in ref 24 , where the dependence of hydrogen binding energy on local curvature was quantified, and its preferential binding on convex sites and spontaneous release by curvature inversion demonstrated by means Density Functional Theory (DFT) calculations and Car-Parrinello (CP) simulations. This study also suggested the possibility of using this effect as the basis for a hydrogen storage device. Curvature control was also proposed to create specific hydrogen decorations for nano-electronic devices25,26. The effect was experimentally confirmed on the naturally corrugated freestanding graphene layer grown by Si evaporation from 6HSiC(0001)27. A number of computational studies, mainly DFT-based, have addressed hydrogen chemisorption and clustering on flat graphene 28 , 29 , 30 , 31 , while the effect of curvature was considered in fullerenes32 and nanotubes33 and in tightly and artificially rippled graphene24,22. Naturally rippled graphene grown by Si evaporation on SiC is a natural platform to experimentally study the effects of curvature on the physic-chemical properties of graphene. On SiC, however, graphene curvature follows a Moiré pattern with a specific supercell periodicity due to the mismatch with the substrate lattice parameters. Another set of DFT studies analyzes specific aspects of the structural and electronic properties of this system, either including the exact 34 or approximate 35 supercell periodicity, portion of the substrate 36 , 37 , intercalants 38 or dopants39. Though the interplay between curvature and hydrogen adhesion is proven, the route towards its exploitation requires further, more specific and quantitative investigations. This study addresses the interplay between curvature, electronic properties, and hydrogen binding on curved graphene, on a model system optimized for reproducing the features of graphene on SiC. The specific aim is to provide an accurate, yet the simplest possible, computational representation of the real system, for direct comparison with experiment and


Page 2 of 25

Page 3 of 25

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


interpretation of measurements. In addition, the relative simplicity of the model allows extensive and systematic calculations. Models with different levels of curvature are evaluated, and their structural and electronic properties investigated. Subsequently, those better corresponding to the experimental ones are selected and hydrogenated in stages, mimicking the process occurring upon exposure to hydrogen. Increasingly hydrogenated systems are analyzed, and their structural and electronic properties related to H-coverage. After the discussion of results, conclusions and future developments are illustrated, with focus on the applications in hydrogen storage and nanoelectronics, respectively. 1. Models and Methods 1.1 Model systems Si evaporation from the (0001) face of SiC produces first a honeycomb carbon lattice, covalently bonded to the substrate. If the process is prolonged, the first completely detached sheet (“free standing”27, or monolayer) is produced, which displays the graphene electronic properties. As a consequence of the small mismatch between SiC and graphene parameter lattices, the monolayer displays a Moiré pattern of corrugation, whose exact periodicity is 6√3×6√3R30 referred to SiC, or 13×13 referred to graphene40 (red boundaries in Fig 1(a,c)). The periodicity of rippling can also be approximately described by a 6×6 supercell (referred to SiC, or 4√3×4√3R30, if referred to graphene, green boundaries in Fig 1(a,c)), thus both supercells are considered in this study. Smaller supercells (specifically the standard unit cell, the minimal rotated by 30 deg, √3×√3R30, and an intermediate one, the 3√3×3√3R30) are considered for comparison. All supercells and their Brillouin Zone (BZ) are reported in Fig 1 (a,b and d). It is to be observed that the Dirac point (“blue” K in Fig 1(d)) is remapped in different special points upon BZ refolding when different supercell are considered.

Fig 1 (a) Supercells considered in this study (in colors), superimposed to the graphene lattice (in black). In (b) The cell periodicity (with respect to graphene and SiC lattice) and the number of atoms included (corresponding colors to (a)). (c) The two supercells with the periodicity of the Moiré lattice (exact in red, and approximate, in green). In the background a representation of the corrugation is reported: light grey areas are protruding areas. (d) The Brillouin zones of the main cells (corresponding colors to (a)). The main symmetry points of the unit and of the 6x6 cell are reported in blue and green, respectively.



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 periodicity along z axis is set at 15.00Å, thus the model systems is a periodic multilayer with negligible inter-sheet interaction. The sheets were increasingly rippled and stretched by uniformly reducing or increasing the lateral cell size (by step of 1%), and subsequently relaxing atomic position. For the 13x13 cell, the formation of ripples was guided by superimposed starting small corrugations, as described in Section 3. The local curvature is measured by the average value of the pyramidalization pseudo-dihedral φ related to a single C site41. Selected systems are then hydrogenated with a protocol emulating the kinetics of adsorption occurring upon exposure to low concentration atomic hydrogen27, and in agreement with binding energy enhancement with curvature24: couples of H atoms are added onto the sites with highest local curvature in dimers or tetramers, on a single side of the sheet. The system is then relaxed and the procedure repeated. The stability of selected hydrogenated systems is then tested by means of simulated annealing with CP dynamics. Selected cases of double side hydrogenations are also considered, for comparison. 1.2 Density Functional Theory calculations and simulations setup DFT calculations were performed expanding the Kohn-Sham orbitals in plane waves and using Vanderbilt ultrasoft pseudopotentials 42 with energy cutoff set at 25 Ry (see SI, S1 for convergence check). van der Waals (vdW) corrections were added, according to the scheme by Grimme43. Perdew-Burke-Ernzerhof (PBE) exchange and correlation functional44 was used for production calculations, but preliminary ones were also performed in Local Density Approximation (LDA)45 for comparison (see SI, S2.1). For electronic structure calculations the Brillouin zone was sampled with Monkhorst-Pack grids with sufficiently dense sampling46. For force and energy calculations, converge checks were performed in different compression conditions with respect to the grid sampling (see SI.1 for details). After noting that cohesive energies changes less or about 1%, we finally chose the Γ point scheme for structure calculations. Gaussian smearing with 0.01Ry spread was used. Local structural minima were searched using the Broyden Fletcher Goldfarb Shanno algorithm47, preceded by coordinated randomization in a few cases, and using standard convergence criteria for electronic self-consistency (10-8au) and forces (10-3au). Molecular dynamics was performed within the CP48 scheme, with the preconditioning scheme for the orbital mass49, and a time step of 0.1205 fs. Equation integration was performed with Verlet algorithm, Nosé thermostat 50 is used to control the temperature. Annealing procedures are applied to selected systems. The temperature is increased by a combination of coordinate randomization and velocity rescaling. Calculations and simulations are performed with Quantum Espresso (QE5.0.151). Inputs creation and post-processing analysis was performed using in-house made software. 2. Results 2.1 Structure and stability of rippled and stretched graphene A model system with a few layers of SiC substrate and hydrogen includes ~2000 atoms. While the simulations of such a large model system is currently in progress for selected cases on high performance computational facilities, here we are interested in optimizing simpler model systems (without including explicitly the SiC substrate), equivalent in terms of structural and electronic properties of the graphene layers, to be used for extensive calculations. To this aim, the 13×13 and 4√3×4√3R30 cells were considered, as they match (exactly or approximately) the Moiré symmetry. These were laterally compressed (stretched) in the range +(-)10% with respect to their


Page 4 of 25

Page 5 of 25

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


relaxed size. The same calculation was also performed on the 3√3×3√3R30 cell, to obtain a rippling with a smaller wavelength and higher level of curvature. Two different functionals (LDA and PBE) with and without vdW interactions were compared for the smaller supercells. The energies of optimized structures are reported in Table S3-S5 in the SI and in Fig 2(a), for the 4√3×4√3R30 (green dots), the 3√3×3√3R30 (magenta dots) and 13x13 cells (red dots). The general behavior of the energy as a function of the contraction/expansion (here measured in terms of C=C distance of the flat structure in the corresponding cell size) displays a minimum at zero contraction (C=C = 1.424 Å), with a harmonic-like behavior around it (irregularities are discussed below). Fig 2(a) also shows that PBE (filled symbols) returns a value of the cohesive energy ~1eV smaller than the LDA (empty symbols) and nearer to the experimental one 52 . vdW interactions bring little change in the cohesive energy (squares vs circles), however they were included in all subsequent calculations because they are likely to play a role in the interaction with hydrogen24. For the smaller cells, and for selected structures, convergence checks on the cohesive energy and on the structure was performed with the respect to the increasing sampling of the BZ (see SI Sec S1). A stabilization of ~0.6% is observed, with negligible changes in the structure, which justify the use of Γ point only to sample the BZ for structure and forces calculations.

Fig 2 (a) Coesive energies per C atom as a function of the strain, measured by the “effective” C-C bondlength (i.e. the bond C-C bond length one would have if the graphene remained flat and regular). Magenta: 3√3×3√3R30 cell, green: 4√3×4√3R30, red=13x13. Filled dots=PBE+vdW, empty dots=LDA+vdW, empty squares=LDA. (b) and (c) top and perspective view of the 2% contracted 4√3×4√3R30° supercell (3x3 repeated cells). In (b) the unit cell boundaries (in green) and the 13x13 cell boundaries (red) are superimposed. Structures colored according to the z coordinate: blue = protruding, red= intruding.



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

We now comment on the structural/energetic features of each kind of supercell. The 13×13 has the periodicity of the graphene grown on SiC. The lattice parameter mismatch between SiC and graphene is such that the freestanding layer is ~0.1-0.2% contracted. We then first performed structure relaxation on a 0.2% contracted cell, starting from different configurations with randomly displaced atoms or with regular displacements with different periodicity (see SI S2.4 for details). The final configurations display average out-of-plane C atoms displacements of ~0.04 Å, which is approximately one order of magnitude smaller than the experimental one measured for the first free standing layer on SiC27. In addition, none of the optimized structures is stabilized in the correct symmetry of the Moiré pattern of graphene on SiC. Increasing the contraction level at ~2% (Fig 2(a), C=C = 1.38Å) an out of plane displacement of ~1Å is reached, which is similar to the experimental one, but the pattern of displacement is still dissimilar. No improvement is obtained at higher contraction (~10%). Our conclusion are (i) a simple contraction of the 13×13 supercell cannot reproduce the experimental corrugation, indicating that the substrate has a major role in determining the Moiré pattern (ii) average out of plane displacements corresponding to experimental ones are obtained with larger levels of compression with respect to the real one, probably to compensate for the absence of the effect of the highly corrugated buffer layer. In light of this, we considered the 4√3×4√3R30 supercells, whose symmetry corresponds approximately to that of the ripples. For values of the contraction up to ~5% we obtain a pattern very similar to the experimental one (see Fig 2(b),(c) for 2% contraction and inset to Fig 3 for 1%), consisting of hills (colored in blue) separated by inter-connected valleys (in red). Noticeably, this pattern is not symmetric with respect to the plane: while this asymmetry is expected in the experiment due to the substrate, it seems to spontaneously persist even without it. The inverted pattern (wells separated by interconnected ripples) also exists with the same energy, in absence of the substrate. By an analysis of the contracted supercells (see SI S2.3) it turns out that the experimental out of plane displacement is obtained for contraction of ~2%. The larger value of the contraction to reproduce the experimental corrugation is necessary to compensate for the absence of the effect of the buffer layer. We conclude that the 4√3×4√3R30 supercell with ~2% contraction is the optimal choice when one wants to reproduce the corrugated structure of a monolayer of graphene on SiC without explicitly including the substrate. The cohesive energy minimum is achieved at a C-C distance of 1.424Å. The behavior, however, soon deviates from harmonicity especially in the contraction region, where deviations appear around 2-3% and at 7-8% contraction. The analysis of the structural properties reveals interesting features. We first consider the strained structures, which display a simpler behavior. Upon stretching, the cell remains flat, but the bond length distribution becomes bimodal. Bimodality starts around 5% in the case of 4√3×4√3R30 cell (Fig 3(a)) and around 1-2% strain in the case of 3√3×3√3R30 (Fig 3(b)) and is due to the formation of a pattern of benzene-like structures (hexagon with “small” bond lengths) separated by longer bonds (insets to, Fig 3(a,b)). This “benzenization” is preceded, at intermediate stretching, by states in which larger hexagons are separated by smaller bonds (6% stretching level in the 4√3×4√3R30 cell, inset of Fig 3(a)). At 10% strain bonds break in the 3√3×3√3R30 and large vacancies form with no defined symmetry. In the contraction region the structural analysis reveals substantially different behavior comparing the 4√3×4√3R30 and 3√3×3√3R30 cells. This is somehow expected, since the different cell size imposes a wavelength of ripples differing by 25%. In the 4√3×4√3R30 cell (Fig 3,(a),(c)) the bond lengths contracts regularly up to -5%, and their distribution width


Page 6 of 25

Page 7 of 25

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


increases. So does the φ dihedrals distributions, coherently with an increase of the average level of curvature. The inspection of structures (insets to Fig 3(a), the complete series reported in the SI, S2.2) reveals the formation of “hills” of the expected symmetry, which increase in amplitude as the contraction increases. Around -5% however, a transition occurs: the form of the hills changes from round to approximately triangular, with more net boundaries. Some of the sites tend to pyramidalize, as shown by the increasing of two peaks at the positive and negative wings of the φ distribution. Correspondingly longer bonds appear, also moving the center of bond length distribution to larger values. Overall these features are the signature of the sp2 to sp3 transition.

Fig 3 Structural properties of stretched and compressed graphene of the 4√3×4√3R30 (a,c) and 3√3×3√3R30 (b,d) cell. The histograms report the C=C distance distributions (a,b), and dihedral angle distribution (c,d) at the different cell compression/strain, from -10% to +10% of the relaxed cell size (for dihedrals, the distributions of strained structures are not reported, since they are flat). The insets report selected strained (one unit supercell) and contracted (3×3 supercells) structures (the level of strain/compression is indicated). In the strained structures the bonds are colored according to their length (shorter in red, longer in blue). In the contracted structures the colors represent the displacement in z direction (out of plane): blue = up, red= down. Dots and error bars in color (green = 4√3×4√3R30, magenta = 3√3×3√3R30 ) report the mean of the distribution and its variance.

As an effect of the smaller periodicity, in the 3√3×3√3R30 cell the transition features appear earlier (at 2-3% contraction) in the bond length and dihedrals distributions (Fig 3, (b) and (d)). There are additional differences: an inspection of the structures reveals at least four different ripples geometries explored by the system (insets of Fig 3(b); SI S2.2 reports the whole series):



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

ripples separated by wells at small contraction (up to -2%), linear wave-like ripples up to -6-7%, which then become wells separated by ripples (-8%) and then hills separated by valleys (-10%) passing again from a wavelike structure (-9%). Correspondingly the bond lengths and dihedral distributions assume variegate shapes, depending on the contraction. The system clearly displays multi stability, manifesting in the change of ripples conformation and geometry. As already observed, multi stabilities are always present in the presence of ripples on isolated graphene, since the structures obtained by exchanging convexities with concavities are degenerate in energy. Because these effects are less pronounced in the larger cell, one can infer that they are generated by the interplay between periodicity and compression level. In fact, it is reasonable that, if convexities and concavities are spatially separate enough (i.e. for ripples with long wavelength) and/or small in amplitude, the two degenerate structures can coexist. Otherwise they mix and tend to become unstable. In support of this picture, a recent work53, shows that at given uni-axial compression, the number of ripples formed increases as the size of the system increase. This indicates the existence of a critical periodicity, and a corresponding critical ripple wavelength, under which the system becomes unstable. A rough estimate from ref 53 would indicate a value between 1 and 2 nm at 5% uni-axial compression. In our systems compression is bi-axial, thus it is reasonable that such instability could appear at smaller wavelength and/or lower compression levels. Indeed, what we observe is that 4√3×4√3R30 cell (~1.6nm cell vector size) has a regular behavior with stable hills and valleys up to 5%, while in the 3√3×3√3R30 cell (~1.2 cell vector size) those instabilities already appear at 1-2% contraction. This indicates a critical ripple wavelength between 1.1 and 1.5 nm for our systems. In summary, the systematic analysis of graphene in the ±10% deformation range reveals a “benzeinzation” transition in the stretching region and a “pyramidalization” transition in the contraction region. The latter is preceded, at smaller contractions, by the formation of ripples with different geometries and amplitudes, which reveals instabilities as the contraction increases and/or the wavelength (i.e. cell size) decrease, coherently with previous literature on the stability of ripples in graphene53,54. 2.2 Electronic properties of rippled and stretched graphene The band structures of the optimized contracted/strained structures of the 4√3×4√3R30 and 3√3×3√3R30 cells were evaluated. The whole set of data is reported in the SI (section S2.2) and the gap values as a function of the strain/contraction are plotted in Fig 4 (a). The plot shows an irregular behavior, which is related to the high variety of geometry and symmetries of the system. In the stretched structures of the 3√3×3√3R30 cell (magenta dots), the gap is seen to increase quite regularly with the strain, except for the 10% stretching, which correspond to a highly defective one (with broken bonds). In the case of 4√3×4√3R30 cell (green dots), it is null up to ~5% and then it starts increasing. An inspection of the global band structures (see SI) indicates that empty bands are elongated towards the higher energies and the filled ones towards the lower energies, which produces an immediate gap opening for the small cell; conversely, this happens later for the larger one. It is quite natural to relate the gap opening with the “benzenization” observed in the stretched structures, since this produces a symmetry breaking. Therefore, we plotted the gap value vs the squared root variance of the bond lengths σ, which is a measure of the splitting of the bond lengths (Fig 4(b)). The plot shows a clear correlation for both the cells, indicating that the gap opening is related to the “benzenization” of the system. The analysis of the contracted structures is more difficult, due to the large variety of geometries and symmetries of the ripples. We first observe that in the 3√3×3√3R30 cell the gap


Page 8 of 25

Page 9 of 25

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


is null up to 7% compression. The inspection of the ripples and band structure in this region reveals very peculiar symmetries: first a regular alternation of wells and hills (up to 2%) and then wavy-like ripples in a one of the basis vector directions. The first symmetry produces only a slight vertical translation in the high energy empty and low energy occupied bands leaving the Dirac point unaltered, while the second produces the horizontal translation of the Dirac point, with, however, quite unaltered band slopes. At ~5% compression, n-type doping is also observed. The gap finally opens at around 8% when the symmetry of ripples changes, forming wells and hills. At larger compressions, n-type doping and rippling combine producing an alternating behavior (See SI S2.2 for details of the band structures).

Fig 4 (a) Band gap as a function of the strain (positive values of the abscissa) and compression (negative values) for the 3√3×3√3R30 (magenta) and 4√3×4√3R30 (green) cells. (b) gap vs the squared root variance of the bond lengths for extended structures (same color code). (c) Squared root variance (bottom) and squared root of the average squared dihedrals (upper) vs the gap for contracted structures (same color code).

As said, in the 4√3×4√3R30 cell the wavy ripples are never observed, and the hills/wells appear already at small contractions. Consequently the gap opening due to this symmetry is immediately visible although initially small, reasonably proportional to the curvature level. In order to quantify this gap-structure relation, we plotted the gap in the contraction region against the variance of the bond length σ and against the squared root of the average squared dihedral φ, related to a global measurement of the curvature. Since in the curved areas the bonds are also stretched, the two quantities are related. In fact, they show similar behavior as a function of the gap (Fig. 4(c)). Data for 4√3×4√3R30 cell (green) clearly show two regimes: at small curvature the gap increases slowly, at larger curvature, more rapidly. The transition occurs at about 5% compression, and is related to the already noted appearing of pyramidalization (see previous section). In the smaller cell (magenta dots), the transition appears to be directly from the gapless regime (with wavy ripples) to the sharp-pyramidalized one, which is also confirmed by the inspection of the structures (see the SI, S2.2.2.). In summary, also in compressed structures, one can distinguish different structural/electronic regimes, which we label “wavy”, “smooth”, and “sharp” regimes. The wavy regime is characterized by wavy smooth ripples and null gap, the smooth by the slowly increasing gap with no change in the bond hybridization, the sharp by the appearance of the pyramidalization towards the sp3 hybridization and rapidly increasing gaps. The appearance of the different regimes depends also on the cell size, which introduces a substantially different periodicity of the ripples in the system, namely ~1.2 nm for the small cell and ~1.7 nm for the larger one.



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 25

2.3 Hydrogen reactivity and stability on rippled graphene Selected rippled structures were then hydrogenated. These are: the 4√3×4√3R30 cell with 2% contraction, chosen as the model system better mimicking graphene on SiC, and the 3√3×3√3R30 cells with 7% contraction, chosen as a model for highly rippled graphene, because it is the un-defected model with the largest curvature level among those studied. Hydrogen atoms were progressively added on the most convex sites of a single side, and the system subsequently fully relaxed at fixed cell. In the case of 3√3×3√3R30 cell a few cases of double side hydrogenation were also considered. The full series of hydrogenated structures is reported in the SI, section S3. The hydrogen binding energies were evaluated as a function of the H coverage in two different ways: (i) the average binding energy per atom (filled dots in Fig 5) and the incremental binding energy per atom, namely the binding energy of the newly added hydrogen at each step, with respect to previously partially hydrogenated structure (empty dots in Fig 5; see also SI Section S3, for numerical values). From inspection of Fig. 5, it is immediately clear that hydrogen on the 7% compressed 3√3×3√3R30 cells (magenta dots) is more stable than on 2% compressed 4√3×4√3R30 cell (green dots) by 1-2 eV, depending on the coverage. This is most likely an effect of compression that increases the average local curvature, and is in agreement with our previous findings, indicating a linear increase of the binding energy with curvature for isolated atom hydrogen adhesion24. Indeed, this result could be considered the extension to the case in which other hydrogen atoms are present in the structure. In particular Fig. 5 shows that in the case of very high compression and corrugation, the chemisorbed hydrogen can also be more stable than molecular hydrogen, as previously suggested24.

Fig 5 Binding energies vs the hydrogen coverage in the 4√3×4√3R30 cell at 2% compression (green dots) and 3√3×3√3R30 cell at 7% compression (magenta dots) and at 2% compression (magenta squares, the two values correspond to two different decoration, before and after a H site hopping, described in Section 3.5). Filled symbols are the average binding energies per H atom at a given coverage (measured as % of hydrogenated sites). Empty dots are incremental binding energies per atom, i.e. measured with respect to the partially hydrogenated structure of the previous step, instead than to naked graphene. Both sets of energies are evaluated with respect to atomic hydrogen (left vertical axis) and to molecular hydrogen (right vertical axis), the two differing by an offset of ~2.23 eV/H atom. Dotted lines are guides for the eye.


Page 11 of 25

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 behavior of the binding energy as a function of the coverage indicates that, on average, bound hydrogen destabilizes as the coverage increases. This might seem in contrast with the previously reported cooperative effect30, which would rather produce stabilization in the presence of already bound hydrogen. This discrepancy is most likely due to the fact that we relax structures but not the cell. It is true that already present hydrogen increases the local convexity because of the pyramidalization, and consequently the bound hydrogen stability. But hydrogen binding has also the effect of increasing the C-C bond length, and consequently the average level of compression – if cell is kept fixed – and consequently, the average energy of the system. Thus, the result plotted in Fig 5 is the combined effect of hydrogen stabilization due to increased convexity and global destabilization due to increased effective compression. The incremental energy behavior (empty symbols) more difficult to interpret because it strongly depends on the decoration of the partially hydrogenate structure at the previous step (see the SI for representation of the structures). Some of the intermediate configuration display empty sites in the hydrogenated clusters, creating a graphane-like local structures, which turns out more stable and with heavenly distributed stress; in the case of 3√3×3√3R30 cells we also tried a double side hydrogenation. Thus overall the behavior of the incremental energy is roughly alternating around the average level, which is, at least in the case of the larger cell, not too far from the average one at large coverage. In agreement with previous literature, this analysis indicates, that local curvature, and specifically the pyramidalization, is the key feature in the process of hydrogenation. However, in addition, our results highlight the interplay between pyramidalization and clusterization of hydrogen on the graphene surface. 2.4 Electronic properties of hydrogenated rippled graphene The electronic bands of hydrogenated structures were evaluated (see SI, section S3). An interesting property to consider is the value of the band gap as a function of the hydrogen coverage. Unfortunately, this task is not as straightforward as it might seem: the band structure turns out very complex and hydrogen-decoration dependent, so that, at a given coverage, the band gap value displayed by different H decoration spans a wide range, as it can be seen from Fig 6. Though a gap increasing trend seems to emerge with increasing coverage, in most cases very small values of the gap persist even in the presence of large coverage. A closer inspection of the bands (see SI, Tables S7 and S8) reveals that small or zero gaps appear in two distinct cases. The first is at low levels of hydrogenation (e.g. in the case of dimers bound in ortho or para conformation, see red inset to Fig 6). In that case the electronic density of states near the Fermi level is in quantitative agreement with results previously obtained theoretically or by Scanning Tunneling Microscopy (STM) experiments27. High electronic density is visible on the hydrogenated sites, but residual graphene-like delocalized areas are present. In this case, the bands show a small gap and a quasi-graphene like linear dispersion near the Fermi level. The second case of vanishingly small gap is due to the appearance of nearly flat mid-gap states, as in the case of ~8% coverage of 4√3×4√3R30 cell (blue insets to Fig 6). The analysis of the electronic density of the states around the Fermi level reveals that these are very localized on the hydrogenated sites, with very little density around the hydrogenated islands. Conversely, the deeper states at -0.3 – -1.5 eV include also some electronic density in the graphene part, but they are separated by a gap of at least 0.5-0.6eV from the corresponding in the conduction band. Given the presence of those localized mid gap state in several of our structures, in those cases we also evaluated the band gaps excluding them (in these cases the bare gap is reported as small dots in Fig 6). In fact, these are not likely to contribute to the STM measurement, because localized and therefore non-conductive. Non-dispersive states are always



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 12 of 25

present, but in some case they are mixed with dispersive ones. This is the case, for instance, of 4√3×4√3R30 cell with 10 H atoms (orange inset), in which one has a real gap of ~0.4 eV formed by localized states mixed up with the dispersive ones. Coherently the corresponding electronic density (only filled bands) shows both localization and residual density in the conducting graphene-like areas.

Fig 6 Electronic band gaps as a function of coverage for 4√3×4√3R30 cell at 2% compression (green large dots) and 3√3×3√3R30 cell at 7% compression (magenta large dots). If localized mid-gap states are present, the gap is evaluated both including (small dots) and excluding (large dots) them. The black dots are experimental data from ref 27. The eye-guiding line in black is a polynomial fit Egap[eV]=3.8 (cov/100)0.6. For 100% coverage the line passes through Egap=3.8 eV, which is approximately the value of the graphane gap. Insets include the band gap plot around the Fermi level (red lines = empty states, blue lines = filled states) and the STM-like representation of the electron density of states in selected energy intervals for selected cases. Namely: enclosed in red, the 4√3×4√3R30 hydrogenated with the para-dimer, DOS integrated ±0.5eV with respect to Fermi level; enclosed in blue: the 4√3×4√3R30 hydrogenated with 8 atoms, DOS integrated ±0.6eV around the Fermi level (upper part in orange) and integrated between Fermi level and -1.6eV (lower part, in red); enclosed in orange, the 4√3×4√3R30 hydrogenated with 10 atoms, DOS integrated between Fermi level and -1.6eV. The arrows connect the insets to the corresponding numerical data for gap/coverage.

Figure 6 also reports (in black) experimental data taken from the experimental measurement of the gap at two different hydrogen coverages27. As it can be seen, they are in agreement with the trend of the gap evaluated excluding the non-dispersive midgap states, when present, confirming our analysis. Based on these data, the increase of band gap as a function of coverage can be estimated ~0.1eV per 1% site occupation, in the intermediate region, although it is better represented by a polynomial with exponent 0.6 (see caption of Fig 6), which also accounts for the saturation towards the 100% coverage (graphane has a theoretically estimated band gap of ~3.7eV3). The large fluctuations around the fitting line, as well as the emergence of the localized mid-gap states, might be explained by the fact that in most cases here we are considering only decoration with specific symmetries, while the real hydrogenation process is likely to have a higher degree of disorder in the hydrogenated sited distribution, being performed


Page 13 of 25

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


at high temperatures. Nevertheless, a correlation between the gap value and the coverage seems to emerge, which could be used to evaluate the chemisorbed H amount by STS measurements. 2.5 Dynamical stability of hydrogenated rippled graphene. We evaluated the dynamical stability of chemisorbed hydrogen at room temperature in one of the most unfavorable cases, namely the 3√3×3√3R30 cell with highest coverage. To this end we performed a CP molecular dynamics simulation, where, after a gradual temperature increase (~1.5 ps long, data not shown), the system is coupled to the thermostat at 300K (see Fig 7 (a)). A hopping event is observed in 0.5 ps on sites highlighted with the blue oval in Fig7 (b) and (c), as shown by the abrupt change in the involved C-H bond length (red and green lines in Fig 7). This event produces a stabilization of about ~1eV in the system which can be seen as corresponding decrease in the EKS energy (dotted line in the central plot) and triggers a macroscopic oscillation in the temperature and energy (see central plot), whose period is 1.41ps, meaning frequency in the THz range. Considering that the characteristic wavelength of this system is ~1-2nm, and the dispersion relations of the graphene phonons, these frequencies might correspond to ZA phonons, i.e. transverse out of plane displacements 55 . In fact, the analysis of the oscillations of the z coordinate of selected atoms near the hopping site (highlighted in yellow in the Fig 7(c)) show a battement at that frequency superimposed to faster oscillations (dotted lines in the top graph). We remark that ZA phonons are precisely the modes that describe the deformation corresponding to rippling in these systems. The macroscopic energy oscillations have certainly a large component due to the exchange with the thermostat, and are not entirely due to the sheet deformation. However, the fact that they are triggered by a hopping event and that they are resonant with ZA phonons is an indication of the strong interplay between rippling, hydrogenation/hopping events and mechanical properties. This relationship could be a possible way for detecting single hydrogen hopping by mean of vibrational analysis.

Fig 7. CP molecular dynamics simulation of the hydrogenated system. (a) Plot of different observables as a function of time. From bottom to top: bond distance corresponding to the two sites involved in the hopping of the proton (highlighted with the blue oval in (b) and (c)); Temperature (solid line scale on the left) and KS energy (doted line, scale on the right) during simulation; z displacement of a selected C atom (highlighted in yellow in (c); the dotted lines in the top plot highlights the battements.



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 25

3. Conclusions In this work, we have reported the structural and electronic properties of corrugated graphene, pristine and hydrogenated. The corrugation was created with the specific intention of mimicking the natural corrugation of graphene grown on SiC. For this reason the 4√3×4√3R30 supercell was used, whose periodicity and symmetry correspond to the pseudo-periodicity of epitaxial monolayer graphene on SiC. This allows to re-create the natural corrugation by slightly compressing the cell bi-axially. (We also verified that the same strategy did not work on the cell with the exact symmetry of graphene on SiC). The 3√3×3√3R30 cell was also considered, to recreate a similar condition but with smaller wavelength and larger local curvature, in order to enhance the effects. Stability, structural and electronic properties of these two model systems subject to a wide range of levels of strain and compression were analyzed. This systematic study lead to a number of side results, in addition to the main ones, regarding the interplay between structural and electronic properties. A “benzenization” transition is observed for large stretching, which changes substantially the conductive properties of graphene. In the compression region, a variety of different symmetries and geometries of ripples are observed, depending on the strain and wavelength of compression, corresponding, in turn, to different electronic properties. These however, can be rationalized relating them to structural order parameters, namely the bond length variance or the average of the pyramidalization angle, both related to the sp2→sp3 transition. The conclusion is that, while at low level of contraction graphene preserves its conductive properties, at high level of contractions pyramidalization appears inducing a transition towards a corrugation induced insulating state. Selected structures, namely those more similar to the naturally corrugated graphene grown on SiC, were hydrogenated following a protocol which mimics the exposure to atomic H: atoms were added progressively to the more reactive sites, namely those more pyramidalized. The binding energy of hydrogen on the surface was seen to decrease to a stable value up to almost half coverage of graphene surface, which is reasonably the largest value reachable with single-side hydrogenation. The system is seen to accept spontaneously atomic hydrogen up to this level. At high compression, the bound hydrogen turns out to be stable also with respect to molecular hydrogen. This is in agreement with our previous results obtained for isolated H atoms added on the corrugated graphene surface, extending them to the case of larger coverage. In addition, the analysis of the binding energy as a function of the coverage indicates a relation between the reactivity of pyramidal sites, the increase of pyramidalization induced by hydrogenation, and the consequent increased mechanical stress. As expected, the electronic band gap increases with the hydrogen coverage. The increase was quantified with a simple polynomial fit, which also superimposes to experimental data. Therefore, this quantitative rule could be used to indirectly evaluate the H coverage by STS measurements. However, this trend emerges only when localized mid-gap states sporadically appearing in calculations are excluded from the gap calculation. We interpreted the emergence of those states, elusive with respect to direct STM measurements because non conductive, as due to specific hydrogen decoration with specific symmetries, whose statistical weight is low in real hydrogenation experiments at room temperature27. Finally, we evaluated the dynamical stability of chemisorbed hydrogen at room temperature, choosing a system with particularly unfavorable hydrogenation, by means of CP molecular dynamics simulations. We observe, a spontaneous hopping event occurring in the ps time scale, by which the potential energy of the system is decreased of ~1eV. As a consequence, energy oscillations in the THz range are excited, with a component of out-of-plane related to ZA


Page 15 of 25

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


flexural phonons (traveling ripples of nm wavelength). This unexpected result is however coherent with our previous findings [24], because it shows the interplay between the curvature dynamical change related to ZA phonons and the hydrogen destabilization. In addition, this opens the possibility of measuring hydrogen hopping on graphene surface by means of vibrational measurements in the THz range. Overall, therefore, this work clarifies several aspects of the interplay between corrugation, electronic properties, vibro-mechanical properties and hydrogenation. These were obtained as “side results” as an effect of the systematic study of pristine and hydrogenated graphene with different levels of average local curvature. Besides this, we optimized a modeling platform for a number of possible developments. We verified that 4√3×4√3R30 supercell can reproduce the main characteristics of the naturally rippled graphene on SiC. Thus this model system can be used for further studies, and with the aim to compare with experimental measurement of graphene on SiC. The quantitative agreement between measured and calculated trend of the band gap vs hydrogen coverage is a validation test for this model system. In addition, we have also reported other “measurable” relations, such as the variation of binding energy vs coverage, or the possibility of detecting hydrogen hopping by means of measurements vibrations in the THz range. We hope that this study will trigger experimental investigation in the same directions. Acknowledgments We thank Dr Chandramathy Surendrean Praveen for precious technical support for QE code, and Dr Walter Rocchia and Dr Stefan Heun for useful discussions. We gratefully acknowledge ﬁnancial support by the EU, 7th FP, Graphene Flagship (contract no. NECT-ICT-604391), the CINECA award “ISCRA C” IsC10_HBG, 2013 and PRACE “Tier0” award Pra07_1544 for resources on FERMI (IBM Blue Gene/Q@CINECA, Bologna Italy), iit “Compunet platform” for providing computational resources and CINECA staff for technical support. Supporting Information Available Supporting information include: (i) convergence checks, (ii) energy data in tabular form, (iii) graphical representation of selected structures and their electronic bands, for pristine and hydrogenated graphene. This information is available free of charge via the Internet at http://pubs.acs.org References 1

Geim, A.K.; Novoselov K.S. The Rise Of Graphene. Nature Mat 2009, 6, 183-191 Elias D.C.; Nair R. R.; Mohiuddin T. M. G.; Morozov S. V.; Blake P.; Halsall M. P.; Ferrari A. C.; Boukhvalov D. W.; Katsnelson M. I.; Geim A. K.; Novoselov K. S. Control of Graphene's Properties By Reversible Hydrogenation: Evidence for Graphane. Science 2009, 323, 610-613 3 Sofo J O.; Chaudhari A S.; Barber G D. Graphane: A Two-Dimensional Hydrocarbon. Phys. Rev. B 2007, 75, 153401 4 Haberer D.; Vyalikh D. V.; Taioli S.; Dora B.; Farjam M.; Fink J.; Marchenko D.; Pichler T.; Ziegler K., Simonucci S., et al. Tunable Band Gap In Hydrogenated Quasi-Free-Standing Graphene. Nano Lett. 2010, 10, 3360–3366 5 Tozzini V.; Pellegrini V. Electronic Structure And Peierls Instability In Graphene Nanoribbons Sculpted In Graphane. Phys. Rev. B 2010, 81, 113404 6 Schmidt M.J.; Loss D. Edge States And Enhanced Spin-Orbit Interaction At Graphene/Graphane Interfaces. Phys. Rev. B 2010, 81, 165439 7 Singh A. K.; Yakobson B. I. Electronics And Magnetism Of Patterned Graphene Nanoroads. Nano Lett. 2009, 9, 1540–1543 2



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

Page 16 of 25

Hernández-Nieves A. D.; Partoens B.; Peeters F. M. Electronic And Magnetic Properties Of Superlattices Of Graphene/Graphane Nanoribbons With Different Edge Hydrogenation. Phys. Rev. B 2010, 82, 165412 9 Huang L. F.; Zheng X. H.; Zhang G. R.; Li L. L.; Zeng Z. Understanding the Band Gap, Magnetism, and Kinetics of Graphene Nanostripes in Graphane. J. Phys. Chem. C 2011, 115, 21088–21097 10 Balog R.; Jørgensen B.; L Nilsson; Andersen M.; Rienks E.; Bianchi M.; Fanetti M.; Lægsgaard E.; Baraldi A.; Lizzit S.; et al. Bandgap Opening In Graphene Induced By Patterned Hydrogen Adsorption. Nature Mat. 2010, 9, 315–319 11 Tozzini V.; Pellegrini V.; Prospects For Hydrogen Storage In Graphene. Phys. Chem. Chem. Phys., 2013, 15, 8089 12 Patchkovskiim S.; Tse J. S.; Yurchenko S. N.; Zhechkov L.; Heine T.; Seifert G. Graphene Nanostructures As Tunable Storage Media For Molecular Hydrogen. Proc. Natl. Acad. Sci. USA, 2005, 102, 10439 –10444 13 Lee S. M.; Lee Y. H. Hydrogen Storage In Single Walled Carbon Nanotubes. Appl. Phys. Lett., 2000, 76, 2877 14 Okamoto Y.; Miyamoto Y. Ab Initio Investigation Of Physisorption Of Molecular Hydrogen On Planar And Curved Garaphene. J. Phys. Chem. B 2001, 105, 3470-3474 15 Cambria I.; Lopez M. J.; Alonso J. A. The Optimum Average Nanopore Size For Hydrogen Storage In Carbon Nanoporous Materials. Carbon, 2007, 45, 2649–2658 16 Dillon A. C.; Jones K. M.; Bekkedahl T. A.; Kiang C. H.; Bethune D. S.; Heben M. J. Storage Of Hydrogen In Single Walled Carbon Nanotubes Nature, 1997, 386, 377 17 Panella B.; Hirscher M.; Roth S. Hydrogen adsorption in different carbon nanostructures. Carbon 2005, 43, 2209 18 Miura Y.; Dino W; Nakanishi H.; First Principles Studies For The Dissociative Adsorption Of H2 On Graphene. J. Appl. Phys., 2003, 93, 3395 19 Parambhath V. B.; Nagar R.; Sethupathi K.; Ramaprabhu S. Investigation Of Spillover Mechanism In Palladium Decorated Hydrogen Exfoliated Functionalized Graphene. J. Phys. Chem. C 2011, 115, 15679–15685 20 Ao Z. M.; Peeters F. M. Electric Field: A Catalyst For Hydrogenation Of Graphene. Appl Phys Lett 2010, 96, 253106 21 Ao Z. M.; Peeters F. M. Electric Field Activated Hydrogen Dissociative Adsorption To Nitrogen-doped Graphene. J. Phys. Chem. C 2010, 114, 14503–14509 22 Boukhvalov D. W.; Katsnelson M. I. Enhancement of Chemical Activity In Corrugated Graphene J. Phys. Chem. C 2009, 113, 14176–14178 23 Ruffieux P.; Gröning O.; Bielmann M.; Mauron P.; Schlapbach. L.; Gröning P. Hydrogen Adsorption On sp2Bonded Carbon: Influence Of The Local Curvature. Phys Rev B 2002, 66, 245416 24 Tozzini V.; Pellegrini V. Reversible Hydrogen Storage By Controlled Buckling Of Graphene Layers. J. Phys. Chem. C 2011, 115, 25523-25528 25 Wang Z. F.; Zhang Y.; Liu F. Formation Of Hydrogenated Graphene Nanoripples By Strain Engineering And Directed Surface Self-Assembly. Phys. Rev. B 2011, 83, 041403(R) 26 Chernozatonskii L. A.; Sorokin P.B. Nanoengineering Structures On Graphene With Adsorbed Hydrogen “Lines”. J. Phys. Chem. C 2010, 114, 3225–3229 27 Goler S.; Coletti C.; Tozzini V.; Piazza V.; Mashoff T.; Beltram F.; Pellegrini V.; Heun S. The Influence Of Graphene Curvature On Hydrogen Adsorption: A Study For Future Hydrogen Storage Devices. J. Phys. Chem. C 2013, 117, 11506−11513 28 Casolo S.; Løvvik O. M.; Martinazzo R.; Tantardini G. F. Understanding Adsorption Of Hydrogen Atoms On Graphene. J. Chem. Phys. 2009, 130, 054704 29 Šljivančanin Ž.; Rauls E.; Hornekær L.; Xu W.; Besenbacher F.; Hammerb B. Extended Atomic Hydrogen Dimer Configurations On The Graphite (0001) Surface. J. Chem. Phys. 2009, 131, 084706 30 Boukhvalov D. W.; Katsnelson M. I.; Lichtenstein A. I. Hydrogen On Graphene: Electronic Structure, Total Energy, Structural Distortions And Magnetism From First-Principles Calculations. Phys. Rev. B 2008, 77, 035427 31 Ferro Y.; Teillet-Billy D.; Rougeau N.; Sidis V.; Morisset S.; Allouche A. Stability And Magnetism Of Hydrogen Dimers On Graphene. Phys. Rev. B 2008, 78, 085417 32 Sha X.; Knippenberg M. T.; Cooper A C.; Pez G.P.; Cheng H. Dynamics of Hydrogen Spillover on Carbon-Based Materials. J. Phys. Chem. C 2008, 112, 17465–17470 33 Bilic A.; Gale J. D. Chemisorption of Molecular Hydrogen on Carbon Nanotubes: A Route to Effective Hydrogen Storage? J. Phys. Chem. C 2008, 112, 12568–12575 34 Bonini N.; Lazzeri M.; Marzari N.; Mauri F. Phonon Anharmonicities In Graphite And Graphene. Phys. Rev. Lett. 2008, 100, 176802 35 Mattausch A.; Pankratov O. Ab Initio Study Of Graphene On SiC. Phys. Rev. Lett. 2007, 99, 076802


Page 17 of 25

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


36

Lee B.; Han S.; Kim Y-S. First-Principles Study Of Preferential Sites Of Hydrogen Incorporated In Epitaxial Graphene On 6h-Sic(0001). Phys. Rev. B 2010, 81, 075432 37 Huang B.; Xiang H.J.; Wei S-H. Controlling Doping In Graphene Through A Sic Substrate: A First-Principles Study. Phys. Rev. B 2011, 83, 161405(R) 38 Deretzis I.; La Magna A. Role Of Covalent And Metallic Intercalation On The Electronic Properties Of Epitaxial Graphene On Sic(0001) Phys. Rev. B 2011, 84, 235426 39 Jayasekera T.; Kong B.D.; Kim K.W.; Buongiorno Nardelli M. Band Engineering And Magnetic Doping Of Epitaxial Grpahene On Sic (0001). Phis. Rev. Lett. 2010, 104, 146801 40 Goler S.; Coletti C.; Piazza V.; Pingue P.; Colangelo F.; Pellegrini V.; Emtsev K.V.; Forti S.; Starke U.; Beltram F.; Heun S. Revealing The Atomic Structure Of The Buffer Layer Between Sic(0001) And Epitaxial Graphene. Carbon 2013, 51, 249–254 41 The dihedral ϕ is related to the out-of-plane protrusion distance (i.e. the distance of a C site from the plane formed by its three first neighbors) by means of the geometric relationship given in24. 42 Vanderbilt D. Soft Self-Consistent Pseudopotential In A Generalized Eigenvalue Formalism, Phys. Rev. B 1990, 41, 7892-7895 43 Grimme S. J. Semiempirical GGA-Type Density Functional Constructed With A Long-Range Dispersion Correction. Comp. Chem. 2006, 27, 1787-1799 44 Perdew J.P.; Burke K.; Ernzerhof M. Generalized Gradient Approximation Made Simple, Phys. Rev. Lett. 1996, 77, 3865-3868 45 Ceperley D.M.; Alder B.J. Ground State of the Electron Gas by a Stochastic Method. Phys. Rev. Lett. 1980, 45, 566; Vosko S.H.; Wilk L.; Nusair M. Accurate Spin-Dependent Electron Liquid Correlation Energies For Local Spin Density Calculations: A Critical Analysis. Can. J. Phys. 1980, 58, 1200-1211 46 Monkhorst H. J.; Pack J. D. Special Points For Brillouin-Zone Integrations. Phys. Rev. B 1976, 13, 5188. MP grids were (10×10×1 for the 3√3×3√3R30 cell, 4√3×4√3R30 cells, 40×40×1 for the unit cell the; band structure was sampled along the main symmetry path with 100 points. 47 Head J. D.; Zerner M.C. A Broyden-Fletcher-Goldfarb-Shanno optimization procedure for molecular geometries Chem phys. lett., 1985,122, 264–270 48 Car R.; Parrinello M. A Unified Approach For Molecular Dynamics And Density Functional Theory. Phys. Rev. Lett. 1985, 55, 2471 49 µ=250a.u. and 5 Ry mass cutoff; Tassone F.; Mauri F.; Car R. Acceleration schemes for ab initio moleculardynamics simulations and electronic-structure calculations. Phys. Rev. B 1994, 50, 10561 50 Nosé, S.; Klein, M. L. Constant Pressure Molecular Dynamics For Molecular Systems. Mol. Phys 1983, 50, 10551076 51 Giannozzi P.; Baroni S.; Bonini N.; Calandra M.; Car R.; Cavazzoni C.; Ceresoli D.; Chiarotti G.L.; Cococciani M.; Debo I.; et al. QUANTUM ESPRESSO: A Modular And Open-Source Software Project For Quantum Simulations Of Materials. J. Phys.:Condens. Matter 2009 21 395502 URL http://www.quantum-espresso.org 52 Kittel C. Introduction to Solid State Physics. (Wiley, 1996). 53 Wang Z. F.; Zhang Y.; Liu F. Formation Of Hydrogenated Graphene Nanoripples By Strain Engineering And Directed Surface Self-Assembly. Phys. Rev. B 2011, 83, 041403(R) 54 de Andres P.L.; Guinea F.; Katsnelson M. I. Density Functional Theory Analysis Of Flexural Modes, Elastic Constants, And Corrugations In Strained Graphene. Phys. Rev. B 2012, 86, 245409 55 Lindsay L.; Broido D.A.; Optimized Tersoff And Brenner Empirical Potential Parameters For Lattice Dynamics And Thermal Transport In Carbon Nanotubes And Graphene. Phys. Rev. B 2010, 81, 205441



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

(a) Supercells considered in this study (in colors), superimposed to the graphene lattice (in black). In (b) The cell periodicity (with respect to graphene and SiC lattice) and the number of atoms included (corresponding colors to (a)). (c) The two supercells with the periodicity of the Moiré lattice (exact in red, and approximate, in green). In the background a representation of the corrugation is reported: light grey areas are protruding areas. (d) The Brillouin zones of the main cells (corresponding colors to (a)). The main symmetry points of the unit and of the 6x6 cell are reported in blue and green, respectively. 508x381mm (72 x 72 DPI)


Page 18 of 25

Page 19 of 25

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


(a) Coesive energies per C atom as a function of the strain, measured by the “effective” C-C bondlength (i.e. the bond C-C bond length one would have if the graphene remained flat and regular). Magenta: 3√3×3√3R30 cell, green: 4√3×4√3R30, red=13x13. Filled dots=PBE+vdW, empty dots=LDA+vdW, empty squares=LDA. (b) and (c) top and perspective view of the 2% contracted 4√3×4√3R30° supercell (3x3 repeated cells). In (b) the unit cell boundaries (in green) and the 13x13 cell boundaries (red) are superimposed. Structures colored according to the z coordinate: blue = protruding, red= intruding. 508x381mm (72 x 72 DPI)



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

Structural properties of stretched and compressed graphene of the 4√3×4√3R30 (a,c) and 3√3×3√3R30 (b,d) cell. The histograms report the C=C distance distributions (a,b), and dihedral angle distribution (c,d) at the different cell compression/strain, from -10% to +10% of the relaxed cell size (for dihedrals, the distributions of strained structures are not reported, since they are flat). The insets report selected strained (one unit supercell) and contracted (3×3 supercells) structures (the level of strain/compression is indicated). In the strained structures the bonds are colored according to their length (shorter in red, longer in blue). In the contracted structures the colors represent the displacement in z direction (out of plane): blue = up, red= down. Dots and error bars in color (green = 4√3×4√3R30, magenta = 3√3×3√3R30 ) report the mean of the distribution and its variance. 508x381mm (72 x 72 DPI)


Page 20 of 25

Page 21 of 25

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


(a) Band gap as a function of the strain (positive values of the abscissa) and compression (negative values) for the 3√3×3√3R30 (magenta) and 4√3×4√3R30 (green) cells. (b) gap vs the squared root variance of the bond lengths for extended structures (same color code). (c) Squared root variance (bottom) and squared root of the average squared dihedrals (upper) vs the gap for contracted structures (same color code). 502x160mm (72 x 72 DPI)



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

Binding energies vs the hydrogen coverage in the 4√3×4√3R30 cell at 2% compression (green dots) and 3√3×3√3R30 cell at 7% compression (magenta dots) and at 2% compression (magenta squares, the two values correspond to two different decoration, before and after a H site hopping, described in Section 3.5). Filled symbols are the average binding energies per H atom at a given coverage (measured as % of hydrogenated sites). Empty dots are incremental binding energies per atom, i.e. measured with respect to the partially hydrogenated structure of the previous step, instead than to naked graphene. Both sets of energies are evaluated with respect to atomic hydrogen (left vertical axis) and to molecular hydrogen (right vertical axis), the two differing by an offset of ~2.23 eV/H atom. Dotted lines are guides for the eye. 189x158mm (150 x 150 DPI)


Page 22 of 25

Page 23 of 25

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


Electronic band gaps as a function of coverage for 4√3×4√3R30 cell at 2% compression (green large dots) and 3√3×3√3R30 cell at 7% compression (magenta large dots). If localized mid-gap states are present, the gap is evaluated both including (small dots) and excluding (large dots) them. The black dots are experimental data from ref 27. The eye-guiding line in black is a polynomial fit Egap[eV]=3.8 (cov/100)0.6. For 100% coverage the line passes through Egap=3.8 eV, which is approximately the value of the graphane gap. Insets include the band gap plot around the Fermi level (red lines = empty states, blue lines = filled states) and the STM-like representation of the electron density of states in selected energy intervals for selected cases. Namely: enclosed in red, the 4√3×4√3R30 hydrogenated with the para-dimer, DOS integrated ±0.5eV with respect to Fermi level; enclosed in blue: the 4√3×4√3R30 hydrogenated with 8 atoms, DOS integrated ±0.6eV around the Fermi level (upper part in orange) and integrated between Fermi level and -1.6eV (lower part, in red); enclosed in orange, the 4√3×4√3R30 hydrogenated with 10 atoms, DOS integrated between Fermi level and -1.6eV. The arrows connect the insets to the corresponding numerical data for gap/coverage. 396x365mm (72 x 72 DPI)



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

CP molecular dynamics simulation of the hydrogenated system. (a) Plot of different observables as a function of time. From bottom to top: bond distance corresponding to the two sites involved in the hopping of the proton (highlighted with the blue oval in (b) and (c)); Temperature (solid line scale on the left) and KS energy (doted line, scale on the right) during simulation; z displacement of a selected C atom (highlighted in yellow in (c); the dotted lines in the top plot highlights the battements. 481x304mm (72 x 72 DPI)


Page 24 of 25

Page 25 of 25

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


Graphical representation of rippled pristine and hydrogenated graphene with their band structures 239x152mm (72 x 72 DPI)


Nano-Scale Corrugations in Graphene: A Density Functional Theory

Recommend Documents