{'fraction_non_alphanumeric': 0.03468282619907169, 'fraction_numerical': 0.0221763795771016, 'mean_word_length': 4.7078734363502575, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Atomistic View of Homogeneous Nucleation of Water into Polymorphic Ices

Water is one of the most abundant substances on Earth, and ice, i.e., solid water, has more than 18 known phases. Normally ice in nature exists only as Ice Ih, Ice Ic, or a stacking disordered mixture of both. Although many theoretical efforts have been devoted to understanding the thermodynamics of different ice phases at ambient temperature and pressure, there still remains many puzzles. We simulated the reversible transitions between water and different ice phases by performing full atom molecular dynamics simulations. Using the enhanced sampling method MetaITS with the two selected X-ray diffraction peak intensities as collective variables, the ternary phase diagrams of liquid water, ice Ih, ice Ic at multiple were obtained. We also present a simple physical model which successfully explains the thermodynamic stability of ice. Our results agree with experiments and leads to a deeper understanding of the ice nucleation mechanism.

Water-to-ice transition is still a mystery 1 which puzzles scientists working in biochemistry 2 , physics 3 and environmental science 4-6 . Understanding the ice interaction with biomacromolecule is vital to functional protein design, such as antifreeze proteins 7 , ice-nucleating proteins 8 or organic polymers 9 . The phase diagram of ice is complicated with more than eighteen solid phases 10 , but only the hexagonal ice (ice Ih, ABAB... stacking 11 ) and the cubic ice (ice Ic, ABCABC… stacking 11 ) are present in nature. At ambient pressure, water freezes into Ih rather than Ic 12 . Ice Ic only exists partially at low temperatures and high pressure 13 . In fact, pure ice Ic without stacking defects was not experimentally accessible until 2020 10 . In most cases, water freezes into stacking disordered ice (ice Isd) consisting of hexagonal and cubic structures [14][15][16][17][18] .
Molecular dynamics (MD) simulation is an efficient tool in seeking the ice stacking mechanism on the molecular level. However, with the limitation imposed by the large conformation space and the timescale 19,20 in crystallization, a direct all-atom simulation of ice nucleation from liquid water is almost impossible 21 . The only one successful case of spontaneous nucleation with all-atom model is the work of Ohmine and co-workers 22 in 2002, which has been difficult to be reproduced. Therefore, many groups used instead enhanced sampling methods to help simulate ice nucleation 23-28 . However, reversible transition between water and ice is still a great challenge. Recently, Niu and co-workers 29 realized the reversible water-ice transition in all-atom MD simulations using the enhanced sampling method MetaITS 30 with a carefully designed set of collective variables (CVs), which includes the X-ray diffraction (XRD) peak intensities 31 . The XRD peak intensity is an efficient CV for the phase transition between the solid and liquid states, and the simulation of reversible transition has been succeeded using the metadynamics with only one XRD peak intensity as CV for the systems of Na, Al 32 , Si 33 , etc. It turned out that the phase transition between water and ice is much more difficult than these covalent crystals. Although a complex combination of seven X-ray diffraction peak intensities and the translation entropy of water molecules were used as CVs, the solid state obtained in Ref. 28 is a mixed phase of stacking disordered ices. An efficient sampling of different ice phases using MD simulation thus remains as a challenge.
In this paper, we performed MD simulation with the all-atom TIP4P/Ice water model. The TIP4P/Ice model 34 provides one of the best predictions for the melting point in comparison to the experiment data, and was widely used in homogeneous nucleation of water 33,35,36 . Rather than the coarse-grained models, like monotonic water (MW) potential 23,25,[37][38][39][40][41] , the all-atom TIP4P/Ice model provides precise dynamics behaviours as well as H-bond connection and dipole properties. Our simulation system contains 1152 water molecules placed in a periodic box of 3.14 nm * 2.72 nm * 4.45 nm, which can match the lattice of ice Ih, Ic or their stacking disordered states (See Section S-Ⅰ for detailed information).
We adopt the enhanced sampling methods MetaITS to facilitate the phase transition between water and ice. MetaITS 30 is a hybrid method that combines metadynamics 42 (MetaD) and integrated tempering sampling 43 (ITS), which increases the sampling efficiency of MD simulation by accelerating the CVs related degrees of freedom (See Section S-Ⅲ for detailed information). The CVs we used in this work are the two intensities of 3D and 2D XRD peaks 33 : Due to the anisotropy of the CV 2 , the plane of sphere packing can only stack at the Z-axis of the simulation system, and up to 12 ice bilayers can be formed at this size of periodic box. Since each bilayer may form a hexagonal or cubic structure, fully iced systems can form the states with different hexagonal/cubic ratios. In the simulation, seven different hexagonal/cubic ratios of structures were observed. The sampled structures as well as the corresponding cubicity (the fraction of cubic ice, using a ternary definition 44 ) are shown Fig. 1(b-i) and Section S-Ⅶ in SI. with different values of cubicity (Fig. 1). We also calculated the ternary phase free energy diagram of the liquid state, ice Ih, and ice Ic (Fig.2b). Similar to the FES in Fig.   2a, this phase diagram also shows 7 ice states as distinguished by their different values of cubicity. To further investigate the thermodynamics properties of the different ice state, we calculated the free energy change for water to ice transition, Δ water→ice , as a function of cubicity from 260K to 280K in Fig. 3a. This figure clearly exhibits that there are 7 stable states with the cubicity C* = 0%, 16.7%, 33.3%, 50%, 66.7%, 83.3% and 100%, which correspond to the seven states in Fig. 1: Ih, Isd(2h10c), Isd(4h8c), Isd(6h6c), Isd(8h4c), Isd(10h2c), and Ic, respectively.
where ΔHcubicity is the free energy difference between bulk hexagonal and cubic ices 11 , ΔHsurface is the free energy difference between ice and liquid water 45 , which is zero pure solid ice, Hfault is the penalty of the mismatched interface, ΔSlayer = R  ln(Ω(C*)) by is the entropy of layer, and Ω( * ) is the number of combinatory way of hexagonal and cubic layers in the system. The main difference between the current and earlier theories is in the means of Ω( * ) calculation. It appears that the current model is in clearly better consistency with the simulation results. We calculated the free energy of ice as a function of cubicity at 270 K using Eq. 3 and show the results in Fig. 3b A number of thermodynamic properties including the ternary phase diagram are calculated from the simulation trajectories, which agree with literature. Our simulations also show that in the simulation system containing a maximum of 12 ice bilayers, seven, but not thirteen (as predicted by earlier theories), stable ice states can exist. We propose an alternative model, which satisfactorily explains the stable structure of ice that may exist in the system, especially in the periodic boundary conditions. We also propose a method to estimate the free energy of various states of ice based on this simple model.
The corresponding calculation provides a qualitatively explanation on the stability of different states of ice in the simulation system. Further studies, including those on the kinetics of ice nucleation and simulation of larger systems etc., are currently underway, which hopefully will lead to further understanding of the fascinating physical properties of ice.