Elasticity of Diamond at High Pressures and Temperatures

Elasticity of Diamond at High Pressures and Temperatures

Maribel Núñez Valdez nunez025@umn.edu Dept. of Chemical Engineering and Materials Science, University of Minnesota, Minneapolis, Minnesota 55455, USA    Koichiro Umemoto Dept. of Geology and Geophysics, University of Minnesota, Minneapolis, Minnesota 55455, USA    Renata M. Wentzcovitch Dept. of Chemical Engineering and Materials Science, University of Minnesota, Minneapolis, Minnesota 55455, USA Minnesota Supercomputing Institute (MSI), University of Minnesota, Minneapolis, Minnesota 55455, USA
July 4, 2019

We combine density functional theory within the local density approximation, the quasiharmonic approximation, and vibrational density of states to calculate single crystal elastic constants, and bulk and shear moduli of diamond at simultaneous high pressures and temperatures in the ranges of GPa and K. Comparison with experimental values at ambient pressure and high temperature shows an excellent agreement for the first time with our first-principles results validating our method. We show that the anisotropy factor of diamond increases to 40% at high pressures and becomes temperature independent.

diamond, first-principles, elasticity, high-pressure, high-temperature
preprint: AIP/123-QED

Diamond has singular physical properties. Its strength and transparency make it remarkably useful for industrial applications in abrasive machining, micro electromechanical devices etc. Kohn, Gluche, and Adamschikr (1999) It is also an anvil material for experimental studies under extreme pressures and temperatures. Pressures greater to that reached at the center of the Earth ( GPa) can be produced in diamond anvil cells (DACs) using single crystal (SC) diamonds as opposing anvils. Applications of diamond to geoscience and planetary science are ideal to study high pressure-temperature behavior of magmas and minerals relevant to the Earth and other planets. For example, in 2004 a high pressure form of (Fe,Mg)SiO, postperovskite, was discovered using DAC techniques Murakami et al. (2004), this mineral exists above the mantle-core boundary of the Earth.

Natural diamond can be SC or polycrystalline (PC). The high quality SC diamond is optically transparent and often used as gems. It also possesses higher strength than the PC type, which is more opaque or translucent with impurities and inclusions. PC diamond, however, is ideal for cutting tools due to its toughness.

Thin films of PC diamond with very low growth rate (m/h) were synthesized and reported in the 1960s using chemical vapor deposition (CVD), but those films were of poor quality. SC diamond synthesis is performed near graphite-diamond phase boundary conditions ( GPa and K) within large-volume presses Shigley, Abbaschian, and Clarke (2002). The main limitations to this production are the size of the press, yielding only diamonds of some millimeters and a very slow growth rate Sumiyaa, Toda, and Satoh (2002). However, in recent years SCs of diamonds larger than 2 cm and more than 10 carats have been successfully fabricated using microwave plasma CVD techniques. The purest CVD diamond was reported to have comparable features to those of high purity natural SC with hardness and toughness that can be tuned by controlling growth and annealing Liang et al. (2009); Meng et al. (2012). In 2003 in Japan, using large-volume presses Irifune et al. (2003), graphite was directly converted into diamond to form nanopolycrystalline diamond (NPD) at simultaneous high pressure and temperature ( GPa and K) within minutes. NPD was found to have high elastic stiffness, very high fracture toughness and similar or higher hardness than that of most single crystal diamonds. These features and the fact that NPD rods can be fabricated up to 1 cm in diameter and length (14.5 carats), make NPD very suitable for larger anvils than the usually available in DACs. Larger sample volume in high pressure apparatuses has open new windows into mineral physics studies of materials that require accurate and simultaneous very high pressures and temperatures such as those found deep in the Earth and other planets Ohfuji et al. (2010).

With these new techniques to fabricate diamonds and their increasing applications, accurate elasticity of diamond at simultaneous high pressures and temperatures acquires new significance. However, it cannot be found in the literature in neither experimental form nor computational modeling. Using Brillouin scattering, elastic moduli of diamond were measured as a function of temperature in the range of 300–1600 K at ambient pressure Zouboulis et al. (1998). A molecular dynamics study attempted the computation of elastic constants between 100 and 1100 K with a modest degree of success Gao et al. (2006). Single crystal elasticity and bulk modulus of diamond were computed as a function of pressure at 300 K in the range of 0–500 GPa Fu et al. (2009), using first-principles methods based on density functional theory (DFT) Hohenberg and Kohn (1964); Kohn and Sham (1965) as implemented in the CRYSTAL03 (CR03) program Saunders, Dovesi, and et al. (2003) and a qusiharmonic Debye model (QHDM) Blanco, Francisco, and Luana (2004). But the lack of experimental and other computational values for such quantities made impossible a comparison for validation of these CRYSTAL03 results.

In this letter, for the first time we address the elasticity of single crystal diamond and of a polycrystalline aggregate at simultaneous high pressures and temperatures. The computational approach is based on DFT Hohenberg and Kohn (1964); Kohn and Sham (1965) within the local density approximation (LDA) Ceperley and Alder (1980). At arbitrary pressures, cubic equilibrium structures of diamond (space group Fd3m, two C atoms/primitive cell) were found using the variable cell-shape damped molecular dynamics approach Wentzcovitch (1991); Wentzcovitch, Martins, and Price (1993) as implemented in the quantum-ESPRESSO (QE) code Giannozzi et al. (2009). The C pseudopotential was of Vanderbilt-type with cutoff radii of 1.3 a.u., and partial core correction. We used a plane-wave kinetic energy cutoff of 50 Ry and 200 Ry for the charge density. The k-point sampling for charge density was determined on a Monkhorst-Pack grid of the Brillouin Zone (BZ) shifted by from the origin. These parameters corresponded to having interatomic forces and stresses smaller than Ry/a.u.

To obtain the independent cubic static elastic constants , , and (in Voigt notation Wallace (1972)) at each pressure, the equilibrium structures were strained, and their internal degrees of freedom re-relaxed. Then elastic constants were extracted using the relationship . Positive and negative strains of 1% magnitude were applied in order to attain accuracy in the limit of zero strain. This method has been applied successfully for over a decade Wentzcovitch et al. (1995); Karki et al. (1997); NúñezValdez, Umemoto, and Wentzcovitch (2010); NúñezValdez, da Silveira, and Wentzcovitch (2011).

Dynamical matrices were calculated using density functional perturbation theory (DFPT) Baroni et al. (2001). At each pressure, a dynamical matrix was obtained on a q-point mesh. Force constants were extracted and interpolated on a mesh to produce VDoS.

Figure 1: (Color online) Temperature dependence of (a) ; (b) ; (c) ; and (d) bulk, , and shear, , moduli. Dashed lines indicate a lesser degree of confidence within the QHA limit.
Reference (GPa) (GPa) (GPa) (GPa) (GPa)
This work LDA 1078.99 140.81 577.37 453.54 531.31
GGA Mounet and Marzari (2005) 1060 125 562 436.67 522.08
MD Gao et al. (2006) 963.6-982.2 148.4-165.8 674.5-694.0 425.76-437.93 550.12-565.99
B3PW Fu et al. (2009) 1097.5 115.5 598.2 442.8 552.74
Exp McSkimin and P. Andreatch (1972) 1079 124 578 442 535
Exp Grimsditch and Ramdas (1975) 1076.4 0.2 125.2 2.3 577.4 1.4 442.3 534.27
Exp Zouboulis et al. (1998) 1080.4 127.0 576.6 444.8 534.32
Table 1: Comparison of elastic mouduli of diamond at ambient conditions (=0 GPa and = 300K)

Then we exploited the dependence of phonon frequencies on compression to determine the thermal contribution to the Helmholtz free energy , within the QHA Wallace (1972),


where is the phonon wave vector, is the normal mode index, is temperature, is the static internal energy at equilibrium volume under isotropic pressure and infinitesimal strain , and are Planck and Boltzmann constants, respectively. Isothermal elastic constants are given by with , the Gibbs energy and . To convert to adiabatic elastic constants, we use , (heat capacity at constant , and entropy). have the static- and phonon- contributions, the latter can be expressed as function of strain and mode Grünesein parameters. Assuming that the angular distribution of phonons is isotropic, isothermal elastic constants can be calculated without performing phonon calculations for strained configurations Wu and Wentzcovitch (2011). This approximation is equivalent to assuming that thermal pressure is isotropic. Even though this assumption is not completely rigorousCarrier, Wentzcovitch, and Tsuchiya (2007), the method is quite accurate compared to experimental uncertainties Carrier, Justo, and Wentzcovitch (2008). Bulk and shear moduli were calculated as Voigt-Reuss-Hill averages Watt (1979) using the adiabatic , i.e., and , where, , and .

The computed elastic constants , , and , and bulk and shear , moduli as a function of temperature are shown in Figure 1. We compared our results with experimental trends in the temperature range of K at ambient pressure and we found excellent agreement. The previous molecular dynamics calculation Gao et al. (2006) was quite unsuccessful at reproducing the experimental values for individually. The quality of our findings with respect to measurements validates our overall computational approach.

Then we obtained the elasticity as a function of pressure between 0 and 500 GPa at several temperatures, (Figure 2). Unfortunately experimental values are only available at ambient conditions for comparison. The only other calculation of elasticity as a function of pressure Fu et al. (2009) reported using CR03 and the QHDM was compared only to measurements at room temperature. Figure 2(a) shows that CR03 results are quite different from ours. There is agreement in the bulk modulus trend [Figure 2(b)], but that probably reflects the fact that the compression curves are similar in both studies.

Since the QHA is valid up to about two thirds of the melting temperature (for diamond K at GPa and it increases in the range GPa Correa, Bonev, and Galli (2006)), dashed-lines in Figures 1 and 2 indicate a lesser degree of confidence for the QHA results Wu (2010). Table 1 compares results on the elasticity of diamond at ambient conditions with the other DFT-based values and experimental measurements. One can see clearly the improvement in the computation of , , and with respect to the other computational methods. Our computed bulk modulus is 2.5% larger than experimental values. This is typical of LDA results that underestimate volume by % and reflects on only.

Figure 2: (Color online) Pressure dependence of (a) , , and ; and (b) bulk, , and shear, , moduli. Dashed lines indicate a lesser degree of confidence within the QHA limit.

Finally, we solved the Christoffel equation (Landau and Lifshitz, 1970) for acoustic velocities in single crystal, , where is elastic wave velocity, is propagation direction, is Kronecker delta, and is density. For a given , there are three solutions, i.e., one P-wave () and two S-waves (, ). Figures 3(a-c) show the variation with crystallographic direction (anisotropy) of , and at 0 and 100 GPa and 300 and 1000 K for comparison. We found similar qualitative behavior independent of pressure and temperature for the three waves. In Figure 3(d) we show the anisotropic factor of diamond , as a function of at several s. One can observe that becomes roughly independent of and . often decreases with Karki et al. (1999), but this is not the case for diamond. As seen in Figures 3(a-c), velocities vary with direction much more at high pressures.

Figure 3: (Color online) Variation of (a) Compressional velocity and (b,c) Shear velocities as a function of direction; and (d) Anisotropy factor of diamond

In conclusion, we have computed the high pressure and temperature elasticity of diamond by first principles using a novel analytical approach. In this new powerful method, vibrational density of states for strained configurations are unnecessary (Wu and Wentzcovitch, 2011). Treatment of strain Grüneisen parameters via isotropic averages reduced greatly the computational cost of the task, which, otherwise, would have been much more demanding. We have shown an outstanding improvement in the calculated elasticity of diamond as a function of temperature at ambient pressure (Figure 1). These results increase the confidence in our approach to calculate high temperature elasticity (Figure 2). In general, accurate elasticity of diamond at high pressures and temperatures is a step forward in improving and opening applications of diamonds in science and technology, particularly in geophysics and planetary sciences where extreme conditions of pressure and temperature are typical. The elastic anisotropy of diamond is shown to increase considerably at high pressures. This result is significant for understanding the performance under pressure of large SC and PC diamonds specimens nowadays synthesized by CVD or as NPDs.

We acknowledge the support of grant NSF/EAR 1047629. Computations were performed at the MSI. The authors also thank Z. Wu for helpful discussions.


  • Kohn, Gluche, and Adamschikr (1999) E. Kohn, P. Gluche,  and M. Adamschikr, Diamond and Rel. Mat. 8, 934 (1999).
  • Murakami et al. (2004) M. Murakami, K. Hirose, K. Kawamura, N. Sata,  and Y. Ohishi, Science 304, 442 (2004).
  • Shigley, Abbaschian, and Clarke (2002) J. E. Shigley, R. Abbaschian,  and C. Clarke, Gems and Gemology 38, 301 (2002).
  • Sumiyaa, Toda, and Satoh (2002) H. Sumiyaa, N. Toda,  and S. Satoh, J. Cryst. Growth 237, 1281 (2002).
  • Liang et al. (2009) Q. Liang, C. S. Yan, Y. Meng, S. K. L. Lai, H. K. Mao,  and R. J. Hemley, Diamond Related Mater. 18, 698 (2009).
  • Meng et al. (2012) Y. F. Meng, C. S. Yan, S. Krasnicki, Q. Liang, J. Lai, Y. Y. H. Shu, A. S. Steele, H. K. Mao,  and R. J. Hemley, Phys. Status Solidi A. 209(1), 101 (2012).
  • Irifune et al. (2003) T. Irifune, A. Kuio, S. Sakamoto, T. Inoue,  and H. Sumiya, Nature 421, 599 (2003).
  • Ohfuji et al. (2010) H. Ohfuji, T. Okada, T. Yagi, H. Sumiya, ,  and T. Irifune, High Pressure Res. 30, 142 (2010).
  • Zouboulis et al. (1998) E. S. Zouboulis, M. Grimsditch, A. K. Ramdas,  and S. Rodriguez, Phys. Rev. B 57, 2889 (1998).
  • Gao et al. (2006) G. Gao, K. V. Workum, J. D. Schall,  and J. A. Harrison, J. Phys. Cond. Matter 18, S1737 (2006).
  • Fu et al. (2009) Z. J. Fu, G. F. Ji, X. R. Chen,  and Q. Q. Gou, Commun. Theor. Phys. 51, 1129 (2009).
  • Hohenberg and Kohn (1964) P. Hohenberg and W. Kohn, Phys. Rev. 136, B864 (1964).
  • Kohn and Sham (1965) W. Kohn and L. J. Sham, Phys. Rev. 140, A1133 (1965).
  • Saunders, Dovesi, and et al. (2003) V. R. Saunders, R. Dovesi,  and C. R. et al.CRYSTAL03 User’s Manual (University of Torino, 2003).
  • Blanco, Francisco, and Luana (2004) M. A. Blanco, E. Francisco,  and V. Luana, Comp. Phys. Commun. 158, 57 (2004).
  • Ceperley and Alder (1980) D. M. Ceperley and B. J. Alder, Phys. Rev. Lett. 45, 566 (1980).
  • Wentzcovitch (1991) R. M. Wentzcovitch, Phys. Rev. B 44, 2358 (1991).
  • Wentzcovitch, Martins, and Price (1993) R. M. Wentzcovitch, J. L. Martins,  and G. D. Price, Phys. Rev. Lett. 70, 3947 (1993).
  • Giannozzi et al. (2009) P. Giannozzi, S. Baroni, N. Bonini, M. Calandra, R. Car, C. Cavazzoni, D. Ceresoli, G. L. Chiarotti, M. Cococcioni, I. Dabo,  and A. D. et al., J. Phys.: Condens. Matter 21, 395502 (2009).
  • Wallace (1972) D. C. Wallace, Thermodynamics of Crystals (John Wiley& Sons, Inc., USA, 1972).
  • Wentzcovitch et al. (1995) R. M. Wentzcovitch, D. A. Hugh-Jones, R. J. Angel,  and G. D. Price, Phys. Chem. Minerals 22, 453 (1995).
  • Karki et al. (1997) B. B. Karki, L. Stixrude, S. Clark, M. C. Warren, G. J. Ackland,  and J. Crain, Am. Mineral. 82, 635 (1997).
  • NúñezValdez, Umemoto, and Wentzcovitch (2010) M. NúñezValdez, K. Umemoto,  and R. M. Wentzcovitch, Geophys. Res. Lett. 37, L14308 (2010).
  • NúñezValdez, da Silveira, and Wentzcovitch (2011) M. NúñezValdez, P. da Silveira,  and R. M. Wentzcovitch, J. Geophys. Res. 116, B12207 (2011).
  • Baroni et al. (2001) S. Baroni, A. D. Corso, S. de Gironcoli,  and P. Gianozzi, Rev. of Mod. Phys. 73(2), 515 (2001).
  • Mounet and Marzari (2005) N. Mounet and N. Marzari, Phys. Rev. B 71, 205214 (2005).
  • McSkimin and P. Andreatch (1972) H. J. McSkimin and J. P. Andreatch, J. Appl. Phys. 43(7), 2944 (1972).
  • Grimsditch and Ramdas (1975) M. H. Grimsditch and A. K. Ramdas, Phys. Rev. B 11(8), 3139 (1975).
  • Wu and Wentzcovitch (2011) Z. Wu and R. M. Wentzcovitch, Phys. Rev. B 83, 184115 (2011).
  • Carrier, Wentzcovitch, and Tsuchiya (2007) P. Carrier, R. M. Wentzcovitch,  and J. Tsuchiya, Phys. Rev. B 87, 064116 (2007).
  • Carrier, Justo, and Wentzcovitch (2008) P. Carrier, J. F. Justo,  and R. M. Wentzcovitch, Phys. Rev. B 78, 144302 (2008).
  • Watt (1979) J. P. Watt, J. Appl. Phys., 50, 6290 (1979).
  • Correa, Bonev, and Galli (2006) A. A. Correa, S. A. Bonev,  and G. Galli, PNAS 103, 1204 (2006).
  • Wu (2010) Z. Wu, Phys. Rev. B 81, 172301 (2010).
  • Landau and Lifshitz (1970) L. D. Landau and E. M. Lifshitz, Theory of Elasticity (Pergamon Press, UK, 1970).
  • Karki et al. (1999) B. B. Karki, R. M. Wentzcovitch, S. de Gironcoli,  and S. Baroni, Science 286, 1705 (1999).
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description