2D ice from first principles: structures and phase transitions
Ji Chen, Georg Schusteritsch, Chris J. Pickard, Christoph G. Salzmann, Angelos Michaelides
22D ice from first principles: structures and phase transitions
Ji Chen,
1, 2
Georg Schusteritsch,
2, 3, 4
Chris J. Pickard,
2, 3, 4
Christoph G. Salzmann, and Angelos Michaelides
1, 2, 5, ∗ London Centre for Nanotechnology,17-19 Gordon Street, London WC1H 0AH, U.K. Thomas Young Centre, University College London,20 Gordon Street, London, WC1H 0AJ, U.K. Department of Physics and Astronomy, University College London,Gower Street, London WC1E 6BT, U.K. Department of Materials Science and Metallurgy, University of Cambridge,27 Charles Babbage Road, Cambridge CB3 0FS, U.K. Department of Chemistry, University College London,20 Gordon Street, London, WC1H 0AJ, U.K.
Abstract
Despite relevance to disparate areas such as cloud microphysics and tribology, major gaps inthe understanding of the structures and phase transitions of low-dimensional water ice remain.Here we report a first principles study of confined 2D ice as a function of pressure. We find thatat ambient pressure hexagonal and pentagonal monolayer structures are the two lowest enthalpyphases identified. Upon mild compression the pentagonal structure becomes the most stable andpersists up to ca. a r X i v : . [ c ond - m a t . m t r l - s c i ] A ug onfined and interfacial water-ice is ubiquitous in nature, playing an important role ina wide range of areas such as rock fracture, friction, and nanofluidics [1–4]. As a resultof a delicate balance of forces (hydrogen bonding, van der Waals, and interaction withthe confining material or substrate) confined and interfacial water forms a rich variety ofstructures [5–10]. Almost every specific system examined has revealed a different structuresuch as a 2D overlayer built from heptagons and pentagons on a platinum surface or thesquare ice observed within layers of graphene [8, 10]. This shows that in contrast to bulkice the phase behavior of 2D ice is much less well understood.From an experimental perspective a full exploration of the phase diagram of 2D ice has notbeen achieved yet. However recent experiments revealed the exciting possibility of exploring2D ice structures at specific conditions by trapping water within layered materials [2, 10, 11].For example by confining water between layers of graphene it is possible to create so callednanocapillaries in which water experiences a pressure in the GPa regime due to the vander Waals forces pulling the sheets together [10]. Using transmission electron microscopy(TEM) square ice structures from a single up to a few layers were observed in such graphenenanocapillaries [10], and in the multilayer regime square ice layers located directly on top ofone another in an AA stacking arrangement. Force field simulations performed as part of thesame study of confined water in graphene layers reproduced the square monolayer ice butfailed to explain the AA stacking [10]. Indeed prior to this recent study there was alreadya long tradition of computer simulation studies of nanoconfined water, mostly involvingclassical force field approaches [12–27]. Such work has been incredibly valuable and providedconsiderable insight into the structures and phase transitions of monolayer and multilayer ice.The 2D ice structures predicted include hexagonal, pentagonal, quasicrystalline, hexatic, andorthogonal phases as well as various amorphous structures. Although the relative stabilityof the structures depends on the particular force field used, for monolayer ice an orthogonalphase has been widely predicted to be stable across a broad pressure regime [13, 14, 19, 23].Given the sensitivity of force field studies of confined ice to the potential used, the aspirationto understand the observation of monolayer square ice, and the unexplained AA stackingorder of square ice layers, a systematic study with an electronic structure method such asdensity functional theory (DFT) is highly desirable. DFT studies of ice do not come withouttheir own sensitivity to the exchange-correlation functional used, however, they have provedto be very useful in predicting and understanding structures of adsorbed water and bulk ice.25–9, 28].Here we report a systematic study of 2D phases of water-ice from first principles basedon an unbiased exploration of monolayer ice structures. The much greater computationalcost of DFT compared to force fields means that we cannot map the entire phase diagramof monolayer confined water. Instead we focus on the phase transitions as a function oflateral pressure and confinement width at zero Kelvin. The stable structures identified atdifferent pressures include an hexagonal structure, a pentagonal Cairo tiling (CT) structure,a flat square structure, and a buckled rhombic structure. The observation of a flat squarestructure is consistent with the recent experimental observation. However, the sequenceof low energy phases identified differs significantly from that predicted in force field basedstudies [13, 14, 19, 21, 23] and a recent DFT report [26]. Interestingly the sequence ofstructures observed depends sensitively on the confinement width used, suggesting thatit should be possible to tune the monolayer ice structures produced in experiments by e.g. varying the confining material. We also propose a double layer square ice phase withinterlayer hydrogen bonds, which further explains the AA stacking order of multi-layer icerecently observed [10].In order to explore ice structures in an unbiased manner, we used the ab initio randomstructure search (AIRSS) technique [29, 30], an approach which has previously predicted newice, 2D, and interfacial structures [31, 32]. Structures from ambient up to a lateral pressureof 10 GPa were considered. The 2D confinement was introduced via a Morse potentialfit to quantum Monte Carlo results for the binding of a water monomer to graphene [33].By tuning the confinement width, we are not only able to study the general properties ofwater under flat and smooth confinement but also to compare with the recent experimentsof ice in graphene confinement [10]. Our electronic structure calculations were carried outusing the Vienna Ab-initio
Simulation Package (VASP) [34] with the DFT+vdW approach[35] in conjunction with the Perdew-Burke-Ernzerhof exchange-correlation functional [36].Tests with other exchange-correlation functionals show that whilst the transition pressuresbetween the various phases depend to some extent on the choice of exchange-correlationfunctional, the overall conclusions do not change. See the supporting information (SI) forthese results as well as further computational details [37].From a preliminary set of calculations we established that the optimal separation betweengraphene sheets in which a monolayer of water is sandwiched is somewhere between 6.0 and3 ) b)d) e)c)
FIG. 1. Monolayer ice structures. The top and side views of the hexagonal (a), the Cairo tiling(CT) (b), the flat square (f-SQ) (c), the rhombic (b-RH) (d) and the buckled square (b-SQ) (e)structures. Red and pink spheres represent oxygen atoms at different heights and white spheresare hydrogen atoms. The green boxes show the primitive unit cells. O. This difference drops to just 2 meV/H O when harmonic zeropoint energy effects are taken into account (Table S2, Fig. S9). Tests with other exchange-correlation functionals generally concur that the energy difference between the two phasesis tiny, with a vanishingly small preference for the hexagonal structure (Table S2) [37].Anharmonic zero point energy differences between the two phases and finite temperatureeffects could also easily exceed the energy difference [38], suggesting that both the hexagonaland the CT phases could be observed at ambient pressures.4 Δ H ( m e V / H O ) hexagonalCairo-tiling (CT)flat-square (f-SQ) buckled-rhombic (b-RH)buckled-square (b-SQ) Å Δ H ( m e V / H O ) Å w=6.5w=6.0 FIG. 2. Enthalpies of the water monolayer phases as a function of lateral pressure under 6.0 ˚A (a)and 6.5 ˚A (b) confinement. ∆H is the relative enthalpy with respect to the CT phase. EnthalpyH = E water +E confinement +P × A × w, where E water is the total energy per water molecule, E confinement is the energy (per water molecule) in the confinement potential, P is the lateral pressure, A is thelateral area per water molecule, w is the width of confinement. The hexagonal structure has p6mm wall-paper group symmetry if only oxygen atoms areconsidered and is built exclusively of six-membered rings (Fig. 1a). Water molecules arethree-fold coordinated, with half of them having one OH bond directed out of the monolayer(a so-called “dangling OH”). The average O-O separation in this monolayer hexagonal phaseis about 2.72 ˚A which is similar to the bulk O-O separation in ice I [39]. However the confinedhexagonal structure identified here is quite flat with the vertical separation between oxygenatoms < p4gm (Fig. 1b).5he unit cell has twelve water molecules, four with dangling OH bonds. One third of thewater molecules are four-fold coordinated and the rest are three fold coordinated. The higheraverage coordination and smaller ring size of the CT phase renders the density of this phasehigher than that of the hexagonal phase (Fig. S4, S5) [37]. Therefore, upon increasingthe pressure the CT phase clearly becomes more stable than the hexagonal phase. Aftersearching for the lowest enthalpy structures at finite pressure we find that the CT structureis clearly the most stable in a broad range of pressures all the way up to ca. O than the hexagonal phase at ambient pressure. In the SI wetrace this difference to the different computational setups [37]. As shown in the SI we areconfident that the hexagonal and CT monolayer ice structures are indeed more stable thanany square ice structure at the low pressure limit.At pressures beyond 2 GPa higher density phases obeying the Bernal-Fowler-Pauling(BFP) ice rules are identified more frequently in the structure searches. This includes aflat square phase (f-SQ) and a buckled rhombic phase (b-RH). Both phases consist of four-fold coordinated water molecules with two donor HBs and two acceptor HBs. The f-SQstructure has a p gm wall-paper group symmetry where the dipoles of the water molecules6 Pressure (GPa)0 100.02 2 C on f i n e m e n t w i d t h ( Å ) FIG. 3. Sketch of the phase diagram of monolayer ice with respect to lateral pressure and confine-ment width. For clarity the pressure axis has an artificial scale with the various pressures indicatedcorresponding to the transition pressures obtained at 6 ˚A confinement. Note that at ambient pres-sure the enthalpies of the hexagonal and Cairo tiled phases are essentially degenerate but by ca. are distributed on two orthogonal antiferroelectric sublattices (Fig. 1c). The HB networkof the b-RH phase is similar to f-SQ but it is buckled and has a higher lateral density(Fig. 1d). The relative stability of the f-SQ and b-RH structures in the 2-4 GPa regimedepends sensitively on the confinement width (Fig. 2). At 6.0 ˚A f-SQ is more stable whileat 6.5 ˚A b-RH has the lower enthalpy. Beyond 4 GPa the b-RH phase is more stable thanany other structure identified. Several other metastable structures belonging to the b-RHfamily with different hydrogen bond ordering have also been observed. However, since amore delicate discussion of hydrogen ordering is beyond the scope of this study, we onlyshow one of the most stable members of the b-RH family at the pressure and confinementconditions considered (Fig. 1d).We also identified a second metastable square phase which we dub “b-SQ” because ofits buckled “basketweave”-like pattern of HBs (Fig. 1e). Its lattice structure resemblesa 2D projection of bulk ice VIII but it is unique in that the two sublayers are hydrogenbonding with each other. The energy of the b-SQ phase is higher than the most stablephases identified, however because of its unique hydrogen bonding arrangement and as it7ight be possible to observe it in systems where the substrate has a square lattice we feelit is worth reporting.Beyond the phase behavior of monolayer ice confined within the 6.0 ˚A to 6.5 ˚A regime, wealso explored a broader range of confinement widths at different pressures. These additionalcalculations support the validity of the conclusions reached but also show that there isscope for altering the relative enthalpies of the various phases by tuning the confinement.The phase diagram for monolayer water with respect to lateral pressure (0-10 GPa) andconfinement width (5-8 ˚A) that emerges from these calculations is shown schematically inFig. 3. At small confinements and low pressures ( < ∼ cm − due todangling OH groups and lower frequency stretching modes < cm − arising from strongHBs. The stretching regions also have quite different total widths for the different phasesand the bending mode of the b-SQ phase is softer than the other phases by ca. cm − .Therefore it should also be possible to discriminate one phase from another with vibrationalspectroscopy.Multi-layer square ice was also observed experimentally with a structure in which oxygenatoms are located directly on top of each other, suggesting an AA stacking square phase [10].However, the force field simulations performed as part of the experimental study obtainedan AB stacking arrangement. The difference was attributed to a possible inaccuracy of the8 ) b)c) f-SQb-SQdl-SQ E n t h a l py ( e V / H O ) FIG. 4. Double layer ice energies and structures. (a) Enthalpy per water molecule for differentdouble layer square ice structures as a function of confinement width at 2 GPa. The solid lines forthe f-SQ and b-SQ phases represent an AA stacking whereas the dashed lines are for AB stacking.(b,c) Structure of the dl-SQ phase. Other AA and AB stacking double layer ice structures areshown in the SI [37]. force field models used. In order to clarify this difference between experiment and theory,we calculated the enthalpy for different confinement widths at 2 GPa for a number of AAand AB stacking structures (Fig. 4a). In agreement with the force field simulations we findAB stacking is preferred for the f-SQ phase. However the other square structure identifiedin this study, b-SQ, is more stable with AA stacking, which agrees with experiment. Themain interaction between the layers in these structures is van der Waals, with the HBnetworks within the layers remaining unperturbed. We find, however, that these stackedvan der Waals bonded structures are considerably less stable than double layer structureswith interlayer hydrogen bonds, the most stable of which we dub “dl-SQ”. The structureof dl-SQ is shown in Fig. 4b,c. It is AA stacked and as can be seen from Fig. 4a it is ca.
80 meV/H O more stable than any of the stacked van der Waals bonded double layerstructures. Thus we conclude that the AA stacked double layer square structure observedin experiments is unlikely to be a van der Waals bonded layered structure but rather onestabilised by interlayer hydrogen bonds.In summary, monolayer ice phases and their phase transitions in confinement have beenstudied with DFT and a random structure search approach. At ambient pressure we havepredicted hexagonal and pentagonal Cairo tiled structure, which are similar in enthalpy9nd more stable than other structures. The CT structure becomes more stable than thehexagonal structure when lateral pressure is applied. Upon increasing the pressure to aboveabout 2 GPa high density square and rhombic phases are observed. Looking forward acomplete description of the temperature dependent phase diagram of 2D water is desirable.Experimentally it would be interesting to explore a broader range of temperatures, waterdensities, and confining materials. From the computational perspective, it would be de-sirable to explore 2D ice at higher temperature with electronic structure methods, such as ab initio molecular dynamics. In addition, in light of the small ring sizes and relativelyshort intermolecular separations, 2D ice also provides further opportunities for investigatingcollective proton quantum dynamics [47–49].
ACKNOWLEDGEMENTS
J.C. and A.M. are supported by the European Research Council under the EuropeanUnion’s Seventh Framework Programme (FP/2007-2013) / ERC Grant Agreement number616121 (HeteroIce project). A.M is also supported by the Royal Society through a RoyalSociety Wolfson Research Merit Award. G.S. is supported by the TOUCAN EPSRC grant(EP/J010863/1). C.G.S is supported by the Royal Society (UF100144). We are also grate-ful for computational resources to the London Centre for Nanotechnology, UCL ResearchComputing, and to the UKCP consortium (EP/ F036884/1) for access to Archer. ∗ [email protected][1] D. Chandler, Nature , 640 (2005).[2] K. Xu, P. Cao, and J. R. Heath, Science , 1188 (2010).[3] J. Carrasco, A. Hodgson, and A. Michaelides, Nature Materials , 667 (2012).[4] Y. Gao, S. Kim, S. Zhou, H.-C. Chiu, D. Nlias, C. Berger, W. deHeer, L. Polloni, R. Sordan,A. Bongiorno, and E. Riedo, Nature Materials , 714 (2015).[5] J. Yang, S. Meng, L. F. Xu, and E. G. Wang, Physical Review Letters , 146102 (2004).[6] J. Carrasco, A. Michaelides, M. Forster, S. Haq, R. Raval, and A. Hodgson, Nature Materials , 427 (2009).
7] M. Forster, R. Raval, A. Hodgson, J. Carrasco, and A. Michaelides, Physical Review Letters , 046103 (2011).[8] S. Nie, P. J. Feibelman, N. C. Bartelt, and K. Th¨urmer, Physical Review Letters , 026102(2010).[9] J. Chen, J. Guo, X. Meng, J. Peng, J. Sheng, L. Xu, Y. Jiang, X.-Z. Li, and E.-G. Wang,Nature Communications , 4056 (2014).[10] G. Algara-Siller, O. Lehtinen, F. C. Wang, R. R. Nair, U. Kaiser, H. A. Wu, A. K. Geim, andI. V. Grigorieva, Nature , 443 (2015).[11] Q. Li, J. Song, F. Besenbacher, and M. Dong, Accounts of Chemical Research , 119 (2015).[12] K. Koga, H. Tanaka, and X. C. Zeng, Nature , 564 (2000).[13] R. Zangi and A. E. Mark, Physical Review Letters , 025502 (2003).[14] K. Koga and H. Tanaka, The Journal of Chemical Physics , 104711 (2005).[15] N. Giovambattista, P. J. Rossky, and P. G. Debenedetti, Physical Review Letters , 050603(2009).[16] S. Han, M. Y. Choi, P. Kumar, and H. E. Stanley, Nature Physics , 685 (2010).[17] J. C. Johnston, N. Kastelowitz, and V. Molinero, The Journal of Chemical Physics ,154516 (2010).[18] N. Kastelowitz, J. C. Johnston, and V. Molinero, The Journal of Chemical Physics ,124511 (2010).[19] J. Bai, C. A. Angell, and X. C. Zeng, Proceedings of the National Academy of Sciences ,5718 (2010).[20] J. Bai and X. C. Zeng, Proceedings of the National Academy of Sciences , 21240 (2012).[21] A. L. Ferguson, N. Giovambattista, P. J. Rossky, A. Z. Panagiotopoulos, and P. G.Debenedetti, The Journal of Chemical Physics , 144501 (2012).[22] N. Giovambattista, P. Rossky, and P. Debenedetti, Annual Review of Physical Chemistry ,179 (2012).[23] W.-H. Zhao, J. Bai, L.-F. Yuan, J. Yang, and X. C. Zeng, Chemical Science , 1757 (2014).[24] W.-H. Zhao, L. Wang, J. Bai, L.-F. Yuan, J. Yang, and X. C. Zeng, Accounts of ChemicalResearch , 2505 (2014).[25] Q. Lu, J. Kim, J. D. Farrell, D. J. Wales, and J. E. Straub, The Journal of Chemical Physics , 18C525 (2014).
26] F. Corsetti, P. Matthews, and E. Artacho, arXiv:1502.03750 [physics] (2015), arXiv:1502.03750.[27] F. Corsetti, J. Zubeltzu, and E. Artacho, arXiv:1506.04668 [cond-mat] (2015), arXiv:1506.04668.[28] G. A. Tribello, B. Slater, and C. G. Salzmann, Journal of the American Chemical Society , 12594 (2006).[29] C. J. Pickard and R. J. Needs, Physical Review Letters , 045504 (2006).[30] C. J. Pickard and R. J. Needs, Journal of Physics: Condensed Matter , 053201 (2011).[31] C. J. Pickard, M. Martinez-Canales, and R. J. Needs, Physical Review Letters , 245701(2013).[32] G. Schusteritsch and C. J. Pickard, Physical Review B , 035424 (2014).[33] J. Ma, A. Michaelides, D. Alf`e, L. Schimka, G. Kresse, and E. Wang, Physical Review B ,033402 (2011).[34] G. Kresse and J. Furthm¨uller, Physical Review B , 11169 (1996).[35] A. Tkatchenko and M. Scheffler, Physical Review Letters , 073005 (2009).[36] J. P. Perdew, K. Burke, and M. Ernzerhof, Physical Review Letters , 3865 (1996).[37] See Supplemental Material for more details on our calculations, additional information on theice phases and various tests.[38] E. A. Engel, B. Monserrat, and R. J. Needs, Physical Review X , 021033 (2015).[39] K. R¨ottger, A. Endriss, J. Ihringer, S. Doyle, and W. F. Kuhs, Acta Crystallographica SectionB , 644 (1994).[40] V. F. Petrenko and R. W. Whitworth, Physics of Ice (Oxford University Press, 2002).[41] C. G. Salzmann, P. G. Radaelli, B. Slater, and J. L. Finney, Physical Chemistry ChemicalPhysics , 18468 (2011).[42] S. Zhang, J. Zhou, Q. Wang, X. Chen, Y. Kawazoe, and P. Jena, Proceedings of the NationalAcademy of Sciences , 2372 (2015).[43] C. Vega, J. L. F. Abascal, M. M. Conde, and J. L. Aragones, Faraday Discussions , 251(2008).[44] H. J. C. Berendsen, J. R. Grigera, and T. P. Straatsma, The Journal of Physical Chemistry , 6269 (1987).[45] J. Zhang, P. Chen, B. Yuan, W. Ji, Z. Cheng, and X. Qiu, Science , 611 (2013).
46] J. Guo, X. Meng, J. Chen, J. Peng, J. Sheng, X.-Z. Li, L. Xu, J.-R. Shi, E. Wang, andY. Jiang, Nature Materials , 184 (2014).[47] X.-Z. Li, M. I. J. Probert, A. Alavi, and A. Michaelides, Physical Review Letters , 066102(2010).[48] C. Drechsel-Grau and D. Marx, Physical Review Letters , 148302 (2014).[49] X. Meng, J. Guo, J. Peng, J. Chen, Z. Wang, J.-R. Shi, X.-Z. Li, E.-G. Wang, and Y. Jiang,Nature Physics , 235 (2015)., 235 (2015).