記住我
The purpose of this Letter is to investigate thermal transport in IB-TBG 2D diamond superstructures and determine the dependence of their thermal conductivity on the 2D diamond fraction of the superstructures based on molecular-dynamics (MD) simulations. As a metric for 2D diamond fraction, we use the fraction, fsp3, of C atoms in the superstructure that participate in interlayer covalent C–C bond formation, NCib: fsp3≡NCib/NCtot, where NCtot is the total number of C atoms in the bilayer. We find that the introduction of interlayer C–C bonds in these bilayer structures leads to an abrupt drop in thermal conductivity compared to that of pristine, non-interlayer-bonded bilayer graphene. However, increased formation of interlayer C–C bonds leads to a monotonic increase in the thermal conductivity of the resulting superstructures with increasing 2D diamond fraction toward the high thermal conductivity of 2D diamond, demonstrating that the thermal conductivity of IB-TBG 2D diamond superstructures can be precisely tuned by controlling the 2D diamond fraction in the superstructures. We also find that the thermal conductivity of interlayer-bonded graphene bilayers with randomly distributed individual interlayer C–C bonds (RD-IBGs) follows a similar trend with that of the IB-TBG superstructures as a function of fsp3; however, at given fsp3, the thermal conductivity of the IB-TBG 2D diamond superstructures exceeds that of RD-IBGs. The results are discussed using effective medium theories and percolation theory and explained on the basis of lattice distortions induced in the bilayer structures as a result of interlayer bonding.
To generate and equilibrate the graphene bilayer structures and compute their thermal conductivity, we carried out MD simulations using the LAMMPS software package.3939. S. Plimpton, J. Comput. Phys. 117, 1–19 (1995). https://doi.org/10.1006/jcph.1995.1039 Our MD simulations account only for the lattice (phonon) contribution to the thermal conductivity; this contribution amounts to ∼99% of the total thermal conductivity of graphene40,4140. D. L. Nika and A. A. Balandin, Rep. Prog. Phys. 80, 036502 (2017). https://doi.org/10.1088/1361-6633/80/3/03650241. S. K. Jaćimovski, M. Bukurov, J. P. Šetrajčić, and D. I. Raković, Superlattices Microstruct. 88, 330–337 (2015). https://doi.org/10.1016/j.spmi.2015.09.027 and is a fairly justified approximation of the thermal conductivity of graphene derivatives. For the description of interatomic interactions, we employed the adaptive interatomic reactive empirical bond-order (AIREBO) potential function,4242. S. J. Stuart, A. B. Tutein, and J. A. Harrison, J. Chem. Phys. 112, 6472–6486 (2000). https://doi.org/10.1063/1.481208 a reliable reactive potential that has been very commonly used in the literature for atomic-scale modeling of graphene-based materials.4343. G. Dhaliwal, P. B. Nair, and C. V. Singh, Carbon 142, 300–310 (2019). https://doi.org/10.1016/j.carbon.2018.10.020 In our simulations, we set the cutoff distance of C–C and Lennard-Jones interactions to 2 and 10.2 Å, respectively, while neglecting any torsional contributions to the potential due to their minimal impact on the computation of thermal transport properties, as reported in other studies.44–4644. X. Zhang and J. Jiang, J. Phys. Chem. C 117, 18441–18447 (2013). https://doi.org/10.1021/jp405156y45. L. Cui, Y. Feng, and X. Zhang, J. Phys. Chem. A 119, 11226–11232 (2015). https://doi.org/10.1021/acs.jpca.5b0799546. J. Chen, Y. Gao, C. Wang, R. Zhang, H. Zhao, and H. Fang, J. Phys. Chem. C 119, 17362–17368 (2015). https://doi.org/10.1021/acs.jpcc.5b02235 We used large-area supercells with a sufficiently thick vacuum layer perpendicular to the graphene plane and periodic boundary conditions in all three directions. The supercell sizes in this study range from 8064 atoms (in pristine graphene bilayers) to 12 096 (C and H) atoms in 2D diamond (diamane) sheets. For the generation of the IB-TBG 2D diamond superstructures, we followed the detailed structure formation algorithm of Ref. 3434. M. Chen, A. R. Muniz, and D. Maroudas, ACS Appl. Mater. Interfaces 10, 28898–28908 (2018). https://doi.org/10.1021/acsami.8b09741. All the interlayer-bonded structures were fully relaxed prior to carrying out MD simulations of thermal transport; this structural relaxation scheme involved conjugate-gradient energy minimization, isothermal-isobaric (NPT) MD simulation, and canonical (NVT) MD simulation in a manner similar to that utilized in a previous study.4747. X. Zhang, Y. Gao, Y. Chen, and M. Hu, Sci. Rep. 6, 22011 (2016). https://doi.org/10.1038/srep22011Figure 1 shows various views of representative, relaxed atomic configurations of interlayer-bonded graphene bilayers. The starting configurations for the formation of these structures are graphene bilayers characterized by either AA (hexagonal) or AB (Bernal) stacking arrangements or TBG with specified twist angle θ and commensurate bilayers resulting in a Moiré superlattice that exhibits AA- and AB-stacked nanodomains; this requires θ<15°, with θ=0° corresponding to AA-stacked bilayers. Figure 1(a) shows an IB-TBG 2D diamond superstructure with AA-stacked 2D nanodiamond clusters at fsp3=22.13% and θ=7.341°; prescribing these three degrees of freedom, namely, twist angle, fsp3, and stacking arrangement of the nanodomain of the Moiré superlattice where interlayer C–C bonds are formed (induced by hydrogenation as shown by the H atoms chemisorbed onto each graphene layer of the bilayer) is required to fully define an IB-TBG superstructure. Another IB-TBG 2D diamond superstructure with AB-stacked nanodiamond clusters with fsp3=11.09% and θ=5.509° is shown in Fig. 1(b). If the chemical functionalization of the graphene layers is not patterned, interlayer C–C bonding does not occur in a periodic arrangement and, instead of 2D diamond superstructures, interlayer-bonded graphene bilayers with randomly distributed individual interlayer bonds (RD-IBGs) are formed. Figure 1(c) shows such a representative RD-IBG configuration for an interlayer-bonded AB-stacked graphene bilayer at fsp3=10%. A fully interlayer-bonded (fsp3=50%) AB-stacked graphene bilayer, corresponding to a 2D cubic diamond layer (diamane), is depicted in Fig. 1(d).For analysis of thermal transport in the interlayer-bonded graphene bilayer structures under consideration, we implemented the equilibrium molecular-dynamics (EMD) approach4848. D. Frenkel and B. Smit, Understanding Molecular Simulation: From Algorithms to Applications ( Elsevier, 2001), Vol. 1. with a time step of 0.1 fs and computed the thermal conductivity of these structures by employing the Green–Kubo (GK) formula.4949. P. K. Schelling, S. R. Phillpot, and P. Keblinski, Phys. Rev. B 65, 144306 (2002). https://doi.org/10.1103/PhysRevB.65.144306 A detailed discussion of EMD simulations in LAMMPS for lattice thermal conductivity computations using the GK formula5050. Z. Fan, L. F. C. Pereira, H.-Q. Wang, J.-C. Zheng, D. Donadio, and A. Harju, Phys. Rev. B 92, 094301 (2015). https://doi.org/10.1103/PhysRevB.92.094301 and comparisons of such EMD predictions for monolayer and bilayer graphene sheets with those of non-equilibrium molecular-dynamics (NEMD) simulations are presented in the supplementary material. An optimum set of EMD simulation parameters (a correlation time of 35 ps and an ensemble size of 100), derived based on thorough convergence tests, was used for accurate heat flux calculations, which resulted in faster convergence of the heat-current autocorrelation function (HCACF). For each bilayer structure considered, three independent EMD simulation runs were performed varying the initial distribution of atomic velocities. The total duration of each EMD simulation run was 5.25 ns, with data for thermal conductivity computation collected after (typically) a couple of ns when higher-amplitude fluctuations in the computed thermal conductivity values had dissipated.Our EMD simulation results for the thermal conductivity of interlayer-bonded graphene bilayers as a function of fsp3 are shown in Fig. 2(a). Graphene bilayer structures with fsp3=0 and fsp3=50% correspond to pristine non-interlayer-bonded bilayer graphene (BG, AA-stacked, and AB-stacked) and 2D diamond sheets (diamanes, 2D AB-stacked cubic diamond, and 2D AA-stacked hexagonal diamond or lonsdaleite), respectively, which are the reference structures for this study. The average EMD-computed lattice thermal conductivity of BG is ∼370 W/mK, in agreement with Ref. 4747. X. Zhang, Y. Gao, Y. Chen, and M. Hu, Sci. Rep. 6, 22011 (2016). https://doi.org/10.1038/srep22011, whereas the corresponding thermal conductivity of 2D diamond sheets ranged between 1250 and 1450 W/mK. These results are in very good agreement with those of a recent theoretical study, which reported a thermal conductivity of 1360 W/mK for AB-stacked diamane according to a machine-learning-optimized potential.5151. S. Chowdhury, V. A. Demin, L. A. Chernozatonskii, and A. G. Kvashnin, Membranes 12, 925 (2022). https://doi.org/10.3390/membranes12100925 The results of Fig. 2(a) reveal an initial abrupt decrease in thermal conductivity with respect to that of pristine non-interlayer-bonded BG with the introduction of interlayer C–C bonds at very low values of fsp3. Nevertheless, the results indicate that, after this initial drop in thermal conductivity, increasing 2D diamond fraction (i.e., increasing fsp3) causes a monotonic increase in the thermal conductivity of IB-TBG 2D diamond superstructures, implying the capability to precisely tune the thermal conductivity of the superstructures by properly controlling fsp3; this trend in the thermal conductivity of the IB-TBG superstructures with fsp3 tends to raise their thermal conductivity toward the very high thermal conductivity of 2D diamond. It should be mentioned that a similar trend was found in a recent study on the interlayer thermal resistance in multilayer graphene and was attributed to the dual role of interlayer crosslinks in such multilayer graphene structures.5252. Y. Chen, J. Wan, Y. Chen, H. Qin, Y. Liu, Q.-X. Pei, and Y.-W. Zhang, Int. J. Therm. Sci. 183, 107871 (2023). https://doi.org/10.1016/j.ijthermalsci.2022.107871 The results also indicate a similar trend, namely, a monotonic increase in thermal conductivity with increasing fsp3, for the RD-IBG bilayers. However, interestingly, it is evident from the results of Fig. 2(a) that the thermal conductivity of the IB-TBG 2D diamond superstructures is consistently higher than that of RD-IBGs at given fsp3. This observation may be explained by noting that at any value of fsp3, an IB-TBG nanodiamond superstructure always contains continuous connected pathways for phonon transport, whereas an RD-IBG bilayer, being a random network of interlayer C–C bonds in bilayer graphene, is statistically more likely to contain random scattering centers that impede phonon transport. Moreover, we notice that bilayer stacking order does not have any significant effect on the thermal conductivity values of interlayer-bonded graphene bilayers. Specifically, the EMD-computed thermal conductivities of AA- and AB-stacked RD-IBGs are practically identical (to within statistical error) at given fsp3. In a similar manner, the thermal conductivities of IB-TBG 2D diamond superstructures with AA- and AB-stacked nanodiamond domains do not show any substantial difference due to stacking order at given fsp3. This finding also is consistent with the results of a recent study, which reported that the stacking patterns and interlayer twist angles play a negligible role in determining the in-plane and interlayer mechanical behaviors of IB-TBG nanodiamond superstructures.5353. S. Liu, Y. Chen, and Y. Liu, J. Appl. Phys. 132, 235104 (2022). https://doi.org/10.1063/5.0128970 Finally, it should be mentioned that we carried out selected tests on RD-IBG bilayer sheets without H passivation and found that they exhibit a monotonic decrease in thermal conductivity with increasing fsp3 in agreement with prior studies.4747. X. Zhang, Y. Gao, Y. Chen, and M. Hu, Sci. Rep. 6, 22011 (2016). https://doi.org/10.1038/srep22011 However, as shown in our prior work,23,37,5423. A. R. Muniz and D. Maroudas, J. Phys. Chem. C 117, 7315–7325 (2013). https://doi.org/10.1021/jp310184c37. A. R. Muniz and D. Maroudas, Phys. Rev. B 86, 075404 (2012). https://doi.org/10.1103/PhysRevB.86.07540454. A. R. Muniz and D. Maroudas, J. Appl. Phys. 111, 043513 (2012). https://doi.org/10.1063/1.3682475 proper passivation of dangling C bonds is necessary to ensure the structural stability of the interlayer-bonded graphene bilayers; moreover, this is also the situation of experimental relevance, wherein chemisorption of H atoms onto graphene layers acts as a trigger for interlayer bonding in few-layer graphene.5555. S. Rajasekaran, F. Abild-Pedersen, H. Ogasawara, A. Nilsson, and S. Kaya, Phys. Rev. Lett. 111, 085503 (2013). https://doi.org/10.1103/PhysRevLett.111.085503IB-TBG 2D diamond superstructures and RD-IBG bilayers can be viewed as graphene–diamond nanocomposite structures with BG considered as the matrix and 2D diamond nanoclusters or, in general, sp3 interlayer C–C bonds considered as the filler. Therefore, the EMD results of Fig. 2(a) can be examined in the context of effective medium theories (models), which give the thermal conductivity of the composite as an effective conductivity, keff, that can be expressed in terms of the thermal conductivities of the matrix and the filler, km and kf, respectively, as a function of the filler volume fraction, ϕ; for the interlayer-bonded graphene bilayer sheets under consideration, km and kf represent the thermal conductivity of BG and 2D diamond, respectively, and ϕ=2fsp3.The outcomes of this analysis are shown in Fig. 2(b), where the main plot and the inset show the EMD results for AA- and AB-stacked bilayers, respectively, or AA- and AB-stacked diamond nanodomains in IB-TBG superstructures, respectively. The effective medium theoretical predictions are shown with solid and dot-dashed lines in the main plot and inset, respectively. The first model employed is the rule of mixtures (ROM), which expresses keff as a linear function of ϕ, The next model used is based on Maxwell's effective medium theory (MEMT),5656. M. L. Levin and M. A. Miller, Sov. Phys. Usp. 24, 904 (1981). https://doi.org/10.1070/PU1981v024n11ABEH004793 which assumes that filler particles do not interact with each other and gives keff as keffkm=1+3ϕ(kf+2km)(kf−km)−ϕ.(2)We have also used a model based on Rayleigh's effective medium theory (REMT),5757. L. Rayleigh, London, Edinburgh, Dublin Philos. Mag. J. Sci. 34, 481–502 (1892). https://doi.org/10.1080/14786449208620364 which accounts for interactions between filler particles and expresses keff as keffkm=1+2ϕc1−ϕ+c2(0.306ϕ4+0.013ϕ8),where c1=(kf−km)(kf+km) and c2=(kf+km)(kf−km).(3)The terms in parentheses in the denominator at the right-hand side of Eq. (3) represent an infinite series in powers of ϕ that has been truncated retaining terms up to O(ϕ8) since higher-order terms make negligible contributions to the series for ϕ≤50%. It is evident from Fig. 2(b) that the effective medium models of Eqs. (1)–(3) fail to fully capture the trends in the dependence of the thermal conductivity of IB-TBG 2D diamond superstructures on nanodiamond fraction (or fsp3), both the abrupt thermal conductivity drop at low fsp3 and the type of monotonic increase in the thermal conductivity with increasing fsp3 at higher fsp3. The models also completely fail to capture the trend in the thermal conductivity of RD-IBGs as a function of fsp3. However, the EMD data for RD-IBGs [open circles in Fig. 2(b)] can be described satisfactorily, beyond the initial sudden drop, by a model based on percolation theory,58,5958. D. Stauffer and A. Aharony, Introduction to Percolation Theory ( Taylor & Francis, 2018).59. S. Torquato, Random Heterogeneous Materials: Microstructure and Macroscopic Properties ( Springer, 2002). which expresses keff according to the power-law scaling relation at ϕ > ϕc, where κ is a scaling factor, ϕc is the percolation threshold, and ν is the critical exponent. Black solid and dot-dashed lines in the main plot and inset of Fig. 2(b), respectively, represent optimal fits of the EMD data for RD-IBGs according to Eq. (4). The resulting values of the fitting parameters κ, ϕc, and ν are 21.40 ± 1.48 W/mK, 5.02 ± 0.57%, and 1.11 ± 0.02, respectively, for RD-IBGs in AA-stacked graphene bilayers and 20.14 ± 1.99 W/mK, 4.80 ± 0.45%, and 1.108 ± 0.03, respectively, for RD-IBGs in AB-stacked graphene bilayers.We attribute the inability of the effective medium models to fully capture the dependence of the thermal conductivity of the various interlayer-bonded bilayer graphene structures on the interlayer bond density (or nanodiamond fraction) to the lattice distortions in the bilayer graphene sheets induced by the interlayer C–C bonding, which are not accounted for in the effective medium theories. These lattice distortions effectively change the shape of the “filler particles” and the interactions between their stress fields, thus affecting phonon transport and the resulting thermal conductivity. Representative lattice distortion effects are shown in Fig. 3. Figure 3(a) depicts the atomic-level stress distribution in an RD-IBG at low interlayer bond density, fsp3=1.25%; the atomic-level stress is computed from the virial stress, and the stress metric used is the second scalar invariant of the stress tensor (von Mises stress). These stress fields decay radially outwards from each C–C bond (red centers of high-stress concentration) in the plane of the bilayer, but they are overlapping, effectively creating a domain wall that impedes the in-plane paths of phonons, increasing phonon scattering and, thus, reducing the in-plane thermal conductivity of the RD-IBG sheet compared to that of unstressed graphene. The effect is analogous in IB-TBG nanodiamond superstructures at low fsp3, the difference being that the interlayer bonds (or small 2D diamond nanoclusters) are regularly instead of randomly distributed. This effect explains the abrupt drop in the thermal conductivities of the interlayer-bonded bilayer sheets at low fsp3 demonstrated in Fig. 2. The effect is stronger in RD-IBGs due to the higher probability for formation of extended domain walls for randomly distributed interlayer C–C bonds compared to those in IB-TBG superstructures, where regularly arranged gaps in such domain walls have a weaker effect on phonon scattering. Figure 3(b) depicts the C–C bond length distribution in an IB-TBG 2D diamond superstructure at a relatively high nanodiamond fraction, fsp3=14.9%, demonstrating interlayer C–C bond stretching that increases radially outwards from the center of each nanodiamond cluster in the superstructure and generates a displacement field that overlaps with that induced in the graphene regions of the superstructure; these overlapping deformation fields effectively create a continuous/connected fast “highway” for in-plane phonon transport in the superstructure, thus resulting in a high thermal conductivity. Again, the effect is analogous in RD-IBG sheets at high fsp3, which results in a continuous percolation network of overlapping randomly arranged deformation fields that accelerates phonon transport and increases thermal conductivity with increasing fsp3. The continuity/connectivity of the in-plane fast phonon highways in the IB-TBG 2D diamond superstructures extending throughout each superstructure, as demonstrated in Fig. 3(b), explains the higher thermal conductivity of the IB-TBG superstructures than that of the RD-IBGs over the same fsp3 range.Our findings demonstrate that the thermal conductivity of these 2D graphene-diamond nanocomposite materials is precisely tunable by controlling the 2D diamond fraction in the nanocomposite (interlayer C–C bond density), which makes these 2D materials particularly appealing for thermal management applications, including their use as heat sinks in nanoelectronic devices.
See the supplementary material for a detailed discussion of EMD computations of lattice thermal conductivity in LAMMPS using the GK formula and comparisons of such EMD predictions for monolayer and bilayer graphene sheets with those of non-equilibrium molecular-dynamics (NEMD) simulations.We acknowledge financial support by the Army Research Laboratory under Grant No. HQ00341520007 and the usage of the facilities of the Massachusetts Green High-Performance Computing Center (MGHPCC).
Conflict of Interest
The authors have no conflicts to disclose.
Author Contributions
Afnan Mostafa: Formal analysis (equal); Investigation (equal); Methodology (equal); Software (equal); Visualization (equal); Writing – original draft (equal). Ashwin Ramasubramaniam: Formal analysis (equal); Funding acquisition (equal); Investigation (equal); Methodology (equal); Project administration (equal); Resources (equal); Software (equal); Supervision (equal); Writing – review & editing (equal). Dimitrios Maroudas: Conceptualization (lead); Data curation (lead); Formal analysis (equal); Funding acquisition (equal); Investigation (equal); Methodology (equal); Project administration (equal); Resources (equal); Software (equal); Supervision (equal); Validation (equal); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal).
The data that support the findings of this study are available from the corresponding author upon reasonable request.
REFERENCES
1. A. K. Geim, Science 324, 1530–1534 (2009). https://doi.org/10.1126/science.1158877,
留言 (0)