Atomistic Simulations of Graphene Growth: From Kinetics to

Mar 1, 2018 - DOI: 10.1021/acs.accounts.7b00592 ... By checking the simulation results, we find that such nonlinearity is caused by lattice mismatch, ...
3 downloads 5 Views 5MB Size
Article Cite This: Acc. Chem. Res. XXXX, XXX, XXX−XXX

pubs.acs.org/accounts

Atomistic Simulations of Graphene Growth: From Kinetics to Mechanism Zongyang Qiu, Pai Li, Zhenyu Li,* and Jinlong Yang Hefei National Laboratory for Physical Sciences at the Microscale, University of Science and Technology of China, Hefei, Anhui 230026, China CONSPECTUS: Epitaxial growth is a promising strategy to produce high-quality graphene samples. At the same time, this method has great flexibility for industrial scale-up. To optimize growth protocols, it is essential to understand the underlying growth mechanisms. This is, however, very challenging, as the growth process is complicated and involves many elementary steps. Experimentally, atomic-scale in situ characterization methods are generally not feasible at the high temperature of graphene growth. Therefore, kinetics is the main experimental information to study growth mechanisms. Theoretically, first-principles calculations routinely provide atomic structures and energetics but have a stringent limit on the accessible spatial and time scales. Such gap between experiment and theory can be bridged by atomistic simulations using first-principles atomic details as input and providing the overall growth kinetics, which can be directly compared with experiment, as output. Typically, system-specific approximations should be applied to make such simulations computationally feasible. By feeding kinetic Monte Carlo (kMC) simulations with first-principles parameters, we can directly simulate the graphene growth process and thus understand the growth mechanisms. Our simulations suggest that the carbon dimer is the dominant feeding species in the epitaxial growth of graphene on both Cu(111) and Cu(100) surfaces, which enables us to understand why the reaction is diffusion limited on Cu(111) but attachment limited on Cu(100). When hydrogen is explicitly considered in the simulation, the central role hydrogen plays in graphene growth is revealed, which solves the long-standing puzzle into why H2 should be fed in the chemical vapor deposition of graphene. The simulation results can be directly compared with the experimental kinetic data, if available. Our kMC simulations reproduce the experimentally observed quintic-like behavior of graphene growth on Ir(111). By checking the simulation results, we find that such nonlinearity is caused by lattice mismatch, and the induced growth front inhomogeneity can be universally used to predict growth behaviors in other heteroepitaxial systems. Notably, although experimental kinetics usually gives useful insight into atomic mechanisms, it can sometimes be misleading. Such pitfalls can be avoided via atomistic simulations, as demonstrated in our study of the graphene etching process. Growth protocols can be designed theoretically with computational kinetic and mechanistic information. By contrasting the different activation energies involved in an atom-exchange-based carbon penetration process for monolayer and bilayer graphene, we propose a three-step strategy to grow high-quality bilayer graphene. Based on first-principles parameters, a kinetic pathway toward the high-density, ordered N doping of epitaxial graphene on Cu(111) using a C5NCl5 precursor is also identified. These studies demonstrate that atomistic simulations can unambiguously produce or reproduce the kinetic information on graphene growth, which is pivotal to understanding the growth mechanism and designing better growth protocols. A similar strategy can be used in growth mechanism studies of other two-dimensional atomic crystals.

1. INTRODUCTION Graphene is an important material that has attracted intense research interest from both academia and industry.1,2 Among the available methods for its synthesis, epitaxial growth is promising for the realization of both high quality and large quantity of the obtained graphene samples.3,4 By using different metal substrates,3,5−7 carbon sources,8 and growth environments (temperature, pressure, etc.),9 different groups have developed very different protocols to grow graphene, which has led to diverse growth behaviors and sample morphologies.10,11 Guided optimization of the growth conditions requires a deep understanding of the growth mechanisms. Unfortunately, characterization of the atomic details of graphene growth in © XXXX American Chemical Society

situ and in operando is still a great challenge, mainly due to the high growth temperature. Kinetic information is thus the main experimental clue to conjecture growth mechanisms. First-principles, semiempirical, and empirical theoretical methods are routinely used in geometry optimizations and energy calculations to compare the relative stabilities of different species on a metal surface. Such results of atomistic resolution certainly provide useful insight for understanding the growth mechanism.12−18 However, the energetics, in principle, is an equilibrium state property and is difficult to be directly Received: November 27, 2017

A

DOI: 10.1021/acs.accounts.7b00592 Acc. Chem. Res. XXXX, XXX, XXX−XXX

Article

Accounts of Chemical Research related to the overall kinetics obtained from experiment. To fill this gap between theory and experiment, explicit simulations of graphene growth are desirable. Due to the complexity, multiscale simulation techniques, from molecular dynamics (MD)19 to kinetic Monte Carlo (kMC),20−22 rate equations,21,23 kinetic Wulff construction,24,25 phase field,26 and even computational fluid dynamics (CFD),27 have been used to study graphene growth. To gain enough mechanistic insight, maintaining atomistic resolution in the simulation is desirable. In this Account, we summarize our recent progress in the atomistic simulation of graphene growth. Based on our simulation results, we can now largely understand the mechanism of graphene growth on different substrates, even for the most complicated and most widely adopted case of chemical vapor deposition (CVD) on a Cu substrate. A direct comparison between the simulation results and experimental data provides an ultimate way to judge mechanistic hypotheses gained from experiment. With kinetic and mechanistic information from simulations, we can also theoretically design growth protocols to tune graphene growth. These results demonstrate that atomistic simulation plays an important role in understanding two-dimensional (2D) van der Waals (vdW) epitaxy.

Figure 1. Schematic diagram of the SOF-kMC model. Solid and hollow circles represent occupied and unoccupied lattice sites. Brown sites have already formed graphene. For simplicity, surface sites for carbon adsorption are represented by a honeycomb lattice. Adapted from ref 20. Copyright 2012 American Chemical Society.

(yellow dashed box); growth front, whose velocity determines the growth rate (red dashed box). If we focus on the last region and use a simulation box moving with the growth front, the socalled standing-on-the-growth-front (SOF) model permits us to use a very small simulation box but still obtain a statistically significant growth rate.20

2. SIMULATION TECHNIQUES The most straightforward atomistic simulation technique is MD. In a classical MD simulation, the atoms move following Newton’s law. To discretize Newton’s equation with a desirable accuracy, the proper time step is about 1 fs. Therefore, the time scale accessible for an MD simulation is typically much shorter than that of graphene growth, which makes MD more suitable for studying specific problems in graphene growth instead of obtaining an overall picture. High temperature and/or high oversaturation may be required in MD simulations to artificially speed up the growth. Alternatively, MD simulations can be used with enhanced sampling techniques, such as metadynamics28 and umbrella sampling,29 to obtain the potential of mean force and kinetic information on key reactions. To go beyond MD, a time scale coarse grain technique such as kMC can be adopted. In a typical kMC simulation, the whole system is represented by a lattice, where different sites can be occupied by different species. The occupation status of all sites determines the configuration of the system. At each step, all events that can change the current system configuration are listed in a table along with their rates, which are typically provided by first-principles calculations within the framework of transition state theory. A specific event is selected with a probability proportional to its rate, which changes the lattice configuration and advances the simulation. However, even with kMC, an explicit simulation of graphene growth may still be computationally unfeasible. Therefore, we developed some system-specific techniques to speed up kMC simulations of graphene growth.

2.2. Time Scale Separation

If there is a large rate gap in the event list, the possibility to select low-rate events, which will determine the overall growth kinetics, becomes extremely low. Most of the simulation time will be wasted on fast events and the system is trapped in local equilibrium states. In this situation, a time scale separation technique can be used.20 If slow events have not been selected for a long time, then the fast events can be temporarily removed from the event list. After a slow event is selected, the full event list is restored. 2.3. Mean Field Approximation for Specific Species

A species that diffuses quickly can reach its equilibrium distribution quickly for any configuration of other species. In an explicit kMC simulation, the diffusion of this species, if abundant, will have a very high probability to be selected and thus consumes a large amount of simulation time. To avoid this situation, we can apply a mean field approximation for this species,22 where its lattice distribution is not bookkept and only its total number is recorded. Fast diffusion events for this species are permanently removed from the event list, and the rates of all other events involving this species are adjusted according to its total number.

3. IDENTIFICATION OF DOMINANT KINETIC PATHWAYS Because it can explicitly simulate the graphene growth process, atomistic modeling provides an independent means to investigate growth mechanisms. For example, dominant kinetic pathways of growth can be identified without any input from experiment.

2.1. Floating Simulation Box

To obtain a statistically significant growth rate, the graphene edge as the growth front must be propagated forward over a long distance, which requires a very large simulation cell. To solve this problem, we divided the whole system into four regions (Figure 1): graphene area, where graphene is already formed (black dashed box); far field, in which a stationary distribution of different surface species can be assumed (blue dashed box); diffusion layer, through which a carbon flux passes

3.1. Hydrogen-Free Growth on a Cu Substrate

Graphene growth on Cu(111) has been found experimentally to be diffusion limited.30 Because the diffusion of carbon atoms on Cu(111) is almost barrierless, this is a surprising observation if graphene growth occurs through the attachment of carbon atoms to graphene edges. One possibility is that carbon atom attachment is not the dominant kinetic pathway of graphene B

DOI: 10.1021/acs.accounts.7b00592 Acc. Chem. Res. XXXX, XXX, XXX−XXX

Article

Accounts of Chemical Research growth. We performed an atomistic simulation study to find the actual feeding species for graphene growth.21 First, we calculated the energy barriers of all possible events on the surface involving carbon species up to trimers from first principles. As expected, monomer attachment to the graphene edges has a high barrier (1.17 eV), which blocks the pathway of graphene growth via carbon monomer attachment. The dimer is energetically much more favorable than the monomer on Cu(111),12 and the formation of a dimer from two monomers is kinetically feasible with a low barrier (0.3 eV). Further increasing the cluster size to form a trimer has a large barrier (1.15 eV). Therefore, many dimers will be formed on the surface. By solving the rate equations, we found that the steadystate concentration of dimer is indeed at least 2 orders of magnitude higher than that of other species. Although the barrier to dimer diffusion (0.49 eV) is still slightly lower than that to attachment at the graphene edge (0.58 eV), the low concentrations of carbon species and graphene islands on the surface may make graphene growth diffusion limited. To determine whether dimer attachment is the dominant growth pathway and whether this pathway will lead to diffusionlimited growth behavior, kMC simulations were performed. First, the dominance of the dimer predicted by the rateequation approach within the mean field approximation was confirmed. More importantly, a carbon concentration gradient around the graphene islands, which is characteristic of diffusionlimited growth, was observed (Figure 2a). We also used kMC

limited. This conclusion is consistent with experimental observations.31 Identification of the dimer as the dominant feeding species can also be used to explain the experimentally observed graphene island morphology transformation at different temperatures.30 Generally, the formation of dendritic or compact islands is determined by the competition between carbon attachment and its diffusion along the edge. Based on such a model, we can estimate the transition temperatures, which reasonably agree with experiment.21 3.2. Hydrogen-Involved Growth on a Cu Substrate

Hydrogen is highly involved in the industry-compatible CVD growth of graphene, where hydrocarbon molecules, such as methane, are used as the carbon source. Our density functional theory (DFT) calculations indicate that the dehydrogenation of methane and other intermediate species is coupled with other elementary steps of graphene growth.32 Therefore, all CxHy species should be considered in H-involved CVD growth on Cu substrates. What makes things more complicated is that, in most CVD experiments, H2 gas is also provided. Because the overall reaction of graphene CVD growth with methane is CH4 → C + 2H2, it is counterintuitive as to why H2, which is a byproduct of the reaction, should be supplied. Both hydrogen-promoted growth33,34 and hydrogen-induced graphene etching35,36 have been observed experimentally, but the underlying mechanism is still under debate.37 Hydrogen can affect the surface species concentration, graphene edge configuration, and carbon attachment/detachment dynamics. Relevant elementary processes form a complicated kinetic network. To understand the mechanisms behind this network, the dominant kinetic pathways should be unambiguously identified.22 DFT calculations indicate that H2 can easily dissociate into H adatoms on Cu(111). In contrast, the dehydrogenation of CH4 is endothermic with a high activation barrier, making the surface sticking coefficient of CH4 extremely low. As a result, the surface concentration of H atoms is typically much higher than those of carbon-containing species, even if the partial pressure of H2 is lower than that of CH4. Our kMC simulations indicate that the surface H atom concentration is almost exclusively determined by the H2 partial pressure. At low H2 pressures, the most abundant carbon-containing species is C2, and CH becomes dominant at high H2 pressure. The graphene edge configuration also depends on the H2 partial pressure. Between metal passivation at low H2 pressures16 and hydrogen saturation at high H2 pressures,38 homogeneous partially hydrogen-saturated configurations exist due to entropic effects. The attachment of different carbon-containing species to the graphene edge is generally easier at metal-passivated sites than at hydrogen-saturated sites. The same trend also stands for carbon detachment. An overall picture of growth can be obtained via kMC simulations. As shown in Figure 3, under a relatively low H2 pressure, CH and C2 are the most commonly attached species. Although the raw attachment frequencies of these two species are comparable, the contribution to graphene growth from C2 is much higher. The reason is that the C2 detachment barrier at metal-passivated edge sites (2.19 eV) is much higher than the CH detachment barrier (1.08 eV). Because the concentration of carbon-containing species in graphene CVD growth is extremely low (