Free Submitted: 14 June 2018 Accepted: 19 July 2018 Published Online: 02 August 2018
J. Chem. Phys. 149, 051101 (2018); https://doi.org/10.1063/1.5044636

The dielectric properties of molecules and nanostructures are usually modified in a complex manner, when assembled into a condensed phase. We propose a first-principles method to compute polarizabilities of sub-entities of solids and liquids, which accounts for multipolar interactions at all orders and is applicable to semiconductors and insulators. The method only requires the evaluation of induced fields in the condensed phase, with no need of multiple calculations for each constituent. As an example, we present results for the molecular polarizabilities of water in a wide pressure and temperature range. We found that at ambient conditions, the dipole-induced-dipole approximation is sufficiently accurate and the Clausius-Mossotti relation may be used, e.g., to obtain molecular polarizabilities from experimental refractive indexes. However with increasing pressure, this approximation becomes unreliable and in the case of ice X the Clausius-Mossotti relation is not valid.
The polarizability of molecules and nano-structures is an important property determining the assembly and behavior of solids and liquids composed of well-defined building blocks (BBs). In addition, molecular polarizabilities play a key role in vibrational spectroscopy, e.g., in the calculation of Raman11. A. Putrino and M. Parrinello, Phys. Rev. Lett. 88, 176401 (2002). https://doi.org/10.1103/physrevlett.88.176401 and sum-frequency generation22. Q. Wan and G. Galli, Phys. Rev. Lett. 115, 246404 (2015). https://doi.org/10.1103/physrevlett.115.246404 spectra, in the determination of van der Waals interactions in solids and liquids,3,43. G. D. Mahan, J. Chem. Phys. 76, 493 (1982). https://doi.org/10.1063/1.4427504. J. Klimeš and A. Michaelides, J. Chem. Phys. 137, 120901 (2012). https://doi.org/10.1063/1.4754130 and in the development of polarizable force fields.
While several methods are available to compute and predict polarizabilities of isolated molecules and nanostructures, their definition and calculation in condensed phases have been challenging and various levels of approximations have been adopted in the literature (e.g., Refs. 5–95. R. P. Feynman, R. B. Leighton, and M. Sands, Feynman Lectures on Physics. Volume 2: The New Millennium Edition: Mainly Electromagnetism and Matter (Basic Books, 2011), Vol. 2.6. C. Kittel, Introduction to Solid State Physics, 8th ed. (Wiley, 2004).7. R. J. Heaton, P. A. Madden, S. J. Clark, and S. Jahn, J. Chem. Phys. 125, 144104 (2006). https://doi.org/10.1063/1.23571518. M. Salanne, R. Vuilleumier, P. A. Madden, C. Simon, P. Turq, and B. Guillot, J. Phys.: Condens. Matter 20, 494207 (2008). https://doi.org/10.1088/0953-8984/20/49/4942079. Q. Wan, L. Spanu, G. A. Galli, and F. Gygi, J. Chem. Theory Comput. 9, 4124 (2013). https://doi.org/10.1021/ct4005307). For example, the Clausius-Mossotti (CM) equation relates the average atomic or molecular polarizability α of a material building block to its electronic dielectric constant ϵ,5,65. R. P. Feynman, R. B. Leighton, and M. Sands, Feynman Lectures on Physics. Volume 2: The New Millennium Edition: Mainly Electromagnetism and Matter (Basic Books, 2011), Vol. 2.6. C. Kittel, Introduction to Solid State Physics, 8th ed. (Wiley, 2004).
 $4πNα3=ϵ∞−1ϵ∞+2,$ (1)
where N is the number density of atoms or molecules. In Eq. (1), if we substitute ϵ with the refractive index of the material n using ϵ = n2, the Lorentz-Lorentz equation is recovered. The validity of the CM relation in condensed phases, such as molecular liquids or assembly of nanostructured solids, depends on the system and general rules to establish its regime of applicability are not available. We note that often times, the variation of polarizabilities from the gas to the condensed phase are neglected. For example, many molecular dynamics simulations of aqueous solutions using force fields assume a fixed molecular polarizability of water,1010. B. Schropp and P. Tavan, J. Phys. Chem. B 112, 6233 (2008). https://doi.org/10.1021/jp0757356 though first-principles electronic structure studies of water at ambient conditions have shown that molecular polarizabilities have rather broad distributions and are not isotropic.7–97. R. J. Heaton, P. A. Madden, S. J. Clark, and S. Jahn, J. Chem. Phys. 125, 144104 (2006). https://doi.org/10.1063/1.23571518. M. Salanne, R. Vuilleumier, P. A. Madden, C. Simon, P. Turq, and B. Guillot, J. Phys.: Condens. Matter 20, 494207 (2008). https://doi.org/10.1088/0953-8984/20/49/4942079. Q. Wan, L. Spanu, G. A. Galli, and F. Gygi, J. Chem. Theory Comput. 9, 4124 (2013). https://doi.org/10.1021/ct4005307 Most electronic structure methods7–97. R. J. Heaton, P. A. Madden, S. J. Clark, and S. Jahn, J. Chem. Phys. 125, 144104 (2006). https://doi.org/10.1063/1.23571518. M. Salanne, R. Vuilleumier, P. A. Madden, C. Simon, P. Turq, and B. Guillot, J. Phys.: Condens. Matter 20, 494207 (2008). https://doi.org/10.1088/0953-8984/20/49/4942079. Q. Wan, L. Spanu, G. A. Galli, and F. Gygi, J. Chem. Theory Comput. 9, 4124 (2013). https://doi.org/10.1021/ct4005307 consider only the dipole-dipole interaction when calculating the variation of polarizabilities upon assembly of molecular fluids or solids, and in many cases, such approximations have remained untested. Recently, substantial progress has been reported in computing polarizabilities of building blocks using maximally localized Wannier functions (MLWFs).1111. N. Marzari, A. A. Mostofi, J. R. Yates, I. Souza, and D. Vanderbilt, Rev. Mod. Phys. 84, 1419 (2012). https://doi.org/10.1103/revmodphys.84.1419 However, the method of Ref. 1111. N. Marzari, A. A. Mostofi, J. R. Yates, I. Souza, and D. Vanderbilt, Rev. Mod. Phys. 84, 1419 (2012). https://doi.org/10.1103/revmodphys.84.1419 requires separate calculations of the dielectric properties of each constituent self-consistently, and it is not based on global induced fields within the condensed system.12,1312. X. Ge and D. Lu, Phys. Rev. B 92, 241107 (2015). https://doi.org/10.1103/physrevb.92.24110713. X. Ge and D. Lu, Phys. Rev. B 96, 075114 (2017). https://doi.org/10.1103/physrevb.96.075114
In this paper, we propose a first-principles method to compute the polarizabilities of building blocks in condensed phases. The method, based solely on electronic structure calculations for the condensed phase, is applicable to semiconductors and insulators (i.e., to any systems with a gap between occupied and virtual states such that maximally localized Wannier functions can be defined1111. N. Marzari, A. A. Mostofi, J. R. Yates, I. Souza, and D. Vanderbilt, Rev. Mod. Phys. 84, 1419 (2012). https://doi.org/10.1103/revmodphys.84.1419). We present results for the molecular polarizabilities of water in a wide pressure-temperature (P-T) range, and we validate the CM relation for water at ambient conditions and the dipole-induced-dipole approximation (DID). We found that the DID becomes increasingly less accurate under pressure and breaks down when covalent bonds are present and oxygen ions are formed within the solid.
We start by summarizing our formulation. A building block (BB) composing of a condensed system (e.g., a molecule in a molecular crystal) is defined by its ionic coordinates and by electronic wave functions spatially localized at the BB site, for example, maximally localized Wannier functions $wi(r→)$1111. N. Marzari, A. A. Mostofi, J. R. Yates, I. Souza, and D. Vanderbilt, Rev. Mod. Phys. 84, 1419 (2012). https://doi.org/10.1103/revmodphys.84.1419 constructed from the Bloch orbitals of the condensed phase. The linearly induced electron polarization density of the BB, in response to a macroscopic field $E→$ is
 $ΔρBB(r→)=2∑iNorbwi*(r→)Δwi(r→)+c.c.,$ (2)
where Norb is the number of localized electronic orbitals (e.g., four doubly-degenerate orbitals for a water molecule with 8 valence electrons) and Δ$wi$ is the variation of the ith Wannier function under the field $E→$. The local field ($E→loc$) acting on the BB is given by two contributions,
 $E→loc=E→+E→env,$ (3)
where $E→env$ denotes the field produced by the environment surrounding the BB, that is, by all the electrons that are not associated, through the definition of the $wi$, to the BB. In most previous studies, $E→env$ was approximated by dipole-induced-dipole (DID) electrostatic interactions.7–97. R. J. Heaton, P. A. Madden, S. J. Clark, and S. Jahn, J. Chem. Phys. 125, 144104 (2006). https://doi.org/10.1063/1.23571518. M. Salanne, R. Vuilleumier, P. A. Madden, C. Simon, P. Turq, and B. Guillot, J. Phys.: Condens. Matter 20, 494207 (2008). https://doi.org/10.1088/0953-8984/20/49/4942079. Q. Wan, L. Spanu, G. A. Galli, and F. Gygi, J. Chem. Theory Comput. 9, 4124 (2013). https://doi.org/10.1021/ct4005307
The polarizability tensor αBB of the BB is defined by the equation
 $μ→BB=αBBE→loc,$ (4)
where $μ→BB$ is the dipole moment computed as
 $μ→BB=−qe∫r→ΔρBB(r→)dr→.$ (5)
In Eq. (5), qe is the elementary charge and ΔρBB is from Eq. (2). To compute αBB, we need to calculate $E→loc$. Since $E→$ is fixed, we only need to determine $E→env$, which is simply
 $E→env=1Norb∑iNorb∫e′→(r→)wi(r→)2dr→ ,$ (6)
where $e′→(r→)$ is the microscopic electric field induced by all the electrons outside the BB.
The microscopic electric field $e′→$ can be evaluated within the random phase approximation (RPA) or by including the variation of the exchange and correlation potential [we denote the latter with density functional theory (DFT)]. Within RPA, $e′→RPA$ is obtained using Gauss’ law
 $∇⋅e′→RPA(r→)=−4πqeΔρ′(r→),$ (7)
where Δρ′ = Δρ − ΔρBB and Δρ is the electron polarization density of the whole system. Historically, it is from the wing part of the inverse dielectric matrix.27,2827. A. Baldereschi, R. Car, and E. Tosatti, Solid State Commun. 32, 757 (1979). https://doi.org/10.1016/0038-1098(79)90749-x28. R. Resta and A. Baldereschi, Phys. Rev. B 23, 6615 (1981). https://doi.org/10.1103/physrevb.23.6615 At the DFT level, the exchange-correlation (xc) potential Vxc also contributes to the microscopic local field,
 $e′→DFT(r→)=e′→RPA(r→)+1qe∇dVxcdρΔρ′(r→).$ (8)
Once Eloc and μ are computed from Eqs. (3) and (5), respectively, αBB is known. Therefore the procedure outlined here to obtain polarizabilities of BB within a condensed system is rather simple. Once the density and single particle wave functions are computed, e.g., by solving the Kohn-Sham equations, the electron polarization density is obtained by performing a single self-consistent calculation for the whole system. Equation (7) or (8) is then solved non-self-consistently, including multipole interactions at all orders. Solvers to obtain linear variations of the electron density exist in most DFT codes, using either density-functional perturbation theory (DFPT)1414. S. Baroni, S. De Gironcoli, A. Dal Corso, and P. Giannozzi, Rev. Mod. Phys. 73, 515 (2001). https://doi.org/10.1103/revmodphys.73.515 or finite fields.1515. I. Souza, J. Íñiguez, and D. Vanderbilt, Phys. Rev. Lett. 89, 117602 (2002). https://doi.org/10.1103/physrevlett.89.117602
We now turn to applying the method outlined above to the study of the molecular polarizabilities of water in a broad P-T range from ambient to supercritical conditions. The electron polarization density was obtained by DFPT, as implemented in the plane-wave pseudopotential code Qbox (http://qboxcode.org/).9,169. Q. Wan, L. Spanu, G. A. Galli, and F. Gygi, J. Chem. Theory Comput. 9, 4124 (2013). https://doi.org/10.1021/ct400530716. F. Gygi, IBM J. Res. Dev. 52, 137 (2008). https://doi.org/10.1147/rd.521.0137 We used Hamann-Schluter-Chiang-Vanderbilt norm-conserving pseudopotentials29,3029. D. R. Hamann, M. Schlüter, and C. Chiang, Phys. Rev. Lett. 43, 1494 (1979). https://doi.org/10.1103/physrevlett.43.149430. D. Vanderbilt, Phys. Rev. B 32, 8412 (1985). https://doi.org/10.1103/physrevb.32.8412 with a plane-wave kinetic energy cutoff of 85 Ry. The MD trajectories of water at ambient conditions were taken from the water PBE400 dataset:3131. W. Dawson and F. Gygi, J. Chem. Phys. 148, 124501 (2018). https://doi.org/10.1063/1.5018116 http://www.quantum-simulation.org, obtained using 64 water molecules in the simulation box. The supercritical water trajectories were from our previous simulations,18,1918. D. Pan, L. Spanu, B. Harrison, D. A. Sverjensky, and G. Galli, Proc. Natl. Acad. Sci. U. S. A. 110, 6646 (2013). https://doi.org/10.1073/pnas.122158111019. D. Pan, Q. Wan, and G. Galli, Nat. Commun. 5, 3919 (2014). https://doi.org/10.1038/ncomms4919 with 128 water molecules in the supercell. The MD trajectories of the Na+-water solution are from Ref. 2020. A. P. Gaiduk and G. Galli, J. Phys. Chem. Lett. 8, 1496 (2017). https://doi.org/10.1021/acs.jpclett.7b00239. At least 60 snapshots from each MD trajectory were employed in our electronic structure calculations. We calculated the polarizability of each water molecule for each of the MD snapshots and obtained averages and standard deviations. For ice VIII and ice X, we used a 96 and 128-molecule supercells, respectively; the accuracy of the cell size was validated using a Monkhorst-Pack k-point mesh of 8 × 8 × 8 with the primitive cells.3232. H. J. Monkhorst and J. D. Pack, Phys. Rev. B 13, 5188 (1976). https://doi.org/10.1103/physrevb.13.5188 The Perdew-Burke-Ernzerhof (PBE) exchange-correlation (xc) functional1717. J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996). https://doi.org/10.1103/physrevlett.77.3865 was used. Although PBE overestimates the molecular polarizability of an isolated water molecule by ∼10%, for water under pressure PBE gives both the static and the electronic dielectric constants in better agreement with experimental values than at ambient conditions, as shown in previous studies.18,1918. D. Pan, L. Spanu, B. Harrison, D. A. Sverjensky, and G. Galli, Proc. Natl. Acad. Sci. U. S. A. 110, 6646 (2013). https://doi.org/10.1073/pnas.122158111019. D. Pan, Q. Wan, and G. Galli, Nat. Commun. 5, 3919 (2014). https://doi.org/10.1038/ncomms4919 Here we used one xc functional to analyze trends of polarizabilities as a function of P and T; however, the method is general and can be used with any functional. In particular, we note that using finite field methods to compute polarizabilities (http://qboxcode.org/),99. Q. Wan, L. Spanu, G. A. Galli, and F. Gygi, J. Chem. Theory Comput. 9, 4124 (2013). https://doi.org/10.1021/ct4005307 calculations with hybrid functionals are readily carried out.
Figure 1 shows that at ambient conditions, the molecular polarizabilities of water given by DFT [Eq. (8)] are anisotropic. The out of plane polarizability is the largest, and the ones in-plane and perpendicular to the water dipole direction are smaller, consistent with the reports of other authors using just DID interactions to compute Eenv.7–97. R. J. Heaton, P. A. Madden, S. J. Clark, and S. Jahn, J. Chem. Phys. 125, 144104 (2006). https://doi.org/10.1063/1.23571518. M. Salanne, R. Vuilleumier, P. A. Madden, C. Simon, P. Turq, and B. Guillot, J. Phys.: Condens. Matter 20, 494207 (2008). https://doi.org/10.1088/0953-8984/20/49/4942079. Q. Wan, L. Spanu, G. A. Galli, and F. Gygi, J. Chem. Theory Comput. 9, 4124 (2013). https://doi.org/10.1021/ct4005307 We found that at high pressures and high temperatures, the anisotropy substantially decreases as shown in Fig. 1. Note that the polarizability of an isolated water molecule is less anisotropic than in the liquid at ambient conditions,9,139. Q. Wan, L. Spanu, G. A. Galli, and F. Gygi, J. Chem. Theory Comput. 9, 4124 (2013). https://doi.org/10.1021/ct400530713. X. Ge and D. Lu, Phys. Rev. B 96, 075114 (2017). https://doi.org/10.1103/physrevb.96.075114 suggesting that the anisotropy is critically related to the formation of hydrogen bonds. Indeed, also in supercritical water, the polarizability components are less dissimilar than at ambient conditions (see Fig. 1).
In Table I, four different methods to compute polarizabilities are compared from conditions varying from ambient to 11 GPa and from 0 to 2000 K. All methods show that with increasing pressure along an isotherm, the average molecular polarizability of water ($α¯mol=13Tr{αmol}$) decreases, while with increasing temperature along an isobar, it increases. Our previous study showed that the average dipole moment of water molecules increases with pressure, but decreases with temperature, so the present results indicate that varying the molecular dipole moments of water becomes more difficult when the values of the moments increase.
TABLE I. The average value of the three diagonal elements of the polarizability tensors of water molecules obtained by four methods: Clausius-Mossotti (CM) relation, dipole-induced-dipole (DID), random phase approximation (RPA), and density functional theory (DFT). The standard deviations of data in the molecular dynamics simulations are shown in parentheses. (Unit: Å3.)
Ambient1.605 (0.007)1.595 (0.082)1.625 (0.082)1.643 (0.083)
1.1 GPa, 1000 K1.638 (0.009)1.626 (0.107)1.668 (0.109)1.681 (0.110)
5.8 GPa, 1000 K1.526 (0.009)1.521 (0.107)1.599 (0.109)1.619 (0.111)
11.4 GPa, 1000 K1.453 (0.009)1.449 (0.109)1.549 (0.114)1.572 (0.116)
5.2 GPa, 2000 K1.639 (0.017)1.628 (0.187)1.712 (0.201)1.730 (0.204)
8.9 GPa, 2000 K1.586 (0.022)1.578 (0.251)1.688 (0.273)1.711 (0.276)
1st solvation shell of Na+1.569 (0.140)1.607 (0.133)1.623 (0.130)
Ice Ih {0001} surface1.610 (0.114)1.607 (0.120)1.622 (0.114)
The polarizabilities obtained by DFT [Eq. (8)] are slightly larger than those from RPA [Eq. (7)] by ∼0.02 Å3. When applying the two methods, we used the same electron polarization density Δρ, which is obtained when the exchange-correlation functional is included. The local electric field mainly comes from the electrostatic interactions, so the DFT and RPA values turn out to be similar.
In order to test the validity of the CM relation, we substituted the electronic dielectric constant ϵ, obtained by DFPT into the CM relation to calculate the average molecular polarizability. It is interesting to see that the CM relation yields nearly identical results as the DID approximation. The standard deviations obtained for the CM relation are smaller than those from the DID approximation by one order of magnitude, as they only arise from the thermal fluctuation of ϵ, not from molecular distributions as shown in Fig. 1. The CM relation holds when the Lorentz relation holds in an isotropic material, i.e., when the field at the center of a fictitious spherical cavity created by molecules inside the cavity vanishes.66. C. Kittel, Introduction to Solid State Physics, 8th ed. (Wiley, 2004). A well-known example for the Lorentz relation is the lattice with cubic symmetry, where only dipole-dipole interactions are considered.66. C. Kittel, Introduction to Solid State Physics, 8th ed. (Wiley, 2004). The agreement between the results obtained using the CM relation and the DID approximation suggests that the Lorentz relation is accurate when we consider only the dipole-dipole interaction for the water systems studied in Table I.
We now turn to comparing molecular polarizabilities of water in various phases. If we substitute the refractive index of 1.333, the experimental value for water at 293 K, and ambient pressure, into the CM relation, we get a molecular polarizability of 1.47 Å3, which is the same as the experimental value for water vapor. At the PBE level of theory, the polarizability of an isolated water molecule is 1.60 Å3, the same as the values obtained at ambient conditions using the CM relation and DID approximation (see Table I), consistent with previous studies.99. Q. Wan, L. Spanu, G. A. Galli, and F. Gygi, J. Chem. Theory Comput. 9, 4124 (2013). https://doi.org/10.1021/ct4005307 Hence, within the DID approximation, the average molecular polarizability of water does not change from gas to liquid phase. However, both the RPA and DFT methods give slightly larger values (∼3% larger than that obtained with CM and DID methods). We note that recently, Ge and Lu reported the molecular polarizabilities of water and ice at ambient conditions calculated using the local dielectric response of orbitals,1313. X. Ge and D. Lu, Phys. Rev. B 96, 075114 (2017). https://doi.org/10.1103/physrevb.96.075114 where the electron polarization density of each molecule Δρmol is evaluated individually. In general, the sum of Δρmol does not equate the total Δρ. In Ref. 1313. X. Ge and D. Lu, Phys. Rev. B 96, 075114 (2017). https://doi.org/10.1103/physrevb.96.075114, the $α¯mol$ of water was found to increase by ∼10% (instead of ∼3%) from gas to the liquid at ambient conditions.
Using the DFT, we also calculated the molecular polarizabilities of water in the first solvation shell of the Na+ ion at ambient conditions:2020. A. P. Gaiduk and G. Galli, J. Phys. Chem. Lett. 8, 1496 (2017). https://doi.org/10.1021/acs.jpclett.7b00239 $α¯mol$ is 1.62 Å3, which is again slightly larger than that obtained by the DID approximation by ∼3%. For water molecules with dangling bonds in the basal surface layer of ice Ih,21,2221. D. Pan, L.-M. Liu, G. A. Tribello, B. Slater, A. Michaelides, and E. Wang, Phys. Rev. Lett. 101, 155703 (2008). https://doi.org/10.1103/physrevlett.101.15570322. M. Watkins, D. Pan, E. G. Wang, A. Michaelides, J. VandeVondele, and B. Slater, Nat. Mater. 10, 794 (2011). https://doi.org/10.1038/nmat3096 the difference in $α¯mol$ given by the DFT and DID approaches is even smaller: 1.62 Å3 vs 1.61 Å3, only ∼1%. Our results suggest that at ambient conditions, the DID approximation works remarkably well.
Table I shows that the polarizabilities obtained by our method and the CM relation or the DID approximation differ with increasing pressure at a fixed temperature. For ice VIII, a high pressure phase consisting of two interpenetrating cubic ice sublattices,2323. V. F. Petrenko and R. W. Whitworth, Physics of Ice (Oxford University Press, 1999). when increasing pressure from 0 to 30 GPa, the difference between DFT values and the DID approximation increases from 6% to 13%, as shown in Fig. 2. These results indicate that interactions higher than dipole-dipole play a bigger role for denser water.
We close by considering the case of extremely dense water, ice X, the highest pressure phase ever determined experimentally.2424. R. J. Hemley, A. P. Jephcoat, H. K. Mao, C. S. Zha, L. W. Finger, and D. E. Cox, Nature 330, 737 (1987). https://doi.org/10.1038/330737a0 In ice X, the oxygen atoms are in a body-centered cubic lattice and the hydrogen atoms sit right between the two nearest O atoms (see Fig. 3). Because the H atom is equidistant to two O atoms, it is no longer possible to define H2O molecules; however, since the four maximally localized Wannier orbitals are still closely localized around O atoms, a new BB can be defined, and the molecular polarizability discussed below refers to the polarizability of the O2− anion.
Figure 3 shows that the CM relation gives the same results as the DID approximation, whereas the molecular polarizabilities given by the DFT and RPA methods are about 20% larger. The reason is that in ice X covalent bonds are present, and indeed the BB identified by our calculation is no longer a water molecule, but rather an anion, for which higher-order interactions play an important role.
For ice X, another interesting finding of our calculation is that the electronic dielectric constant ϵ has a minimum at around 250 GPa, and accordingly the bandgap increases up to 150 GPa and then decreases slowly, as shown in Fig. 4. The inverse correlation between the electronic dielectric constant and the bandgap of ice X is consistent with the Penn model,25,2625. D. R. Penn, Phys. Rev. 128, 2093 (1962). https://doi.org/10.1103/physrev.128.209326. W.-L. Yim, H. Shi, Y. Liang, R. J. Hemley, and J. S. Tse, in Correlations in Condensed Matter Under Extreme Conditions (Springer, 2017), Chap. 9. and differs from what we found in ice VII/VIII and water up to 30 GPa in our previous study.1919. D. Pan, Q. Wan, and G. Galli, Nat. Commun. 5, 3919 (2014). https://doi.org/10.1038/ncomms4919 Generally, ϵ increases when both the molecular polarizability and material density become larger. With increasing pressure, the molecular polarizability of ice X decreases as shown in Fig. 3, whereas the material density increases due to volume shrinking, hence the molecular polarizability and the material density of ice X are two competing factors determining ϵ; this is also the reason why the variation of ϵ is not very pronounced (see Fig. 4). From 50 to 250 GPa, the molecular polarizability dominates the change of ϵ, but above 250 GPa, the rate of its decrease becomes slower and thus the material density becomes a more important factor. As a result, ϵ increases slowly as shown in Fig. 4.
In order to predict the properties of solids and liquids composed of well-defined building blocks, it is important to determine the variation of the dielectric properties of the isolated molecular or nano-scale constituents upon assembly. Hence the ability to compute dipole moments and polarizabilities of building blocks in condensed phases is critical. In this paper, we proposed a first-principles method to compute polarizabilities of sub-entities in condensed phases, which includes multipole interactions at all orders and is applicable to semiconductors and insulators. The method only requires a single self-consistent calculation for the entire condensed system, as opposed to multiple calculations for each building block, and it is readily applicable within and beyond the RPA.
As an example, we presented results for the molecular polarizabilities of liquid water in a wide pressure and temperature range. We found that at ambient conditions, the dipole-induced-dipole approximation is sufficiently accurate and the Clausius-Mossotti relation may be used, e.g., to obtain molecular polarizabilities from experimental refractive indexes. However with increasing pressure, this approximation becomes unreliable. For example, in ice VIII, the contribution of multipoles beyond the dipole is ∼13% at 30 GPa, and in ice X, the difference between all multipoles and the DID contribution is about ∼20% at 350 GPa, indicating that when hydrogen bonds are replaced by covalent bonds, higher-order interactions cannot be ignored. In the case of ice X, the CM relation is not valid, though the Lorentz relation still holds under the DID approximation. We also found that the bandgap of ice X has a maximum, while the electronic dielectric constant of ice X has a minimum, as a function of pressure. Finally we note that the knowledge of the polarizabilities of sub-entities under pressure may help to design polarizable force fields suitable for extreme P-T conditions. The method presented here can be used to study the local dielectric response of a wide range of semiconductors and insulators.
We thank Deyu Lu, He Ma, and Ikutaro Hamada for their helpful discussions. D.P. acknowledges support from Hong Kong Research Grants Council (Project No. ECS-26305017), the National Natural Science Foundation of China (Project No. 11774072), the Alfred P. Sloan Foundation through the Deep Carbon Observatory, and the Croucher Foundation through the Croucher Innovation Grant. M.G. and G.G. were supported by MICCoM, as part of the Computational Materials Sciences Program funded by the U.S. Department of Energy, Office of Science, Basic Energy Sciences, Materials Sciences and Engineering Division. This research used resources of the Research Computing Center at the University of Chicago, the Argonne Leadership Computing Facility, which is a DOE Office of Science User Facility supported under Contract No. DE-AC02-06CH11357.
1. 1. A. Putrino and M. Parrinello, Phys. Rev. Lett. 88, 176401 (2002). https://doi.org/10.1103/physrevlett.88.176401, Google ScholarCrossref
2. 2. Q. Wan and G. Galli, Phys. Rev. Lett. 115, 246404 (2015). https://doi.org/10.1103/physrevlett.115.246404, Google ScholarCrossref
3. 3. G. D. Mahan, J. Chem. Phys. 76, 493 (1982). https://doi.org/10.1063/1.442750, Google ScholarScitation, ISI
4. 4. J. Klimeš and A. Michaelides, J. Chem. Phys. 137, 120901 (2012). https://doi.org/10.1063/1.4754130, Google ScholarScitation, ISI
5. 5. R. P. Feynman, R. B. Leighton, and M. Sands, Feynman Lectures on Physics. Volume 2: The New Millennium Edition: Mainly Electromagnetism and Matter (Basic Books, 2011), Vol. 2. Google Scholar
6. 6. C. Kittel, Introduction to Solid State Physics, 8th ed. (Wiley, 2004). Google Scholar
7. 7. R. J. Heaton, P. A. Madden, S. J. Clark, and S. Jahn, J. Chem. Phys. 125, 144104 (2006). https://doi.org/10.1063/1.2357151, Google ScholarScitation, ISI
8. 8. M. Salanne, R. Vuilleumier, P. A. Madden, C. Simon, P. Turq, and B. Guillot, J. Phys.: Condens. Matter 20, 494207 (2008). https://doi.org/10.1088/0953-8984/20/49/494207, Google ScholarCrossref, ISI
9. 9. Q. Wan, L. Spanu, G. A. Galli, and F. Gygi, J. Chem. Theory Comput. 9, 4124 (2013). https://doi.org/10.1021/ct4005307, Google ScholarCrossref, ISI
10. 10. B. Schropp and P. Tavan, J. Phys. Chem. B 112, 6233 (2008). https://doi.org/10.1021/jp0757356, Google ScholarCrossref
11. 11. N. Marzari, A. A. Mostofi, J. R. Yates, I. Souza, and D. Vanderbilt, Rev. Mod. Phys. 84, 1419 (2012). https://doi.org/10.1103/revmodphys.84.1419, Google ScholarCrossref, ISI
12. 12. X. Ge and D. Lu, Phys. Rev. B 92, 241107 (2015). https://doi.org/10.1103/physrevb.92.241107, Google ScholarCrossref
13. 13. X. Ge and D. Lu, Phys. Rev. B 96, 075114 (2017). https://doi.org/10.1103/physrevb.96.075114, Google ScholarCrossref
14. 14. S. Baroni, S. De Gironcoli, A. Dal Corso, and P. Giannozzi, Rev. Mod. Phys. 73, 515 (2001). https://doi.org/10.1103/revmodphys.73.515, Google ScholarCrossref, ISI
15. 15. I. Souza, J. Íñiguez, and D. Vanderbilt, Phys. Rev. Lett. 89, 117602 (2002). https://doi.org/10.1103/physrevlett.89.117602, Google ScholarCrossref, ISI
16. 16. F. Gygi, IBM J. Res. Dev. 52, 137 (2008). https://doi.org/10.1147/rd.521.0137, Google ScholarCrossref, ISI
17. 17. J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996). https://doi.org/10.1103/physrevlett.77.3865, Google ScholarCrossref, ISI
18. 18. D. Pan, L. Spanu, B. Harrison, D. A. Sverjensky, and G. Galli, Proc. Natl. Acad. Sci. U. S. A. 110, 6646 (2013). https://doi.org/10.1073/pnas.1221581110, Google ScholarCrossref
19. 19. D. Pan, Q. Wan, and G. Galli, Nat. Commun. 5, 3919 (2014). https://doi.org/10.1038/ncomms4919, Google ScholarCrossref
20. 20. A. P. Gaiduk and G. Galli, J. Phys. Chem. Lett. 8, 1496 (2017). https://doi.org/10.1021/acs.jpclett.7b00239, Google ScholarCrossref, ISI
21. 21. D. Pan, L.-M. Liu, G. A. Tribello, B. Slater, A. Michaelides, and E. Wang, Phys. Rev. Lett. 101, 155703 (2008). https://doi.org/10.1103/physrevlett.101.155703, Google ScholarCrossref
22. 22. M. Watkins, D. Pan, E. G. Wang, A. Michaelides, J. VandeVondele, and B. Slater, Nat. Mater. 10, 794 (2011). https://doi.org/10.1038/nmat3096, Google ScholarCrossref
23. 23. V. F. Petrenko and R. W. Whitworth, Physics of Ice (Oxford University Press, 1999). Google Scholar
24. 24. R. J. Hemley, A. P. Jephcoat, H. K. Mao, C. S. Zha, L. W. Finger, and D. E. Cox, Nature 330, 737 (1987). https://doi.org/10.1038/330737a0, Google ScholarCrossref, ISI
25. 25. D. R. Penn, Phys. Rev. 128, 2093 (1962). https://doi.org/10.1103/physrev.128.2093, Google ScholarCrossref
26. 26. W.-L. Yim, H. Shi, Y. Liang, R. J. Hemley, and J. S. Tse, in Correlations in Condensed Matter Under Extreme Conditions (Springer, 2017), Chap. 9. Google ScholarCrossref
27. 27. A. Baldereschi, R. Car, and E. Tosatti, Solid State Commun. 32, 757 (1979). https://doi.org/10.1016/0038-1098(79)90749-x, Google ScholarCrossref
28. 28. R. Resta and A. Baldereschi, Phys. Rev. B 23, 6615 (1981). https://doi.org/10.1103/physrevb.23.6615, Google ScholarCrossref
29. 29. D. R. Hamann, M. Schlüter, and C. Chiang, Phys. Rev. Lett. 43, 1494 (1979). https://doi.org/10.1103/physrevlett.43.1494, Google ScholarCrossref, ISI
30. 30. D. Vanderbilt, Phys. Rev. B 32, 8412 (1985). https://doi.org/10.1103/physrevb.32.8412, Google ScholarCrossref
31. 31. W. Dawson and F. Gygi, J. Chem. Phys. 148, 124501 (2018). https://doi.org/10.1063/1.5018116, Google ScholarScitation, ISI
32. 32. H. J. Monkhorst and J. D. Pack, Phys. Rev. B 13, 5188 (1976). https://doi.org/10.1103/physrevb.13.5188, Google ScholarCrossref, ISI