^{1}, Naoki Karasawa

^{2}and William A. Goddard III

^{1,a)}

### Abstract

As assemblies of graphene sheets, carbon nanotubes, and fullerenes become components of new nanotechnologies, it is important to be able to predict the structures and properties of these systems. A problem has been that the level of quantum mechanics practical for such systems (density functional theory at the PBE level) cannot describe the London dispersion forces responsible for interaction of the graphene planes (thus graphite falls apart into graphene sheets). To provide a basis for describing these London interactions, we derive the quantum mechanics based force field for carbon (QMFF-Cx) by fitting to results from density functional theory calculations at the M06-2X level, which demonstrates accuracies for a broad class of molecules at short and medium range intermolecular distances. We carried out calculations on the dehydrogenated coronene (C24) dimer, emphasizing two geometries: parallel-displaced X (close to the observed structure in graphite crystal) and PD-Y (the lowest energy transition state for sliding graphene sheets with respect to each other). A third, eclipsed geometry is calculated to be much higher in energy. The QMFF-Cx force field leads to accurate predictions of available experimental mechanical and thermodynamics data of graphite (lattice vibrations, elastic constants, Poisson ratios, lattice modes, phonon dispersion curves, specific heat, and thermal expansion). This validates the use of M06-2X as a practical method for development of new first principles based generations of QMFF force fields.

This research was partially supported by the Functional Engineered Nano Architects (FENA) via the Microelectronics Advanced Research Corporation (MARCO) with the prime award (2009–NT-2048) at UCLA (PI Kang Wang). T.A.P. thanks the National Science Foundation and the U.S. Department of Energy (CSGF) for graduate fellowships. W.A.G. acknowledges the WCU program at the Korea Advanced Institute of Science and Technology (KAIST) for additional support.

I. INTRODUCTION

II. EXPERIMENTAL DATA FOR GRAPHITE

A. Crystal structure

B. Elastic constants

C. Lattice vibrations

D. Thermal properties

III. COMPUTATIONAL DETAILS

A. QM calculations on the DHC dimer

B. The FF

1. vdW terms

2. Valence terms

C. Fitting procedure

D. Thermodynamic properties

E. Thermal expansion

IV. RESULTS

A. Structure and energies of DHC dimer

B. QMFF-Cx parameters

1. Comparison to experiment

2. Comparison with other two-body vdW parameters

3. Graphene-graphene interactions

C. Mechanical properties of graphite

1. C44 shear stiffness

2. C13 stiffness

3. Lattice modes

D. Phonon dispersion curve

E. Thermodynamic properties of graphite

1. Surface energy

2. Specific heat capacity

3. Lattice parameters and thermal expansion

F. Rhombohedral graphite

V. DISCUSSION

VI. CONCLUDING REMARKS

### Key Topics

- Graphite
- 80.0
- Lattice constants
- 22.0
- Elasticity
- 21.0
- Graphene
- 18.0
- Heat capacity
- 16.0

## Figures

(a) DHC unit used in this study. Since coronene has armchairlike edges, there are no dangling bonds: Each empty -orbital forms a bond with its neighbor. All carbon atoms are equivalent, making this molecule a reliable model for bulk graphite. (b) The optimized PD-X configuration from level of DFT. This is the most stable of the three configurations tested and is most similar to bulk graphite . (c) The optimized PD-Y configuration, which is 0.66 kcal/mol higher than PD-X and represents a saddle point for sliding two DHC units. Our optimized geometry is similar to QM studies on the CD. (d) The high energy eclipsed structure is 6.05 kcal/mol higher than PD-X.

(a) DHC unit used in this study. Since coronene has armchairlike edges, there are no dangling bonds: Each empty -orbital forms a bond with its neighbor. All carbon atoms are equivalent, making this molecule a reliable model for bulk graphite. (b) The optimized PD-X configuration from level of DFT. This is the most stable of the three configurations tested and is most similar to bulk graphite . (c) The optimized PD-Y configuration, which is 0.66 kcal/mol higher than PD-X and represents a saddle point for sliding two DHC units. Our optimized geometry is similar to QM studies on the CD. (d) The high energy eclipsed structure is 6.05 kcal/mol higher than PD-X.

Comparison of basis set for the counterpoise corrected binding energy of the PD-X configuration of DHC at the MO6-2X level DFT theory. Data are smooth using cubic splines (dashed lines) for presentation purposes. Single point energies in the various basis set were calculated using the minimum energy structure and the counterpoise corrections are then estimated. The basis set was selected as a good compromise of speed (1056 basis functions) and accuracy.

Comparison of basis set for the counterpoise corrected binding energy of the PD-X configuration of DHC at the MO6-2X level DFT theory. Data are smooth using cubic splines (dashed lines) for presentation purposes. Single point energies in the various basis set were calculated using the minimum energy structure and the counterpoise corrections are then estimated. The basis set was selected as a good compromise of speed (1056 basis functions) and accuracy.

Comparison of the QMFF-Cx and M06-2X vdW curves for the DHC (a) PD-X and (b) eclipsed geometries. The bottom of the PD-X curve was used in the fitting for all the QMFF-Cx potentials, while the PD-Y energy was used to fit the X6 and Morse potentials. The high energy, eclipsed geometry is well described only for the X6S (stretched X6) potential, which reproduces the QM value distance (3.565 Å) and energy (−10.68 kcal/mol). Inset: examination of equilibrium positions (3.4–3.8 Å) for the both geometries.

Comparison of the QMFF-Cx and M06-2X vdW curves for the DHC (a) PD-X and (b) eclipsed geometries. The bottom of the PD-X curve was used in the fitting for all the QMFF-Cx potentials, while the PD-Y energy was used to fit the X6 and Morse potentials. The high energy, eclipsed geometry is well described only for the X6S (stretched X6) potential, which reproduces the QM value distance (3.565 Å) and energy (−10.68 kcal/mol). Inset: examination of equilibrium positions (3.4–3.8 Å) for the both geometries.

PECs for the four potentials considered in this study. The parameters were all determined from the interaction of the DHC dimer from M06-2X DFT theory. The X6 and Morse potentials lead to predictions in excellent agreement with the QM for the PD-X and PD-Y geometrics, but are not able to reproduce the energetics of the high energy DHC eclipsed structure. The X6S potential has an additional parameter relating to the inner wall curvature parameter of 26.9 (compared to 20.7 for X6) and is able to reproduce the high energy structure.

PECs for the four potentials considered in this study. The parameters were all determined from the interaction of the DHC dimer from M06-2X DFT theory. The X6 and Morse potentials lead to predictions in excellent agreement with the QM for the PD-X and PD-Y geometrics, but are not able to reproduce the energetics of the high energy DHC eclipsed structure. The X6S potential has an additional parameter relating to the inner wall curvature parameter of 26.9 (compared to 20.7 for X6) and is able to reproduce the high energy structure.

(a) PES and (b) contour plot for sliding a periodic graphene sheet (96 atoms) over another in 0.1 Å displacements, using the X6 potential. The interplanar distance was optimized at each displacement. All energies are referenced to the PD-X structure, the global minima. The PD-Y structure is a low energy minima . The PD-X PD-Y barrier is 0.0156 kcal/mol C and a low energy pathway for sliding is obtained by tracing the edges of the hexagon unit. The AA stacked eclipsed graphite structure is high energy .

(a) PES and (b) contour plot for sliding a periodic graphene sheet (96 atoms) over another in 0.1 Å displacements, using the X6 potential. The interplanar distance was optimized at each displacement. All energies are referenced to the PD-X structure, the global minima. The PD-Y structure is a low energy minima . The PD-X PD-Y barrier is 0.0156 kcal/mol C and a low energy pathway for sliding is obtained by tracing the edges of the hexagon unit. The AA stacked eclipsed graphite structure is high energy .

Phonon dispersion curve for all vibrational modes of (hexagonal) graphite at 0 K using the X6 potential.

Phonon dispersion curve for all vibrational modes of (hexagonal) graphite at 0 K using the X6 potential.

Phonon dispersion curves for the low frequency modes of (hexagonal) graphite at 300 K. Solid lines from theory and symbols from experimental data (Refs. 18, 19, 42, and 77).

Phonon dispersion curves for the low frequency modes of (hexagonal) graphite at 300 K. Solid lines from theory and symbols from experimental data (Refs. 18, 19, 42, and 77).

Specific heat of hexagonal graphite as computed with the QM-FF X6 potential. Low temperature results obtained from the thin plate approximation, other results obtained from the uniform grid method. Experimental results from different sources are indicated, as reported in Refs. 27 and 46. Values for rhombohedral graphite are not plotted since the lines would be essentially superimposed on the hexagonal graphite lines.

Specific heat of hexagonal graphite as computed with the QM-FF X6 potential. Low temperature results obtained from the thin plate approximation, other results obtained from the uniform grid method. Experimental results from different sources are indicated, as reported in Refs. 27 and 46. Values for rhombohedral graphite are not plotted since the lines would be essentially superimposed on the hexagonal graphite lines.

(a) In-plane lattice parameter of graphite as a function of temperature, calculated with the QM-FF X6 potential. The experimental results (Ref. 52) (red squares) are compared the calculated values (black triangles). The solid black line is the least-squares line to the FF using cubic spline regression. (b) Out of plane lattice parameter of graphite as a function of temperature, calculated with the QM-FF X6 potential. The experimental results (Ref. 52) (red squares) are compared the calculated values (black triangles). The solid black line is the least-squares line by cubic spline regression.

(a) In-plane lattice parameter of graphite as a function of temperature, calculated with the QM-FF X6 potential. The experimental results (Ref. 52) (red squares) are compared the calculated values (black triangles). The solid black line is the least-squares line to the FF using cubic spline regression. (b) Out of plane lattice parameter of graphite as a function of temperature, calculated with the QM-FF X6 potential. The experimental results (Ref. 52) (red squares) are compared the calculated values (black triangles). The solid black line is the least-squares line by cubic spline regression.

## Tables

Comparison of interaction energies (kcal/mol) for PD-X, PD-Y, and eclipsed DHC dimer structures between QM (counterpoise corrected and uncorrected ) and QMFF-Cx. These can be compared to similar numbers on the CD from DFT and MP2 calculations. PD-X is topologically similar to the graphite structure, while PD-Y is a saddle point corresponding to the barrier between two adjacent PD-X minima. The interplane displacements for the QM structures are quoted in brackets.

Comparison of interaction energies (kcal/mol) for PD-X, PD-Y, and eclipsed DHC dimer structures between QM (counterpoise corrected and uncorrected ) and QMFF-Cx. These can be compared to similar numbers on the CD from DFT and MP2 calculations. PD-X is topologically similar to the graphite structure, while PD-Y is a saddle point corresponding to the barrier between two adjacent PD-X minima. The interplane displacements for the QM structures are quoted in brackets.

Optimized FF parameters for graphite at 0 K.

Optimized FF parameters for graphite at 0 K.

Calculated properties for graphite at 0 K (numbers in parenthesis indicate values at 300 K). The , and , experimental lattice modes are used in the FF optimization.

Calculated properties for graphite at 0 K (numbers in parenthesis indicate values at 300 K). The , and , experimental lattice modes are used in the FF optimization.

Comparison of QMFF-Cx X6 with published vdW parameters for carbon. In each, the valence parameters are adjusted to match the correct experimental in-plane (a) lattice parameter. Errors relative to the experimental value are indicated in brackets for properties sensitive to the vdW.

Comparison of QMFF-Cx X6 with published vdW parameters for carbon. In each, the valence parameters are adjusted to match the correct experimental in-plane (a) lattice parameter. Errors relative to the experimental value are indicated in brackets for properties sensitive to the vdW.

Phonon frequencies of graphite and derivatives at the high-symmetry points , , , and in . The X6 (0 K) is compared to results from *ab initio* DFT studies and experiment (300 K).

Phonon frequencies of graphite and derivatives at the high-symmetry points , , , and in . The X6 (0 K) is compared to results from *ab initio* DFT studies and experiment (300 K).

Article metrics loading...

Full text loading...

Commenting has been disabled for this content