Many-body contributions to polymorphism and polyhexaticity in a water monolayer

TL;DR

Confined water monolayers exhibit polymorphism and hexatic ordering whose stability is governed by hydrogen-bond cooperativity. The authors extend the Franzese–Stanley model by introducing three-body ($J_ heta$) and five-body ($J_\sigma$) HB-MBIs and map the $P$–$T$ phase diagram using $NPT$ Monte Carlo; they find low-density square ice (LDi) and high-density triangular ice (HDi) separated from the liquid by hexatic phases (LDh, HDh), with a line of density maxima (TMD). The three-body term promotes crystallization, while five-body cooperativity stabilizes hexatic fluctuations, producing polyamorphism and a nonmonotonic specific-heat signature that converges toward a liquid–liquid critical point (LLCP) only under suitable MBIs balance. These results demonstrate that HB-MBIs are central to the phase behavior of confined water, with wide implications for nanofluidics, interfacial science, biology, and pharmaceutical technologies where water under confinement governs function.

Abstract

Nanoconfined water plays a crucial role in nanofluidics, biology, and cutting-edge technologies. The process of melting water monolayers and quasi-two-dimensional confined water involves, as an intermediate stage, the hexatic phase--a state that lies between solid and liquid and is characterized by quasi-long-range orientational order and short-range translational order. However, the influence of hydrogen bond (HB) cooperativity in this process has not been thoroughly investigated. This gap hampers our understanding of the phase behavior of confined water and limits the accuracy of our models. To address this, we extend the water model developed by Franzese and Stanley, which explicitly includes many-body interactions (MBIs) of HBs. We distinguish the contributions of three-body and five-body HB-MBIs. Our Monte Carlo calculations in the isobaric-isothermal ensemble produces a detailed pressure-temperature phase diagram, revealing polymorphism and polyhexaticity: low-density square ice and high-density triangular ice are separated from the liquid phase by distinct hexatic phases. Three-body interactions notably promote crystallization and can destabilize the low-density hexatic phase, while cooperative five-body interactions help restore it, thus modifying the thermodynamic landscape. These findings demonstrate that HB-MBIs are key in determining the phase behavior of confined water, influencing phenomena such as the non-monotonic specific heat, maximum density lines, and the accessibility of the liquid-liquid critical point. Beyond advancing theoretical understanding, these results have wide-ranging implications for nanofluidics, interfacial science, and applications in biology, food technology, and pharmaceutics, where controlling water under confinement is essential.

Figures (7)

• Figure 1: Polymorphism and polyhexaticism in a water monolayer.$P-T$ phase diagram with parameters $J/4\epsilon=0.5$, $J_{\sigma} /4\epsilon= 0.1$ and for: (a)$J_{\theta}/4\epsilon= 0.08$, (b)$J_{\theta}/4\epsilon= 0.10$, and (c)$J_{\theta}/4\epsilon= 0.12$. Other parameters are detailed in Table \ref{['Table1']}. We identify the liquid (white region), LD hexatic (LDh, pale blue), HD hexatic (HDh, pale red), LD ice (LDi, blue), and HD ice (HDi, red) phases, as well as the line of temperatures of maximum density (TMD, line with open circles) in the liquid phase. The gas phase at higher temperatures is not shown. The boundaries of the hexatic phases (marked by lines with open squares and triangles) are determined through structural analysis, as detailed in SI (./supplementary.pdf). Letters in (a) denote the state points considered in text.
• Figure 2: Structural changes in a water monolayer at different state-points. Radial distribution function $g(r)$ for the model with $J_{\theta}/4\epsilon =0.08$ evaluated at: (Top)$Pv_0/4\epsilon=0.05$ and (A) $Tk_B/4\epsilon=0.15$ in the liquid phase, (B) $Tk_B/4\epsilon=0.11$ in the LDh phase and (C) $Tk_B/4\epsilon=0.06$ in the LDi phase; (Bottom) at (D) $Pv_0/4\epsilon=0.2$ and $Tk_B/4\epsilon=0.18$ in the liquid phase, (E) $Pv_0/4\epsilon=0.5$ and $Tk_B/4\epsilon=0.17$ in the HDh phase and (F) $Pv_0/4\epsilon=0.8$ and $Tk_B/4\epsilon=0.12$ in the HDi phase. Insets: Configurations of the molecules' center of masses corresponding to the state points (Top) A, B, C, and (Bottom) D, E, F. The letters refer to the state-points marked in Fig.\ref{['fig:PhaseDiagram']}a.
• Figure 3: Specific heat decomposition of a water monolayer at low pressure. For $J_{\theta}/4\epsilon=0.08$ at $Pv_0/4\epsilon =0.05$, we calculate the specific heat $C_P$ (dashed line with gray triangles) from enthalpy ($\mathcal{H} +PV$) fluctuations and decompose it into its components: i) the van der Waals component $\mathcal{H}_{\epsilon} +PV$ (yellow open squares); ii) the two-body HB component $\mathcal{H}_{J}$ (blue squares); iii) the three-body component $\mathcal{H}_{3}$ (green circles); iv) the first coordination shell (up to) five-body component $\mathcal{H}_{5}$ (red open circles). The lines are guides for the eye. The star marks the minimum in $C_P$. The shaded regions mark the LD ice (blue), LD hexatic (pale blue), the liquid (white), and the gas phase (yellow).
• Figure 4: Locus of maxima of specific heat in the $P$-$T$ phase diagram. The line of $C_P$ maxima (green with full squares) and the maximum in the two-body HB component of $C_P$(green open squares with gray line) converges at $Pv_0/4\epsilon \simeq 0.25$ and $Tk_B/4\epsilon \simeq 0.115$ (purple circle), regardless of $J_{\theta}$. This point resembles a critical point at the end of a line of first-order phase transitions. The first-order phase-transition line coincides with the $C_P$-maximum line at high pressure. The TMD (blue) line is shown as a reference. Model parameters are as in Fig.\ref{['fig:PhaseDiagram']}: $J_{\theta}4\epsilon=0.08$ in (a), $0.10$ in (b), and $0.12$ in (c). The lines serve as guides to the eye.
• Figure 5: Correlation functions $g_{\alpha}(r)$ and $g_{\vec{G}}(r)$ (solid lines) and their corresponding fitting functions (dashed lines), computed at different state points. The state point labels refer to those indicated in Fig. 1a of the main text: A at $(Tk_B/4\epsilon=0.15, Pv_0/4\epsilon=0.05)$, C at $(Tk_B/4\epsilon=0.06, Pv_0/4\epsilon=0.05)$, D at $(Tk_B/4\epsilon=0.18,Pv_0/4\epsilon=0.2)$, F at $(Tk_B/4\epsilon=0.12, Pv_0/4\epsilon=0.8)$. a) 4-fold bond-orientational correlation function in A (red) and C (gray); b) 6-fold bond-orientational correlation function in D (red) and F (gray); c) translational square correlation function in A (red) and C (gray); d) translational triangular correlation function in D (red) and F (gray).