PCCP View Article Online
Published on 30 July 2015. Downloaded by Stockholms Universitet on 26/08/2015 15:49:10.
PAPER
Cite this: Phys. Chem. Chem. Phys., 2015, 17, 22832
View Journal | View Issue
Gas-phase dynamics in graphene growth by chemical vapour deposition† Gan Li,a Sheng-Hong Huang*b and Zhenyu Li*a Chemical vapour deposition on a Cu substrate is becoming a very important approach to obtain high quality graphene samples. Previous studies of graphene growth on Cu mainly focus on surface processes. However, recent experiments suggest that gas-phase dynamics also plays an important role in graphene growth. In this article, gas-phase processes are systematically studied using computational
Received 20th April 2015, Accepted 30th July 2015
fluid dynamics. Our simulations clearly show that graphene growth is limited by mass transport under
DOI: 10.1039/c5cp02301g
rate at different positions in the tube furnace and the concentration of different gas phase species are
ambient pressures while it is limited by surface reactions under low pressures. The carbon deposition calculated. Our results confirm that the previously realized graphene thickness control by changing the
www.rsc.org/pccp
position of the Cu foil is a result of gas-phase methane decomposition reactions.
Introduction Graphene is a two-dimensional material composed of carbon atoms with sp2 hybridization. It has a unique electronic structure and excellent physical properties1 which lead to a wide range of applications.2,3 Therefore, the production of highquality graphene films at a large scale is very desirable. Several approaches have been used to prepare graphene samples, including micro-mechanical exfoliation of graphite, epitaxial growth on different substrates, and reduction of graphene oxide.4 Among them, graphene growth on metal surfaces via chemical vapour deposition (CVD) is a promising way for largescale, inexpensive, high-quality graphene production.5,6 Cu is currently the most widely used substrate material due to its low carbon solubility and the resulting flexibility for graphene layer number control.7 Mechanisms of graphene growth on a Cu substrate have been intensively studied theoretically.8 The structure of graphene edges is determined by surface passivation9 and also by partial pressure of hydrogen.10 Active species for graphene growth is suggested to be of more than one carbon atom11 and may contain H atoms.12,13 The complexity of both the graphene edge structure and active species composition can lead to complex growth behaviour.14 These previous studies mainly focused on processes that occurred on the Cu surface. a
Hefei National Laboratory for Physical Science at Microscale, University of Science and Technology of China, Hefei, Anhui 230026, China. E-mail: [email protected] b Department of Modern Mechanics, University of Science and Technology of China, Hefei, Anhui 230026, China. E-mail: [email protected] † Electronic supplementary information (ESI) available: More computational details and results at intermediate pressure (2666 Pa). See DOI: 10.1039/ c5cp02301g
22832 | Phys. Chem. Chem. Phys., 2015, 17, 22832--22836
The gas phase effects have been considered only from a thermodynamic point of view.12 To the best of our knowledge, a systematic study on gas phase dynamics in graphene growth has not been reported yet. On the experimental side, however, several recent studies have suggested that gas phase dynamics is important in graphene growth. Important evidence comes from strong sensitivity of the growth behaviour to pressure. At low pressures (LP) or in ultrahigh vacuum, uniform single-layer graphene is usually obtained due to low carbon solubility in Cu.15 Whereas if graphene is grown at atmospheric pressures (AP), bilayer or multilayer graphene can be obtained.16 A possible explanation is that the growth kinetics is determined by different processes involving gas phase dynamics at LP and AP. Another experiment which has clearly demonstrated the importance of gas-phase dynamics is the dependence of graphene thickness on the position of the Cu foil in the tube furnace.17 In that experiment, it is found that graphene grown at downstream positions becomes thicker than that at upstream positions. The differences in different positions are expected to mainly come from the gas phase. For example, the concentration of methyl radicals in the gas phase, which is active for graphene growth,18 can be different. Motivated by these results, we perform a systematic computational fluid dynamics (CFD) study on gas-phase dynamics during CVD growth of graphene in a horizontal tube furnace. The surface deposition rate of carbon on a Cu surface is calculated to understand the experimental observations for CVD growth of graphene with methane as the precursor gas. Our results clearly show that graphene growth is mainly controlled by mass transport at atmospheric pressure and by surface reactions at low pressure. Additional simulations also
This journal is © the Owner Societies 2015
View Article Online
Paper
PCCP
coefficients of each species in the gas mixture are given by the following formula19 1 Xi Di;m ¼ P Xj Dij
(2)
j; jai
where Xi is the mole fraction of gas i in the mixture. Fig. 1
Schematic of the CVD reactor model used in this study.
Published on 30 July 2015. Downloaded by Stockholms Universitet on 26/08/2015 15:49:10.
Reaction model confirm that methane decomposition is initiated in the gas phase.
Computational details CVD reactor
Atomic details of graphene growth on a Cu substrate are still elusive, involving many elementary processes. In this study, we only consider the overall reaction CH4 - Chsi + 2H2, where Chsi means solid state carbon. The effective reaction constant ks can be written as E
ks ¼ AeRT
Fig. 1 shows a schematic of a general CVD reactor, including a gas-inlet pipe (8 cm in length), a tube reactor (100 cm in length), a furnace (46 cm in length), and a gas-outlet pipe (8 cm in length). The tube reactor has a diameter of 4 cm and is connected with the gas-inlet and gas-outlet pipes. Both pipes have a diameter of 1 cm. The furnace is located in the middle part of the tube reactor, providing required heating power. A Cu foil is placed in the furnace zone of the tube reactor and graphene is primarily grown on the Cu foil when gas is flowing through it. Fluid dynamics model and numerical methods The fluid dynamics model adopted here is based on 3D, laminar, compressible flow with surface chemical reaction, governed by equations of continuity and momentum/energy conservation. The flow field in the CVD reactor can be obtained by solving these equations. The CFD code FLUENT19 is used in this study, with the following numerical setup: (1) the pressurebased solver is chosen for low speed compressible flow, in which the SIMPLC (Semi-Implicit Method for Pressure-Linked Equations Consistent) algorithm is used for pressure velocity coupling; (2) the second-order upwind difference is used to discretize all convective terms due to its relatively high resolution, and the central difference is applied to discretize the diffusion term for its low numerical diffusivity; and (3) the steady state solution is obtained using an under-relaxation algorithm for momentum, energy, as well as species variables. Gas mixture Physical properties of Ar, H2, and CH4, such as heat capacity, viscosity, and thermal conductivity, are taken from the literature.20 We consider their mixture as a compressible ideal gas. The binary diffusion coefficients at temperature T and under pressure P are estimated using the following equation,20 Dij ¼
0:0026T 3=2 PMij 1=2 sij 2 OD
(1)
where Mij = 2[(1/Mi) + (1/Mj)]1. Mi and Mj are molecular weights si þ sj , with si and sj being characterof gas species i and j. sij ¼ 2 istic Lennard-Jones lengths of species i and j, respectively. OD is the dimensionless diffusion collision integral (see ESI†). The diffusion
This journal is © the Owner Societies 2015
(3)
where A is the pre-exponential factor (consistent units) and E is the activation energy. The effective activation energy is set to 2.6 eV as observed in experiments.21 The pre-exponential factor is set to 1 109 s1, to make the calculated surface deposition rate agree with experimental results21 at 1270 K within the order of magnitude. The surface deposition rate of carbon is calculated using RD = MCksCs, where MC is the molecular weight of carbon and Cs is the molar concentration of methane at the Cu surface. Computational domain, mesh arrangement, and boundary conditions The computational domain considered in this study includes a 1/2 symmetric CVD reactor. The computational mesh arrangement and boundary conditions are given in Fig. S1 (ESI†). As shown in Fig. S1 (ESI†), a structured rectangle mesh is applied for boundary surfaces and a prismatic 3D mesh is applied for the whole computational domain. To resolve the boundary layer flow near the wall and the Cu foil, a boundary layer mesh geometrically adaptive to the wall is generated. After several mesh convergence checks, about 500 thousand cells are adopted in the present investigation. The boundary conditions are set as the following: at the gas inlet, the mass flow rate, the species mass fraction, and the total temperature in the CVD reactor are directly specified. The flow velocity, pressure, and density are calculated automatically based on these specified values. At the gas outlet, a static pressure is specified and other conditions are extrapolated from the interior of the domain. At Cu foil surfaces and the furnace surface, same constant temperature is specified. Besides, for the Cu foil surface, the surface reaction model is activated. Finally, to study the sensitivity of the graphene growth behaviour to the pressure, different operating pressure conditions are specified for the whole computational domain. The main boundary parameters for different cases are listed in Table 1.
Results and discussion Fig. 2(a) shows contours of temperature on the symmetry plane at AP and LP when the furnace temperature is 1270 K. The temperature distribution is uniform above the reaction wall,
Phys. Chem. Chem. Phys., 2015, 17, 22832--22836 | 22833
View Article Online
PCCP Table 1
Paper Different boundary conditions used in this study
Table 2
P/Pa
T/K
Flow rate CH4 : H2 : Ar/sccm
T (K)
AP:101 325 LP:83
1230–1330 1230–1330
2 : 50 : 450 15 : 7 : 0
1230
Published on 30 July 2015. Downloaded by Stockholms Universitet on 26/08/2015 15:49:10.
1250
which is desirable for graphene growth. Notice that the temperature and the velocity near the boundaries of the furnace zone are not uniform at AP, which is caused by thermal buoyancy due to the relatively large density change under high pressure conditions. The distribution of velocity above the reaction wall is uniform, suggesting a laminar flow region close to the reaction wall in the furnace zone. Fig. 2(c) shows that there is a methane concentration gradient above the reaction wall at AP, with the lowest concentration on the surface. At LP, the concentration of methane is distributed uniformly. The concentration at the surface is almost the same as that in the bulk. Further calculations at different heating temperatures (1230–1330 K) give a similar velocity field, temperature field, and concentration of methane. Therefore, the thickness of the boundary layer (d) is determined mainly by pressure in the temperature range considered here. According to Fick’s law of diffusion, the mass transport coefficient of methane hg can be estimated using the formula hg = Dg/d. The values of the mass transport coefficient hg and the surface reaction constant ks are shown in Table 2 at different pressures and temperatures. In all
Fig. 2 Contours of (a) temperature, (b) velocity, and (c) molar concentration of methane distribution on the symmetry plane at different pressures. The furnace temperature is 1270 K.
22834 | Phys. Chem. Chem. Phys., 2015, 17, 22832--22836
1270 1290 1310 1330
Values of DCH4, d, hg and ks at AP and LP
DCH4 (m2 s1) AP LP AP LP AP LP AP LP AP LP AP LP
2.58 0.83 2.65 0.853 2.73 0.876 2.80 0.9 2.87 0.924 2.95 0.947
1004 04
10
1004 04
10
1004 04
10
d (m)
hg (m s1)
ks (m s1)
0.015 0.02 0.015 0.02 0.015 0.02 0.015 0.02 0.015 0.02 0.015 0.02
0.0172 41.5 0.0177 42.65 0.0182 43.8 0.0187 45 0.0191 46.2 0.0197 47.35
0.022 0.0328 0.048 0.069 0.099 0.14
cases, hg at LP is several orders of magnitude larger than that at AP. According to the Grove model,22 the rate of methane diffusion through the boundary layer is MCH4hg(Cg Cs), where MCH4 is the molecular weight of methane, and Cg and Cs are the molar concentrations of methane in the gas phase and on the Cu foil, respectively. The rate of methane being consumed on the Cu foil is MCH4ksCs. In the steady state, these two rates must be equal, so we have Cs = hgCg/(ks + hg), which leads to the surface deposition rate of carbon RD = MC[kshg/(ks + hg)]Cg. Since hg is sensitive to pressure, the surface deposition rate may be controlled either by ks or by hg at different pressures. Fig. 3 shows the surface deposition rate at different heating temperatures. Under LP, it increases exponentially with heating temperature, while it increases only moderately with temperature at AP. According to eqn (1), we have hg p T3/2. As shown in Fig. 3, the deposition rate at AP and LP can be fitted well by the 3/2 power law and the exponent law. In the latter case, the effective activation energy obtained by fitting the simulated deposition rate as an exponential function is 2.57 eV, which agrees well with the input experimental value of 2.6 eV. According to Fig. 3, it is concluded that the surface deposition rate is mainly determined by hg at AP, while at LP, it is mainly determined by ks. That is to say, under low pressure, graphene growth by CVD on a Cu substrate is mainly controlled by reaction dynamics on the foil surface reaction due to a rapid
Fig. 3
Calculated surface carbon deposition rate at different temperatures.
This journal is © the Owner Societies 2015
View Article Online
Published on 30 July 2015. Downloaded by Stockholms Universitet on 26/08/2015 15:49:10.
Paper
mass transport process caused by large diffusivity in a low density environment; under atmospheric pressure, graphene growth by CVD on a Cu substrate is mainly controlled by a mass transport process due to low diffusivity in a high density environment. We now shift our attention to the effects of Cu foil position on the deposition rate. The Cu foil is placed at different positions in the tube. A three-dimensional model (Fig. S2, ESI†) with seven Cu foil positions equally distributed (6 cm interval) in a tube is established. The Cu reaction surface length is 1.5 cm. The temperature of the furnace zone is fixed at 1270 K. Firstly, the case with a single Cu foil activated at seven different positions independently is studied, which requires seven simulations. Then, the case with seven Cu foils activated simultaneously is computed for comparison. The results at AP are shown in Fig. 4(a). It is observed that, in both cases, the surface deposition rate dramatically decreases from positions 1 to 5 or 6. Then, an increase of deposition rate downstream is observed. Variations in the case with seven Cu foils activated simultaneously are more remarkable. This is because the boundary layer heights (see Fig. S3, ESI†) at AP are very different at different positions due to thermal buoyancy effects as shown in Fig. 2 and Fig. S4 (ESI†) (first increase then gradually decrease). The deposition rate is mainly determined by the local boundary layer height which is inversely proportional to the mass transport coefficient. For the case with seven Cu foils activated simultaneously, the methane concentration decreases along the flow direction due to the consumption on sequential foils, resulting in more remarkable change in the deposition rate. However, as shown in Fig. 4(b), the situation at LP is different. The deposition rate remains almost unchanged in both cases, since it is mainly controlled by reaction dynamics on the Cu foil surface determined by the input parameter of wall temperature. As a result, the variations of the boundary layer height and the methane concentration at different positions have very small effects on the deposition rate. Notice that, the methane concentration decreases along the flow direction due to the consumption by the overall reaction on upstream Cu foils, when seven Cu foils are activated simultaneously. Therefore, it is reasonable to observe a small decrease of the surface deposition rate in Fig. 4(b) in the simultaneous activation case. Interestingly, thicker graphene samples are obtained downstream in the recent experiment under medium pressure (MP).17 As a test, we have also performed simulations using the same pressure of 2666 Pa. The variations of the calculated deposition rate are situated between those at AP and LP (Fig. S5, ESI†). Therefore, new dynamics should enter into the growth mechanisms to compensate the lower deposition rate downstream observed in this study. One thing missed here is gas phase decomposition of methane. To study the gas-phase decomposition reactions, we built a 2D axisymmetric grid. The tube reactor was assumed to be axisymmetric (Fig. S6, ESI†), the copper foil was not considered. The heating zone is set after 36 cm from the reactor entrance and 43 cm long. The flow rate was fixed at 50 sccm and 30 sccm
This journal is © the Owner Societies 2015
PCCP
Fig. 4 Variation of the calculated surface average deposition rate with the position of the Cu surface in the furnace zone under (a) AP conditions and (b) LP conditions. The sign ‘simul’ represents that the reaction wall is placed at seven different positions simultaneously in the tube and the sign ‘indep’ represents that the reaction wall is placed at seven different positions independently.
for methane and hydrogen respectively. The pressure was fixed at 2666 Pa. The chemical mechanism was implemented with the associated kinetic laws, taken from the literature.23,24 The simulation results show that methane is decomposed in the tube reactor. The generated unsaturated species, such as
Fig. 5 Mole fraction of methane and unsaturated species along the reactor axis under 2666 Pa, the black triangles mark the positions of Cu foils in the tube reactor in the experiment.17
Phys. Chem. Chem. Phys., 2015, 17, 22832--22836 | 22835
View Article Online
PCCP
Published on 30 July 2015. Downloaded by Stockholms Universitet on 26/08/2015 15:49:10.
CH3, C2H3 and C2H5, are expected to be important for graphene growth.18 Without considering the surface reaction, we find that the concentration of unsaturated species becomes significantly higher along the tube in the heating zone as can be seen from Fig. 5. These results explain the experimental observation that graphene grown at downstream positions is thicker than those at upstream positions, which confirm the experimental conjecture that gas-phase reactions are critical for CVD growth of graphene.17
Conclusions In summary, we have simulated the flow field and calculated the surface deposition rate of CVD growth of graphene on Cu foils via CFD. We have demonstrated numerically that when the reaction chamber is at low pressure, graphene growth is mainly controlled by a surface reaction process; in the reaction chamber at atmospheric pressure, graphene growth is mainly controlled by a mass transport process. Without considering the gas-phase methane decomposition reactions, simulations predict that the surface deposition rate decreases downstream. At the same time, gas phase simulations give higher concentrations of unsaturated species downstream along the reactor tube. Our simulation results thus confirm that gas-phase reactions are important integrants of the mechanism of CVD synthesis of graphene.
Acknowledgements This work was partially supported by MOST (2011CB921404 and 2014CB932700), NSFC (21173202, 21222304, and 21421063), CAS (XSB01020300), CUSF, and by USTC-SCC, SCCAS, Tianjin, and Shanghai Supercomputer Centers.
Notes and references 1 A. C. Neto, F. Guinea, N. Peres, K. S. Novoselov and A. K. Geim, Rev. Mod. Phys., 2009, 81, 109. 2 M. J. Allen, V. C. Tung and R. B. Kaner, Chem. Rev., 2009, 110, 132–145. 3 D. A. Brownson, D. K. Kampouris and C. E. Banks, J. Power Sources, 2011, 196, 4873–4885. 4 K. S. Novoselov, V. Fal, L. Colombo, P. Gellert, M. Schwab and K. Kim, Nature, 2012, 490, 192–200.
22836 | Phys. Chem. Chem. Phys., 2015, 17, 22832--22836
Paper
5 S. Bae, H. Kim, Y. Lee, X. Xu, J.-S. Park, Y. Zheng, J. Balakrishnan, T. Lei, H. R. Kim and Y. I. Song, Nat. Nanotechnol., 2010, 5, 574–578. 6 Y. Zhang, L. Zhang and C. Zhou, Acc. Chem. Res., 2013, 46, 2329–2339. 7 X. Li, W. Cai, J. An, S. Kim, J. Nah, D. Yang, R. Piner, A. Velamakanni, I. Jung and E. Tutuc, Science, 2009, 324, 1312–1314. 8 P. Wu, W. Zhang, Z. Li and J. Yang, Small, 2014, 10, 2136–2150. 9 H. Shu, X. Chen, X. Tao and F. Ding, ACS Nano, 2012, 6, 3243–3250. 10 X. Zhang, L. Wang, J. Xin, B. I. Yakobson and F. Ding, J. Am. Chem. Soc., 2014, 136, 3040–3047. 11 P. Wu, Y. Zhang, P. Cui, Z. Li, J. Yang and Z. Zhang, Phys. Rev. Lett., 2015, 114, 216102. 12 W. Zhang, P. Wu, Z. Li and J. Yang, J. Phys. Chem. C, 2011, 115, 17782–17787. 13 H. Shu, X.-M. Tao and F. Ding, Nanoscale, 2015, 7, 1627–1634. 14 P. Wu, H. Jiang, W. Zhang, Z. Li, Z. Hou and J. Yang, J. Am. Chem. Soc., 2012, 134, 6045–6051. 15 X. Li, W. Cai, L. Colombo and R. S. Ruoff, Nano Lett., 2009, 9, 4268–4272. 16 S. Bhaviripudi, X. Jia, M. S. Dresselhaus and J. Kong, Nano Lett., 2010, 10, 4128–4133. 17 Z. Li, W. Zhang, X. Fan, P. Wu, C. Zeng, Z. Li, X. Zhai, J. Yang and J. Hou, J. Phys. Chem. C, 2012, 116, 10557–10562. 18 Y. Carissan and W. Klopper, ChemPhysChem, 2006, 7, 1770–1778. 19 A. Fluent, Fluent Inc, 2006. 20 B. E. Poling, J. M. Prausnitz, O. C. John Paul and R. C. Reid, The properties of gases and liquids, McGraw-Hill, New York, 2001, vol. 5. 21 H. Kim, C. Mattevi, M. R. Calvo, J. C. Oberg, L. Artiglia, S. Agnoli, C. F. Hirjibehedin, M. Chhowalla and E. Saiz, ACS Nano, 2012, 6, 3614–3623. 22 A. S. Grove, Ind. Eng. Chem., 1966, 58, 48–55. 23 A. Holmen, O. Olsvik and O. Rokstad, Fuel Process. Technol., 1995, 42, 249–267. 24 P. Trinsoutrot, C. Rabot, H. Vergnes, A. Delamoreanu, A. Zenasni and B. Caussat, Chem. Vap. Deposition, 2014, 20, 51–58.
This journal is © the Owner Societies 2015
Published on 30 July 2015. Downloaded by Stockholms Universitet on 26/08/2015 15:49:10.
PAPER
Cite this: Phys. Chem. Chem. Phys., 2015, 17, 22832
View Journal | View Issue
Gas-phase dynamics in graphene growth by chemical vapour deposition† Gan Li,a Sheng-Hong Huang*b and Zhenyu Li*a Chemical vapour deposition on a Cu substrate is becoming a very important approach to obtain high quality graphene samples. Previous studies of graphene growth on Cu mainly focus on surface processes. However, recent experiments suggest that gas-phase dynamics also plays an important role in graphene growth. In this article, gas-phase processes are systematically studied using computational
Received 20th April 2015, Accepted 30th July 2015
fluid dynamics. Our simulations clearly show that graphene growth is limited by mass transport under
DOI: 10.1039/c5cp02301g
rate at different positions in the tube furnace and the concentration of different gas phase species are
ambient pressures while it is limited by surface reactions under low pressures. The carbon deposition calculated. Our results confirm that the previously realized graphene thickness control by changing the
www.rsc.org/pccp
position of the Cu foil is a result of gas-phase methane decomposition reactions.
Introduction Graphene is a two-dimensional material composed of carbon atoms with sp2 hybridization. It has a unique electronic structure and excellent physical properties1 which lead to a wide range of applications.2,3 Therefore, the production of highquality graphene films at a large scale is very desirable. Several approaches have been used to prepare graphene samples, including micro-mechanical exfoliation of graphite, epitaxial growth on different substrates, and reduction of graphene oxide.4 Among them, graphene growth on metal surfaces via chemical vapour deposition (CVD) is a promising way for largescale, inexpensive, high-quality graphene production.5,6 Cu is currently the most widely used substrate material due to its low carbon solubility and the resulting flexibility for graphene layer number control.7 Mechanisms of graphene growth on a Cu substrate have been intensively studied theoretically.8 The structure of graphene edges is determined by surface passivation9 and also by partial pressure of hydrogen.10 Active species for graphene growth is suggested to be of more than one carbon atom11 and may contain H atoms.12,13 The complexity of both the graphene edge structure and active species composition can lead to complex growth behaviour.14 These previous studies mainly focused on processes that occurred on the Cu surface. a
Hefei National Laboratory for Physical Science at Microscale, University of Science and Technology of China, Hefei, Anhui 230026, China. E-mail: [email protected] b Department of Modern Mechanics, University of Science and Technology of China, Hefei, Anhui 230026, China. E-mail: [email protected] † Electronic supplementary information (ESI) available: More computational details and results at intermediate pressure (2666 Pa). See DOI: 10.1039/ c5cp02301g
22832 | Phys. Chem. Chem. Phys., 2015, 17, 22832--22836
The gas phase effects have been considered only from a thermodynamic point of view.12 To the best of our knowledge, a systematic study on gas phase dynamics in graphene growth has not been reported yet. On the experimental side, however, several recent studies have suggested that gas phase dynamics is important in graphene growth. Important evidence comes from strong sensitivity of the growth behaviour to pressure. At low pressures (LP) or in ultrahigh vacuum, uniform single-layer graphene is usually obtained due to low carbon solubility in Cu.15 Whereas if graphene is grown at atmospheric pressures (AP), bilayer or multilayer graphene can be obtained.16 A possible explanation is that the growth kinetics is determined by different processes involving gas phase dynamics at LP and AP. Another experiment which has clearly demonstrated the importance of gas-phase dynamics is the dependence of graphene thickness on the position of the Cu foil in the tube furnace.17 In that experiment, it is found that graphene grown at downstream positions becomes thicker than that at upstream positions. The differences in different positions are expected to mainly come from the gas phase. For example, the concentration of methyl radicals in the gas phase, which is active for graphene growth,18 can be different. Motivated by these results, we perform a systematic computational fluid dynamics (CFD) study on gas-phase dynamics during CVD growth of graphene in a horizontal tube furnace. The surface deposition rate of carbon on a Cu surface is calculated to understand the experimental observations for CVD growth of graphene with methane as the precursor gas. Our results clearly show that graphene growth is mainly controlled by mass transport at atmospheric pressure and by surface reactions at low pressure. Additional simulations also
This journal is © the Owner Societies 2015
View Article Online
Paper
PCCP
coefficients of each species in the gas mixture are given by the following formula19 1 Xi Di;m ¼ P Xj Dij
(2)
j; jai
where Xi is the mole fraction of gas i in the mixture. Fig. 1
Schematic of the CVD reactor model used in this study.
Published on 30 July 2015. Downloaded by Stockholms Universitet on 26/08/2015 15:49:10.
Reaction model confirm that methane decomposition is initiated in the gas phase.
Computational details CVD reactor
Atomic details of graphene growth on a Cu substrate are still elusive, involving many elementary processes. In this study, we only consider the overall reaction CH4 - Chsi + 2H2, where Chsi means solid state carbon. The effective reaction constant ks can be written as E
ks ¼ AeRT
Fig. 1 shows a schematic of a general CVD reactor, including a gas-inlet pipe (8 cm in length), a tube reactor (100 cm in length), a furnace (46 cm in length), and a gas-outlet pipe (8 cm in length). The tube reactor has a diameter of 4 cm and is connected with the gas-inlet and gas-outlet pipes. Both pipes have a diameter of 1 cm. The furnace is located in the middle part of the tube reactor, providing required heating power. A Cu foil is placed in the furnace zone of the tube reactor and graphene is primarily grown on the Cu foil when gas is flowing through it. Fluid dynamics model and numerical methods The fluid dynamics model adopted here is based on 3D, laminar, compressible flow with surface chemical reaction, governed by equations of continuity and momentum/energy conservation. The flow field in the CVD reactor can be obtained by solving these equations. The CFD code FLUENT19 is used in this study, with the following numerical setup: (1) the pressurebased solver is chosen for low speed compressible flow, in which the SIMPLC (Semi-Implicit Method for Pressure-Linked Equations Consistent) algorithm is used for pressure velocity coupling; (2) the second-order upwind difference is used to discretize all convective terms due to its relatively high resolution, and the central difference is applied to discretize the diffusion term for its low numerical diffusivity; and (3) the steady state solution is obtained using an under-relaxation algorithm for momentum, energy, as well as species variables. Gas mixture Physical properties of Ar, H2, and CH4, such as heat capacity, viscosity, and thermal conductivity, are taken from the literature.20 We consider their mixture as a compressible ideal gas. The binary diffusion coefficients at temperature T and under pressure P are estimated using the following equation,20 Dij ¼
0:0026T 3=2 PMij 1=2 sij 2 OD
(1)
where Mij = 2[(1/Mi) + (1/Mj)]1. Mi and Mj are molecular weights si þ sj , with si and sj being characterof gas species i and j. sij ¼ 2 istic Lennard-Jones lengths of species i and j, respectively. OD is the dimensionless diffusion collision integral (see ESI†). The diffusion
This journal is © the Owner Societies 2015
(3)
where A is the pre-exponential factor (consistent units) and E is the activation energy. The effective activation energy is set to 2.6 eV as observed in experiments.21 The pre-exponential factor is set to 1 109 s1, to make the calculated surface deposition rate agree with experimental results21 at 1270 K within the order of magnitude. The surface deposition rate of carbon is calculated using RD = MCksCs, where MC is the molecular weight of carbon and Cs is the molar concentration of methane at the Cu surface. Computational domain, mesh arrangement, and boundary conditions The computational domain considered in this study includes a 1/2 symmetric CVD reactor. The computational mesh arrangement and boundary conditions are given in Fig. S1 (ESI†). As shown in Fig. S1 (ESI†), a structured rectangle mesh is applied for boundary surfaces and a prismatic 3D mesh is applied for the whole computational domain. To resolve the boundary layer flow near the wall and the Cu foil, a boundary layer mesh geometrically adaptive to the wall is generated. After several mesh convergence checks, about 500 thousand cells are adopted in the present investigation. The boundary conditions are set as the following: at the gas inlet, the mass flow rate, the species mass fraction, and the total temperature in the CVD reactor are directly specified. The flow velocity, pressure, and density are calculated automatically based on these specified values. At the gas outlet, a static pressure is specified and other conditions are extrapolated from the interior of the domain. At Cu foil surfaces and the furnace surface, same constant temperature is specified. Besides, for the Cu foil surface, the surface reaction model is activated. Finally, to study the sensitivity of the graphene growth behaviour to the pressure, different operating pressure conditions are specified for the whole computational domain. The main boundary parameters for different cases are listed in Table 1.
Results and discussion Fig. 2(a) shows contours of temperature on the symmetry plane at AP and LP when the furnace temperature is 1270 K. The temperature distribution is uniform above the reaction wall,
Phys. Chem. Chem. Phys., 2015, 17, 22832--22836 | 22833
View Article Online
PCCP Table 1
Paper Different boundary conditions used in this study
Table 2
P/Pa
T/K
Flow rate CH4 : H2 : Ar/sccm
T (K)
AP:101 325 LP:83
1230–1330 1230–1330
2 : 50 : 450 15 : 7 : 0
1230
Published on 30 July 2015. Downloaded by Stockholms Universitet on 26/08/2015 15:49:10.
1250
which is desirable for graphene growth. Notice that the temperature and the velocity near the boundaries of the furnace zone are not uniform at AP, which is caused by thermal buoyancy due to the relatively large density change under high pressure conditions. The distribution of velocity above the reaction wall is uniform, suggesting a laminar flow region close to the reaction wall in the furnace zone. Fig. 2(c) shows that there is a methane concentration gradient above the reaction wall at AP, with the lowest concentration on the surface. At LP, the concentration of methane is distributed uniformly. The concentration at the surface is almost the same as that in the bulk. Further calculations at different heating temperatures (1230–1330 K) give a similar velocity field, temperature field, and concentration of methane. Therefore, the thickness of the boundary layer (d) is determined mainly by pressure in the temperature range considered here. According to Fick’s law of diffusion, the mass transport coefficient of methane hg can be estimated using the formula hg = Dg/d. The values of the mass transport coefficient hg and the surface reaction constant ks are shown in Table 2 at different pressures and temperatures. In all
Fig. 2 Contours of (a) temperature, (b) velocity, and (c) molar concentration of methane distribution on the symmetry plane at different pressures. The furnace temperature is 1270 K.
22834 | Phys. Chem. Chem. Phys., 2015, 17, 22832--22836
1270 1290 1310 1330
Values of DCH4, d, hg and ks at AP and LP
DCH4 (m2 s1) AP LP AP LP AP LP AP LP AP LP AP LP
2.58 0.83 2.65 0.853 2.73 0.876 2.80 0.9 2.87 0.924 2.95 0.947
1004 04
10
1004 04
10
1004 04
10
d (m)
hg (m s1)
ks (m s1)
0.015 0.02 0.015 0.02 0.015 0.02 0.015 0.02 0.015 0.02 0.015 0.02
0.0172 41.5 0.0177 42.65 0.0182 43.8 0.0187 45 0.0191 46.2 0.0197 47.35
0.022 0.0328 0.048 0.069 0.099 0.14
cases, hg at LP is several orders of magnitude larger than that at AP. According to the Grove model,22 the rate of methane diffusion through the boundary layer is MCH4hg(Cg Cs), where MCH4 is the molecular weight of methane, and Cg and Cs are the molar concentrations of methane in the gas phase and on the Cu foil, respectively. The rate of methane being consumed on the Cu foil is MCH4ksCs. In the steady state, these two rates must be equal, so we have Cs = hgCg/(ks + hg), which leads to the surface deposition rate of carbon RD = MC[kshg/(ks + hg)]Cg. Since hg is sensitive to pressure, the surface deposition rate may be controlled either by ks or by hg at different pressures. Fig. 3 shows the surface deposition rate at different heating temperatures. Under LP, it increases exponentially with heating temperature, while it increases only moderately with temperature at AP. According to eqn (1), we have hg p T3/2. As shown in Fig. 3, the deposition rate at AP and LP can be fitted well by the 3/2 power law and the exponent law. In the latter case, the effective activation energy obtained by fitting the simulated deposition rate as an exponential function is 2.57 eV, which agrees well with the input experimental value of 2.6 eV. According to Fig. 3, it is concluded that the surface deposition rate is mainly determined by hg at AP, while at LP, it is mainly determined by ks. That is to say, under low pressure, graphene growth by CVD on a Cu substrate is mainly controlled by reaction dynamics on the foil surface reaction due to a rapid
Fig. 3
Calculated surface carbon deposition rate at different temperatures.
This journal is © the Owner Societies 2015
View Article Online
Published on 30 July 2015. Downloaded by Stockholms Universitet on 26/08/2015 15:49:10.
Paper
mass transport process caused by large diffusivity in a low density environment; under atmospheric pressure, graphene growth by CVD on a Cu substrate is mainly controlled by a mass transport process due to low diffusivity in a high density environment. We now shift our attention to the effects of Cu foil position on the deposition rate. The Cu foil is placed at different positions in the tube. A three-dimensional model (Fig. S2, ESI†) with seven Cu foil positions equally distributed (6 cm interval) in a tube is established. The Cu reaction surface length is 1.5 cm. The temperature of the furnace zone is fixed at 1270 K. Firstly, the case with a single Cu foil activated at seven different positions independently is studied, which requires seven simulations. Then, the case with seven Cu foils activated simultaneously is computed for comparison. The results at AP are shown in Fig. 4(a). It is observed that, in both cases, the surface deposition rate dramatically decreases from positions 1 to 5 or 6. Then, an increase of deposition rate downstream is observed. Variations in the case with seven Cu foils activated simultaneously are more remarkable. This is because the boundary layer heights (see Fig. S3, ESI†) at AP are very different at different positions due to thermal buoyancy effects as shown in Fig. 2 and Fig. S4 (ESI†) (first increase then gradually decrease). The deposition rate is mainly determined by the local boundary layer height which is inversely proportional to the mass transport coefficient. For the case with seven Cu foils activated simultaneously, the methane concentration decreases along the flow direction due to the consumption on sequential foils, resulting in more remarkable change in the deposition rate. However, as shown in Fig. 4(b), the situation at LP is different. The deposition rate remains almost unchanged in both cases, since it is mainly controlled by reaction dynamics on the Cu foil surface determined by the input parameter of wall temperature. As a result, the variations of the boundary layer height and the methane concentration at different positions have very small effects on the deposition rate. Notice that, the methane concentration decreases along the flow direction due to the consumption by the overall reaction on upstream Cu foils, when seven Cu foils are activated simultaneously. Therefore, it is reasonable to observe a small decrease of the surface deposition rate in Fig. 4(b) in the simultaneous activation case. Interestingly, thicker graphene samples are obtained downstream in the recent experiment under medium pressure (MP).17 As a test, we have also performed simulations using the same pressure of 2666 Pa. The variations of the calculated deposition rate are situated between those at AP and LP (Fig. S5, ESI†). Therefore, new dynamics should enter into the growth mechanisms to compensate the lower deposition rate downstream observed in this study. One thing missed here is gas phase decomposition of methane. To study the gas-phase decomposition reactions, we built a 2D axisymmetric grid. The tube reactor was assumed to be axisymmetric (Fig. S6, ESI†), the copper foil was not considered. The heating zone is set after 36 cm from the reactor entrance and 43 cm long. The flow rate was fixed at 50 sccm and 30 sccm
This journal is © the Owner Societies 2015
PCCP
Fig. 4 Variation of the calculated surface average deposition rate with the position of the Cu surface in the furnace zone under (a) AP conditions and (b) LP conditions. The sign ‘simul’ represents that the reaction wall is placed at seven different positions simultaneously in the tube and the sign ‘indep’ represents that the reaction wall is placed at seven different positions independently.
for methane and hydrogen respectively. The pressure was fixed at 2666 Pa. The chemical mechanism was implemented with the associated kinetic laws, taken from the literature.23,24 The simulation results show that methane is decomposed in the tube reactor. The generated unsaturated species, such as
Fig. 5 Mole fraction of methane and unsaturated species along the reactor axis under 2666 Pa, the black triangles mark the positions of Cu foils in the tube reactor in the experiment.17
Phys. Chem. Chem. Phys., 2015, 17, 22832--22836 | 22835
View Article Online
PCCP
Published on 30 July 2015. Downloaded by Stockholms Universitet on 26/08/2015 15:49:10.
CH3, C2H3 and C2H5, are expected to be important for graphene growth.18 Without considering the surface reaction, we find that the concentration of unsaturated species becomes significantly higher along the tube in the heating zone as can be seen from Fig. 5. These results explain the experimental observation that graphene grown at downstream positions is thicker than those at upstream positions, which confirm the experimental conjecture that gas-phase reactions are critical for CVD growth of graphene.17
Conclusions In summary, we have simulated the flow field and calculated the surface deposition rate of CVD growth of graphene on Cu foils via CFD. We have demonstrated numerically that when the reaction chamber is at low pressure, graphene growth is mainly controlled by a surface reaction process; in the reaction chamber at atmospheric pressure, graphene growth is mainly controlled by a mass transport process. Without considering the gas-phase methane decomposition reactions, simulations predict that the surface deposition rate decreases downstream. At the same time, gas phase simulations give higher concentrations of unsaturated species downstream along the reactor tube. Our simulation results thus confirm that gas-phase reactions are important integrants of the mechanism of CVD synthesis of graphene.
Acknowledgements This work was partially supported by MOST (2011CB921404 and 2014CB932700), NSFC (21173202, 21222304, and 21421063), CAS (XSB01020300), CUSF, and by USTC-SCC, SCCAS, Tianjin, and Shanghai Supercomputer Centers.
Notes and references 1 A. C. Neto, F. Guinea, N. Peres, K. S. Novoselov and A. K. Geim, Rev. Mod. Phys., 2009, 81, 109. 2 M. J. Allen, V. C. Tung and R. B. Kaner, Chem. Rev., 2009, 110, 132–145. 3 D. A. Brownson, D. K. Kampouris and C. E. Banks, J. Power Sources, 2011, 196, 4873–4885. 4 K. S. Novoselov, V. Fal, L. Colombo, P. Gellert, M. Schwab and K. Kim, Nature, 2012, 490, 192–200.
22836 | Phys. Chem. Chem. Phys., 2015, 17, 22832--22836
Paper
5 S. Bae, H. Kim, Y. Lee, X. Xu, J.-S. Park, Y. Zheng, J. Balakrishnan, T. Lei, H. R. Kim and Y. I. Song, Nat. Nanotechnol., 2010, 5, 574–578. 6 Y. Zhang, L. Zhang and C. Zhou, Acc. Chem. Res., 2013, 46, 2329–2339. 7 X. Li, W. Cai, J. An, S. Kim, J. Nah, D. Yang, R. Piner, A. Velamakanni, I. Jung and E. Tutuc, Science, 2009, 324, 1312–1314. 8 P. Wu, W. Zhang, Z. Li and J. Yang, Small, 2014, 10, 2136–2150. 9 H. Shu, X. Chen, X. Tao and F. Ding, ACS Nano, 2012, 6, 3243–3250. 10 X. Zhang, L. Wang, J. Xin, B. I. Yakobson and F. Ding, J. Am. Chem. Soc., 2014, 136, 3040–3047. 11 P. Wu, Y. Zhang, P. Cui, Z. Li, J. Yang and Z. Zhang, Phys. Rev. Lett., 2015, 114, 216102. 12 W. Zhang, P. Wu, Z. Li and J. Yang, J. Phys. Chem. C, 2011, 115, 17782–17787. 13 H. Shu, X.-M. Tao and F. Ding, Nanoscale, 2015, 7, 1627–1634. 14 P. Wu, H. Jiang, W. Zhang, Z. Li, Z. Hou and J. Yang, J. Am. Chem. Soc., 2012, 134, 6045–6051. 15 X. Li, W. Cai, L. Colombo and R. S. Ruoff, Nano Lett., 2009, 9, 4268–4272. 16 S. Bhaviripudi, X. Jia, M. S. Dresselhaus and J. Kong, Nano Lett., 2010, 10, 4128–4133. 17 Z. Li, W. Zhang, X. Fan, P. Wu, C. Zeng, Z. Li, X. Zhai, J. Yang and J. Hou, J. Phys. Chem. C, 2012, 116, 10557–10562. 18 Y. Carissan and W. Klopper, ChemPhysChem, 2006, 7, 1770–1778. 19 A. Fluent, Fluent Inc, 2006. 20 B. E. Poling, J. M. Prausnitz, O. C. John Paul and R. C. Reid, The properties of gases and liquids, McGraw-Hill, New York, 2001, vol. 5. 21 H. Kim, C. Mattevi, M. R. Calvo, J. C. Oberg, L. Artiglia, S. Agnoli, C. F. Hirjibehedin, M. Chhowalla and E. Saiz, ACS Nano, 2012, 6, 3614–3623. 22 A. S. Grove, Ind. Eng. Chem., 1966, 58, 48–55. 23 A. Holmen, O. Olsvik and O. Rokstad, Fuel Process. Technol., 1995, 42, 249–267. 24 P. Trinsoutrot, C. Rabot, H. Vergnes, A. Delamoreanu, A. Zenasni and B. Caussat, Chem. Vap. Deposition, 2014, 20, 51–58.
This journal is © the Owner Societies 2015
