Freezing and Melting Transitions under Mesoscalic Confinement

Subscriber access provided by UNIV OF NEBRASKA - LINCOLN

Article

Freezing and Melting Transitions under Mesoscalic Confinement: Application of the Kossel-Stranski Crystal-Growth Model. Daria Kondrashova, and Rustem Valiullin J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/jp510467t • Publication Date (Web): 27 Jan 2015 Downloaded from http://pubs.acs.org on February 1, 2015

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.

The Journal of Physical Chemistry C 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 37

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

The Journal of Physical Chemistry

76x66mm (96 x 96 DPI)

ACS Paragon Plus Environment


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

Freezing and melting transitions under mesoscalic confinement: Application of the Kossel-Stranski crystal-growth model. D. Kondrashova and R. Valiullin∗ University of Leipzig, Institute for Experimental Physics I, Leipzig, Germany E-mail: [email protected]

∗

To whom correspondence should be addressed

1 ACS Paragon Plus Environment

Page 2 of 37

Page 3 of 37

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


Abstract The Kossel-Stranski model of crystal growth has been adopted to study freezing and melting of fluids in mesoporous materials. The model is found to exhibit the key features observed in the experiments, including shifted solid-liquid and liquid-solid transition temperatures, irreversibility between freezing and melting, and strong impact of the pore geometry. By first analyzing fluids confined to cylindrical pores, we obtain several important insights into the transition mechanisms. In particular, we establish the conditions for the occurrences of equilibrium and metastable transitions and derive exact analytical equations interrelating the transition temperatures, confinement size and the interaction parameters of the Kossel crystal. Variation of the channel diameter along the channel axis, mimicking disorder in real porous materials, is shown to strongly affect both freezing and melting. The model predicts that the freezing transition in disordered materials is governed by the pore blocking mechanism. The melting transition is found to result as an interplay of two different transition mechanisms.

Introduction It is well documented that freezing and melting transitions for fluids confined to porous solids are found to occur at temperatures deviating from the bulk, equilibrium transition temperatures T0 . 1,2 Using a broad spectrum of experimental techniques, including differential scanning calorimetry, 3,4 nuclear magnetic resonance cryoporometry 5–8 and relaxometry, 9–11 X-ray 12–14 and neutron 15–17 scattering techniques, various aspects of these phenomena have thoroughly been addressed. In particular, these studies focused on establishing correlations between the confinement size, pore geometry, surface wetting and transition temperatures and on understanding irreversibility between freezing and melting. In earlier works, the experiments were mostly performed using disordered pore networks, such as porous sol-gel glasses 18,19 or controlled pore glasses. 20–22 Complex pore structure of these materials made, however, direct comparison of the experimental results to theoretical pre2 ACS Paragon Plus Environment


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 37

dictions questionable. The recent advent of mesoporous materials with well-ordered pores has notably contributed to advancements in the fundamental knowledge about confined fluids. They provided important options to either completely decouple pure confinement effects from those arousing from complexities of the pore geometries 13,23–31 or to control the degree of their coupling using materials with ink-bottle pore geometries. 32,33 Notably, the deeper understanding of fluid behavior on a single-pore level allowed more elaborate analyses of the experimental data obtained for disordered materials. Despite this progress brought about by the availability of the ordered materials, some aspects of thermodynamics of fluids confined in disordered porous materials have still remained scarcely studied. 34–39 The shifts of the liquid-solid and solid-liquid transition temperatures to lower values observed in porous materials as compared to that in the bulk state are typically associated with the excess surface free energy of the confined crystals. Due to the common phenomenon of ’interfacial melting’, i.e. the formation of a thin, pre-molten layer between two solids causing a partial reduction of the excess interfacial energy, 40,41 the crystals formed in the pore interiors are found to be surrounded by one-to-two monolayers thick liquid-like layers adjacent to the pore walls. Thermodynamics of the freezing and melting processes of such systems can be explored by analyzing their total free energy as a function of the thickness of the liquid layer. 35,42–44 Under certain assumptions made about the interactions and the pore geometry and by considering only the leading terms in the free energy functional, the transition temperatures can be shown to conform to a modified Gibbs-Thomson law

T0 − T = K/(d − τ ),

(1)

where d is the pore size, τ is the thickness of the non-frozen surface layer, and the proportionality factor K depends on the thermodynamical properties of the confined substances, the geometry of the pore space, and the phase transition mechanisms. Corrections to this equation can further be obtained by taking into account additional mechanisms contributing to the phase equilibria. 3 ACS Paragon Plus Environment

Page 5 of 37

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


As mentioned above, for predicting the freezing or melting transition temperature by this approach, i.e. for quantifying K, one has to assume a particular phase transition pathway. Thus, for freezing to occur in ideal cylindrical pores (the most simple pore geometry frequently used in the literature), a crystal nucleus has to be formed first. The highest temperature Tn , at which a spherically-shaped crystal seed can be nucleated, can be estimated from the condition that, at this temperature, the channel is able to accommodate a nucleus of the critical size. 1 The nucleation process is effective at sufficiently low temperatures only, i.e., it requires a significant degree of undercooling of the confined fluid. Hence, as soon as a stable crystal seed of the critical size is nucleated, it will immediately initiate freezing in the whole channel. The nucleation-related supercooling can be, however, avoided by providing seeds of the frozen phase at the pore openings. The latter corresponds to experimental situations in which porous materials containing capillary-condensed liquids are brought into contact with the bulk frozen liquid. Under this condition, freezing occurs by invasion of the frozen phase from the outer boundaries of the porous monoliths. The onset of this process is, obviously, found to occur at temperatures higher than Tn . Because this temperature depends on the confinement size, geometric disorder of the pore space can lead to quite complex percolation patterns, thus giving rise, once again, to the formation of metastable liquid phases in the pores. Analogous phenomena may be noted for the melting transition. Similarly to the previously considered case of the liquid-solid transition, melting is triggered by nucleation of seeds of the liquid phase. In contrast, however, the existence of liquid-like layers at the pore walls may facilitate the nucleation process. Because nucleation requires the overcoming of barriers in the free energy (crystals in the pore interior are, at these temperatures, in a metastable state), the formation of any liquid phase domain over the whole channel cross-section will cause melting of the frozen phase in the whole pore. Structural disorder may as well alter the melting process. It may, in particular, lead to situations where different phenomena may act in a concerted way and may, thus, result in complex patterns of phase coexistence.



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

Already this compendious overview illustrated how the disorder-induced interplay between different transition mechanisms, possibly further complicated by metastability effects, may give rise to different transition pathways. The direct application of theoretical predictions based on ideal pore systems to disordered materials may, therefore, lead to faulty conclusions. This may have important implications for methods of structural characterization, based on the application of Eq. 1 for the determination of pore sizes, 45 or for the prediction of the solid-liquid equilibria in industrially or environmentally-relevant processes. In order to gain deeper insights into freezing and melting in disordered porous materials, the application of microscopic approaches, free of the a priori assumptions inherent to macroscopic approaches, becomes therefore crucial. In the most simple way, this can be done using lattice models, in which the molecules are represented by the lattice sites and the interaction energies between the nearest neighbors are considered to depend on whether a particular site belongs to crystal or liquid. Despite their simplicity, these approaches have proven to be extremely useful. 46–48 In particular, they contributed significantly to an understanding of crystal growth processes in bulk systems. Among them, the Kossel-Stranski model with the simplest cubic symmetry was frequently used. 49 In this work, we employ it for an in-depth study of the generic properties of the freezing and melting transitions of fluids confined to pores with different pore geometries. 37 In addition, we also address the effect of the interactions of the confined substances with the pore walls, which may notably complicate the freezing and melting phenomena. 36,38,50 In what follows, we first elaborate on the lattice model used and its implementation for confined fluids. Thereafter, we apply it to examine different phenomena of phase transitions in pores with cylindrical pore geometry. The results obtained in this section will then be used to discuss fluid behavior in pore spaces with structural disorder.


Page 6 of 37

Page 7 of 37

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


Lattice model Let us consider an arbitrary cubic lattice configuration Ψj , in which different lattice sites can be in either solid (pore walls), fluid (liquid) or crystal (frozen liquid) states. The solid sites are used to model the pore geometry and remain in this fixed state. The probability P for the occurrence of a certain configuration is { [ P (Ψj ) = Q−1 exp −

∑ i

∑

1 β(−Nii ϕii − Nis ϕis + Ni Fi ) + (− βNik ϕik ) 2 i̸=k

]} ,

(2)

where β = 1/kT is the inverse temperature, the subscripts i and k can assume two values, namely c (crystal) or f (fluid), the subscript s stands for solid, Nik is the number of the ik bonds, Ni is the number of sites of type i, Q is the partition function, and -ϕik is the potential energy of the nearest neighbor interaction between i-th and k-th sites. Fc and Ff in eq 2 are the internal energies of the crystal and fluid sites, respectively. Taking account of only the nearest neighbor interactions, the number of sites Ni can be expressed as 1 1 3Ni = Nii + Nik + Nis , 2 2

(3)

where k ̸= i. With Nii given by eq 3, eq 2 becomes { P (Ψj ) = Q−1 exp −β

∑

} Ni (Fi − 3ϕii ) − ωcf Ncf − ωcs Ncs − ωf s Nf s

,

(4)

i

where the following definitions have been used: ωcf = β[−ϕcf +(ϕcc +ϕf f )/2], ωcs = β[−ϕcs + ϕcc /2], and ωf s = β[−ϕf s + ϕf f /2]. The chemical potential µi , associated with a site i, is µi = Fi − 3ϕii . The total number of sites is fixed, Nc + Nf = C. Thus, P (Ψj ) = Q−1 exp {βNc ∆µ − βCµf − ωcf Ncf − ωcs Ncs − ωf s Nf s } ,

(5)

where ∆µ = µf − µc . The probability P (Ψl ) to find a configuration Ψl which contains one



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 37

extra crystal site as compared to the configuration Ψj (with the bond numbers adjusted correspondingly and which are marked by an apostrophe) is { } ′ ′ P (Ψl ) = Q−1 exp β(Nc + 1)∆µ − βCµf − ωcf Ncf − ωcs Ncs − ωf s Nf′ s .

(6)

According to the principle of microscopic reversibility, the transition probabilities p(Ψj → Ψl ) and p(Ψl → Ψj ) fulfill the condition P (Ψj )p(Ψj → Ψl ) = P (Ψl )p(Ψl → Ψj ). This leads to p(Ψj → Ψl ) = exp {β∆µ − ωcf ∆Ncf − ωcs ∆Ncs − ωf s ∆Nf s } , p(Ψl → Ψj )

(7)

′ ′ where ∆Ncf = Ncf − Ncf , ∆Ncs = Ncs − Ncs , and ∆Nf s = Nf′ s − Nf s are the changes in the

bond numbers. In what follows, we will confine ourselves to considering only fluids wetting the pore walls. In addition, we will assume that the interaction energy of the fluid sites with the walls is larger than that of the crystal sites. In this way, we model the occurrence of the interfacial pre-melting as discussed in the introduction. Deviations from these conditions will be touched upon in a separate subsection. Furthermore, for the sake of simplicity we will assume that ϕf f = ϕcf . Under these conditions, ωcf = β[(ϕcc − ϕf f )/2] = βL/6, where L = 3(ϕcc − ϕf f ) is the heat of fusion. Finally, taking into account that ∆Nf s = −∆Ncs and defining ωc ≡ β[−ϕcs + ϕf s + (ϕcc − ϕf f )/2], eq 7 for the ratio between the ’creation’ (p+ ) and ’annihilation’ (p− ) probabilities of a crystal site results as { } 1 p+ = exp β∆µ − β∆Ncf L − ωc ∆Ncs . p− 6

(8)

The difference in the chemical potential ∆µ upon converting a site from the fluid to the crystal state at a given temperature T is ∆µ = L(T0 − T )/T0 . By noting that ∆Ncf = 2(3 − Nii ), where Nii is the resulting number of neighbors of the same type i after the site


Page 9 of 37

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


conversion, the exponent on the right hand side of eq 8 can be defined as ( −β∆G = βL

T Nii − 3 T0

) − ωc ∆Ncs

(9)

Finally, for the calculations, the probabilities p+ and p− can be taken as { } 1 1 p = exp ∓ β∆G 2 2 ±

(10)

to comply with eq 8. In the present work, we analysed the freezing and melting processes of the Kossel crystal confined to straight channel-like pores. Two different types of channels were considered. In the first case, referred to as ideal pores, the channel diameter was kept constant along the channel axis (z-axis). To obtain such channels, all sites satisfying a2 + b2 < r2 , where r is the channel radius and a and b are the site indices in the xy-plane (the pore center is taken at a = b = 0), were assigned to the pore sites, while all remaining ones formed the solid sites. Note that, for relatively small pore sizes, the actual shape of the pore boundary appears to be corrugated. In the second case, referred to as disordered pores, the channels were composed of short sections of a length l. Within each section the channel diameter was held constant (as for ideal pores). The channel diameter, however, varied from section to section. The radii of two adjacent pore sections were taken to be completely uncorrelated and were chosen to follow a Gaussian distribution with a width w centered around r0 . To adopt the model to real situations, two cut-offs, rmin and rmax , were applied to the Gaussian distribution. At the channel openings, direct contact to the sites resembling the bulk phase, sometime referred to as crystallizer, was provided (see Figure 1b). The system evolution was studied using dynamic Monte Carlo simulations. 51 In our simulations, one Monte Carlo step consisted of two subsequent runs. In the first run, a fluid site, randomly chosen from the pool of the pore and interfacial (see the next paragraph) sites, was attempted to be converted to the crystal site. The conversion was accepted with 8 ACS Paragon Plus Environment


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

a

bulk

solid

crystal/fluid

x

b

y

z

z

x y

Figure 1: Schematics showing (a) the xy-plane cross-section of a disordered pore and (ii) a three-dimensional perspective of an ideal pore. (a) corresponds to a temperature just slightly below T0 , at which all sites, except the bulk ones, are in the fluid state. (b) shows a typical configuration attained at sufficiently low temperatures with the crystal sites forming a frozen core in the pore interior. the probability p+ as given by eq 10. This procedure was repeated Nf times. In the second run, the same procedure was applied Nc times to the crystal sites attempting to convert them into the fluid ones with the probability p− . The Monte Carlo steps were repeated t times, with t denoting the Monte Carlo time of our simulations. The results presented were typically averaged over 200 pores (otherwise, it is explicitly indicated in the text). Throughout the simulations, the bulk sites were kept in either the crystal state (for temperatures below T0 ) or the fluid state (for temperatures above T0 ). In this way, any uncontrollable homogenous nucleation delays for freezing in the bulk phase were intentionally removed. Notably, this closely resembles the conditions of typical thermoporometry experiments, in which the cooling-heating cycles are performed in the presence of the excess bulk fluid. In addition, before the experiments, the bulk phase is first intentionally crystallized at low temperatures and then the temperature is raised to a temperature just slightly below than T0 . 7,8 Finally, all sites in the interfacial layer between the crystallizer and the porous


Page 10 of 37

Page 11 of 37

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


material were allowed to change their states (see Figure 1b). In what follows, the results will be presented in the form of the relative (equilibrium) fractions f of the pore sites obtained in the fluid state as a function of temperature. In this way, the results reported will correspond to the raw-data measured by NMR cryoporometry. 7,8 To provide a closer parallelism to the experimental situations, each lattice site may be associated with a molecule, which is either part of the crystal or liquid phases. Therefore, we fixed the site length lm to 0.5 nm, which is a typical molecular diameter of most low-molecular-weight organic liquids. Accordingly, the pore diameters will be indicated in nanometres. The bulk transition temperature is arbitrarily chosen to be T0 = 278 K. Two representative values for the heat of fusion, L = 10 kJ/mol and L = 20 kJ/mol, will be considered. These two values describe entropy changes typical for organic liquids. With these two parameters, L and T0 , and for the cubic crystal structure, the Jackson α-factor (α = Lη1 /ZRT0 , where R is the universal gas constant, η1 is the number of bonds in one layer, and Z is the crystal coordination number) 49 results as α = 2.9 and α = 5.8 for L = 10 kJ/mol and L = 20 kJ/mol, respectively. Notably, α > 2 holds for situations in which a layer-by-layer growth of the crystal is preferred over surface roughening. Thus, for α = 2.9 stronger fluctuations of the crystal-fluid are expected than for α = 5.8. To model the condition of interfacial pre-melting, ωc in eq 9 was set to ωc = 3βL.

Melting and freezing in ideal pores Melting-freezing hysteresis Figure 2 shows a typical result obtained for the most simple case of an ideal channel open at both ends. The freezing and melting transitions are found to be relatively sharp and irreversible. The shape of the hysteresis loop formed between freezing and melting closely resembles those found in experiments with nanoporous materials having cylindrical pores, like SBA-15. 13,25 At low temperatures, below the freezing temperature, f has non-zero values, 10 ACS Paragon Plus Environment


associated with the formation of a liquid film between the channel walls and the frozen core in the channel interior. The thickness of this film is found to be temperature-dependent. Notably, this finding is in a good agreement with the literature data providing convincing evidences for the existence of such layers possessing in high-mobility, liquid-like states and with thicknesses varying with temperature. 11,41,52–54

1,0

Liquid fraction, f

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 12 of 37

0,8

0,6

0,4

254

256

258

260

262

264

266

Temperature, K

Figure 2: Fractions of the lattice sites in the fluid state within the channel with a diameter of d = 6 nm and a length of l = 10 nm as a function of temperature on the cooling (triangles) and heating (circles) branches. The full and open symbols refer to the channels open at both ends and at one end only, respectively. At temperatures below T0 = 278 K, direct contact to the frozen bulk phase at the channel openings was provided. As revealed by the data of Figure 2, both the freezing and melting transition temperatures can be identified with high accuracy. Their suppression with respect to the bulk, ∆T = T0 − T , is shown in Figure 3 as a function of the channel diameter. The suppression is found to be in a qualitative agreement with the Gibbs-Thomson law, predicting proportionality between ∆T and the reciprocal channel diameter d−1 . Before any qualitative comparison with eq 1 will be made, we shall first explore in more detail the mechanisms controlling freezing and melting. Indeed, as it has been discussed in the introduction, the numerical constants in eq 1 are determined not only by the geometry of the pore space, but also by the 11 ACS Paragon Plus Environment

Page 13 of 37

assumptions made about the transition mechanisms. It is, therefore, necessary to establish first the thermodynamic conditions for either of the transitions.

40 T, K

heating

Temperature supression

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


cooling

30

20

10

0

4

6

8

10

12

Channel diameter d, nm

Figure 3: Suppressions of the transition temperatures of freezing (circles) and melting (triangles) in cylindrical pores open at both ends as a function of the channel diameter for two values of the heat of fusion L = 10 kJ/mol (full symbols) and L = 20 kJ/mol (open symbols).

Equilibrium transition Let us first identify which of the transitions, if there is any, is an equilibrium one and what is the equilibrium transition temperature in ideal cylindrical channels. In what follows, this temperature will be referred to as the pore equilibrium temperature, T0,p . To obtain T0,p , we have intentionally prepared an arbitrary phase configuration (typically resembling a large, continuous domain composed of the crystal sites surrounded by the bath of the sites in the liquid state) in the channel and followed its evolution at a given temperature T . At equilibrium, i.e., if T = T0,p , the thus-prepared phase configuration should first relax to a configuration providing the lowest energy state (e.g., by forming a liquid layer between the frozen core and the pore walls or/and rounding the crystal shape) and then, at longer times,



fluctuate around this latter configuration. If, however, T is higher or lower than T0,p , then there will be a long-time kinetics leading to complete pore melting or freezing, respectively. a

b 1,0

T = 259 K

T = 260,65 K

T = 258.5 K T = 258 K

0,8

0,7 T = 257.5 K 0,6

0,5

T = 261,15 K

0,9

0,9

Liquid fraction, f

1,0

Liquid fraction, f

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 37

0,8

T = 260,0 K

0,7

0,6

T = 259,65 K

0,5 T = 259,15 K

T = 257 K 0,4

0,4 0,3 0,0

3

2,0x10

3

4,0x10

3

6,0x10

3

8,0x10

4

1,0x10

4

0,0

1,2x10

4

2,0x10

Time, MC steps

4

4,0x10

4

6,0x10

4

8,0x10

5

1,0x10

Time, MC steps

Figure 4: Freezing and melting kinetics in cylindrical pores with d = 6 nm obtained by preparing an initial configuration with 30% of the sites in the crystal (one domain) and 70% in the liquid states and by tracing the system evolutions at temperatures indicated in the figure (L = 10 kJ/mol in a and L = 20 kJ/mol in b). At the equilibrium transition temperature, the magnitudes of f fluctuate around their mean values as shown by the solid lines. Figure 4 shows the equilibration kinetics obtained with a system initially containing about 30% of the sites in the crystal state (forming a continuous crystal domain). For providing the initial configuration for which the lowest-energy configuration can be attained most efficiently, it was prepared by letting the crystal phase grow from the crystallizer at low temperatures. The process was stopped when the crystal phase in the channel occupied about 30% of the total channel volume. If the thus prepared system is allowed to evolve at a certain temperature, the thickness of the fluid layers between the wall and the frozen core attains quickly its equilibrium value. Further evolution on a longer time scale, which we are interested in, is determined by the attainment of the global minimum in the free energy. As revealed by the data of Figure 4, there is a well-defined temperature T0,p at which the phase composition in the pore, as quantified by f , fluctuates around some mean value, while for T > T0,p and T < T0,p occurrence of, respectively, complete melting and complete freezing is observed. Most notably, the thus-obtained temperatures T0,p coincide, within the limit of the computational accuracy, with the freezing transition temperatures shown in Figure 3. 13 ACS Paragon Plus Environment

Page 15 of 37

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


This finding allows to conclude that, if the nucleation barrier is removed, freezing in ideal cylindrical pores occurs as the equilibrium transition. The mechanism of the transition is the axial advancement of the solid-fluid interface and is reminiscent of the desorption transition from cylindrical channels. In the same line it may be argued that, whenever there will be a mechanism for the formation of a solid-fluid interface across the channel, melting as well will occur at T0,p . The latter situation can be encountered, e.g., in channels closed at one end. Because a liquid-like layer is formed at the closed end, it will facilitate the melting process for temperatures just below T0,p . As exemplified by the data of Figure 2, for such systems indeed no hysteresis is observed (note that a small hysteresis in Figure 2 has purely kinetic origin associated with slow kinetics of freezing, which will be discussed later, and can be eliminated by providing sufficiently long simulation times) and the reversible transition occurs at T0,p . Once again, this is in the complete analogy to the reversible gas-liquid and liquid-gas transitions in capped cylindrical channels as predicted by Cohan. 55–57 This finding, namely identifying freezing and melting in ideal pores as the equilibrium and metastable transitions, respectively, is in line with a part of the theoretical studies made earlier. 42–44 There were many subsequent works devoted to validations of these predictions using porous materials with well-defined, cylindrical pore geometries (like MCM-41 or SBA15). Some results obtained experimentally have shown deviations from these theoretical predictions (see, e.g, related discussions in Ref. 31 In particular, it was found that the melting transition for water in MCM-41 and SBA-15 materials with different pore diameters was reasonably captured by eq 1. At the same time, freezing exhibited much more scattering, questioning the applicability eq 1 to the freezing transition. 13,24 Most importantly, fit of eq 1 to the melting data (taking account of one monolayer of the non-frozen water and using independently obtained thermodynamical parameters for water) revealed equilibrium character of this transition. We would like to note, however, that, in this particular case, one has to be aware of the fact that water is a substance having a relatively low Jackson factor



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 37

α ≈ 2.0 (Ref., 49 page 302). This value is very close to the critical one, when the layer-by-layer crystal growth is not preferred over dendritic one, i.e., the solid-fluid interface is subject to strong fluctuations. This scenario is further supported by the findings that, upon freezing in small pores, ice exhibits a very defective structure deviating from bulk hexagonal ice. 38 As it will be shown in the following sections, such interface fluctuations tend to shift the melting temperatures towards the equilibrium ones. Thus, in sufficiently small pores melting may occur close to the equilibrium transition temperature and the coincidence of the melting and equilibrium transition temperatures in Refs. 13,24 may be due to thermodynamic fluctuations. Although the effect of these fluctuations is barely understood and is rarely discussed in the literature, such lowering of the melting temperature, which may ultimately reach T0,p , has already been pointed out by Petrov and Furó. 8

Freezing temperature By resolving the freezing mechanism in ideal pores, the freezing temperature Tf of the Kossel crystal in cylindrical pores can readily be established. Using eq 9, the condition of equilibrium, ∆G = 0, is given by (

Nii T − 3 T0

)

( −

ωc βL

) ∆Ncs = 0.

(11)

Here, the bars denote the ensemble average. For the bulk Kossel crystal, ∆Ncs = 0, hence Nii = 3. Thus, changing the state of a site having exactly three neighbours in the crystal state does not change the overall energy of the whole system. In the literature, this is known as the concept of a repeatable or kink site. In turn, if Nii is known (for bulk systems Nii = Z/2, where Z is the total number of bonds), the equilibrium transition temperature can be found using eq 11. Nii for the confined Kossel crystals can be evaluated as follows. If the crystallizer is supplied at the pore openings, the crystal will grow in axial direction of the channel at


Page 17 of 37

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


temperatures just slightly below T0,p . In radial direction, the crystal is confined by the channel walls. Because the liquid is considered to wet the pore walls, a non-frozen layer will be found between the wall and the crystal core in the channels. Hence, one may distinguish between two different types of repeatable sites facilitating the axial growth of the crystal. For any site which has no direct contact with the sites belonging to the non-frozen layer, Nii,m = 3 like for the bulk system. The sites in contact with the non-frozen layer may, however, form on average only Z = 5 bonds with sites in the crystal state. Hence, for these sites Nii,p = 5/2. The mean Nii can, thus, be found as

Nii = pm Nii,m + pp Nii,p ,

(12)

where pm and pp are the relative fractions of the two types of the repeatable sites. For cylindrical symmetry, pm = 1 − 4/Nd and pp = 4/Nd , where Nd is the number of sites along the diametric line of the crystal. Finally, Nii results as Nii = 3 − 2/Nd . By noting that Nd = dc /lm , where dc is the crystal diameter and lm is the site dimension, eq 11 with ∆Ncs = 0 (because of a very low probability to have the crystal-solid bonds) yields ∆Tf = T0 − T0,p =

2 lm T0 . 3 dc

(13)

Figure 5 shows the data of Figure 3 plotted versus the inverse pore and crystal diameters. In the latter case, the mean crystal diameters for temperatures in the vicinity of Tf were obtained by direct analysis of the simulation data. The data of Figure 5b are found to be in perfect agreement with the predictions of eq 13. This confirms that the classical GibbsThomson law corrected for the thickness of the non-frozen layers, as given by eq 1, can accurately describe the freezing transition in ideal pores. Notably, this type of correction is a common practice for the determination of actual pore sizes. 13,24,27,54



a

b 0,14

0,11

0,11

f

T /T

0

0,14

0 f

T /T

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 18 of 37

0,07

0,04

0,07

0,04

0,00

0,00 0,1

0,2

0,3

0,1

-1

0,2

0,3

0,4 -1

Inverse pore diameter 1/d, nm

Inverse crystal diameter 1/d, nm c

Figure 5: Freezing temperature suppression ∆T normalized by the equilibrium transition temperature T0 as a function of the inverse (a) pore and (b) crystal diameters for L = 10 kJ/mol (full symbols) and L = 20 kJ/mol (opens symbols). The solid lines show ∆Tf /T0 resulting from eq 13 with lm = 0.5 nm.

Melting temperature As shown in Figure 2, for ideal channels with two open ends, at which the crystallizer is supplied, there is a significant difference between the freezing and melting temperatures. Melting in this case occurs at substantially higher temperatures due to the lack of the (radial) crystal-liquid interface allowing axial shrinking of the crystal at temperatures above T0,p . Hence, melting starts only at temperatures at which the formation of the liquid domains forced by fluctuations becomes possible. Fluctuations play an essential role in this process and become increasingly important with decreasing pore diameters. This is, in particular, revealed by the data for the melting temperature suppressions obtained in simulations which are shown in Figure 6. It is not straightforward, however, to describe the effect of fluctuations quantitatively. For the purpose of this work, we estimate therefore only the highest melting temperature at which there will be no barriers in the free energy for the process of radial shrinking. Formally, this process may be considered to occur at equilibrium for the outermost layers of the crystal. In contrast to the flat interfaces analyzed in the preceding section, here we deal, however, with convex cylindrical interfaces. With the last remark, the melting temperature, which in what follows will be denoted as T0,m , can be established using eq 11. Nii can most easily be calculated by considering wide 17 ACS Paragon Plus Environment

Page 19 of 37

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


channels allowing to neglect the finite size of the crystal lattice. Let us consider the outermost layer of the frozen core in a channel. In one cross-section, it contains approximately N = πdc /lm sites. The layer adjacent to it contains approximately N ′ = π(dc − 2lm )/lm sites. Hence, N − N ′ sites in the outermost layer will form only 5/2 bonds, while the rest will form 3 bonds. Using eq 12 with all relevant quantities, Nii results as (5 + N ′ /N )/2. Finally, the melting temperature suppression is readily found to be two times smaller than ∆Tf :

∆Tm = T0 − T0,m =

1 lm T0 . 3 dc

(14)

By comparing the prediction of eq 14 with the simulation data in Figure 6b, one finds reasonable agreement for large pores, but notable deviations for smaller ones. In particular, the most dramatic deviation is found for narrow channels and the crystal with the lower heat of fusion, L = 10 kJ/mol. Presumably, this is caused by the increasing role of fluctuations with decreasing crystal diameter and decreasing L. Indeed, fluctuations of the crystal-fluid interface with amplitudes comparable to the pore diameter lead to the formation of liquid bridges triggering melting by axial shrinking of the crystal (recall that the equilibrium transition temperature T0,p is notably lower than T0,m ). Therefore, for one and the same temperature, fluctuations will be more effective for smaller pores. At the same time, according to eq 9, the energy penalty for fluctuations leading to increase of the crystal-liquid interfacial area is directly proportional to L. Thus, for crystals with lower bond strength stronger deviations from eq 14 can be expected. At this point, it is important to make a remark about water, which is often used as a test liquid for freezing and melting experiments in porous materials. This choice is partially determined by the fact that all thermodynamical parameters for bulk water are known independently and with high accuracies. However, as it has already been mentioned before, water has a sufficiently low Jackson factor α ≈ 2.0. For the Kossel crystal with T0 = 278.0 K, as considered in the preceding sections, this corresponds to the heat of fusion L = 7 kJ/mol. This reveals that the effect of thermodynamic fluctuations upon melting will be 18 ACS Paragon Plus Environment


even stronger than for the case of L = 10 kJ/mol as exemplified in Figure 6b. b

a

0,07

0,07

m

m

T /T

0

0,11

0

0,11

T /T

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 20 of 37

0,04

0,04

0,00

0,00 0,1

0,2

0,3

0,1

-1

0,2

0,3

0,4 -1

Inverse pore diameter 1/d, nm


Figure 6: Melting temperature suppression ∆T normalized by the equilibrium transition temperature T0 as a function of the inverse (a) pore and (b) crystal diameters for L = 10 kJ/mol (full symbols) and L = 20 kJ/mol (opens symbols). The dashed and solid lines show ∆Tm /T0 resulting from eqs 13 and 14 with lm = 0.5 nm, respectively.

Transition kinetics It is interesting to note that Figure 4 indicates that, for small, identical deviations ±(T −T0,p ), the kinetics of (equilibrium) freezing transitions are slower than the kinetics of (equilibrium) melting transitions. Remarkably, this effect becomes stronger with decreasing L, although the location of the equilibrium transition temperature T0,p is not affected by L. This observation may be explained by referring to the fact that, for the values of L and T0 used in our simulations resulting in Jackson α-factors larger than 2, the addition/removal of one layer to/from the crystal surface requires overcoming of nucleation barriers. The latter are higher for systems with larger bond energies. In full agreement, both the freezing and melting rates are found to be notably slower in Figure 4b for L = 20 kJ/mol than in Figure 4a for L = 10 kJ/mol. To explain the difference between the freezing and melting kinetics for identical L, one has to take into account the finite length of the systems considered and the occurrence of a liquid layer in contact with the crystal. Thus, the formation of a new layer (in the cross-sectional plane of the channel) will require nucleation of a critical cluster on the crystal facet, while the removal can be facilitated by the already existing liquid sites. In the case 19 ACS Paragon Plus Environment

Page 21 of 37

of weak bonds (low L), a rougher surface of the crystal will help to overcome the nucleation barriers, making the difference between the freezing and melting rates more modest. This line of reasoning is also supported by the simulation data.

Effect of the pore wall interaction

20

10

0

Temperature shift T -T, K

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


0

-10

-20

-0,5

0,0

0,5

(

-

cs

1,0

)/(

fs

-

cc

1,5

2,0

)

fc

Figure 7: Shift of the freezing (filled squares) and melting (filled circles) temperatures T0 − T as a function of the ratio b = (ϕcs − ϕf s )/(ϕcc − ϕf c ) in the channels with d = 6 nm. The simulations were performed with L = 20 kJ/mol. For T < T0 (upper half of the figure) the bulk phase at the pore openings was kept in the frozen state, for T > T0 (lower half of the figure) in the liquid state. So far, only liquids wetting the pore walls and having much better affinity to the pore walls as compared to their crystalline phases have been considered. This was intentionally provided by setting ωc in eq 9 to 3βL, i.e. by setting a sufficiently high energy penalty for having a crystal site in contact with a solid one. Because ωc uniquely combines all interactions parameters, this allows to compile a global phase diagram as a function of ωc . This diagram may further be used to analyse the effect of the variation of different interaction terms. It is, however, more convenient to introduce a new parameter b = (ϕcs − ϕf s )/(ϕcc − ϕf c ), 20 ACS Paragon Plus Environment


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

so that ωc = ( 12 − b) βL . In particular, all previously considered situations (ωc = 3βL) 3 corresponded to b = −17/2. To build the phase diagram, we have performed dynamic Monte Carlo simulations for different values of ωc in eq 9 and studied how the freezing and melting temperatures are altered by its variation. A typical result is shown in Figure 7. It demonstrates the shifts ∆Tf and ∆Tm of the transition temperatures as a function of b as obtained for the channels with d = 6 nm. For b < 0 (ϕcs < ϕf s ), i.e. for liquids having better affinities to the pore walls as compared to their solid phases (the occurrence of the interfacial pre-melting), ∆Tf and ∆Tm are found to be independent of b. This behavior emerges due to the short-range character of the interactions in our model, where only the first-neighbour interactions between adjacent sites are considered. For these interaction strengths and for temperatures not significantly lower than the transition temperature, one layer of the sites located between the pore walls and the frozen core always remains in the liquid state. The transition temperatures are, thus, determined by eqs 13 and 14 with the confinement size corrected for the thickness of the liquid layer. For 0 < b < 1, i.e. for interaction strengths ϕf s < ϕcs < ϕcc resembling better affinity of the crystal phase to the pore walls than the fluid one, ∆Tf is found to decrease with increasing b. In this regime, the last term on the right side hand of eq 9 (or eq 11) cannot be neglected (∆Ncs ̸= 0). Hence, ∆Tf is altered due to lowering of the excess surface energy of the confined crystal. In particular, for b = 1 the excess surface energy becomes zero and, therefore, the transition occurs at T0 (∆Tf = 0). The melting transition is found to be dramatically affected by the occurrence of the favorable crystal-solid bonds, i.e. by the diminishment of the liquid layers. Indeed, as it has been discussed earlier, these layers facilitate nucleation of the liquid bridges. Only for relatively small values of b, the composition fluctuations in the layers adjacent to the pore walls may trigger the melting process at temperatures T > T0 . For larger b, nucleation of the liquid bridges is strongly suppressed and melting is, therefore, postponed to the bulk transition temperature T0 . In this case, the liquid phase supplied at


Page 22 of 37

Page 23 of 37

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


the pore openings invades into the channels. The situations with b > 1 (ϕcs > ϕcc ), i.e. when the pore walls may facilitate heterogeneous nucleation of the crystal phase, can be discussed in the same way by noting the symmetry between freezing and melting. They are found to coincide upon two subsequent reflections over the lines b = 1 and ∆T = 0. In this case, the above discussion is identically applicable upon interchanging the crystal and fluid sites. It is important to mention that, in our model, the crystalline structures of the pore walls and of the crystal growing in the pore space are identical because they both are determined by the lattice. The mismatch between the crystalline structures in real system may lead to more complex behavior, which is not directly captured by our model. 58

Melting and freezing transitions in disordered pores One-dimensional channels with disorder Disorder in the pore structure may strongly affect the freezing and melting processes. This is caused by the temperature-dependent character of these transitions in confined spaces and by the interplay of two different transition mechanisms contributing to the melting transition. To address the effect of geometric disorder, we made use of a most simple model of disordered porous materials, considering it to be composed of channel-like pores with the variation of the channel diameters along their axes. It turns out that already this one-dimensional model shows all essential features for rationalizing the complex transition phenomena in materials with three-dimensional pore spaces. In more detail, we composed long channels by joining short sections, represented by ideal channels of length l. Diameters of neighbouring channel sections were uncorrelated, following a discretised Gaussian distribution. In what follows, we present the results for the channels with an average pore diameter d = 8 nm and a distribution width σ = 2 nm. Cut-offs at d = 4 nm and d = 14 nm were used. The length l of the individual 22 ACS Paragon Plus Environment


sections was 2.5 nm (we have proven that choosing longer section lengths resulted in similar results, but required sufficiently longer simulation times). Notably, any particular detail of the disorder (as modeled here by assuming missing correlation between the diameters of neighbouring pores and by using relatively short section lengths), have only quantitative, but not qualitative consequences. a

b 1,0

1,0

0,8

0,8

Liquid fraction, f

Liquid fraction, f

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 24 of 37

0,6

0,4

0,6

0,4

0,2

0,2 240

250

260

270

240

Temperature (K)

250

260

270

Temperature (K)

Figure 8: (a) Freezing and melting transitions in disordered channels with total lengths of 100 nm (open circles) and 750 nm (filled triangles). The simulations were performed with L = 10 kJ/mol. For comparison, the filled squares show the results obtained with L = 20 kJ/mol and for a channel length of 750 nm. In both cases, the formation of liquid layers was ensured (b = −17/2). The results shown for the short pores were averaged over 70 different channel configurations, for the long pores over 10 configurations. The dashed line approaches the change of the volume fraction occupied by the non-frozen surface layer. (b) Freezing and melting data obtained with L = 20 kJ/mol and for a channel length of 750 nm replotted from the figure (a). The solid and dashed vertical lines show, respectively, the melting and freezing transition temperatures obtained in an ideal channel with the pore diameter d = 8 nm and with L = 20 kJ/mol. Figure 8 shows the results obtained for the channels with disorder. First of all, we note the occurrence of a relatively wide hysteresis loops formed by the freezing and melting branches. The freezing transition is found to strongly depend on the channel length, while the melting one shows no length dependence. Finally, the transition enthalpy has a minor, but still notable effect on both freezing and melting. As revealed by the data of Figure 8b, the hysteresis widths obtained for disordered channels is found to notably exceed that obtained in ideal pores with similar pore diameters. This observation may be rationalized by recalling the equilibrium character of the freezing 23 ACS Paragon Plus Environment

Page 25 of 37

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


transition occurring via the invasion of the crystal from the pore openings. Because Tf is determined by the channel diameter, the growth of the invading crystal phase will be controlled by the channel sections with the smallest diameters (the pore-blocking effect). Notably, this scenario also naturally emerges in the macroscopic theories assuming that freezing occurs at equilibrium. 42–44 Due to the channel-like pore geometry, the narrowest channel sections located (statistically) close to the pore openings may effectively postpone freezing in the whole channel to substantially low temperatures. Indeed, the onset of freezing is found at temperatures corresponding to freezing in ideal channels with d = 4 nm. The longer the channels are, the higher fractions of the sites will remain in the liquid state until these temperatures. For the shorter pores with the total lengths comparable to an average separation between the narrowest pores, pore blocking will not be the dominant mechanism and the temperature dependence of the freezing transition will be determined by the pore size distribution. The difference between these two situations is nicely illustrated in Figure 8a. In contrast to freezing, melting in disordered channels may involve two distinct mechanisms. If a crystal in a pore section of diameter d will melt via nucleation of a liquid bridge at the temperature T ≈ T0,m (d) as given by eq 14, it may trigger axial melting in the neighbouring sections for which T0,p ≤ T0,m (d). This phenomenon, referred to as advanced melting, has already been discussed qualitatively in the literature and has been shown to impact the accuracy of the pore size determination. 59 Within our model, this effect may be considered on a quantitative basis. With eqs 13 and 14, the condition of the advanced melting means that melting in any channel section of diameter d = dc + τ (where τ is the thickness of the liquid layer) will be accompanied by successive melting in the neighbouring sections with diameters up to 2dc + τ . This condition is visualized in Figure 9. Thus, at relatively low temperatures, melting will predominantly occur as a concerted action of these two mechanisms. At higher temperatures, equilibrium melting of large pores will dominate the overall process due to the occurrence of a large number of the liquid domains. It is evident that, in this scenario, the channel length does not play any decisive role.



0,10

axial melting 0,08

T

0

radial melting 0,06

T

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 26 of 37

0,04

0,02

0,00 0,10

0,15

0,20

0,25

0,30

0,35

0,40

-1


Figure 9: Diagram showing the interplay between two different melting mechanisms in disordered pores. The lines are the equilibrium and metastable melting transition temperatures versus inverse crystal size dc . Melting of the frozen liquid in the channels with the diameters up to d = dc + τ by nucleating liquid bridges will give rise to the formation of continuous domains with the molten liquid due to axial melting of the adjacent sections with the channel diameters up to d = 2dc + τ . Finally, Figure 8 demonstrates the effect of the transition enthalpy L on both melting and freezing in disordered pores. Notably, eqs 13 and 14 predict that, in ideal pores, the transition temperatures should not depend on L. These equations, however, do not include the effect of the liquid-crystal interface fluctuations, which increase with decreasing L. 49 This has already earlier been demonstrated for the melting transition in ideal pores. The same phenomenon can also contribute to a more effective freezing in disordered channels. Indeed, fluctuations may facilitate the propagation of the solid-liquid interface through a narrow channel section. This may clearly occur if the channel lengths are sufficiently short. Thus, for fluids with higher L both transition temperatures in disordered channels may be expected to be shifted towards higher temperatures. This prediction is found to be in nice agreement with the simulation results shown in Figure 3.


Page 27 of 37

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


Structural information accessible Measurements of solid-liquid transitions are often used as a means for structural characterization of nanoporous solids. It is therefore instructive to discuss in more detail which type of structural information can be obtained for disordered porous materials based on the transition mechanisms established. Let us first consider the melting transition, which is most frequently used to assess the pore size distributions. With the relative fraction f of the liquid phase in the sample measured upon heating, the crystal size distribution (relative number of the channel sections Nc with crystal diameter x) can be found as dNc 1 d(f − fm ) 1 d(f − fm ) dT ∝ 2 = 2 , dx x dx x dT dx

(15)

where fm is the relative volume fraction of the liquid layer. The thus obtained distribution can be considered to approach the pore size distribution with d = x + τ , where τ is the thickness of the liquid layers. The derivative dT /dx in eq 15 is readily obtained from eq 13 or eq 14, resulting in dNc 1 d(f − fm ) ∝ 4 . dx x dT

(16)

At relatively low temperatures the melting process is controlled by radial melting of the narrowest channels. Eq 16 combined with eq 14 can therefore be used to assess the pore sizes corresponding to smallest ones in the real pore size distribution. Although eq 16 can correctly predict the range of the pore sizes, as it is demonstrated in Figure 10a, the information about their relative fractions will be corrupted by the simultaneous occurrence of axial melting in the neighbouring pores. In the same spirit, the range of the largest pore sizes can be obtained by the combined use of eq 16 and eq 13. This procedure applied to the freezing data of Figure 8 as well yields the correct range of the largest pore sizes (see Figure 10a). Notably, this latter procedure is only valid for materials with relatively broad PSDs, i.e. having heavy tails with d > 2dmin , 26 ACS Paragon Plus Environment


where dmin represents the typical diameter of the smallest pore size. Otherwise, the frozen liquid in materials with narrow pore size distributions can melt in one step, simultaneously with the onset of radial melting. a

b 1,0

1,0

0,8

0,8

Frequency

Frequency

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 28 of 37

0,6

0,4

0,2

0,6

0,4

0,2

0,0

0,0 4

6

8

10

12

14

4

6


8

10

12

14


Figure 10: Pore size distributions derived from the melting (a) and the freezing (b) transitions shown in Figure 8 for L = 10 kJ/mol using eq 16. The solid line is the real pore size distribution used in the simulations. The filled circles and triangles in (b) represent the data obtained for the long and short channels, respectively. The fractions fm in eq 16 were approached by the dashed line shown in Figure 8. The freezing transition also contains important structural information. Because freezing under the condition of crystal phase growth starting from the pore openings is controlled by the pore blocking mechanism, eq 16 and eq 13 can be applied to the freezing data to obtain the neck size distribution. This is especially valid for long channels, resembling monolithic porous materials. With decreasing the channel length, more larger pores become accessible by this procedure. These features are nicely illustrated in Figure 10b.

Conclusions In the present work, we have studied the freezing and melting phenomena in porous materials using the Kossel-Stranski crystal growth model. According to this model, molecules in bulk systems arrive at the crystal surface at a certain, specified rate and leave it at a rate depending on the number of adjacent crystal sites. We have modified the model by including additional sites representing the porous material as a solid matrix. The binding energies 27 ACS Paragon Plus Environment

Page 29 of 37

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


between two crystal sites and between the crystal and sold sites have been considered to be different. Thus, the crystal growth could occur only in the pore spaces and the detachment rate from the crystal surface was determined by not only the adjacent crystal sites, but also the solid ones. With only this modification, a number of phenomena inherent to solid-liquid equilibria and phase transitions in confined spaces have naturally emerged. The flexibility of the model in terms of the possibility to vary different interactions and the pore space geometry, permitted us to establish some generic features in the freezing and melting behavior. Their deeper understanding, obtained on this microscopic route, may serve for a more proficient analysis of the experimental results obtained for porous materials with complex pore structures. This may not always be accessible in the frame of macroscopic thermodynamics requiring certain a priori assumptions. By addressing foremost ideal, channel-like pores, we examined the status of the GibbsThomson law for the Kossel crystal by deriving its analogues for two different transitions. The first one was the freezing transition for which the crystal nucleation barriers were intentionally removed by providing the crystalline phase at the pore openings, as often done in thermoporometry experiments. In this case, we show unequivocally that the freezing transition occurs at equilibrium, i.e. without change in the overall energy. This finding is in agreement with the theoretical models put forward earlier. 42–44 The second transition concerned melting occurring by nucleation of liquid bridges within the channels, which is intrinsically metastable. These two transitions were considered because they are most frequently encountered in the experiments. By comparing the analytical results with the data obtained using Monte Carlo simulations we identified conditions for the applicability of the Gibbs-Thomson law. As one of the significant results of this work, we found that the Jackson α-factor, a parameter quantifying the entropy change during the transition, plays an important role for the melting and freezing transitions in confined spaces. Thus, for materials with large entropy changes, the Gibbs-Thomson law, corrected for the occurrence of non-frozen layers adjacent



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

to the pore walls, is found to be in a satisfactory agreement with the simulation data. For materials with low α-factors, however, the fluctuations are found to lead to notable deviations from this law. In particular, it is found that for the melting transition these deviations are notably stronger than for freezing. This finding suggests that the data analysis of the freezing and melting experiments in porous solids should be performed taking account of the Jackson factor. This might especially be important for water having a sufficiently small α-factor of about 2. 49 By varying the binding energies, we explored how the interaction between the crystal, the liquid, and the pore wall affects the transition temperature. For pore walls favoring liquid over crystal sites, the model exhibited the emergence of a liquid layer minimizing the overall energy of the system, resembling the phenomenon of interfacial pre-melting. 11,41,52–54 Its occurrence is shown to lead to the correction of the Gibbs-Thomson equation for the thickness of this layer. Notably, this correction is widely used for the pore size determination. 13,24,27,54 Under these conditions, we found that the transition temperatures do not depend on the interaction strength. The same is found to be valid in the opposite case, when the crystalsolid bonds lower the total energy, giving rise to an epitaxial growth of the crystal on the pore walls. In this case, both melting and freezing appeared to occur at temperatures higher than the bulk transition temperature, but irrespective of the substrate-crystal binding energies. For intermediate interaction energies, we find a transition between these two limits. Transition mechanisms established using ordered pores allowed for a more robust analysis of the phase transformations occurring in disordered materials. We showed that structural disorder renders freezing to occur under strong pore blocking control, i.e. crystal growth is restrained by the smallest pores in the porous material. This finding is as well found to be in agreement with theoretical predictions considering freezing as equilibrium transitions. Exactly these smallest pores were found to determine the onset of the melting process because of smaller nucleation barriers in these pores. Importantly, the thus formed domains of the liquid phase permitted advanced melting in the neighbouring pores, giving rise to a concerted


Page 30 of 37

Page 31 of 37

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


action of two different melting mechanisms. With the increasing number of liquid domains, melting was found to be controlled by the largest pores, postponing melting to temperatures of their equilibrium transitions. Thus, depending on the particular geometry of the pore structure and the pore sizes in disordered materials or ordered materials with defects, the melting transition may exhibit features of metastable or of equilibrium transition, or of their mixture. To assess the underlying mechanisms in these situations, more elaborate experiments probing scanning freezing and freezing behavior are required. 37,39 What type of information is delivered by this type of the experiments may be efficiently educed by using our model. Such studies are currently under progress.

Acknowledgement The authors thank DFG (the German Science Foundation), in particular in the frame of the research project FOR-877 ”From local constraints to macroscopic motion”, for the financial support.

References (1) Christenson, H. K. Confinement effects on freezing and melting. J. Phys.-Condes. Matter 2001, 13, R95–R133. (2) Alba-Simionesco, C.; Coasne, B.; Dosseh, G.; Dudziak, G.; Gubbins, K. E.; Radhakrishnan, R.; Sliwinska-Bartkowiak, M. Effects of confinement on freezing and melting. J. Phys.: Condens. Matter 2006, 18, R15–R68. (3) Brun, M.; Lallemand, A.; Quinson, J. F.; Eyraud, C. New Method for Simultaneous Determination of Size and Shape of Pores - Thermoporometry. Thermochim Acta 1977, 21, 59–88. (4) Landry, M. R. Thermoporometry by differential scanning calorimetry: experimental considerations and applications. Thermochim Acta 2005, 433, 27–50. 30 ACS Paragon Plus Environment


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

(5) Overloop, K.; Vangerven, L. Freezing Phenomena in Adsorbed Water as Studied by NMR. J. Magn. Reson. A 1993, 101, 179–187. (6) Strange, J. H.; Rahman, M.; Smith, E. G. Characterization of Porous Solids by Nmr. Phys. Rev. Lett. 1993, 71, 3589–3591. (7) Mitchell, J.; Webber, J. B. W.; Strange, J. H. Nuclear magnetic resonance cryoporometry. Phys. Rep. 2008, 461, 1–36. (8) Petrov, O. V.; Furó, I. NMR cryoporometry: Principles, applications and potential. Prog. Nucl. Magn. Reson. Spectrosc. 2009, 54, 97–122. (9) Rennie, G. K.; Clifford, J. Melting of Ice in Porous Solids. J. Chem. Soc.-Faraday Trans. I 1977, 73, 680–689. (10) Hills, B. P.; Lefloch, G. NMR-studies of non-freezing water in randomly packed-beds of porous particles. Mol. Phys. 1994, 82, 751–763. (11) Valiullin, R.; Furo, I. The morphology of coexisting liquid and frozen phases in porous materials as revealed by exchange of nuclear spin magnetization followed by H-1 nuclear magnetic resonance. J. Chem. Phys. 2002, 117, 2307–2316. (12) Handa, Y. P.; Zakrzewski, M.; Fairbridge, C. Effect of restricted geometries on the structure and thermodynamic properties of ice. J. Phys. Chem. 1992, 96, 8594–8599. (13) Morishige, K.; Kawano, K. Freezing and melting of water in a single cylindrical pore: The pore-size dependence of freezing and melting behavior. J. Chem. Phys. 1999, 110, 4867–4872. (14) Koga, K.; Tanaka, H.; Zeng, X. C. First-order transition in confined water between high-density liquid and low-density amorphous phases. Nature 2000, 408, 564–567. (15) Li, J. C.; Ross, D. K.; Benham, M. J. Small-angle neutron-scattering studies of water and ice in porous Vycor glass. J. Appl. Crystallogr. 1991, 24, 794–802. 31 ACS Paragon Plus Environment

Page 32 of 37

Page 33 of 37

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


(16) Takamuku, T.; Yamagami, M.; Wakita, H.; Masuda, Y.; Yamaguchi, T. Thermal property, structure, and dynamics of supercooled water in porous silica by calorimetry, neutron scattering, and NMR relaxation. J. Phys. Chem. B 1997, 101, 5730–5739. (17) Baker, J. M.; Dore, J. C.; Behrens, P. Nucleation of ice in confined geometry. J. Phys. Chem. B 1997, 101, 6226–6229. (18) Awschalom, D. D.; Warnock, J. Supercooled liquids and solids in porous glass. Phys. Rev. B 1987, 35, 6779–6785. (19) Ishikiriyama, K.; Todoki, M. Pore Size Distribution Measurements of Silica Gels by Means of Differential Scanning Calorimetry: II. Thermoporosimetry. J. Colloid Interface Sci. 1995, 171, 103–111. (20) Jackson, C. L.; McKenna, G. B. The melting behavior of organic materials confined in porous solids. J. Chem. Phys. 1990, 93, 9002–9011. (21) Molz, E.; Wong, A. P. Y.; Chan, M. H. W.; Beamish, J. R. Freezing and melting of fluids in porous glasses. Phys. Rev. B 1993, 48, 5741 LP – 5750. (22) Unruh, K. M.; Huber, T. E.; Huber, C. A. Melting and Freezing Behavior of Indium Metal in Porous Glasses. Phys. Rev. B 1993, 48, 9021–9027. (23) Morishige, K.; Kawano, K. Freezing and melting of methyl chloride in a single cylindrical pore: Anomalous pore-size dependence of phase-transition temperature. J. Phys. Chem. B 1999, 103, 7906–7910. (24) Schreiber, A.; Ketelsen, I.; Findenegg, G. H. Melting and freezing of water in ordered mesoporous silica materials. Phys. Chem. Chem. Phys. 2001, 3, 1185–1195. (25) Morishige, K.; Iwasaki, H. X-ray study of freezing and melting of water confined within SBA-15. Langmuir 2003, 19, 2808–2811.



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

(26) Morishige, K.; Uematsu, H.; Tateishi, N. Comparative study of liquid-solid phase transition of nitrogen in controlled pore glass and SBA-15: Accurate pore size of conventional mesoporous materials. J. Phys. Chem. B 2004, 108, 7241–7246. (27) Jahnert, S.; Chavez, F. V.; Schaumann, G. E.; Schreiber, A.; Schonhoff, M.; Findenegg, G. H. Melting and freezing of water in cylindrical silica nanopores. Phys. Chem. Chem. Phys. 2008, 10, 6039–6051. (28) Schaefer, C.; Hofmann, T.; Wallacher, D.; Huber, P.; Knorr, K. Melting and Freezing of Argon in a Granular Packing of Linear Mesopore Arrays. Phys. Rev. Lett. 2008, 100, 175701–4. (29) Findenegg, G. H.; Jhnert, S.; Akcakayiran, D.; Schreiber, A. Freezing and Melting of Water Confined in Silica Nanopores. ChemPhysChem 2008, 9, 2651–2659. (30) Deschamps, J.; Audonnet, F.; Brodie-Linder, N.; Schoeffel, M.; Alba-Simionesco, C. A thermodynamic limit of the melting/freezing processes of water under strongly hydrophobic nanoscopic confinement. Phys. Chem. Chem. Phys. 2010, 12, 1440–1443. (31) Petrov, O.; Furo, I. A study of freezing-melting hysteresis of water in different porous materials. Part II: surfactant-templated silicas. Phys. Chem. Chem. Phys. 2011, 13, 16358–16365. (32) Morishige, K.; Yasunaga, H.; Denoyel, R.; Wernert, V. Pore-Blocking-Controlled Freezing of Water in Cagelike Pores of KIT-5. J. Phys. Chem. C 2007, 111, 9488–9495. (33) Khokhlov, A.; Valiullin, R.; Kärger, J.; Steinbach, F.; Feldhoff, A. Freezing and melting transitions of liquids in mesopores with ink-bottle geometry. New J. Phys. 2007, 9, 272. (34) Sliwinska-Bartkowiak, M.; Gras, J.; Sikorski, R.; Radhakrishnan, R.; Gelb, L.; Gubbins, K. E. Phase transitions in pores: Experimental and simulation studies of melting and freezing. Langmuir 1999, 15, 6060–6069. 33 ACS Paragon Plus Environment

Page 34 of 37

Page 35 of 37

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


(35) Wallacher, D.; Knorr, K. Melting and freezing of Ar in nanopores. Phys. Rev. B 2001, 63, 104202. (36) Xu, Q.; Sharp, I. D.; Yuan, C. W.; Yi, D. O.; Liao, C. Y.; Glaeser, A. M.; Minor, A. M.; Beeman, J. W.; Ridgway, M. C.; Kluth, P. et al. Large Melting-Point Hysteresis of Ge Nanocrystals Embedded in SiO[sub 2]. Phys. Rev. Lett. 2006, 97, 155701–4. (37) Kondrashova, D.; Reichenbach, C.; Valiullin, R. Probing pore connectivity in random porous materials by scanning freezing and melting experiments. Langmuir 2010, 26, 6380–6385. (38) Jelassi, J.; Castricum, H. L.; Bellissent-Funel, M. C.; Dore, J.; Webber, J. B. W.; SridiDorbez, R. Studies of water and ice in hydrophilic and hydrophobic mesoporous silicas: pore characterisation and phase transformations. Phys. Chem. Chem. Phys. 2010, 12, 2838–2849. (39) Petrov, O.; Furo, I. A study of freezing-melting hysteresis of water in different porous materials. Part I: Porous silica glasses. Microporous Mesoporous Mat. 2011, 138, 221– 227. (40) Dash, J. G. History of the search for continuous melting. Rev. Mod. Phys. 1999, 71, 1737–1743. (41) Dash, J. G.; Rempel, A. W.; Wettlaufer, J. S. The physics of premelted ice and its geophysical consequences. Rev. Mod. Phys. 2006, 78, 695–741. (42) Vanfleet, R. R.; Mochel, J. M. Thermodynamics of Melting and Freezing in Small Particles. Surf. Sci. 1995, 341, 40–50. (43) Denoyel, R.; Pellenq, R. J. M. Simple phenomenological models for phase transitions in a confined geometry. 1: Melting and solidification in a cylindrical pore. Langmuir 2002, 18, 2710–2716. 34 ACS Paragon Plus Environment


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

(44) Petrov, O.; Furo, I. Curvature-dependent metastability of the solid phase and the freezing-melting hysteresis in pores. Phys. Rev. E 2006, 73, 011608. (45) Kondrashova, D.; Valiullin, R. Improving structural analysis of disordered mesoporous materials using NMR cryoporometry. Microporous Mesoporous Mat. 2013, 178, 15–19. (46) Bennema, P.; Vandereerden, J. P. Crystal-growth from solution - Development in computer-simulation. J. Cryst. Growth 1977, 42, 201–213. (47) Nenow, D. Surface premelting. Prog. Cryst. Growth Charact. Mater. 1984, 9, 185–225. (48) Conrad, E. H. Surface roughening, melting, and faceting. Prog. Surf. Sci. 1992, 39, 65–116. (49) Jackson, K. A. Kinetic Processes: Crystal Growth, Diffusion, and Phase Transitions in Materials; Wiley-VCH: Weinheim, 2004. (50) Liu, Z.; Bando, Y.; Mitome, M.; Zhan, J. Unusual Freezing and Melting of Gallium Encapsulated in Carbon Nanotubes. Phys. Rev. Lett. 2004, 93, 095504–4. (51) Binder, K. Applications of Monte Carlo methods to statistical physics. Rep. Prog. Phys. 1997, 60, 487–559. (52) Cahn, J. W.; Dash, J. G.; Fu, H. Y. Theory of ice premelting in monosized powders. J. Cryst. Growth 1992, 123, 101–108. (53) Engemann, S.; Reichert, H.; Dosch, H.; Bilgram, J.; Honkim,; auml,; ki, V.; Snigirev, A. Interfacial Melting of Ice in Contact with SiO2. Phys. Rev. Lett. 2004, 92, 205701. (54) Petrov, O.; Vargas-Florencia, D.; Furo, I. Surface Melting of Octamethylcyclotetrasiloxane Confined in Controlled Pore Glasses: Curvature Effects Observed by 1H NMR. J. Phys. Chem. B 2007, 111, 1574–1581.


Page 36 of 37

Page 37 of 37

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


(55) Cohan, L. H. Sorption hysteresis and the vapor pressure of concave surfaces. J. Am. Chem. Soc. 1938, 60, 433–435. (56) Schneider, D.; Valiullin, R.; Monson, P. A. Filling Dynamics of Closed End Nanocapillaries. Langmuir 2014, 30, 1290–1294. (57) Zeng, Y. H.; Phadungbut, P.; Do, D. D.; Nicholson, D. Anatomy of Adsorption in Open-End and Closed-End Slit Mesopores: Adsorption, Desorption, and Equilibrium Branches of Hysteresis Loop. J. Phys. Chem. C 2014, 118, 25496–25504. (58) Sun, G.; Giovambattista, N.; Wang, E.; Xu, L. Effects of surface structure and solvophilicity on the crystallization of confined liquids. Soft Matter 2013, 9, 11374– 11382. (59) Hitchcock, I.; Holt, E. M.; Lowe, J. P.; Rigby, S. P. Studies of freezing-melting hysteresis in cryoporometry scanning loop experiments using NMR diffusometry and relaxometry. Chem. Eng. Sci. 2011, 66, 582–592.


Freezing and Melting Transitions under Mesoscalic Confinement

Recommend Documents