The structural origin of anomalous properties of liquid water

Water is the most common liquid in nature, with unusual properties that could be linked to the peculiar hydrogen-bonding network holding the molecules together. Here, Nilsson and Pettersson review recent progress in searching the connections between local configurations and thermodynamic responses of water.

W ater is the most important liquid for our existence and plays an essential role in physics, chemistry, biology and geoscience. What makes water unique is not only its importance but also the anomalous behaviour of many of its macroscopic properties. The ability to form up to four hydrogen bonds (H-bonds), in addition to the non-directional interactions seen in simple liquids, leads to many unusual properties such as increased density on melting, decreased viscosity under pressure, density maximum at 4°C, high surface tension and many more (see, for example, http://www.lsbu.ac.uk/water/index.html). If water would not behave in this unusual way it is most questionable if life could have developed on planet Earth 1 . Figure 1 shows the temperature dependence of the isobaric heat capacity (C P ) and the isothermal compressibility (k T ) that initially decrease with decreasing temperature but, just above the ambient temperature regime, there is a deviation and both the heat capacity and compressibility begin to increase on further cooling, resulting in a minimum at 308 K (35°C) and 319 K (46°C), respectively [2][3][4][5][6] . These properties are related to fluctuations in the liquid where k T derives from volume, or equivalently density, fluctuations, whereas C P is related to entropy fluctuations 2,7 . Similarly, the density (r) of a simple liquid increases with decreasing temperature and this is true also for water when it is hot. However, as it enters the ambient regime, the rate of increase decreases and a maximum in the density is seen at 277 K (4°C) after which the density instead decreases on further cooling [8][9][10] . The density variation can be derived from the thermal expansion coefficient (a P ), dependent on the cross-correlation between fluctuations in density and entropy, which thus becomes negative for water below the density maximum 2,7 . As these thermodynamic properties measure the response to external perturbations, such as temperature or pressure, they are referred to as thermodynamic response functions 11 . Their unusual behaviour for water is denoted anomalous where the anomaly becomes much stronger on cooling. Indeed, when fitted to a power law each of these properties seems to diverge towards a temperature of 228 K (ref. 4).
To gain fundamental understanding of the origin of these anomalies we need to address the instantaneous local structure of the liquid at various thermodynamic state points and establish how this structure couples to the dynamics of the molecular motion. The starting question is: on a larger scale liquid water is homogeneous, but is it also locally homogeneous or can it be heterogeneous (see Box 1)? Not in terms of a static structural picture, but induced by fluctuations on some length and time scales between specific classes of local structures. Many different plausible explanations exist for the unusual properties of water where maybe both homogeneous and heterogeneous models could be viable and sophisticated structural and dynamical experimental data are needed to determine their validity.
The other major question is: how large are distortions in the H-bonding network? In a locally homogeneous model the distortions should be around a near-tetrahedral H-bond arrangement, whereas in a heterogeneous model there would be distortions within each class of configurations in addition to the distinction in local structure between the classes. The more popular heterogeneous models build on fluctuations between two main classes of contrasting structures with notations such as tetrahedral and distorted 12 , symmetrical and asymmetrical 13,14 , locally favoured and normal 15 , and low-density liquid (LDL) and high-density liquid (HDL) 7,[16][17][18] to mention the most recent proposals; these refer to the same two general structural classes where the notation reflects which specific properties are brought forward by the various experimental 7,11,12,16-18 and modelling [13][14][15] techniques.

Box 1 | A simplified picture of heterogeneous-homogeneous fluctuations
Fluctuations in a liquid could appear on various length and time scales. If the fluctuations are correlated over distances that extend over more than the molecular size, the fluctuations could give rise to different local structural regions that we denote patches. Here, in the figure below, we envision a heterogeneous picture with two local structural environments (blue and yellow), while structures representing a structural average in between these two extremes are indicated as green. The patches are not static but will rearrange into each other through structural fluctuations on some time-scale that is longer than the hydrogen bond lifetime. We also have the situation where the liquid only encompasses a homogeneous single structural environment with normal thermal fluctuations that is denoted a simple liquid.
In an illustrative picture we can liken the dynamics of the fluctuations in the heterogeneous case to that of a pendulum moving between two turning points. The speed is zero at the turning points and maximum at the point in the middle. Here the black curve shows such a scenario as time spent on various points along the fluctuating trajectory where the molecules of the liquid will spend most of the time at the extreme points and little time as they are converting in between. The timescales for the collective motion is expected to be a minimum of a few picoseconds but varying with temperature. This leads to an extreme heterogeneous situation where two well-defined structural classes develop on a temperature-dependent length and time scale. The length scale should correspond as a minimum to a factor 2 to 3 of the molecular dimension allowing for the existence of several shells. The boundary in between contains few molecules, as the time spent in this region is very low. Another much debated question is whether such heterogeneous fluctuating structures could develop into metastable macroscopic phases on extreme supercooling 19,20 . Would this also lead to the existence of a second critical point 19 ? A liquid-liquid critical point (LLCP) has never been observed for a one-component fluid, but the apparent power-law divergence of compressibility and heat capacity on supercooling water has been proposed to be an indication of an LLCP 19 .
Here we will review developments in the last few years around the question whether pure bulk water is structurally heterogeneous or homogeneous and how the answer may relate to its unique anomalous properties. At the end we propose a unified picture that can explain many observations in both the ambient and supercooled regimes.
Electronic structure and vibrational spectroscopy Transitions between different electronic or vibrational states give specific spectral signatures that can shed light on possible unique structural environments. The goal of this section is to establish trends in these spectral features based on experimental observations of their temperature dependence and perturbations due to added salt.
Local structure probed by X-ray spectroscopies. Here we first review electronic spectroscopies involving core level excitations to and from the valence levels. Figure 2a shows the temperature dependence in X-ray absorption spectroscopy (XAS), which probes the unoccupied states and was recently obtained also for slightly supercooled water 21 . There is general agreement that the pre-(535 eV) and main-edge peaks (537-538 eV) fingerprint distorted H-bonds, whereas the post edge (540-541 eV) is associated with strong H-bonds and is further enhanced for tetrahedral H-bond structures 13,[22][23][24][25] . This is consistent with the temperature dependence where the post edge, which dominates in hexagonal ice, loses intensity to the pre-edge and main edge on heating as seen in Fig. 2a (refs 21,22). The spectral features do not undergo any major shift or broaden with increasing temperature, but instead there are changes in intensity at an almost fixed energy. It has also recently been shown that the main-edge intensity becomes enhanced on formation of high-density forms of ice such as high-density amorphous (HDA) ice 26 and various crystalline high-pressure forms such as ice II, VI, VII and VIII 27 . The effect on adding NaCl salt is a decrease of the post-edge spectral feature and increase in the pre-edge and main-edge features similar to the changes with temperature 22,28 ; adding NaCl salt has furthermore been shown to affect water structure similar to increased pressure, that is, reducing the tetrahedrality 29 (see next section). Figure 2b shows that the occupied non-bonding 1b 1 lone-pair in X-ray emission spectroscopy (XES) is split into two components where the intensities show dependence on temperature 30,31 and on adding NaCl 32 . Here, again there is no broadening or major shift in the two spectral components, only intensity changes. The origin of the split has been discussed as due to two different instantaneous structural environments giving different electronic emission energies in the spectroscopic process 31,33,34 or giving different O-H stretch dynamics in the core ionized intermediate state 35,36 . Through resonant excitations corresponding to various parts of the XAS spectrum, a direct connection has been shown between the two spectroscopies where the low-energy 1b 1 corresponds to the post edge, and that at high energy correlates with the pre-and main edges fully consistent with the temperature and NaCl concentration dependences 31,33 . Curve fitting of the XAS and XES spectra at ambient temperature gives around 20-30% locally tetrahedral molecules 25,31 , which is consistent with estimations in two-state thermodynamic models 15,37 .
Structure and dynamics probed by vibrational spectroscopies. We next describe trends from vibrational spectroscopy in both the OH stretch region (Fig. 2c) 38,39 and in the low-energy vibrational modes (Fig. 2d) 40 . Important information about H-bond dynamics has been obtained using two-dimensional spectroscopy 41,42 , but here we discuss linear spectroscopy to connect the trends with those obtained from X-ray spectroscopies. Often the focus has been on the OH (OD) stretch spectrum in HDO (HDO is water with one proton (H) substituted by deuterium (D)) in D 2 O (H 2 O), which effectively decouples the oscillator from the environment and makes the OH (OD) group a local probe of the H-bonding 41 . In this case only a broad spectral feature without fine structure is observed and the temperature dependence shows a redistribution of intensity from one side towards the other, opening for many different interpretations in terms of both homogeneous and heterogeneous distributions 14,39,41 . However, when the coupling to the surrounding liquid is turned on in neat H 2 O, the spectrum becomes different due to coupling of resonances through H-bonds with neighbouring molecules 41,43,44 .
The OH stretch vibration in H 2 O has been studied for many years 45 but recently Sun 38,39 compared measurements over an extended range of temperatures for pure water and a range of concentrations of NaCl solution at fixed temperature and found consistent behaviour of the well-resolved spectral feature on the low-frequency side of the OH spectrum close to the broad spectral distribution in hexagonal ice 46 . We note that this  40 , involving low-energy vibrations in the H-bonding network, there are clearly identified signatures of two components in the low-frequency region with the same temperature dependence as the other spectroscopies (Fig. 2d). A further indication of heterogeneities in water is given by recent pump-probe vibrational spectroscopy, which finds a significantly weaker coupling to other chromophores when exciting on the blue side of the spectrum than on the red 47 . Consistent with this weaker interaction, vibrationally resolved XES spectra 48 , using resonant excitation at the pre-edge in XAS, result in an OH-stretch frequency very near that of gas phase water supporting the interpretation of the pre-edge as due to molecules with broken or weakened donated H-bonds 24,25 . In this resonantly excited case, only the high-energy spectral component was observed in the 1b 1 valence band region 16 .
Connecting X-ray and vibrational spectroscopy information.
Thus, the X-ray and vibrational spectroscopies show similar trends with spectral features at fixed energy and only changes in intensity with temperature and NaCl concentration. These observations are consistent with two distinct, different structural classes where the population in each varies with temperature and salt concentration.
As it is well known that increasing temperature 49 and NaCl concentration both remove the tetrahedral coordination 29 , it is clear that the increase in the post-edge shoulder, the low energy 1b 1 component, the 3,200 cm À 1 OH stretch and the 225 cm À 1 low-energy vibrational regions point to an environment connected to tetrahedral structures. These features are also fully aligned with the corresponding features of hexagonal ice. If the tetrahedral structures with strong H-bonds are collective and involve several other tetrahedral water molecules, as discussed in the next section, it could potentially explain the appearance of the 3,200 cm À 1 stretch vibration (not seen in HDO diluted water) as due to coupling to neighbouring tetrahedral water molecules 41 . On the contrary, the pre-edge and main-edge XAS resonances, the high-energy 1b 1 component, the 3,400 cm À 1 OH stretch and the 180 cm À 1 low-energy vibrational regions can be related to non-tetrahedral or more distorted structures. These features completely dominate when tetrahedral order has been removed.  Magnitude of distortions in disordered structures. An important question is how to determine the magnitude of distortions that could give rise to the pre-edge feature in XAS, as it is clearly related to the blue part of the OH stretch spectrum and to the high-energy component in the split 1b 1 peak in XES. One way to investigate such distortions around the first shell is through theoretical energy decomposition of H-bonds in ab initio molecular dynamics (MD) simulations, as recently done by Khaliullin and Kühne 13,50 . This reveals a significant asymmetry in terms of H-bond strength, both for the two donor and for the two acceptor bonds. Figure 2e shows the asymmetry parameter Y for acceptor and donor bonds, where Y ¼ 0 when the two H-bond strengths are equal and Y ¼ 1 when the two bonds are extremely asymmetric with one nearly broken and the other very strong 13 .
In the simulations an asymmetry of around Y ¼ 0.5 is predominant, which in terms of geometry corresponds to small deviations in the H-bond length, and thus a predominance of tetrahedral structures. However, there is a fraction of molecules in the upper right corner of Fig. 2e with very high asymmetry that Khaliullin and Kühne demonstrated contribute high pre-edge peak intensity in the computed XAS spectrum of water 25 . In this particular ab initio MD simulation related to Fig. 2e, the fraction of highly asymmetrical species was small but the simulation is overstructured compared with experiment 13,50 , indicating that in real water this fraction should be significantly larger.
Simulations and X-ray scattering The most direct way to obtain structural information in liquids is from the pair-distribution functions (PDFs) that measure correlations at various distances. Let us first visit a simulation using the TIP4P/2005 force field 51 to inspect how the O-O PDF (Fig. 3a) appears for different structural environments distinguished based on the local structure index (LSI) 52 . What characterizes tetrahedral structures is a well-defined separation of the first shell at 2.8 Å and the second shell at 4.5 Å (high LSI), whereas structures that are less tetrahedral contain interstitials that fill in the region between 3 and 4 Å (low LSI) 53 . The interstitials are typically seen for high-pressure ices due to inwards collapse of the second shell. Very similar results using an LSI classification have been obtained by Car and colleagues 54 based on advanced ab initio MD simulation techniques. From neutron-scattering experiments similar O-O PDFs have been derived by extracting data at various temperatures and pressures, and extrapolating to pure water phases, to correspond to extremely tetrahedral environments, or with a high degree of interstitials that were denoted LDL and HDL phases, respectively 17,18 . This notation comes from that the interstitial structures typically appear when applying pressure to generate a high-density phase in either water or ice. The same occurs, as discussed in the previous section, when NaCl is added to water where the 4.5 Å correlation disappears and instead the number of interstitials increases. The tetrahedral structures, on the other hand, create more open space and are therefore denoted a low-density phase. We note that for the low-LSI PDF there is only a first shell peak at around 2.8 Å but the rest is rather unstructured. This shows that the interstitial-related environment is much more disordered, but the first peak is also much lower and broader, indicating that also the first coordination shell is much more distorted, which is in line with the discussion of the distorted structures in the preceding section.
Pair-distribution functions. The ability to perform X-ray-scattering experiments on water, both in the low- 11 49). This shows that as the first shell peak height decreases, correlations are shifted from 2.8 to the 3.0-3.3 Å region, meaning that interstitials build up not only from the collapse of the second shell but also from distortions of the first shell. Figure 3b shows the temperature dependence of the positions of the first two shells as extracted from the O-O PDFs 49 . As expected, the radial distance of the first shell increases linearly with increasing temperature, as the liquid expands with increasing temperature. The second shell follows the same trend with a similar increase in distance up to a temperature around 320 K, but then there is a dramatic increase and the second peak in the PDF becomes more undefined. This change occurs around the same temperature as the minimum in the isothermal compressibility, shown in Fig. 1. There are some even more striking changes in the liquid if we look at shells at longer distances. Figure 3c shows that the new high-quality data clearly resolve even up to eight shells with the farthest at around 17 Å and they show very distinct differences in their temperature dependence. Previous work has indicated that the feature around 11 Å can be affiliated with local regions of tetrahedral structures 53,57 . Figure 3d shows the temperature dependence of the peak height of the correlation at 11 Å. At high temperatures the height is near zero but then the 11 Å correlation starts to increase around 340-320 K as the liquid is cooled. This shows that it is only around the compressibility minimum that H-bonding becomes decisive and collective tetrahedral fluctuations begin to result in well-defined regions that have a radial extent around 11 Å. This explains the observation in Fig. 3b in terms of the second shell, as it is at 320 K that the tetrahedral structures become well defined as collective regions, whereas at higher temperatures these are much less well defined.
Small-angle X-ray scattering. Such a spatial extent of heterogeneities has earlier been proposed based on small-angle X-ray scattering measurements of water where a temperaturedependent enhancement of the structure factor at low momentum transfer Q has been reported 16,53,55 . Such an enhancement is typically seen for nanoparticles and aggregates in solutions allowing determination of sizes and shapes, but this should also be applicable to density fluctuations, as these are significantly slower than the attosecond X-ray photon-scattering time and an average spatial extent of B10 Å at ambient conditions was proposed based on the measured correlation length B3.1 Å (ref. 16). This interpretation was criticized by Clark et al. 58,59 and Soper et al. 60 who instead proposed that the data were fully consistent with a homogeneous description. More insight has recently been obtained through work by Overduin and Patey 57 and Wikfeldt et al. 61,62 , where TIP4P/2005 force-field simulations show trends qualitatively consistent with the experimental data and provide information regarding the spatial distribution and extent of the structures in the liquid consistent with the originally proposed spatial extent 16 . It is noteworthy that the concept of correlation length is related to the asymptotic decay of the correlation function and the size range of regions is nearly a factor 4 larger 11 , which is consistent with the development of the experimentally observed 11 Å correlation. What is interesting is that the enhancement at low Q appears around 320 K on cooling, which correlates well with the changes in Fig. 3d. This is even more clearly seen in the comparison of the total and normal contributions to the isothermal compressibility as a function of temperature in the work by Overduin and Patey 57 . At high temperatures the normal component is nearly equal to the total and both decrease with decreasing temperature as expected for a simple liquid (Fig. 1a). However, close to the minimum around 320 K the normal component begins to deviate from the total with another anomalous component becoming increasingly significant as the liquid is cooled further 57 .
Anomalous properties and long-range correlations. All these data indicate that, at the compressibility minimum around 320 K, water becomes anomalous not only in terms of the thermodynamic response functions but also in the structure with local regions of tetrahedral structures B11 Å in average radial extent appearing as fluctuations in the distorted/interstitial structurally dominated liquid; these furthermore grow to even larger sizes on further cooling. Here it is important to realize that many-body effects are known to give rise to H-bond cooperativity effects where water molecules prefer to bind with molecules that are in a similar surrounding, leading to local regions 57,63 . At temperatures above 320 K, where structural fluctuations are rare and water behaves as a simple liquid, the structure is more homogeneous.

Supercooled water
In the supercooled regime the nomenclature HDL and LDL has some historical connection with the LLCP hypothesis that water could exist as HDL and LDL macroscopic phases 19 . This links strongly to the glassy state of water in terms of HDA and low-density amorphous ices, and the transformation between them 7,64 . As the structural characteristics of LDL and HDL phases are related to, respectively, a well-defined peak at 4.5 Å and the presence of interstitials around 3-4 Å (ref. 17), we will, for simplicity, use this notation also for local structures appearing as fluctuations with O-O PDFs similar to those of the macroscopic phases, although the local density as such is not a well-defined property.
Validity of a two-state model. The very recent works of Holten and Anisimov 37 , and Russo and Tanaka 15 building complete thermodynamic equations-of-state for water show the validity of a two-state picture with fluctuations between HDL and LDL local structures appearing in the transition between the two phases. These equations-of-state reproduce the temperature dependence of the thermodynamic response functions shown in Fig. 1 and indicate that an LLCP at positive pressure is necessary to describe the data. Here we will use some of this modelling together with new experimental data in the deeply supercooled regime to demonstrate consistency with the LLCP and discuss trends with respect to simulations.
Two liquid phases. There has recently been an intense debate regarding the validity of the two-state model with a liquid-liquid transition (LLT) and a LLCP at positive pressure 20,[65][66][67] . This was based on simulations, as no experiments can yet sample this region of the metastable water phase diagram due to rapid ice crystallization. Figure 4a shows the free-energy surface from Palmer et al. 20 using the ST2 water model 68 for which coexisting metastable basins of HDL and LDL are found at the thermodynamic conditions corresponding to an LLT. There is a free-energy barrier between the metastable LDL and the thermodynamically stable cubic ice phase, indicating a LLT and LLCP in the ST2 model. By modifying the H-bond angular flexibility in ST2 it is furthermore possible to shift the LLCP to become located in the stable part of the phase diagram with respect to ice formation, confirming the validity of a LLCP in atomistic simulations 69 . Based on the ST2 model, Holten et al. 70 estimated the fraction of LDL local structures as a function of temperature. The resulting curves indicate the existence of an LLT as shown in Fig. 4b. The derivative becomes infinite at the LLT and LLCP, and becomes smaller as the pressure is lowered into the one-phase region where the change in population very far from the LLCP approaches a linear relationship. This clearly demonstrates that the fluctuations between LDL and HDL (see Box 1) extend far out from the critical point. Similar work by Russo and Tanaka 15 using the TIP4P/2005 force field and an order parameter to find the fraction of the LDL component (there denoted locally favoured structure) showed almost identical trends as for the ST2 model with only a shift in temperature and pressure. Analysis of the TIP4P/2005 data shows that the derivative becomes infinite at a point that seems to indicate the existence of LLT and LLCP 15 but it remains to be determined whether the LLCP is real or only virtual 37 in the TIP4P/2005 model.
Probing water in 'No-Man's Land'. It has been difficult to experimentally probe water deep in the supercooled regime due to rapid ice crystallization. Recently Sellberg et al. 71 conducted an experiment at the Linac Coherent Light Source X-ray laser, where they could very rapidly cool micron-size water droplets and probe the liquid structure with X-ray scattering as shown in Fig. 4c. The lowest temperature that was achieved while maintaining the liquid state was 227 K, which is 5 K below the previous limit of homogeneous ice nucleation 72 defining the upper onset of what has become denoted 'no-man's land'. Figure 4d shows the development of the tetrahedrality of the liquid, based on the height (g 2 ) of the second O-O PDF shell at 4.5 Å, as function of temperature with an accelerated change towards a tetrahedral or LDL dominated liquid at the lowest temperatures. The change in the temperature-dependent g 2 slope is larger and occurs at slightly higher temperatures than the change in TIP4P/2005 water and is much larger than for the SPC/E model.
As it is the LDL component that gives the g 2 parameter as the excess above one in the second shell (see Fig. 3a where HDL is an almost flat line), they are directly related. The TIP4P/2005 model gives a similar trend as ST2 in Fig. 4b with the slope of the change becoming infinite at the LLCP at 1,350 bar 15 . In the experiment we would expect a near discontinuity in g 2 (as indicated with the schematic LLT line in Fig. 4d) if the approach were to an LLT or if an LLCP were located at negative 73 or ambient pressure. It is noteworthy also from Fig. 4b that in this case the LDL fraction would remain nearly constant coming from high temperatures until very close to the discontinuity and thus also g 2 versus temperature would remain constant. In contrast, the measured g 2 in Fig. 4d exhibits at near-ambient pressure a significant and increasing slope more reminiscent of going through the one-phase region in Fig. 4b at pressures lower than the critical. Noting that the measured change becomes more accelerated at low temperatures than in the TIP4P/2005 simulation, we expect the LLCP of real water, if it exists, to be at much lower pressure than in the simulation. From this we can conclude that a rough estimate for the LLCP should lie somewhere between ambient pressure and the TIP4P/2005 LLCP type of behaviour, maybe in the range 500-1,500 bar. This would be within the range of estimates based on measurements on emulsified water by Mishima 74 and decompression-induced melting of ice IV 75 but new experiments are needed to determine whether an LLCP really exists or whether it is only virtual. In the latter case it would cause fluctuations as if there would be an LLCP but with ice crystallization occurring too rapidly for criticality to fully develop 76 .

Proposed unified picture
Here we want to use the latest developments to propose a simple interpretation giving what we consider a more unified qualitative picture in Fig. 5a of water incorporating a broad temperature range and extending the discussion from the supercooled regime also to ambient temperatures, as the anomalies persist up to B320 K. There could exist two separate liquid phases, HDL and LDL, with a coexistence line in the P-T diagram deep in the supercooled regime and at elevated pressure 7,19 . This LLT line ends with decreasing pressure in an LLCP and its extension into the one-phase region corresponds to the Widom line 77 . At the Widom line the density fluctuations would reach a maximum consistent with equal population of molecules in HDL and LDL 62 . It is worth noting that HDL is on the ambient-temperature side of the Widom line, whereas LDL is on the low-temperature side. This explains why the high-density, H-bond-distorted species dominate at ambient conditions. The origin of the anomalous properties of water is the increase in structural fluctuations, as water is cooled down and approaches the Widom line, leading to fluctuations into tetrahedral patches growing in size as directional H-bonding becomes relatively more dominant 16,78 . Figure 5b schematically illustrates the temperature dependence of a thermodynamic response function (k T or C P ) at various pressures relative to the Widom line. The response function shows a maximum at the Widom line with the height and width dependent on the distance to the critical point. The fluctuating funnel-like anomalous region (Fig. 5a) around the Widom line becomes wider the further away from the LLCP we are in pressure. The closer we are, the narrower the region becomes; however, in this region the anomalous behaviour becomes more strongly enhanced, leading to steeper rise in the response functions 79 .
We can relate to the fundamental discussion about a heterogeneous contra homogeneous distribution of structures (Box 1) by inspecting more closely the funnel-like anomalous region beyond the critical point. Outside we have the macroscopic phases of HDL and LDL, each a simple liquid that would correspond to the homogeneous description with thermal fluctuations and distortions around only one class of structural motifs. Changes in temperature do not affect the populations consistent with Fig. 4b, where at 2 kbar and above 260 K the tetrahedral LDL fraction is low and near constant. Inside the anomalous region we are in the region of heterogeneous fluctuations; however, as we get further away from the LLCP the size of the fluctuations becomes smaller. Eventually, we will approach a limit where the sizes reach molecular dimensions as in an ideal mixture with many intermediate structures resulting in an almost flat line in the populations as in the top figure of Box 1. This would correspond to an almost linear dependence in the tetrahedral LDL population with temperature as in the 0 bar case in Fig. 4a or further extrapolated to negative pressures. As MD force-field simulations give an LLCP at higher pressures than what is discussed here for real water, the resulting structures at ambient pressure would be closer to the ideal mixture limit. We should also note that this is a simplified picture and the boundary at the funnel-like anomalous region is not sharp, and there will be an extension of small contributions of LDL/HDL local structures but with more of an ideal mixture-like dispersion into the simple liquid regions. This is why tetrahedral structures are seen in the spectroscopic data also at temperatures above the here-indicated anomalous region but without resulting in spatially separated regions.
What is unique with water is that at ambient pressure the location of the LLCP is such that the anomalous region with fluctuations extends up to around 320 K (47°C). This means that water is anomalous at temperatures where life is sustained and where most processes of importance in nature and to our society occur. It would be interesting in the future to gain insight into whether this is a pure coincidence or has significant implications for understanding biology.

Future perspective
Coming back to Box 1, a burning question is to determine the degree of heterogeneity in water and the sharpness of the boundaries between the fluctuating regions. This will rely partly on theoretical developments to describe the observed spectroscopic information. Would a large boundary region with more average structure between the two regions show a clear spectral contribution or would it correspond to a projection onto both spectral components? We note here that, in the inherent structure (quenched to 0 K) in simulations, the distribution of the LSI order parameter is highly bimodal with a minimum between HDL-and LDL-like character 54,62,80 ; however, the distributions are not completely separated, indicating molecules in environments of mixed character. In the real structure in simulations (including temperature) the distribution is on the other hand continuous, indicating a more homogeneous liquid with more smeared-out heterogeneities (moving more towards the right in Box 1). It will be important to further investigate the difference in inherent and real structure using models that better include cooperativity effects and bring the LLCP towards lower pressure than TIP4P/ 2005 or ST2. Having the LLCP at lower pressure would shift the funnel in Fig. 5a and make the onset of fluctuations at ambient pressure more distinct, which would improve, for example, the MD isothermal compressibility, which rises too slowly 81 and the g 2 parameter in Fig. 4d.
A challenge that requires more experimental developments is to investigate the 'no-man's land' region more deeply and determine whether or not an LLCP and an LLT exist and whether there is a free-energy barrier to rapid ice formation. There could be a virtual LLCP defined as a point of instability in the P-T phase diagram giving rise to fluctuations in its neighbourhood as indicated in Fig. 5a, but where too rapid ice nucleation prevents real criticality in terms of diverging correlation length 76 . Here, the experiments recently demonstrated at Linac Coherent Light Source show an interesting avenue if they can be conducted also under pressure, a most challenging task. Another way is to come up from below 'no-man's land' starting from the glassy state. It has been seen that starting from low-density amorphous or HDA ice leads to different temperatures of the glass transition 82 . The question is whether there is a real liquid state with translational motion at these low temperatures.
There has been a large interest in recent years in H-bond dynamics through the development of infrared pump-probe experiments and two-dimensional spectroscopies 41,42,83,84 . It is important to understand which part of the dynamics can be related to the anomalous contribution from the HDL-LDL fluctuations that causes the enhancement of k T at low temperatures. These fluctuations are naturally related to more collective motions of many water molecules to form and dissolve local LDL regions, which is beyond the simple H-bond dynamics probed through the local OH stretch mode. There are indeed indications from three-dimensional infrared spectroscopy that subassemblies of water molecules at ambient conditions have different dynamics on timescales longer than 0.5 ps 84 . Related to this, a relaxation time of the order of 1.5 ps has been derived from high-resolution inelastic scattering measurements of ambient water 85 . However, it is not yet known whether this dynamics can be related to the anomalous fluctuations connected to the negative slope of the temperature dependence of k T shown in Fig. 1. In particular, it is expected that the dynamics will change abruptly on passing the Widom line in Fig. 5a where the more rigidly H-bonded structure in an LDL-dominated liquid could slow down the dynamics by several orders of magnitude and maybe give rise to the speculated fragile-to-strong liquid transition 86 , that is, from non-Arrhenius to Arrhenius-like temperature dependence of the viscosity. Here, it is also important to develop techniques that probe the true equilibrium dynamics in a probe-probe scheme instead of pump-probe.
Water represents a severe challenge to theoretical simulations due to the delicate balance between different counteracting interactions of similar magnitude and also to the emergent character of the resulting properties, which requires a large number of molecules to develop. It is clear that highly accurate computational approaches are required to reliably capture this delicate balance, but these need also to be applicable to large-enough systems and long-enough trajectories to allow fluctuations on relevant length and time scales to develop. Classical force-field simulations can easily be applied to very large boxes for rather extended simulations, but with questionable accuracy. Recent force-field developments with parameterization of up to three-body interactions at the CCSD(T) level of quantum chemical accuracy hold promise to remedy this 87 , but whether this approach is sufficient to capture all aspects of the water phase diagram remains to be seen. In another approach, VandeVondele and colleagues 88 have exploited algorithmic and computer developments to perform Monte Carlo simulations of 64 water molecules at the fully ab initio periodic quantum chemical MP2 level; however, considering the experimental observation of an eighth radial coordination shell at B17 Å, significantly larger simulation boxes will be necessary, which will be a challenge to this approach. Simulations based on Density Functional Theory (DFT) 89,90 scale better with size of the system but a box size of 30-40 Å still represents a severe challenge. Quantum effects are furthermore important 91,92 and have been proposed to be decisive in determining the difference in entropic contribution from HDL and LDL 11 . Including quantum effects adds additional cost to simulations, but recent developments in terms of coloured noise have reduced the expense of including this important effect 93 . An additional challenge to simulators is to predict or extract reliably experimental observables from the simulation where, in particular, techniques to reproduce the X-ray spectroscopies are currently being under strong development both within DFT 89,94,95 and traditional quantum chemistry 96,97 .