PHYSICAL REVIEW B
96
, 104202 (2017)
Xray absorption of liquid water by advanced
ab initio
methods
Zhaoru Sun,
1
Mohan Chen,
1
Lixin Zheng,
1
Jianping Wang,
1
Biswajit Santra,
2
Huaze Shen,
3
Limei Xu,
3
Wei Kang,
4
Michael L. Klein,
1,5,6
and Xifan Wu
1,6,*
1
Department of Physics, Temple University, Philadelphia, Pennsylvania 19122, USA
2
Department of Chemistry, Princeton University, Princeton, New Jersey 08544, USA
3
International Center for Quantum Materials, School of Physics, Peking University, Beijing 100871, China
4
College of Engineering, Peking University, Beijing 100871, China
5
Department of Chemistry, Temple University, Philadelphia, Pennsylvania 19122, USA
6
Institute for Computational Molecular Science, Temple University, Philadelphia, Pennsylvania 19122, USA
(Received 17 April 2017; revised manuscript received 6 July 2017; published 11 September 2017)Oxygen
K
edge xray absorption spectra of liquid water are computed based on conﬁgurations from advanced
ab initio
molecular dynamics simulations, as well as an electron excitation theory from the GW method.One the one hand, the molecular structures of liquid water are accurately predicted by including both vander Waals interactions and a hybrid functional (PBE0). On the other hand, the dynamic screening effects onelectron excitation are approximately described by the recently developed enhanced static Coulombhole andscreenedexchange approximation of W. Kang and M. S. Hybertsen [Phys. Rev. B
82
, 195108 (2010)]. The
resulting spectra of liquid water are in better quantitative agreement with the experimental spectra due to thesoftened hydrogen bonds and the slightly broadened spectra srcinating from the better screening model.DOI: 10.1103/PhysRevB.96.104202
I. INTRODUCTION
Water is arguably one of the most important materials onEarth and needs to be thoroughly understood [1]. However,
the understanding of liquid water is itself a challenge in manyaspects. Unlike other liquids, water shows a lot of anomaliessuch as the density maximum at 4
◦
C and the isobaric heatcapacity minimum at 35
◦
C, among many other things [2,3].
Understanding the microscopic structures and dynamics of liquid water, in particular its hydrogenbond (HB) network,is the key to understanding these anomalies [4–7]. Recently,
highresolution oxygen
K
edge corelevel spectroscopy, suchasxrayabsorptionspectra(XAS)andxrayRamanscattering,has emerged as a powerful experimental technique to probethe electronic structure and infer the molecular structure of water and ice, as well as aqueous solutions [6,8–22]. Excited
from the oxygen 1
s
core level, the electron excitation probesthe unoccupied electronic states, which are antibonding alongthe covalent OH bonds and particularly sensitive to the HBs.Therefore, the XAS technique serves as a local probe for theHB structures of liquid water and ice.The experimental XAS of water have three distinct featuresas a function of increasing excitation energies: a preedge startsfrom the absorption threshold at 533 to 536 eV with a peakcentered at 535 eV, a main edge spans from 537 to 539 eV, anda postedge exists from 539 eV and beyond [11,18,20–22]. Ex
perimentally,thepreedgefeatureispresentinbothliquidwaterandicebutismoreintenseintheformer.Therelativeintensitiesbetween the main edge and postedge of liquid water andcrystalline ice Ih are substantially different [6,8–11,20–24].
As one of the most qualitative differences, the intensity of themain edge is higher than that of the postedge in the XAS of liquid water, while the opposite trend is true for the spectra inthe ice.
*
Corresponding author: xifanwu@temple.edu
The unambiguous assignments of the XAS features to theunderlying HB structures are prerequisites for the physical interpretationoftheexperimentalspectra,whichcanbeachievedby ﬁrstprinciples methods including both the modeling of the molecular structure of liquid water and the electronholeexcitation process. With snapshots of the represented molecular conﬁgurations from an equilibrated molecular dynamicstrajectory, the XAS can be computed with the knowledge of electronic structures of the excited core hole. In this regard,variousapproximations[17,25,26]forexcitedelectronicstates
havebeenproposedwithintheframeworkofdensityfunctionaltheory (DFT) [27,28]. In the seminal work of Prendergast
and Galli [17], the proposed excited electron and corehole
approximation yielded XAS in close agreement with experimental measurements. More rigorously, the XAS of liquidwatercanbecomputedbysolvingtheBetheSalpeterequation(BSE)describingtheelectronhole interaction[15,29,30].The
BSEapproachinvolvescalculationsoftheselfenergyoperatorand quasiparticles, which are, in general, computationallyexpensive for liquid water. An approximate way of solvingthe BSE with fewer computational resources was introducedby Chen
et al.
[9] based on a model electronscreeningfunction in the static Coulombhole and screenedexchange(COHSEX) approximation [31], which was evaluated usingthe maximally localized Wannier functions as basis [32],
reducing the computational cost signiﬁcantly. In the abovework, the molecular srcins of the spectra features as well asthe intensity difference between the main edge and postedgein water and ice have been validated [9].
Despite these recent developments, two uncertainties stillremaininthetheoreticalXASofliquidwaterobtainedwiththeapproximate solutions of the BSE. The ﬁrst uncertainty comesfrom the drawback in modeling the liquidwater structure byemploying the generalized gradient approximation (GGA) inthe framework of DFT [33–43]. It is known that GGA predicts
overstructured liquid water [40–43]. Thus, a signiﬁcantly
elevated temperature is often adopted to generate a softer HB
24699950/2017/96(10)/104202(9) 1042021 ©2017 American Physical Society
ZHAORU SUN
et al.
PHYSICAL REVIEW B
96
, 104202 (2017)
structureclosertotheexperimentalmeasurement[41–43].The
lack of physics in accurately describing the water structure isdue to the neglected van der Waals (vdW) interactions and thespurious selfinteraction error [44] in the GGA functional.
Speciﬁcally, by including the vdW interactions, the waterpopulation in the interstitial region between the ﬁrst andsecond coordination shells of water molecules is increasedto better match experiment [43]. Furthermore, by mitigatingthe selfinteraction error through the hybrid functional, thedirectional HB strength between water molecules is weakenedso it is closer to the experiment; as a result, the protonsare less easily donated to neighboring water molecules [43].
However, the effects of this improved water structure on thetheoretical XAS have not been elucidated. Second, in theseries of computational works adopting the static COHSEXapproximation for the electronhole excitation of liquid water,a homogeneous electronic screening model was ﬁrst adopted[9] and then extended by using the HybertsenLouie ansatz
[45] to account for the inhomogeneous screening effects from
themolecularenvironment[10].However,somediscrepanciesstill exist, and it is not yet clear to what extent the dynamicscreening effect will affect the quasiparticle wave functions(QWs) and the computed XAS. For example, it was observedthat the width of the theoretical XAS by static COHSEX isslightly narrower than the experimental data [10].
In an effort to address the above issues, we adopted asystematic way to study the XAS of liquid water at ambientconditions.Speciﬁcally,weusedmoreadvanced
abinitio
modeling of molecular structures and electronic excitations. Wegeneratedliquidwatertrajectoriesfromthe
abinitio
moleculardynamics (AIMD) [46] simulations by employing a hierarchyof exchangecorrelation (XC) functionals of Perdew, Burke,and Ernzerhof (PBE) [47], PBE with the vdW interactions in
the form of Tkatchenko and Schefﬂer (PBE
+
vdW) [48], and
the hybrid functional PBE0 [49,50] with vdW interactions
(PBE0
+
vdW). By utilizing the enhanced static COHSEXmethod to treat the excitations [51], we ﬁnd that the XAScomputed from the snapshot generated by the PBE0
+
vdWfunctional agree well with the experiment among the three XCfunctionals studied. The vdW interactions soften the waterstructures by increasing the population of water molecules inthe interstitial region, while the hybrid functional mitigatesthe selfinteraction error and weakens the HB strength so itis closer to experiment [43]. Both structural corrections and
improved excitation theory are crucial in giving rise to anoverall improvement of the three edges of XAS. In particular,the postedge feature in the highenergy region of XAS is inslightly better agreement with experiment because of betterscreening modeling using the enhanced static COHSEX. Inaddition, we also compare a set of XAS computed fromdifferent excitation theories to the experiment in order toshow the importance of selfconsistently diagonalized QWsin capturing the qualitative features of XAS.
II. METHODS
We performed AIMD simulations to generate liquidwater trajectories using a modiﬁed version of the
QUANTUMESPRESSO
package [52]. We simulated 128 water molecules ina periodic cubic cell with a cell length of 15.68 ˚A using theCarParrinello molecular dynamics (CPMD) [46] within the
canonical (
NVT
) ensemble. We employed normconservingpseudopotentials in the form of Troullier and Martins [53] and
set the kineticenergy cutoff of the electronic wave functionsas 71 Ry. We used a hierarchy of XC functionals, includingPBE, PBE
+
vdW, and PBE0
+
vdW, as mentioned. The hybridfunctional PBE0 with a mixing of 25% exact exchange wasevaluated in a linearscaling manner by taking advantageof maximally localized Wannier functions [32]. The ionictemperatures were controlled through NoséHoover chainthermostats with a chain length of 4 for each ion [54–56].
All AIMD simulations were performed at 330 K, where anincrease of 30 K has been found to mimic the nuclear quantumeffectinstructuralquantitiessuchastheoxygenoxygenradialdistribution function in DFTbased simulations of liquid water[57]. A time step of 4.0 a.u. and a ﬁctitious electron mass of
300 a.u. were chosen.We calculated the xray absorption cross section usingFermi’s golden rule:
σ
(
ω
)
=
4
π
2
α
0
¯
hω
f

M
if

2
δ
(
ω
if
−
ω
)
,
(1)where
α
0
is the ﬁnestructure constant and ¯
hω
is theabsorbed photon energy matching the energy difference¯
hω
if
=
E
f
−
E
i
.
E
f
and
E
i
are the eigenvalues of the ﬁnaland initial states, respectively.
M
if
are the transition matrixelements between the initial state

φ
i
and ﬁnal state

φ
f
,which can be evaluated within the electricdipole approximation as
M
if
∼
φ
i

x

φ
f
, averaging over the three Cartesiandirections.Wetakethe1
s
atomiccorewavefunctionfromDFTcalculations as the initial state

φ
i
. For the ﬁnal state

φ
f
, weapply a selfconsistent diagonalization procedure within theenhanced static COHSEX approach [51], and the details areas follows.The ﬁnal state

φ
f
is obtained by utilizing the enhancedCOHSEX approach, which has been implemented withinthe framework of the CPMD [46] scheme. Speciﬁcally, the
COHSEXmethodhasbeenimplementedintheCPMDmodulewithin the
QUANTUM ESPRESSO
package [52]. For each set of
input wave functions

φ
f
, we ﬁx the ion positions and dampthe wave functions of the system in a selfconsistent way,as explained below. Note that the formulas we describe hereare suitable only for the excitation theory we adopt and areindependent of the other CPMD simulations we performedwith vdW and PBE0 functionals.First, the Lagrangian from the CarParrinello approach is
L
=
µ
2
i
˙
ψ
i

˙
ψ
i
+
12
I
M
I
˙
R
2
I
−
E
tot
(
R
,
{
ψ
}
)
+
λ
ij
(
ψ
i

ψ
j
−
δ
ij
)
,
(2)where
ν
is a ﬁctitious mass of electrons and
ψ
i
is the orbitalof state
i
.
M
I
is the mass for atom
I
that is located at
R
I
.
E
tot
(
R
,
{
ψ
}
) is the total energy calculated from ﬁrstprinciplesmethods. The last part is the orthogonality constraint imposedon orbitals by the Lagrangian multiplier to be
λ
ij
. Note thatthe initial ion positions are ﬁxed, and only wave functions aregenerated. Therefore, we need to damp the wave functionstowards the ground state, which is realized via the equations
1042022
XRAY ABSORPTION OF LIQUID WATER BY ADVANCED . .. PHYSICAL REVIEW B
96
, 104202 (2017)
of motion of electrons in planewave basis set,
µ
¨
ψ
i
= −
H
(
R
,
{
ψ
}
)
ψ
i
+
j
λ
ij
ψ
j
.
(3)Here,
H
(
R
,
{
ψ
}
) is the Hamiltonian matrix of the system.Furthermore, the Hamiltonian part can be evaluated in realspace as
H
ψ
(
r
)
=
−
12
∇
2
(
r
)
+
V
ext
(
r
,
R
)
+
V
H
(
r
)
ψ
(
r
) (4)
+
d
r
′
(
r
,
r
′
,E
)
ψ
(
r
′
)
,
(5)where
V
ext
(
r
,
R
) is the external potential and
V
H
(
r
) is theHartree potential. In particular,
(
r
,
r
′
,E
) is the selfenergyoperator that is nonlocal in real space and depends on theselfenergy
E
. In the static COHSEX approximation, theselfenergy operator can be approximated as
staticCOHSEX
(
r
,
r
′
)
=
staticCOH
(
r
,
r
′
)
+
staticSEX
(
r
,
r
′
)
.
(6)The ﬁrst part is
staticCOH
(
r
,
r
′
)
=
12
δ
(
r
−
r
′
)
W
p
(
r
,
r
′
;
E
=
0) (7)
=
12
δ
(
r
−
r
′
)(
W
−
v
)
,
(8)where
W
p
is the Coulomb hole,
W
is the screened Coulombinteraction, and
v
is the bare Coulomb interaction. TheHybertsenLouie ansatz [45] proposed that
W
generallyfollows the local charge density and has the form
W
(
r
,
r
′
;
E
=
0)
=
12
{
W
[
r
−
r
′
;
ρ
(
r
′
)]
+
W
[
r
′
−
r
;
ρ
(
r
)]
}
,
(9)where
W
can be written as
W
[
r
′
−
r
;
ρ
(
r
)]
=
12
π
3
ǫ
−
1
[
q
;
ρ
(
r
)]
v
(
q
)
e
i
q
·
(
r
−
r
′
)
d
q
.
(10)Here, we use the Bechstedt model [58] for the dielectric
function, 0
ǫ
[
q
,ρ
(
r
)]
=
1
+
(
ǫ
0
−
1)
−
1
+
α
(
q/q
TF
)
2
+
q
4
43
k
2
F
q
2TF
−
1
,
(11)where
q
TF
is the ThomasFermi wave vector,
k
F
is the Fermiwave vector,
ǫ
0
is taken from experiment, and
α
is ﬁxed bymatching the Bechstedt model to the
q
2
dependence of thePenn model [59]. Next, we can transform
W
to
W
[
r
′
−
r
;
ρ
(
r
)]
=
v
(
r
−
r
′
)
ǫ
0
−
1
a
(
x
1
−
x
2
)

r
′
−
r
×
e
i
√
x
1

r
′
−
r

x
1
−
e
i
√
x
2

r
′
−
r

x
2
.
(12)Here,
x
1
,
2
=
(
−
b
±√
b
2
−
4
ac
)
/
2
a
,
a
=
(
43
k
2
F
q
2TF
)
−
1
,
b
=
α/q
2TF
, and
c
=
ǫ
0
/
(
ǫ
0
−
1). The second part is
staticSEX
(
r
′
,
r
)
= −
occ
i
ψ
i
(
r
)
ψ
∗
i
(
r
′
)
W
(
r
,
r
′
;
E
=
0)
.
(13)Numerically, it has been shown that most of the error fromusing the static COHSEX approximation, when compared tothe GW approximation, comes from the shortwavelength partof the assumed adiabatic accumulation of the Coulomb hole
W
p
, namely, the COH part, while the SEX term in the staticCOHSEXapproximationyieldsvaluesrelativelyclosetothoseof the GW calculations [51]. Therefore, the enhanced staticCOHSEX was proposed to introduce a universal function
f
to approximately include the dynamics screening in thesrcinal static model in the COHSEX formula. Speciﬁcally,the enhanced static COH term can be evaluated as
newCOH
(
r
,
r
′
)
=
δ
(
r
−
r
′
)2
W
p
(
q
;
E
=
0)
f
(
q/k
f
)
e
−
i
q
·
r
d
q
,
(14)where
q
isaplanewaveand
k
f
istheFermivector.Thescalingfunction is
f
(
x
)
=
1
+
a
1
x
+
a
2
x
2
+
a
3
x
3
+
a
4
x
4
+
a
5
x
5
+
a
6
x
6
1
+
b
1
x
+
b
2
x
2
+
b
3
x
3
+
b
4
x
4
+
b
5
x
5
+
b
6
x
6
,
(15)where
a
1
=
1
.
9085,
a
2
= −
0
.
542572,
a
3
= −
2
.
45811,
a
4
=
3
.
08067,
a
5
= −
1
.
806,
a
6
=
0
.
410031,
b
1
=
2
.
01317,
b
2
=−
1
.
55088,
b
3
=
1
.
58466,
b
4
=
0
.
368325,
b
5
= −
1
.
68927,and
b
6
=
0
.
599225.Above,wedescribedtheprocedurestocomputeXASbasedon an excited water molecule in a snapshot. We then excitedevery water molecule in the 128molecule supercell in orderto sample the different local environments of the disorderedliquidwater structure. In our case, we found converged XASafter randomly exciting 64 water molecules in the snapshot.Wenotethatinapreviousstudy,10individualanduncorrelatedsnapshots of 32 water molecules were chosen, and onlysmalldifferenceswereobservedbetweenthesesnapshots[17].
Therefore, we consider one snapshot of a large cell to besufﬁcient to yield meaningful results and take a representativesnapshotfromeachAIMDtrajectoryreﬂectingtheequilibratedstructure of liquid water to compute XAS.Due to the different local environments in liquid water foreach excited oxygen atom, we adopted the corehole energyshift of each excitation by following Ref. [60], which is astandard approach to compute the corelevel shifts. We usedGaussian broadening of 0.4 eV for all spectra, which wasused in previously calculated XAS of liquid water [9,10].
The computed XAS were aligned to the onset of the preedge(535eV)andthennormalizedtothesameareaofexperimentaldata ranging from 533 to 546 eV.To analyze the realspace locations of QWs in terms of the three edges in XAS of liquid water, we also deﬁne theonedimensional density of the QW as
ρ
i
(
r
)
=

ψ
i
(
r,θ,φ
)
2

dθdφ
(16)along the radial direction, and the origin is taken to bethe position of the excited oxygen with a core hole. Inthe equation,
i
represents the index of the conduction bandwith excitation energy
ε
i
. In order to present the differentlocalizations of the QWs, the index
i
is chosen so that
ε
i
∈
[534
.
5 eV
,
535
.
5 eV] for QWs representing the preedge
1042023
ZHAORU SUN
et al.
PHYSICAL REVIEW B
96
, 104202 (2017)
0.000.010.02 preedge0123g
OO
(r)0.000.010.02 mainedge0123g
OO
(r)0.000.010.02 postedge0 1 2 3 4 5 6 7
Distance (Å)
0123g
OO
(r)
(a)(b)(c)
FIG. 1. DensitydistributionsoftheQW(greendashedline)of(a)the preedge, (b) main edge, and (c) postedge as a function of oxygenoxygen distance computed from the snapshot of the PBE0
+
vdWtrajectory using the enhanced static COHSEX method. The
g
OO
(
r
)(red line) from the PBE0
+
vdW trajectory is shown for comparison.The insets show the representative QWs of the three edges around theexcited water molecule. Water molecules residing within the secondcoordination shell of the excited oxygen are shown. Red, white, andyellow spheres represent oxygen, hydrogen, and oxygen atoms witha core hole, respectively. QWs with opposite signs are depicted inblue and green.
of XAS,
ε
i
∈
[537 eV
,
539 eV] for QWs representing themainedge of XAS, and
ε
i
∈
[540 eV
,
542 eV] for QWsrepresenting the postedge of XAS.
III. RESULTS AND DISCUSSIONA. HB structure probed by the preedge, main edge,and postedge of XAS
The QWs can be strongly perturbed by the local liquidwater structures. The preedge, main edge, and postedge of XAS are found to have distinguishable molecular signaturesthat relate to different spatial regions of the HB network[8,9,18]; the preedge has 4
a
1
character, while the main edgeand postedge have
b
2
character, both of which srcinate fromthe molecular excitations in the gas phase [8]. In order to
quantitatively study the spatial regions in terms of differentXAS edges, we present the density distributions of QWs as afunction of oxygenoxygen distance in Fig. 1. The oxygen
oxygen radial distribution function
g
OO
(
r
) and an excitedoxygen with QWs distributed within the HB network (insets)are also shown for comparison. Overall, Figs. 1(a), 1(b),
and 1(c) show that the density distributions of QWs becomemore delocalized from the preedge to the main edge to thepostedge.ThedensitydistributionofQWsofthepreedgeillustratedinFig. 1(a) has the highest peak at 1.7 ˚A and is mostly localizedwithin 2.75 ˚A, the latter of which is the ﬁrst peak positionof
g
OO
(
r
). The inset in Fig. 1(a) shows that the QW of thepreedge resembles the ﬁrst excited state of a water molecule inthe gas phase with 4
a
1
symmetry. In this regard, our result isconsistent with a previous assignment [9] of the preedge to abound exciton state,wherethe electronorbitalwasfound tobemostly localized within the ﬁrst coordination shell. Therefore,the preedge features can be largely affected by the shortrangestructures, such as broken HBs and covalent bond strength,aroundtheexcitedoxygens. WiththeweakenedHBsdescribedby the hybrid DFT functional (PBE0) and vdW interactions,that is, the shortrange HB network, therefore, the computedpreedge of XAS is expected to be improved in both energiesand intensities.As shown in Fig. 1(b), the density distribution of the QWof the main edge is more delocalized than that of the preedge.The inset shows that the QW of the main edge can be foundnot only on the excited molecule itself but also on its ﬁrstand secondshell neighbors. In addition, a clear
b
2
charactercan be identiﬁed for a typical QW of the main edge. Theabove is consistent with the fact that the main edge wasfound to srcinate from the second excited state of a watermolecule in the gas phase. By comparing the localization of the mainedge density of the QW to that of the preedge, it canbe seen that the former is more localized between the ﬁrst andsecondcoordinationshellsoftheliquidwaterstructure.Inthisregard, it is expected that the mainedge feature of XAS willbe sensitive to the intermediaterange order ofthe liquidwaterstructure, i.e., water molecules in the interstitial region.In contrast to the density distributions of QWs of thepreedge and main edge, the QW of the postedge shown inFig. 1(c) is the most delocalized one but still preserves
b
2
character. The strong delocalization can be clearly seen by theincreased density of the QW as a function of the distanceaway from the excited water molecule. Hence, the watermolecules in longrange order are critical in determining thepostedge features. Using a small simulation cell containing 32water molecules fails to yield a postedge intensity of XAS incloseagreementwithexperimentaldata[9].Thedelocalizationfeatureisconsistentwiththefactthatthepostedgeisaresonantexciton state.
B. XAS calculated from PBE, PBE
+
vdW,and PBE0
+
vdW AIMD trajectories
Among the three levels of XC functionals investigated,the XAS computed from the snapshot obtained with PBEshow the least agreement with the experimental spectra inFig. 2(a). Four major discrepancies can be identiﬁed. First,the intensity of the computed preedge is underestimatedcompared to experiment. Second, the theoretical main edgeis centered around 538.5 eV, showing a large blueshift (around1 eV) compared to the experimental value at 537.5 eV.Third, a signiﬁcantly overestimated postedge intensity leadsto the fact that the main edge and post edge have almostthe same intensities, which contradicts the experimental fact
1042024
XRAY ABSORPTION OF LIQUID WATER BY ADVANCED . .. PHYSICAL REVIEW B
96
, 104202 (2017)
(e) (f)
2 3 4 5 6
r (Å)
00.511.522.53
g
O  O
( r )
Exp.PBE2 3 4 5 6
r (Å)
00.511.522.53
g
O  O
( r )
Exp.PBE+vdW534 536 538 540 542 544 546
Energy (eV)
02468101214
I n t e n s i t y ( A r b . U n i t s )
Exp.PBE534 536 538 540 542 544 546
Energy (eV)
02468101214
I n t e n s i t y ( A r b . U n i t s )
Exp.PBE+vdW534 536 538 540 542 544 546
Energy (eV)
02468101214
I n t e n s i t y ( A r b . U n i t s )
Exp.PBE0+vdW
(a) (b) (c)(d)
2 3 4 5 6
r (Å)
00.511.522.53
g
O  O
( r )
Exp.PBE0+vdW
FIG. 2. Computed XAS and
g
OO
(
r
) from three levels of exchangecorrelation functionals used in the AIMD simulations, namely, PBE,PBE
+
vdW, and PBE0
+
vdW. A representative snapshot consisting of 128 water molecules from each equilibrated AIMD trajectory was usedfor spectra calculation. The enhanced static COHSEX method is adopted as the excitation theory. The experimental (Exp.) data of XAS [21]and
g
OO
(
r
) [61] are also shown for comparison.
that the main edge is more prominent than the postedgein liquid water. Fourth, the width of the computed XAS isslightly narrower than the experimental XAS, especially in thehighenergy region. Experimentally, the XAS of crystallineice Ih with an intact HB network show a more prominentpeak of the postedge than that of the main edge. Therefore,the discrepancies imply that the HB network of liquid waterfrom the PBE AIMD trajectory is overstructured. Indeed,
g
OO
(
r
) computed from the PBE AIMD trajectory signiﬁcantlydeviates from experimental measurement [61], as shown in
Fig. 2(d). Speciﬁcally, the ﬁrst and second peaks of
g
OO
(
r
)from simulations are signiﬁcantly overestimated, and the ﬁrstminimum is largely underestimated. Consistently, the averagenumber of HBs per water molecule is found to be 3.76 inthe PBE trajectory according to the popular HB deﬁnitionproposed by Luzar and Chandler [62], which is the highestamongallthefunctionalsstudiedherein.Alltheaboveindicatethat the HB network of liquid water is overstructured from thePBE AIMD trajectory. Not surprisingly, the theoretical XASpredicted by an overstructured HB network from the PBEAIMD trajectory yield icelike spectra with a relatively moreprominent postedge feature.As shown in Fig. 2(b), the XAS computed from the
PBE
+
vdW trajectory are largely improved with regard to theexperimentalspectracomparedtothoseobtainedfromthePBEtrajectory. The improvement can be seen by a higher preedgeintensity, a shift of the main edge towards lower energy, anda lower postedge intensity. The better agreement is attributedto the improved description of the HB network by includingthe vdW interactions in the AIMD simulation. An explicitaccount of vdW forces strengthens the attractive interactionsamong water molecules, which signiﬁcantly increase thepopulation of water molecules in the interstitial region. Theincreased population of water molecules brings
g
OO
(
r
) incloser agreement with experiment within the ﬁrst and secondcoordination shells. Moreover, the increased population of water molecules in the interstitial region weakens the HBsamong the water molecules in the ﬁrst coordination shell,resulting in a reduced ﬁrst peak in
g
OO
(
r
). The averagenumber of HBs per molecule is found to be 3.56, whichis about 5% smaller than that of PBE. Hence, an excitedoxygen atom experiences a more disordered environment bythe surrounding water molecules. First of all, the preedgeintensity is increased due to the breaking of more HBs. Inorder to verify this, we selected two excited water molecules,one with broken HBs and the other one with four intact HBs,and plotted their density distributions of QWs in Fig. 3(a).The QW of the preedge from the excited molecule with thebroken HBs is more localized with enhanced
p
character dueto the more disordered shortrange molecular environment.Therefore, larger amplitudes of transition matrices
M
ij
areobtainedaccordingtoFermi’sgoldenrule,asshowninEq.(1).
Asaresult,thepreedgeintensityisincreasedwithmorebrokenHBs. Second, the mainedge intensity is also enhanced byan increased population of water molecules in the interstitialregionduetovdWinteractions.Inordertofurtherexplainhowthe water molecules in the interstitial region affect the mainedge features, we also selected two representative excitedwater molecules and calculated their density distributions of QWs; one of the excited water molecules has four intact HBs,and the other one has extra neighboring water moleculeswithin the interstitial region in addition to the four intactHBs. Figure 3(b) suggests that the density of the mainedgeQW of the excited molecule with extra neighboring water
1042025