Nucleation and Growth of Hexagonal Ice by Dynamical Density

Nov 16, 2016 - Synopsis. By considering the enhanced and ordered hydrogen bonding in crystal, a dynamic density functional theory has been proposed to...
1 downloads 14 Views 747KB Size
Subscriber access provided by University of Otago Library

Article

Nucleation and Growth of Hexagonal Ice by Dynamical Density Functional Theory Yi Ye, Nanying Ning, Ming Tian, Liqun Zhang, and Jianguo Mi Cryst. Growth Des., Just Accepted Manuscript • DOI: 10.1021/acs.cgd.6b01289 • Publication Date (Web): 16 Nov 2016 Downloaded from http://pubs.acs.org on November 16, 2016

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Crystal Growth & Design is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 17

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

Nucleation and Growth of Hexagonal Ice by Dynamical Density Functional Theory Yi Yea,b, Nanying Ningb, Ming Tiana,b*, Liqun Zhanga,b, and Jianguo Mia* a

State Key Laboratory of Organic-Inorganic Composites, Beijing University of Chemical Technology, Beijing 100029, China b

Key Laboratory of Beijing City on Preparation and Processing of Novel Polymer Materials, Beijing University of Chemical Technology, Beijing, China

We propose a dynamic density functional theory (DDFT) to describe the transition process from disordered liquid water to ordered ice crystal. With the special consideration of the enhanced and ordered hydrogen bonding in crystal and the hydrodynamic interactions induced by the flow of water molecules, the theory is able to reproduce the process of nucleation and growth of ice crystal. Some theoretical predictions are in general agreement with the available computational and experimental data, indicating that the theory provides a reasonable description of the dynamic process of ice crystal growth.

Corresponding authors: [email protected]; [email protected].

ACS Paragon Plus Environment

Crystal Growth & Design

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 17

1. Introduction The freezing of water to ice is arguably the most important liquid-to-solid transition. Although a rich and complex phase diagram of ice has been well-established experimentally, understanding of the subtle balance of intermolecular interactions which give rise to this richness is incomplete. The common ground in any field of ice research is the urge to understand the ordering mechanisms of phase transition, interfacial structure and dynamics.1-3 In this regard, the molecular-level understanding of ice crystallization plays an important role in a number of fields. Over the past decade and a half, a multitude of experiments4,5 and molecular simulations6-8 for the liquid-to-solid phase transitions have been successfully performed to understand the relative contribution of hydrogen bonding and van der Waals dispersion forces to the cohesive properties of water. However, our understanding of the freezing is far from complete. In experiments, the homogeneous nucleation needs to avoid heterogeneous nucleation with small droplets (10−100 µm). In molecular simulations, the quiescent period is too long to be probed in direct simulation processes at the experimental temperatures. From a theoretical perspective, the varied densities of crystal structure, with the molecules fixed on well-defined lattices, afford an excellent opportunity to quantify and assess the contribution of these intermolecular interactions.9,10 Understanding how ice grows has recently been identified as one of the top open questions in the development of crystallization theory.11 A recent phase field crystal model has been gaining widespread recognition, which captures the essential concepts of crystal nucleation and growth.12,13 In spite of these successes, the mesoscale model overlooks most of the relevant atomic scale physics that leads to defect interactions and coarse-grained boundary nucleation and migration. For a more realistic molecular force field, the model lacks a generalized free-energy formulation that allows the study of important phase transformation. In liquid water, the hydrogen-bond association between individual water molecules yields a disordered hydrogen-bond network ACS Paragon Plus Environment

Page 3 of 17

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

and a large number of possible network configurations.14 In ice crystal, however, the hydrogen bonding transfers to an ordered network. It is obvious that the phase field crystal model cannot describe the hydrogen-bond induced crystallization morphology and process. As a consequence, it is still an open challenge to develop a reliable crystallization theory to reproduce the freezing of ‘real’ water into ice. In this work, we propose a DDFT approach based on the Smoluchowski equation of motion for obtaining insight into the fundamental factors that control the nucleation and growth of ice crystal. From the thermodynamic point of view, the hexagonal ice has the most stable structure at moderate environment; hence we assume that the crystal has the hexagonal morphology. The density functional is constructed with the fundamental measure theory for hard-spherical repulsion and the direct correlation function derived from atomic interactions for dispersive and hydrogen-bond attractions. In particular, the effect of hydrodynamic interactions induced by water molecule flow on the freezing process is also considered in the density functional. The theory is expected to provide the required molecular resolution with the potential to yield structural insights, dynamics, and estimates of thermodynamic quantities.

2. Theory and equations By constituting the fundamental, nonlinear, deterministic equation for the time-evolution of dimensionless density ρ (r, t ) under the liquid-flow mediated hydrodynamic interactions, the DDFT equation can be generalized as the following form15

∂ρ ( r, t ) ∂t

 δ F  ρ ( r, t )  δ F (ρ )  = D∇  ρ ( r, t ) ∇ + ∫ dr′ρ (2) (r, r′, t )ω(r − r′) ⋅∇  δρ ( r, t ) δρ (r′, t )  

(1)

where t is the time evolution, D is the diffusion constant of 1.18 × 10 −9 m2/s,16 ρ (r) denotes the density distribution with configuration r,

ρ (2) (r, r′, t ) stands for the two-body density with

ρ (2) (r, r′, t ) = ρ (r, t ) ρ (r′, t ) g ( r − r′ ) and g (r ) = 1 , ω (r − r ′) is the Rotne-Prager hydrodynamic tensor,17

ACS Paragon Plus Environment

Crystal Growth & Design

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 17

F [ ρ (r, t )] accounts for the intrinsic Helmholtz free-energy of the system. We begin by separating the total free energy into four components15,18

F  ρ ( r, t )  = kBT ∫ dr ρ ( r, t ) ln [ ρ (r, t )] − 1 + F rep [ ρ (r, t )] + F att [ ρ (r, t )] + ∫ drρ (r, t )V ( r, t )

(2)

representing the contributions of ideal reference, hard-sphere repulsion, attractive interaction, and external potential, respectively. For the hard-sphere repulsion, the fundamental measure theory19 is postulated to have the form

F rep [ ρ ] = ∫ dr Φ hs [nm (r )]

(3)

and Φ hs [ nm (r ) ] is the free-energy density modified with the tensor version of White Bear II

Φhs [ nm (r)] = −n0 ln(1− n3 ) +

n1n2 − nV1 ⋅ nV 2 3[n ⋅ n ⋅ n − n n ⋅ n − tr(n3 ) + n2tr(nT2 )] ϕ1(n3 ) + ϕ2 (n3 ) V 2 T V 2 2 V 2 V 2 2 T 1− n3 16π (1− n3 )

(4)

where nm (r) is the weighted density, n0 , n1, n2 , n3 are scalar weighted densities, nV 1 , nV 2 are vector weighted densities, and nT is tensor weighted density. The attractive contribution due to the van der Waals interactions and hydrogen bonding in disordered water can be integrated into F att [ ρ ( r , t ) ] = −

1 d r d r ′′ρ ( r , t ) ρ ( r ′′, t ) c OattO ( r - r ′′ , ρ 0 ) 2∫ ∫

(5)

att ( r -r′′ , ρ0 ) represents the attractive part of the direct correlation function of oxygen−oxygen where cOO

12 6 under the force field u(r ) = 4ε OO [(dOO / rOO ) − (dOO / rOO ) ] + 1/ 4πε 0 ∑∑ qi q j / r , and has been accurately i

j

determined20 by solving the integral equation21 with the TIP4P/ice force field.22 In liquid water, the number of hydrogen bonds per molecule is about 3.2 in the vicinity of room temperature.23,24 In crystalline ice, however, each molecule is hydrogen bonded to exactly four others in an approximately tetrahedral arrangement.25 Therefore, roughly a fifth of the increase of association strength contributes to the phase transition. Such contribution acts as the external force V (r, t ) = −

1 ′ ∑ H (r ′, t )Cassoc OO (r − r , ρ 0 ) 5

ACS Paragon Plus Environment

(6)

Page 5 of 17

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

where H (r ′, t ) is the Heaviside function, H (r ′, t ) = 0 for ρ (r′, t ) < ρp / 2 , and ρp is the peak of local assoc density in the crystal cell, otherwise H (r ′, t ) = 1 , COO (r − r′, ρ0 ) is a vector version of the direct assoc assoc ′ ′ ]cOO ( r − r′ , ρ0 ) , arising directly from the ordered correlation function, Cassoc OO (r − r , ρ 0 ) = δ [(r − r ) − r

hydrogen bonding. In the theoretical calculations, the crystal lattice parameters should be determined in advance. The density profile in a crystal cell can be accurately parametrized by a sum of Gaussian density distributions around the lattice sites

ρ s ( r ) = (α / π )

3/ 2

∑ e α( −

r −R )

2

(7)

where α is a measure of the width of the Gaussian peaks, R is the crystal lattice vector in real-space, R (i, j , k ) = ia + jb + kc + c′ (i, j, and k are the integer pillar’s indices in the directions parallel and

perpendicular to the array’s edge, c ′ is a nonideal distortion of the c-axis). In order to simplify the process for determining of lattice parameters, we apply the reported angles of the staggered molecules, which are ∠(a, b) = 60° , ∠(a, c) = ∠(b, c) = 90° , and ∠(c, c′) = 0° or 180° .26 Apart from these temperature-independent angles, other crystal parameters in terms of temperature variation are calculated by minimizing the free-energy of crystal cell. The hydrogen-bond length is obtained from the radial distribution function of oxygen–hydrogen.20 Accordingly, these lattice parameters can be written as: 2

a = a(1/ 2, 3 / 2,0) , b = a(1/ 2, − 3 / 2,0) , c = ma(0,0,1) , ma = 4 c′ + a 2 + 4 c′ , here

a 2 + 4 c′

2

corresponds to the distance between the closest oxygen atoms, which is the sum of the intramolecular hydrogen–oxygen bond length (0.95 Å) and intermolecular hydrogen-bond length (1.79 Å). At given temperature, the free-energy of ice crystal is interpolated as a function of the lattice constants a and α .27 The interpolation is performed with a quadratic spline. A representative free-energy surface is shown in Figure 1 for a unit cell. The values are calculated at the intersections of solid lines, and the surface shows the interpolated free-energy. The array of points used in the free-energy interpolation is spaced sufficiently fine ACS Paragon Plus Environment

Crystal Growth & Design

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 17

that the interpolation is smooth, and a well-defined minimum exists. After determination of the minimum free-energy, the crystal structure is then recalculated through a full minimization of the grand potential

Ω [ ρ (r )] = F [ ρ (r ) ] − ∫ drρ (r ) µ , where µ is the chemical potential of the bulk ice with µ = δ F [ ρ b ] / δρb . The iteration procedure is repeated until the average fractional difference over any three-dimensional grid point of ( ρold − ρnew ) / ρold is less than 10−3 . As a result, the finial lattice vectors of lengths are:

a = b = 2.55 Å, c = 3.60 Å, c′ = 0.45 Å, and the average density of ice crystal is 0.92 g/cm3.

Figure 1. Interpolated free-energy surface for a unit crystal at 273.15 K. Free-energy values at the intersections of the black curves are calculation results. Surface coloring is a guide to the energy minimum.

3. Results and discussion DDFT allows for a deterministic expression of the time development of the density profile. At time zero, we run a dynamic process starting from a stable fluid phase with thermal fluctuations, which are externally imposed as an external potential and then switched off after a unit cell has been formulated.28 Here we choose the Brownian timescale τ B = σ 2 / D of 8.5 ×10−11 s as the unit timescale.15 Figure 2 illustrates the crystal formation process at the degree of supercooling ∆T = 25 K, corresponding to 248.15 K. Three typical snapshots along the path of nucleation are depicted in the inset of Figure 2. The free-energy profile and snapshots clearly show the metastable water phase and the thermodynamically stable ice phase. At ACS Paragon Plus Environment

Page 7 of 17

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

initial stage ( t = 0.1τ B ), disordered water molecules are reorganized to formulate compact nucleus. At the growing state ( t = 0.3τ B ), a unit cell has been developed surrounded with immature crystal cells, and the increase of free-energy could be caused by the rearrangement of crystal cells. The saddle point of the free-energy at the transition state ( t = 0.8τ B ) is referred to the critical nucleus. In view of the height of the nucleation free-energy barrier ∆Ω [ ρ (r )] = Ω [ ρ (r) ] − Ω b for the cluster, one sees that the value is 192.4

kBT . In this case, the critical nucleus radius Rc is 2.2 nm. The two values are in good agreement with the simulation ones of 175 kBT and 2.3 nm, respectively.29 When the ordered structure has been reformulated, a quick relaxation of the density field in the neighborhood of the imposed nucleation array allows the nucleus begins to grow with sharp increase of radius.

Figure 2. Time-evoluted free-energy and two-dimensional cuts of density distribution in (0001) crystal plane at different moments of ice crystal nucleation at 248.15 K.

Employing the forward method, we study homogeneous nucleation of ice crystal from supercooled liquid water in the region of 258.15−233.15 K. The nucleation free-energy barriers and critical radii at different temperatures are summarized in Figure 3. In general, free-energy barrier provides an estimate of the concentration of critical clusters in the metastable fluid as ρ0 exp(−∆Ω / kBT ) , and ρ0 is the density of the liquid water.29 In this regard, bigger energy barrier corresponds to lower concentration of critical clusters. As ACS Paragon Plus Environment

Crystal Growth & Design

the energy barrier increases, the volume of water necessary for nucleating a critical cluster increases. Once the volume is larger than the one of the whole hydrosphere, it is virtually impossible to observe homogeneous nucleation. It is clear that both the energy barrier and the critical radius decrease as the temperature decreases, and finally reach a plateau, indicating that a further decline of temperature has only trivial effect on the nucleation behavior. 4.0 500 3.5

Critical nucleation size Nucleation free-energy barrier

3.0

400

300 2.5

200

2.0

100

1.5

15

∆Ω/kBT

Rc/nm

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 17

20

25

30

35

40

∆T/K

Figure 3. The nucleation free-energy barriers and the critical radii for ice crystal nucleation with different degrees of supercooling. The spots are the simulation data.29

Within a unit cell of growing crystal, the calculated density contour is plotted in Figure 4(a), which shows more pronounced features on the staggered grids, and the positions and heights of the nearest neighbor peaks are essentially invariant. Moreover, each water molecule displays aspherical density contour, which can be attributed to the overlapping of hydrogen-bond association. The lower clouds of intermediate density connecting neighboring water molecules report the direction of hydrogen-bond network and the three-dimensional topology of hexagonal column ice. A two-dimensional cut of density contour in the (0001) crystal plane through one layer of the lattice sites is displayed in Figure 4(b). It shows the peak of local density in the range of 4 to 5 times the density of bulk crystal, which is also close to the simulation value of 4.5 times.30

ACS Paragon Plus Environment

Page 9 of 17

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

Figure 4. (a) The contour plot of density in a unit cell by the rendered gradient colors with heavier color indicating higher packing. (b) Two-dimensional cut of local density distribution in (0001) crystal plane.

Information on the crystal–liquid interface plays a central role in understanding the initial stage of crystallization. In Figure 5, the considerable variations in three-dimensional structure for the basal interface can be seen clearly. Hence in much of the broad interfacial region, molecules have some liquid-like and some crystal-like properties. The density oscillations imply that ordered arrangements of water molecules exist in the broad interfacial region. As evident upon careful inspection of the profiles, the peaks in density profile decline gradually from the crystal side to the liquid side. One sees that the basal interfaces have the peaks with about four interfacial layers in both (0001) and (1 100 ) crystal planes. In particular, the predicted interfacial tensions of these two crystal planes for basal orientations of ice are γ (0001) = 28.6 mN/m and

γ (1 100) = 30.8 mN/m, which are similar to the experiment value γ = 29.1 ± 0.8 mN/m.31

ACS Paragon Plus Environment

Crystal Growth & Design

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 17

Figure 5. Two-dimensional cuts of density profiles at both (0001) –liquid and (1 100 )–liquid interfaces at 273.15 K. The crystal and liquid sides are given to display the whole interfaces.

Based on the properties of Fourier transform, the density distribution in real-space can be easily transformed to wave vector resolved information in the reciprocal space. The positions and intensity variations of the bright facula are characteristics of the crystalline species and each one of them originates from a particular set of crystal lattice planes. This is similar to the principle of x-ray scattering,32 which has been widely used to determine the configuration of molecules. As a result, we can read the information from the Fourier transform density distributions as from the X-ray scattering intensity maps. Qualitatively, the amorphous water molecules have X-ray patterns which show only broad diffuse halos, as shown in Figure 6(a). We then observe an increase of the intensity as time evolution (Figure 6(b)), and find more and more Bragg peaks (Figures 6(c) and 6(d)). The Bragg peaks are expected to display approximately a six-fold symmetry in a crystalline array, corresponding to the ordered water molecules in ice crystal.

ACS Paragon Plus Environment

Page 11 of 17

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

Figure 6. Two-dimensional cuts of density profiles of the growing crystal in the reciprocal space. The four subfigures correspond to the states at t = 0.1τ B , 0 .3τ B , 0.5τ B , 0.8τ B , respectively.

According to the calculated grand potentials at various temperatures, the homogeneous nucleation rates

J [cm-3s-1] can be quantitatively determined through the following equation33 J = A exp( −∆Ω cri / k BT )

(8)

where ∆Ωcri is the free-energy barrier for nucleation, A is the kinetic pre-factor, and can be defined as 1/ 2

A=

24 ρ0 Dλ (2 Rc )1/6  ∆µ    l2 6π kB  T 

1/2

 ∆µ  = A0 Rc1/6    T 

(9)

Here ∆µ is the chemical potential difference between liquid and solid, A0 is a kinetic constant in liquid with the value of ln A0 = 114 .32 Figure 7 shows the calculated nucleation rates of homogeneous ice. According to the conception of “impossible nucleation rate”, 29 we can assure that the water cannot transform to ice with homogeneous nucleation process when the temperature is higher than 248 K. This phenomenon refers to the fact that, in a situation with little heterogeneous nuclei, moderately supercooled water or drink will shows the liquid properties in bottles, and transforms to ice when poured into cup since some heterogeneous nuclei have been added or formed in the system. By combining the figure with Figure 3, we can predict that the nucleation rate cannot infinitely increase, since the free-energy barrier and nucleation size would reach a plateau at low enough temperatures. In the high supercooling range ( ∆T ≥ 25K ), it is shown that the theoretical results are 2−4 orders of magnitude below the experimental ones.34 The reason for the discrepancy lies in the structure and morphology of ice crystal. In the theoretical framework, the crystal is assumed as the perfect structure with ordered stacking, whereas the crystals in experiment stack disorderly.

ACS Paragon Plus Environment

Page 12 of 17

0

-3 -1

-50

Theoretical Calculations 32 Experiments by Kramer et al. -100

36

-3 -1

log10(J/(m s ))

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

log10(J/(m s ))

Crystal Growth & Design

-150

-200 15

20

25

∆T/K

∆T/K

37

8

6

30

35

40

Figure 7. Predicted nuclei number as a function of the degree of supercooling. The spots represent the corresponding experimental data. The horizontal dashed line is called “impossible nucleation rate”.29 Above the line, the homogeneous nucleation becomes possible. Insets: the same nuclei number as a function of temperature.

4. Conclusion In summary, we have presented a new DDFT model to describe the ice crystal nucleation and growth. With the given angles of different molecules in the staggered arrangements, the cell parameters have been determined through a full minimization of the grand potential. By considering the enhanced and ordered hydrogen bonding as well as the hydrodynamic interactions, we have explicitly modeled the process of disorder-to-order phase transition and crystal growth. The results indicate that, indeed, it is impossible for homogeneous ice nucleation at temperature lower than 248 K. Compared to the experimental measurements, the theoretical calculations underestimate nucleation rates merely by 2–4 orders of magnitude, showing a verified high predictive capability. This work is expected to provide a general theoretical model for understanding ice crystal formation at microscopic scale. Further investigation is to optimize the lattice parameters without any prerequisite upon the combination of charge density distribution.

Acknowledgements ACS Paragon Plus Environment

Page 13 of 17

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

This work is supported by the National Natural Science Foundation of China (Nos. 21476007 and 51525301), and by Chemcloudcomputing of Beijing University of Chemical Technology.

References (1) Pirzadeh, P.; Kusalik, P. G. On Understanding Stacking Fault Formation in Ice. J. Am. Chem. Soc. 2010, 133, 704−707.

(2) Pamuk, B.; Soler, J. M.; Ramírez, R.; Herrero, C. P.; Stephens, P. W.; Allen, P. B.; Fernández-Serra, M.-V. Anomalous Nuclear Quantum Effects in Ice. Phys. Rev. Lett. 2012, 108, 193003. (3) Bonnet, N.; Marzari, N. Static Dielectric Permittivity of Ice from First Principles. Phys. Rev. Lett. 2014, 113, 245501.

(4) Kuhs, W. F.; Sippel, C.; Falenty, A.; Hansen, T. C. Extent and Relevance of Stacking Disorder in “Ice Ic”. P. Natl. Acad. Sci. USA 2012, 109, 21259−21264. (5) Mishima, O. Volume of Supercooled Water under Pressure and the Liquid-liquid Critical Point. J. Chem. Phys. 2010, 133, 144503.

(6) Engel, E. A.; Monserrat, B.; Needs, R. J. Anharmonic Nuclear Motion and the Relative Stability of Hexagonal and Cubic Ice. Phys. Rev. X 2015, 5, 021033. (7) Pfrommer, B. G.; Mauri, F.; Louie, S. G. NMR Chemical Shifts of Ice and Liquid Water: the Effects of Condensation. J. Am. Chem. Soc. 2000, 122, 123−129. (8) Kienle, D. F.; Kuhl, T. L. Density and Phase State of a Confined Nonpolar Fluid. Phys. Rev. Lett. 2016, 117, 036101.

(9) Archer, A. J.; Evans, R. Dynamical Density Functional Theory and its Application to Spinodal Decomposition. J. Chem. Phys. 2004, 121, 4246−4254.

ACS Paragon Plus Environment

Crystal Growth & Design

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 17

(10) Zimmermann, U.; Smallenburg, F.; Löwen, H. Flow of Colloidal Solids and Fluids through Constrictions: Dynamical Density Functional Theory versus Simulation. J. Phys.: Cond. Mat. 2016, 28, 244019. (11) Bartels-Rausch, T.; Bergeron, V.; Cartwright, J. H. E.; Escribano, R.; Finney, J. L.; Grothe, H.; Gutiérrez, P. J.; Haapala, J.; Kuhs, W. F.; Pettersson, J. B. C.; Price, S. D.; Ignacio Sainz-Díaz, C.; Stokes, D. J.; Strazzulla, G.; Thomson, E. S., Trinks, H.; Uras-Aytemiz, N. Ice Structures, Patterns, and Processes: A View Across the Icefields. Rev. Mod. Phys. 2012, 84, 885. (12) Elder, K. R.; Katakowski, M.; Haataja, M.; Grant, M. Modeling Elasticity in Crystal Growth. Phys. Rev. Lett. 2002, 88, 245701.

(13) Berry, J.; Grant, M. Modeling Multiple Time Scales during Glass Formation with Phase-Field Crystals. Phys. Rev. Lett. 2011, 106, 175702.

(14) Matsumoto, M.; Saito, S.; Ohmine, I. Molecular Dynamics Simulation of the Ice Nucleation and Growth Process Leading to Water Freezing. Nature 2002, 416, 409−413. (15) Rex, M.; Löwen, H. Dynamical Density Functional Theory with Hydrodynamic Interactions and Colloids in Unstable Traps. Phys. Rev. Lett. 2008, 101, 148302. (16) Krynicki, K.; Green, C. D.; Sawyer, D. W. Pressure and Temperature Dependence of Self-diffusion in Water. Faraday Discussions of the Chemical Society 1978, 66, 199−208. (17) Rotne, J.; Prager, S. Variational Treatment of Hydrodynamic Interaction in Polymers. J. Chem. Phys. 1969, 50, 4831−4837. (18) van Teeffelen, S.; Backofen, R.; Voigt, A.; Löwen, H. Derivation of the Phase-Field-Crystal Model for Colloidal Solidification. Phys. Rev. E 2009, 79, 051404.

ACS Paragon Plus Environment

Page 15 of 17

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

(19) Oettel, M.; Görig, S.; Härtel, A.; Löwen, H.; Radu, M.; Schilling, T. Free Energies, Vacancy Concentrations, and Density Distribution Anisotropies in Hard-Sphere Crystals: A Combined Density Functional and Simulation Study. Phys. Rev. E 2010, 82, 051404. (20) Wang, X.; Wang, S.; Xu, Q.; Mi, J. Thermodynamics of Ice Nucleation in Liquid Water. J. Phys. Chem. B 2015, 119, 1660−1668.

(21) Schweizer, K. S.; Curro, J. G. Integral-Equation Theory of the Structure of Polymer Melts. Phys. Rev. Lett. 1987, 58, 246.

(22) Abascal, J. L. F.; Sanz, E.; Fernández, R. G.; Vega, C. A Potential Model for the Study of Ices and Amorphous Water: TIP4P/Ice. J. Chem. Phys. 2005, 122, 234511. (23) Smith, J. D.; Cappa, C. D.; Wilson, K. R.; Messer, B. M.; Cohen, R. C.; Saykally, R. J. Energetics of Hydrogen Bond Network Rearrangements in Liquid Water. Science 2004, 306, 851−853. (24) Kumar, R.; Schmidt, J. R.; Skinner, J. L. Hydrogen Bonding Definitions and Dynamics in Liquid Water. J. Chem. Phys. 2007, 126, 204107.

(25) Prendergast, D.; Galli, G. X-ray Absorption Spectra of Water from First Principles Calculations. Phys. Rev. Lett. 2006, 96, 215502.

(26) Malkin, T. L.; Murray, B. J.; Brukhno, A. V.; Anwarc, J.; Salzmannd, C. G. Structure of Ice Crystallized from Supercooled Water. P. Natl. Acad. Sci. USA 2012, 109, 1041−1045. (27) Howard, M. P.; Milner, S. T. Self-Consistent Real Space Free Energy Calculations for Polyethylene and Isotactic Polypropylene Crystals and Interfaces. Macromolecules, 2014, 47, 8335−8350. (28) Neuhaus, T.; Schmiedeberg, M.; Löwen, H. Compatibility Waves Drive Crystal Growth on Patterned Substrates. New Journal of Physics 2013, 15, 073013.

ACS Paragon Plus Environment

Crystal Growth & Design

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 17

(29) Sanz, E.; Vega, C.; Espinosa, J. R.; Caballero-Bernal, R.; Abascal, J. L. F.; Valeriani, C. Homogeneous Ice Nucleation at Moderate Supercooling from Molecular Simulation. J. Am. Chem. Soc. 2013, 135, 15008−15017. (30) Santra, B.; Klimeš, J.; Alfè, D.; Tkatchenko, A.; Slater, B.; Michaelides, A.; Car, R.; Scheffler. M. Hydrogen Bonds and van der Waals Forces in Ice at Ambient and High Pressures. Phys. Rev. Lett. 2011, 107, 185701.

(31) Hardy, S. C. A Grain Boundary Groove Measurement of the Surface Tension between Ice and Water. Philos. Mag. 1977, 35, 471−484.

(32) Vos, W. L.; Megens, M.; van Kats, C. M.; Bösecke, P. X-ray Diffraction of Photonic Colloidal Single Crystals. Langmuir 1997, 13, 6004−6008. (33) Li, T.; Donadio, D. Russo, G., Gallicd, G. Homogeneous Ice Nucleation from Supercooled Water. Phys. Chem. Chem. Phys. 2011, 13, 19807−19813.

(34) Krämer, B.; Hübner, O.; Vortisch, H.; Wöste, L.; Leisner, T.; Schwell, M.; Rühl, E.; Baumgärtel, H. Homogeneous Nucleation Rates of Supercooled Water Measured in Single Levitated Microdroplets. J. Chem. Phys. 1999, 111, 6521−6527.

ACS Paragon Plus Environment

Page 17 of 17

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

For Table of Contents Use Only

Nucleation and Growth of Hexagonal Ice by Dynamical Density Functional Theory Yi Ye, Nanying Ning, Ming Tian, Liqun Zhang, and Jianguo Mi

By considering the enhanced and ordered hydrogen bonding in crystal, a dynamic density functional theory has been proposed to describe the transition dynamics from disordered liquid water to ordered ice crystal. Verified by computational and experimental data, the theory is expected to provide a general method for understanding ice crystal formation at microscopic scale.

ACS Paragon Plus Environment