crystal nucleation in liquids has countless practical consequences 
 in science and technology , and it also affects our everyday experience . 
 one obvious example is the formation of ice , which influences global 
 phenomena such as climate change , as well as 
 processes happening at the nanoscale , such as intracellular freezing . on the other hand , controlling nucleation of molecular crystals 
 from solutions is of great importance to pharmaceuticals , particularly 
 in the context of drug design and production , as the early stages 
 of crystallization impact the crystal polymorph obtained .
even the multibillion - dollar oil industry is affected by the nucleation 
 of hydrocarbon clathrates , which can form inside pipelines , endangering 
 extraction .
finally , crystal nucleation is involved 
 in many processes spontaneously occurring in living beings , from the 
 growth of the beautiful nautilus shells to the dreadful formation in our own brains of amyloid fibrils , 
 which are thought to be responsible for many neurodegenerative disorders 
 such as alzheimer s disease .
each 
 of the above scenarios starts from a liquid below its melting 
 temperature .
this supercooled liquid(12 ) is doomed , according to thermodynamics , to face a first - order 
 phase transition , leading to a crystal . before this 
 can happen , however , a sufficiently large cluster of crystalline atoms 
 ( or molecules or particles ) must form within the liquid , such that 
 the free energy cost of creating an interface between the liquid and 
 the crystalline phase will be overcome by the free energy gain of 
 having a certain volume of crystal .
this event stands at the heart 
 of crystal nucleation , and how this process has been , is , and will 
 be modeled by means of computer simulations is the subject of this 
 review .
the past few decades have witnessed an impressive body 
 of experimental 
 work devoted to crystal nucleation .
for instance , thanks to novel 
 techniques such as transmission electron microscopy at very low temperatures 
 ( cryo - tem ) , we are now able to peek in real time into the early stages 
 of crystallization
. a substantial effort 
 has also been made to understand which materials , in the form of impurities 
 within the liquid phase , can either promote or inhibit nucleation 
 events , a common scenario known as heterogeneous 
 nucleation . however , our understanding of crystal nucleation is far 
 from being complete .
this is because the molecular ( or atomistic ) 
 details of the process are largely unknown because of the very small 
 length scale involved ( nanometers ) , which is exceptionally challenging 
 to probe in real time even with state - of - the - art measurements .
hence , 
 there is a need for computer simulations , and particularly molecular 
 dynamics ( md ) simulations , where the temporal evolution of the liquid 
 into the crystal is more or less faithfully reproduced .
unfortunately , 
 crystal nucleation is a rare event that can occur on time scales of 
 seconds , far beyond the reach of any conventional md framework . in 
 addition , a number of approximations within the computational models , 
 algorithms , and theoretical framework used have been severely questioned 
 for several decades .
although the rush for computational methods able 
 to overcome this time - scale problem is now more competitive than ever , 
 we are almost always forced to base our conclusions on the ancient 
 grounds of classical nucleation theory ( cnt ) , a powerful theoretical 
 tool that nonetheless dates back 90 years to volmer and weber .
in fact , these are exciting times for 
 the crystal nucleation community , 
 as demonstrated by the many reviews covering several aspects of this 
 diverse field .
this particular review is focused almost exclusively on md simulations 
 of crystal nucleation of supercooled liquids and supersaturated solutions . 

we take into account several systems , from colloidal liquids to natural 
 gas hydrates , highlighting long - standing issues as well as recent 
 advances .
although we review a substantial fraction of the theoretical 
 efforts in the field , mainly from the past decade , our goal is not 
 to discuss in detail every contribution .
instead , we try to pinpoint 
 the most pressing issues that still prevent us from furthering our 
 understanding of nucleation .
we 
 introduce the theoretical framework of cnt ( section 1.1 ) , the state - of - the - art experimental techniques 
 ( section 1.2 ) , and 
 the md - based simulation methods ( section 1.3 ) that in the past few decades have provided 
 insight into nucleation . in section 2 , we put such computational approaches into context , 
 describing both achievements and open questions concerning the molecular 
 details of nucleation for different types of systems , namely , colloids 
 ( section 2.1 ) , lennard - jones 
 ( lj ) liquids ( section 2.2 ) , atomic liquids ( section 2.3 ) , water ( section 2.4 ) , nucleation from solution ( section 2.5 ) , and natural gas hydrates ( section 2.6 ) . in the third 
 and last part of the article ( section 3 ) , we highlight future perspectives and open challenges 
 in the field .
almost 
 every computer simulation of crystal nucleation in liquids invokes 
 some elements of classical nucleation 
 theory ( cnt ) .
this theory has been discussed in great detail elsewhere , and we describe it here for the sake of completeness and also to 
 introduce various terms used throughout the review .
cnt was formulated 90 years ago 
 through the contributions of volmer and weber , farkas , becker and dring , and zeldovich , on 
 the basis of the pioneering ideas of none other than gibbs himself .
cnt was created to describe the condensation 
 of supersaturated vapors into the liquid phase , but most of the concepts 
 can also be applied to the crystallization of supercooled liquids 
 and supersaturated solutions . according to cnt ,
clusters of crystalline 
 atoms ( or particles or molecules ) of any size are treated as macroscopic 
 objects , that is , homogeneous chunks of crystalline phase separated 
 from the surrounding liquid by a vanishingly thin interface .
this 
 apparently trivial assumption is known as the capillarity approximation , 
 which encompasses most of the strengths and weaknesses of the theory . 
 according to the capillarity approximation , the interplay between 
 the interfacial free energy , , and the difference in 
 free energy between 
 the liquid and the crystal , , fully describes the thermodynamics of 
 crystal nucleation . in three dimensions , the free energy of formation , , for a spherical crystalline nucleus of 
 radius r
can thus be written as the sum of a surface 
 term and a volume term1this function , sketched in figure 1 , displays a maximum corresponding 
 to the so - called critical nucleus size n*2where   is the number density of the crystalline 
 phase .
the critical nucleus size represents the number of atoms that 
 must be included in the crystalline cluster for the free energy difference , , to match the free energy cost due to the 
 formation of the solid  liquid interface .
clusters of crystalline 
 atoms occur within the supercooled liquid by spontaneous , infrequent 
 fluctuations , which eventually lead the system to overcome the free 
 energy barrier for nucleation3triggering the actual crystal growth ( see figure 1 ) .
sketch of the free energy 
 difference , , as a function of the crystalline nucleus 
 size n. a free energy barrier for nucleation , , must 
 be overcome to proceed from the ( metastable ) 
 supercooled liquid state to the thermodynamically stable crystalline 
 phase through homogeneous nucleation ( purple ) .
heterogeneous nucleation 
 ( green ) can be characterized by a lower free energy barrier , , and a smaller critical nucleus size , nhet * , whereas 
 in the case of spinodal decomposition ( orange ) , the supercooled 
 liquid is unstable with respect to the crystalline phase , and the 
 transformation to the crystal proceeds in a barrierless fashion .
the 
 three snapshots depict a crystalline cluster nucleating within the 
 supercooled liquid phase ( homogeneous nucleation ) or as a result of 
 the presence of a foreign impurity ( heterogeneous nucleation ) , as 
 well as the simultaneous occurrence of multiple crystalline clusters 
 in the unstable liquid .
this scenario is often labeled as spinodal 
 decomposition , although the existence of a genuine spinodal decomposition 
 from the supercooled liquid to the crystalline phase has been debated 
 ( see text ) .
the kinetics of crystal 
 nucleation is typically addressed by assuming 
 that no correlation exists between successive events increasing or 
 decreasing the number of constituents of the crystalline nucleus . 
 in other words , the time evolution of the nucleus size
is presumed 
 to be a markov process , in which atoms in the liquid either order 
 themselves one by one in a crystalline fashion or dissolve one by 
 one into the liquid phase .
in addition , we state that every crystalline 
 nucleus lucky enough to overcome the critical size n * quickly grows to macroscopic dimensions on a time scale much smaller 
 than the long time required for that fortunate fluctuation to come 
 about .
if these conditions are met , the 
 nucleation rate , that is , the probability per unit time per unit volume 
 of forming a critical nucleus does not depend on time , leading to 
 the following formulation of the so - called steady - state nucleation 
 rate 4where kb is the 
 boltzmann constant and   is a prefactor that we discuss 
 later .
the 
 steady - state nucleation rate is the central quantity in the description 
 of crystallization kinetics , as much as the notion of critical nucleus 
 size captures most of the thermodynamics of nucleation .
all 
 quantities specified up to now depend on pressure and most 
 notably temperature . in most cases , the interfacial free energy , ,
is assumed to be linearly 
 dependent on 
 temperature , whereas the free energy difference between the liquid 
 and solid  phases , , is proportional to the supercooling ,   ( or the supersaturation ) .
several approximations 
 exist to treat the temperature dependence of ( 36 ) and , which can 
 vary substantially for different supercooled liquids . in any case
, it follows from eq 3 that the free energy barrier for nucleation ,   decreases with supercooling . in other words , 

the farther one is from the melting temperature , the larger the thermodynamic driving force 
 for nucleation is .
interestingly , in the case of supercooled 
 liquids , kinetics goes 
 the other way , as the dynamics of the liquid slow down with supercooling , 
 thus hindering the occurrence of nucleation events .
in fact , although 
 a conclusive expression for the prefactor the latter is still lacking , it is usually written within cnt 
 as5where   is the number of possible 
 nucleation sites 
 per unit volume ,   is the zeldovich 
 factor ( accounting for the fact that several postcritical 
 clusters might 
 still shrink without growing into the crystalline phase ) , and   is a kinetic prefactor .
the latter should represent the attachment 
 rate , that is , the frequency with which the particles in 
 the liquid phase 
 reach the cluster rearranging themselves in a crystalline fashion . 

however , in a dense supercooled liquid ,   also 
 quantifies the ease with which the 
 system explores configurational space , effectively regulating the 
 amplitude of the fluctuations possibly leading to the formation of 
 a crystalline nucleus . in short
, we can safely say that   involves the atomic or molecular 
 mobility 
 of the liquid phase , more often than not quantified in terms of the 
 self - diffusion coefficient , which obviously 
 decreases with supercooling .
thus , for a supercooled liquid , the competing 
 trends of   and   lead , in the case of diffusion - limited 
 nucleation , to a maximum in the nucleation 
 rate , as depicted in figure 2 .
the same arguments apply when dealing with processes such 
 as the solidification of metallic alloys . in the case of nucleation from solutions ,   and  
however , 
 the dependence of the kinetic prefactor on supersaturation is much 
 weaker than the temperature dependence of   characteristic 
 of supercooled liquids . 
 as a result
, there is usually no maximum in the nucleation rate as 
 a function of supersaturation for nucleation from solutions .
illustration of how certain quantities from cnt vary as 
 a function 
 of supercooling , t , for supercooled liquids . 

the free energy difference between the liquid and the solid phase , the interfacial free energy , and the kinetic prefactor   are reported as functions of t in a generic case of diffusion - limited nucleation , characterized 
 by a maximum in the steady - state nucleation rate .  
is zero at the melting temperature , and   is vanishingly 
 small at the glass transition 
 temperature .
although   is supposed to play a minor role 
 compared 
 to the exponential term in eq 4 , the kinetic prefactor has been repeatedly blamed for the 
 quantitative disagreement between experimental measurements and computed 
 crystal nucleation rates .
atomistic simulations 
 could , in principle , help to clarify the temperature dependence as 
 well as the microscopic origin of   and also 
 of the thermodynamic ingredients 
 involved in the formulation of cnt .
however , quantities such as   are not only infamously 
 difficult to converge 
 within decent levels of accuracy but can even 
 be ill - defined in many situations .
for instance , it remains to be 
 seen whether , which , in principle , 
 refers to a planar 
 interface under equilibrium conditions , can be safely defined when 
 dealing with small crystalline clusters of irregular shapes .
in fact , 
 the early stages of the nucleation process often involve crystalline 
 nuclei whose size and morphology fluctuate on a time scale shorter 
 than the structural relaxation time of the surrounding liquid . moreover , 
 the dimensions of such nuclei can be of the same order as the diffuse 
 interface between the liquid and the solid phases , thus rendering 
 the notion of a well - defined   value quite dangerous . 
 as an example , joswiak 
 et al .
recently showed that , for liquid 
 water droplets ,   can strongly depend on 
 the curvature of 
 the droplet .
the mismatch between the macroscopic interfacial free 
 energy and its curvature - dependent value can spectacularly affect 
 water - droplet nucleation , as reported by atomistic simulations of 
 droplets characterized by radii on the order of 0.51.5 
 nm .
some other quantities , such as the size of the critical cluster , 
 depend in many cases rather strongly on the degree of supercooling . 

this is the case , for example , for the critical nucleus size n * , which can easily span 2 orders of magnitude in just 
 10 c of supercooling .
given the old 
 age of cnt , it is no surprise that substantial efforts have been devoted 
 to extend and/or improve its original theoretical framework .
the most 
 relevant modifications possibly concern the issue of two - step nucleation . 

many excellent works have reviewed this subject extensively ( see , 
 e.g. , refs ( 18 , 24 , 52 , and 53 ) ) , so that 
 we provide only the essential concepts here . in the original 
 formulation of cnt , the system has to overcome a single free energy 
 barrier , corresponding to a crystalline nucleus of a certain critical 
 size , as depicted in figure 3 . when dealing with crystal nucleation from the melt , it is 
 rather common to consider the number of crystalline particles within 
 the largest connected cluster , n , as the natural 
 reaction coordinate describing the whole nucleation process . in many 
 cases ,
the melt is dense enough that local density fluctuations are 
 indeed not particularly relevant and the slow degree of freedom is 
 in fact the crystalline ordering of the particles within the liquid 
 network .
however , one can easily imagine that , in the case of crystal 
 nucleation of molecules in solutions , for example , the situation can 
 be quite different .
specifically , in a realistically supersaturated 
 solution , a consistent fluctuation of the solute density ( concentration ) 
 could be required just to bring a number n of solute molecules close enough to form a connected cluster . assuming 
 that the molecules involved in such a density fluctuation will also 
 order themselves in a crystalline fashion on exactly the same time 
 scale is rather counterintuitive .
schematic comparison of one - step versus 
 two - step nucleation for 
 a generic supersaturated solution .
( a ) sketch of the free energy difference 
 gn , n as a function of the number of solute molecules 
 in the largest  connected
 cluster ( they can be ordered 
 in a crystalline fashion or not ) ( n ) and of the number of crystalline molecules within the largest connected 
 cluster ( n ) .
the one - step mechanism predicted by 
 cnt ( purple ) is characterized by a single free energy barrier for 
 nucleation , gn , n,one - step*. in contrast , the two - step nucleation 
 requires a free energy barrier , gn,two - step * , to be overcome through a local density 
 fluctuation of the solution , leading to a dense , but not crystalline - like , 
 precursor .
the latter can be unstable ( green ) or stable ( orange ) with 
 respect to the liquid phase , being characterized by a higher ( green ) 
 or lower ( orange ) free energy basin .
once this dense precursor has 
 been obtained , the second step consists of climbing a second free 
 energy barrier , gn , two - step * , corresponding 
 to the ordering of the solute molecules within the precursor from 
 a disordered state to the crystalline phase .
( b ) one - step ( purple ) 
 and two - step ( green and orange ) nucleation mechanisms visualized in 
 the density ( n)ordering 
 ( n ) plane .
the one - step mechanism proceeds along 
 the diagonal , as both n and n increase at the same time in such a way that a single 
 free energy barrier has to be overcome . in this scenario ,
the supersaturated 
 solution transforms continuously into the crystalline phase . on the 
 other hand , within a two - step nucleation scenario ,
the system has 
 to experience a favorable density fluctuation along n first , forming a disordered precursor that , in 
 a second step , orders itself in a crystalline fashion , moving along 
 the ( n ) coordinate and ultimately leading to the 
 crystal .
in fact , the formation of crystals 
 from molecules in solution often 
 occurs according to a two - step nucleation mechanism that has no place 
 in the original formulation of cnt . in the prototypical scenario
depicted 
 in figure 3 , a first 
 free energy barrier , gn,two - step * , has to be overcome by means of a density 
 fluctuation of the solute , such that a cluster of connected molecules 
 of size n * is formed .
this object does not yet have any 
 sort of crystalline order , and depending on the system under consideration , 
 it can be either unstable or stable with respect to the supersaturated 
 solution ( see figure 3 ) .
subsequently , the system has to climb a second free energy barrier , 
 gn , two - step * , to order the molecules within the dense 
 cluster in a crystalline - like fashion .
a variety of different nucleation 
 scenarios have been loosely labeled as two - step , from crystal nucleation 
 in colloids ( see section 2.1 ) or lennard - jones liquids ( see section 2.2 ) to the formation of crystals of urea 
 or nacl ( see section 2.5 ) , not to mention biomineralization ( see , e.g. , refs ( 18 and 53 ) ) and protein crystallization ( see , e.g. , refs ( 54 and 55 ) ) . in all of these cases , cnt as it is formulated is simply 
 not capable 
 of dealing with two - step nucleation .
this is why , in the past few 
 decades , a number of extensions and/or modifications of cnt have been 
 proposed and indeed successfully applied to account for the existence 
 of a two - step mechanism . here , we mention the phenomenological theory 
 of pan et al . , who wrote an expression 
 for the nucleation rate assuming a free energy profile similar to 
 the one sketched in figure 3 , where dense metastable states are involved as intermediates 
 on the path toward the final crystalline structure .
the emergence 
 of so - called prenucleation clusters ( pncs ) , namely , stable states 
 within supersaturated solutions , which are known to play a very important 
 role in the crystallization of biominerals , for example , was also 
 recently fit into the framework of cnt by hu et al .
they proposed a modified expression for the excess free 
 energy of the nucleus that takes into account the shape , size and 
 free energy of the pncs as well as the possibility for the pncs to 
 be either metastable or stable with respect to the solution .
a comprehensive 
 review of the subject is offered by the work of gebauer et al .
it is worth noticing that these extensions of 
 cnt are mostly quite recent , as they were triggered by overwhelming 
 experimental evidence for two - step nucleation mechanisms .
cnt is 
 also the tool of the trade for heterogeneous crystal nucleation , that 
 is , nucleation that occurs on account of the presence of a foreign 
 phase ( see figure 1 ) .
in fact , nucleation in liquids occurs heterogeneously more often 
 than not , as in some cases , the presence of foreign substances in 
 contact with the liquid can significantly lower the free energy barrier .
a typical 
 example is given by the formation 
 of ice : as we shall see in sections 2.4.1 and 2.4.2 , it 
 is surprisingly difficult to freeze pure water , which invariably takes 
 advantage of a diverse portfolio of impurities , from clay minerals 
 to bacterial fragments , to facilitate 
 the formation of ice nuclei .
heterogeneous nucleation is customarily 
 formulated within the cnt framework in terms of geometric arguments .
specifically6where f( )  
 1 is the shape factor , a quantity that accounts for 
 the fact that three different interfacial free energies must  be  balanced : , , and .
for instance , considering a supercooled 
 liquid nucleating on top of an ideal planar surface offered by the 
 foreign phase , we obtain the so - called young s relationwhere  is the contact angle , namely , 
 a measure of the extent to which the crystalline nucleus wets the foreign surface .
thus , the contact angle determines whether 
 and how much it could be easier for a critical nucleus to form in 
 an heterogeneous fashion , as for 0   <  , 
 the volume - to - surface energy ratio   is larger for 
 the spherical cap nucleating 
 on the foreign surface than for the sphere nucleating in the liquid . 

this simple formulation is clearly only a rough approximation of what 
 happens in reality . at first
, the contact angle is basically a macroscopic 
 quantity , of which the microscopic equivalent is in most cases ill - defined 
 on the typical length scales involved in the heterogeneous nucleation 
 process .
in addition , in most cases , 
 the nucleus will not be shaped like a spherical cap , and to make things 
 more complicated , many different nucleation sites with different morphologies 
 typically exist on the same impurity . finally , the kinetic prefactor , , becomes even more obscure in 
 heterogeneous 
 nucleation , as it is plausible that the foreign phase will affect 
 the dynamical properties of the supercooled liquid . moving toward strong 
 supercooling , several things can happen to the 
 supercooled liquid phase . whether one can avoid the glass transition 
 largely depends on the specific liquid under consideration and on 
 the cooling rate ( see , e.g. , ref ( 58 ) ) . assuming that the system can be cooled sufficiently 
 slowly , hence avoiding both the glass transition and crystal nucleation , 
 one can , in principle , enter a supercooled regime in which the liquid 
 becomes unstable with respect to the crystalline phase .
this region 
 of the phase diagram is known as the spinodal region , where the tiniest perturbation , for example , of the local density 
 or the degree of ordering leads the system toward the crystalline 
 phase without paying anything in terms of free energy ( see figure 1 ) .
in fact , below 
 a certain critical temperature , , the free energy barrier for nucleation 
 is zero , and the liquid transforms spontaneously into the crystal 
 on very short time scales .
the same picture holds for molecules in 
 solution , as nicely discussed by gebauer et al . , and it can not , by definition , be described by conventional 
 cnt , according to which a small   value 
 persists even at the strongest supercoolings . although spinodal regimes have been observed 
 in a variety of scenarios , the existence 
 of a proper spinodal decomposition from the supercooled liquid to 
 the crystalline phase has been debated ( see , e.g. , ref ( 61 ) ) .
enhanced - sampling md 
 simulations , which we discuss in section 2.2 , have suggested 
 that barrierless crystal nucleation is possible at very strong supercooling , 
 whereas other works claim that this is not the case ( see , e.g. , ref ( 63 ) ) .
here , we simply note 
 that , at strong supercooling  not necessarily within the presumed 
 spinodal regime  a number of assumptions on which cnt relies 
 become , if not erroneous , ill - defined .
the list is long , and in fact , 
 a number of nucleation theories able 
 to at least take into account the emergence of a spinodal decomposition 
 exist , although they have mostly been formulated for condensation 
 problems . in any case
, the capillarity approximation is most likely 
 to fail at strong supercoolings , as the size of the critical nucleus 
 becomes exceedingly small , down to losing its meaning in the event 
 of a proper spinodal decomposition .
moreover , we shall see , for instance , 
 in section 2.2 , that 
 the shape of the crystalline clusters is anything but spherical at 
 strong supercooling and , at the same time , the kinetic prefactor assumes 
 a role of great importance .
in fact , nucleation at strong supercooling 
 might very well be dominated by , as the mobility of the supercooled 
 liquid 
 is what really matters when the free energy barrier for nucleation 
 approaches vanishingly small values .
strong supercooling is important 
 because this is the regime in which most computational studies have 
 been performed .
large values of t imply high 
 nucleation rates and smaller critical nuclei , although as one moves 
 away from , most of the assumptions of cnt 
 are progressively 
 invalidated . at this point ,
given the substantial approximations 
 of cnt and especially its old age , the 
 reader might 
 be waiting for us to introduce the much more elegant , accurate , and 
 comprehensive theories that experiments and simulations surely embrace 
 today . sadly , this is not the case .
many of them , such as dynamical nucleation theory , mean - field kinetic nucleation theory , and coupled flux theory , are mainly limited to condensation 
 problems , and some others have only rarely been applied , for example , 
 to crystallization in glasses , such as 
 diffuse interface theory .
several improvements 
 on cnt have been proposed , targeting specific aspects such as the 
 shape of the crystalline nuclei or the 
 finite size of the nonsharp crystal  liquid interface .
nucleation theories largely unrelated to cnt 
 can also be found , such as classical density functional theory ( cdft ) ( classical , not to be confused with the celebrated 
 quantum mechanical framework of hohenberg and kohn ) .
a fairly complete inventory of nucleation theories , together 
 with an excellent review of nucleation in condensed matter , can be 
 found elsewhere . here
, we do not discuss 
 the details of any of these approaches , as indeed none of them has 
 been consistently used to model crystal nucleation in liquids .
this 
 is because cnt , despite having many shortcomings , is a simple yet 
 powerful theory that is able to capture at least qualitatively the 
 thermodynamics and kinetics of nucleation for very different systems , 
 from liquid metals to organic crystals .
it has been extended to include 
 heterogeneous nucleation , and it is fairly easy to modify it to take 
 into consideration multicomponent systems such as binary mixtures 
 as well .
several different 
 experimental approaches have been employed to understand the thermodynamics 
 and kinetics of crystal nucleation in liquids .
although this review 
 discusses theory and simulations almost exclusively , we present in 
 this section a concise overview of the state - of - the - art experimental 
 techniques to highlight their capabilities as well as their limitations .
a schematic synopsis focusing on both spatial and temporal resolutions 
 is sketched in figure 4 , and an inventory of notable applications is reported in table 1 . as already stated , 
 nucleation is a dynamical process usually occurring on very small 
 time and length scales ( nanoseconds and nanometers , respectively ) . 

thus , obtaining the necessary spatial and temporal resolutions is 
 a tough technical challenge .
overview of some of the experimental methods 
 that have been applied 
 to characterize nucleation .
ranges of the spatial and temporal resolutions 
 typical of each approach are reported on the x and y axes , respectively .
for instance , colloids offer a playground where 
 simple microscopy can image the particles involved in the nucleation 
 events , which occur on such long time scales ( seconds ) that a full 
 characterization in time of the process has been achieved .
specifically , confocal microscopy has led to three - dimensional imaging 
 of colloidal systems , unraveling invaluable information about the 
 critical nucleus size , for example . in a 
 similar fashion ,
sleutel et al . achieved molecular resolution 
 of the formation of two - dimensional glucose isomerase crystals by 
 means of atomic force microscopy .
this 
 particular investigation featured actual movies showing both crystal 
 growth and the dissolution of precritical clusters , as well as providing 
 information about the influence of the substrate .
in addition , cryo - tem 
 techniques have recently provided two - dimensional snapshots of nucleation 
 events at very low temperatures . in selected cases , where the time 
 scales involved are again on the order of seconds , dynamical details 

have been obtained , as in the cases of caco3 , metal phosphate , and magnetite .
however , more often than not , crystal 
 nucleation in liquids takes 
 place within time windows too small ( nanoseconds ) to allow for a sequence 
 of snapshots to be taken with high - spatial - resolution instruments . 
 in these cases ,
microscopic insights can not be obtained , and much 
 more macroscopic measurements have to be performed . in this 
 context ,
several experimental approaches aim at examining 
 a large number of independent nucleation events for a whole set of 
 rather small configurations of the system , basically performing an 
 ensemble average .
for example , in droplet experiments , nucleation 
 is characterized as a function of time or temperature .
freezing is 
 identified for each nucleation event within the ensemble of available 
 configurations by techniques such as femtosecond x - ray scattering , optical microscopy , and powder x - ray diffraction . from these data , the nucleation rate is often reconstructed by measuring 
 either metastable zone widths or induction times ( several examples are listed in , e.g. , refs ( 111115 ) ) , thus providing a solid connection to theoretical frameworks such 
 as cnt ( see section 1.1 ) . an essential technical detail within this class of measurements 

is that the volume available for each nucleation event has to be as 
 small as possible to reduce the occurrence of multiple nucleation 
 events within the same configuration .
high - throughput devices such 
 as the lab - on - a - chip can significantly 
 improve the statistics of the nucleation events , thus enhancing the 
 capabilities of these approaches .
another line of action focuses 
 on the study of large , macroscopic 
 systems .
freezing is detected by techniques such as differential scanning 
 calorimetry , fourier transform infrared spectroscopy ( ftirs ) , and analytical ultracentrifugation or by some flavor of chamber experiments . in this case , the frozen fraction of the overall system and/or the 
 nucleation temperatures can be obtained , and in some cases , nucleation 
 rates have been extracted ( see table 1 ) .
finally , experimental methods that can detect 
 nucleation and the 
 formation of the crystal ( predominantly by means of optical microscopy ) 
 but do not provide any microscopic detail have helped to shed light 
 on issues such as the role of the solvent or impurities .
this is usually 
 possible by examining the amount of crystalline phase obtained along 
 with its structure .
even though there are a large number of 
 powerful experimental techniques 
 and new ones emerging ( e.g. , ultrafast x - ray ) , it is still incredibly challenging to obtain microscopic - level 
 insight into nucleation from experiments .
as we shall see now , md 
 simulations provide a powerful complement to experiments .
when dealing 
 with crystal nucleation in liquids , atomistic simulations should provide 
 a detailed picture of the formation of the critical nucleus . the simplest 
 way to achieve
this is by so - called brute - force md simulations , which 
 involve cooling the system to below the freezing temperature and then 
 following its time evolution until nucleation is observed .
brute - force 
 simulations are the antagonist of enhanced - sampling simulations , where 
 specific computational techniques are used to alter the dynamics of 
 the system so as to observe nucleation on a much shorter time scale . 

monte carlo ( mc ) techniques , although typically coupled with enhanced 
 sampling techniques , can be used to recover , but the calculation of   requires 
 other methods , such as kinetic 
 monte carlo ( kmc ) .
the natural choice 
 to simulate nucleation events is instead md simulations , which directly 
 provide the temporal evolution of the system .
isobaric 
 ensemble ( npt ) , where p ( usually 
 ambient pressure ) and t <   are kept constant by means of a barostat 
 and a thermostat , respectively . such computational tweaking is a double - edged 
 sword .
in fact , nucleation and most notably crystal growth are exothermic 
 processes , and within the length scale 
 probed by conventional atomistic simulations ( 110 nm ) , it is necessary to keep the system at constant temperature . 

on the other hand , in this way , dynamical and structural effects in 
 both the liquid and the crystalline phases due to the heat developed 
 during the nucleation events are basically neglected . although the actual extent of these effects is not yet clear , forcing 
 the sampling of the canonical ensemble is expected to be especially 
 dangerous when dealing with very small systems affected by substantial 
 finite - size effects . more importantly
small coupling constants and clever approaches 
 ( e.g. , stochastic thermostats ) can be 
 employed to limit the effects of the thermostats , but in general , 
 care must be taken . the same reasoning applies for p and barostats as well .
a density change of the system is usually 
 associated with nucleation , the crystalline 
 phase being more ( or less , in the case of , e.g. , water ) dense than 
 the liquid parent phase .
three conditions must be fulfilled 
 to extract   from brute - force 
 md simulations:(1)the system must be allowed to evolve 
 in time until spontaneous fluctuations lead to a nucleation event.(2)the system size must be 
 significantly 
 larger than the critical nucleus.(3)significant statistics of nucleation 
 events must be collected .
the system must be allowed to evolve 
 in time until spontaneous fluctuations lead to a nucleation event .
the most daunting obstacle is probably the first one because of the 
 so - called time - scale problem . in most cases ,
nucleation 
 is a rare event , meaning that it usually occurs on a very long time 
 scale ; precisely how long depends strongly on t .
a rough estimate of the number of simulation steps required to 
 observe a nucleation event within a molecular dynamics run is reported 
 in figure 5 . under 
 the fairly optimistic assumption that classical md simulations can 
 cope with up to 10 molecules on a time scale of 
 nano-/microseconds
, there is only a very narrow set of conditions 
 for which brute - force classical md simulations could be used to investigate 
 nucleation , usually only at strong supercooling .
time scales typical 
 of first - principles simulations , also reported in figure 5 assuming up to 10 molecules , indicate that unbiased ab initio simulations of 
 nucleation events are unfeasible .
nucleation rate   as a 
 function of the simulation time needed 
 within an md simulation to observe a single nucleation event .
the 
 blue shaded region highlights the approximate simulation times currently 
 affordable by classical md simulations ; clearly , only very fast nucleation 
 processes can be simulated with brute - force md . for homogeneous ice 
 nucleation ,   and   can typically be observed for t = 30 k and 
 t = 80 k , respectively . 
 in the derivations of classical and ab initio simulation times , 10 and 10 molecules , respectively ,
were considered , 
 together with the number density of a generic supercooled liquid ,   = 0.01 molecules .
the number 
 of atoms ( or molecules ) in the system defines the time scale accessible 
 to the simulation and , thus , the severity of the time - scale problem . 

the reason large simulation boxes , significantly larger than the size 
 of the critical nucleus , are needed is because periodic boundary conditions 
 will strongly affect nucleation ( and growth ) if even the precritical 
 nuclei are allowed to interact with themselves .
this issue worsens at mild 
 supercooling , where the critical nucleus size rapidly increases toward 
 dimensions not accessible by md simulations .
third , it is not 
 sufficient to collect information on just one 
 nucleation event .
nucleation is a stochastic event following a poisson 
 distribution ( at least ideally ; see section 1.1 ) , and so to obtain the nucleation rate , 
 one needs to accumulate decent statistics . taking
these issues 
 into consideration , various approaches for 
 obtaining   have emerged . 

one approach , known as the 
 yasuoka  matsumoto method , involves 
 simulating a very large system , so that different nucleation events 
 can be observed within a single run . in this case ,
large simulation 
 boxes are needed to collect sufficient statistics and to avoid spurious 
 interactions between different nuclei .
another family of methods involves 
 running many different simulations using much smaller systems , which 
 is usually computationally cheaper .
once a collection of nucleation 
 events has been obtained , several methods for extracting   can be employed . 

the simplest ones ( mean 
 lifetime and survival probability methods ) involve the fitting of the nucleation times to poisson 
 statistics .
a more in - depth technique , the so - called mean first - passage 
 method , allows for a detailed analysis 
 of the nucleus population but requires a probability distribution 
 in terms of nucleus size .
the literature offers a notable number 
 of works in which brute - force 
 md simulations have been successfully applied .
most of them rely on 
 one approach to circumvent the above - mentioned issues , particularly 
 the time - scale problem . as we shall see in section 2 , to simulate nucleation events , one almost 
 always has to either choose a very simple system or increase the level 
 of approximation sometimes dramatically , for instance , by coarse - graining 
 the interatomic potential used . in 
 the previous section
, we introduced the time - scale problem , the main 
 reason brute - force md simulations are generally not feasible when 
 studying crystal nucleation . enhanced sampling methods
alter how the 
 system explores its configurational space , so that nucleation events 
 can be observed within a reasonable amount of computational time . 
 broadly speaking
, one can distinguish between free - energy methods 
 and path - sampling methods , both of which have been extensively discussed 
 elsewhere ( see , e.g. , refs ( 156 and 162164 ) ) .
thus , only the briefest of introductions is needed here . of the many enhanced sampling methods , only a handful
this is because information 
 is needed about both the thermodynamics of the system ( the free energy 
 barrier for nucleation ) and 
 the kinetics of the nucleation process 
 ( the kinetic prefactor ) .
when dealing with crystal nucleation 
 in supercooled liquids , free - energy - based methods are rather common , 
 such as umbrella sampling ( us ) and metadynamics . in both cases , and indeed in almost all enhanced sampling methods 
 currently available , the free energy surface of the actual system 
 is coarse - grained by means of one or more order parameters or collective 
 variables .
the choice of the order parameter is not trivial and can 
 have dramatic consequences .
an external bias is then applied to the 
 system , leading to a modified sampling of the configurational space 
 that allows for the reconstruction of the free energy profile with 
 respect to the chosen order parameter and , thus , for the computation 
 of the free energy barrier .
however , there is a price to be paid : upon introduction 
 of an extra term into the system hamiltonian , the actual dynamics 
 of the system is to some extent hampered , and much of the insight 
 into the nucleation mechanism is lost .
one 
 needs complementary methods , usually 
 aimed at estimating the probability for the system on top of the nucleation 
 barrier  in the space of the selected order parameter  to 
 get back to the liquid phase or to evolve into the crystal .
most frequently , 
 such methods are based on some flavor of transition state theory , such as the bennett  chandler formulation , and require a massive set of md or kmc simulations to be performed . on the other hand , the ever - growing family of path - sampling methods 
 can provide direct access to the kinetics of the nucleation process . 

these approaches again rely on the definition of an order parameter , 
 but instead of applying an external bias potential , an importance 
 sampling is performed so as to enhance the naturally occurring fluctuations 
 of the system . within the majority of the path - sampling approaches 
 currently used , including transition interface sampling ( tis ) and forward flux sampling ( ffs ) , the 
 ensemble of paths connecting the liquid and the crystal is divided 
 into a series of interfaces according to different values of the order 
 parameter . by sampling the probability with which the system crosses 
 each of these interfaces , a cumulative probability directly related 
 to the nucleation rate can be extracted .
other path - sampling techniques 
 such as transition path sampling ( tps ) rely 
 instead on the sampling of the full ensemble of the reactive trajectories . 
 in both cases , by means of additional simulations involving , for example , 
 committor analysis distribution and 
 thermodynamic integration , one can subsequently 
 extract the size of the critical nucleus and the free energy difference 
 between the solid and liquid phases , respectively .
many different 
 path - sampling methods are available , but to our knowledge , only tps , 
 tis , and most prominently ffs have allowed for estimates of crystal 
 nucleation rates . under certain conditions , path - sampling methods 
 do not alter the dynamics of the system , allowing for invaluable insight 
 into the nucleation mechanism .
however , they are particularly sensitive 
 to the slow dynamics of strongly supercooled systems , which hinder 
 the sampling of the paths and makes them exceptionally expensive computationally . 
 although the past few decades have taught us that enhanced sampling 
 techniques are effective in tackling crystal nucleation of colloids 
 ( see section 2.1 ) , 
 lennard - jones melts ( see section 2.2 ) , and other atomic liquids ( see section 2.3 ) , only recently have these 
 techniques been applied to more complex systems .
one challenging 
 scenario for simulations of nucleation is provided 
 by the formation of crystals from solutions characterized by very 
 low solute concentration .
although this occurrence is often encountered 
 in real systems of practical interest , it is clearly extremely difficult 
 for md simulations , even if aided by conventional enhanced sampling 
 techniques , to deal with just a few solute molecules dissolved within 
 1010 solvent molecules . in these cases , 

the diffusion of the solute plays a role of great relevance , and the 
 interaction between solvent and solute can enter the nucleation mechanism 
 itself .
thus , obtaining information about the thermodynamics , let 
 alone the kinetics , of nucleation at very low solute concentrations 
 is presently a formidable task .
however , efforts have been devoted 
 to further our understanding of solute migration and solute - nuclei 
 association , for example , as demonstrated by the pioneering works 
 of gavezzotti and co - workers and more recently 
 by kawska and co - workers . in the latter work , 

the authors illustrate an approach that relies on the modeling of 
 the subsequent growth step , where solute particles ( often ions ) are 
 progressively added to the ( crystalline or not ) cluster . after each 
 of these growth steps , a structural optimization of the cluster and 
 the solvent by means of md simulations
although this 
 method can not provide quantitative results in terms of the thermodynamics 
 and/or the kinetics of nucleation , it can , in principle , provide valuable 
 insight into the very early stages of crystal nucleation when dealing 
 with solutions characterized by very low solute concentrations . on a final note , we mention seeded md simulations .
this technique 
 relies on simulations in which a crystalline nucleus of a certain 
 size is inserted into the system at the beginning of the simulation . 

although useful information about critical nucleus size can be obtained 
 in this way , the method does not usually allow for a direct calculation of the 
 nucleation rate . however , seeded md simulations are one of the very 
 few methods by which it is currently possible to investigate solute 
 precipitate nucleation ( see , e.g. , knott et al . ) . in this case
, the exceedingly low attachment rate of 
 the solute often prevents both free - energy- and transition - path - sampling - enhanced 
 sampling methods from being applied effectively . as we shall 
 see in the next few sections , the daunting computational 
 costs , together with the delicate choice of order parameter and the 
 underlying framework of cnt ,
still make enhanced - sampling simulations 
 of crystal nucleation in liquids an intimidating challenge .
we have chosen to review 
 different classes of systems , which we 
 present in order of increasing complexity .
we start in sections 2.1 and 2.2 with colloids 
 and lennard - jones liquids , respectively .
these systems are described 
 by simple interatomic potentials that allow large - scale md simulations , 
 and thus with them , many aspects of cnt can be investigated and nucleation 
 rates calculated . in some cases , the latter can can be directly compared 
 to experimental results . as such , colloids and
lennard - jones liquids 
 represent a sort of benchmark for md simulations of crystal nucleation 
 in liquids , although we shall see that our understanding of crystal 
 nucleation is far from satisfactory even within these relatively easy 
 playgrounds . in section 2.3 ,
we discuss selected atomic liquids of technological interest 
 such as liquid metals , supercooled liquid silicon , and phase - change 
 materials for which nucleation occurs on very small time scales . as 
 the first example of a molecular system ,
we review the body of computational 
 work devoted to unraveling both the homogeneous ( section 2.4.1 ) and heterogeneous ( section 2.4.2 ) formation 
 of ice , offering a historical perspective guiding the reader through 
 the many advances that have furthered our understanding of ice nucleation 
 in the past decades .
next , we present an overview of nucleation 
 from solution ( section 2.5 ) , where simulations 
 have to deal with solute and solvent .
we take into account systems 
 of great practical relevance such as urea molecular crystals , highlighting 
 the complexity of the nucleation mechanism , which is very different 
 from what cnt predicts .
finally , section 2.6 is devoted to the formation of gas hydrates . as a general rule ,
increasing the complexity of the system raises 
 more questions about the validity of the assumptions underpinning 
 cnt .
the reader will surely notice that simulations have revealed 
 many drawbacks of cnt along the way and that reaching decent agreement 
 for the nucleation rate   between experiments 
 and simulations still 
 remains a formidable task .
one reason for this is 
 the simplicity of the interatomic potential customarily used to model 
 them : the only interaction a hard - sphere particle experiences comes 
 from elastic collisions with other particles . because there is no 
 attractive force between particles , a hard - sphere system is entirely 
 driven by entropy . as a consequence ,
the phase diagram is very simple 
 and can be entirely described with one single parameter , the volume 
 fraction . only two different phases are possible : a fluid 
 and a crystal . at volume fractions
 < 0.494 , the system 
 is in its fluid state ; at 0.492 <  < 0.545 , the system 
 will be a mixture of fluid and crystalline states ; and at  
 > 0.545 , the thermodynamically most stable phase is the crystal . 

the 
 transformation from fluid to crystal occurs through a first - order 
 phase transition . despite their simplicity , 

colloids made of polymers are commonly used for this purpose , the 
 most prominent example being poly(methyl methacrylate ) ( pmma ) spheres 
 coated with a layer of poly(12-hydroxystearic acid ) .
after the spheres 
 have been synthesized , they are dissolved in a mixture of cis - decaline and tetraline , which enables the use of a wide 
 range of powerful optical techniques to investigate nucleation .
the possibility of using these large hard spheres in nucleation 
 experiments has two major advantages : first , a particle size larger 
 than the wavelength used in microscopy experiments makes it possible 
 to track the particle trajectories in real space . in additional , nucleation 
 occurs in a matter of seconds , which allows experimentalists to follow 
 the complete nucleation process in detail . compared to other systems
, 
 it is therefore possible to observe the critical nucleus directly , 
 for example , by confocal microscopy ( see section 1.2 ) , which is of crucial importance for understanding 
 nucleation .
these qualities of hard - sphere systems make them ideal 
 candidates to advance our understanding of nucleation . as such ,
it 
 is not surprising that the freezing of hard spheres is better characterized 
 than any other nucleation scenario , and in fact , a number of excellent 
 reviews in this field already exist .
our aim here is thus not to give a detailed overview of the field 
 but to highlight some of the milestones and key discoveries and connect 
 them to other nucleation studies . to keep the discussion reasonably 
 brief , we limit the latter to neutral and perfectly spherical hard - sphere 
 systems . however , we note that a sizable amount of work has been devoted 
 to a diverse range of colloidal systems , such as nonspherical particles , charged particles , and mixtures of different colloidal particles , to name just a few .
readers interested in the state of the 
 art in about 2000 are referred to other reviews . in the
early 2000s , two major advances in the field were made , one 
 on the theoretical side and the other experimentally . in 2001 , auer 
 and frenkel computed absolute nucleation 
 rates of a hard - sphere system using kmc simulations .
they did so by 
 calculating pcrit , the probability of 
 forming a critical nucleus spontaneously , and , the 
 kinetic prefactor .
they found that the experimental and 
 theoretical nucleation rates disagreed by several orders of magnitude . 

this was surprising , because simulations did really well in describing 
 all sorts of properties of hard spheres before .
it was worrisome because 
 only very few sound approximations were made by auer and frenkel to 
 obtain their nucleation rates .
the authors suggestion , that the problem 
 lay in experiments or , more precisely , in the interpretation of experiments , 
 showed a possible way to resolve the discrepancy . in the same year , 
 gasser et al . conducted ground - breaking 
 experiments , imaging the nucleation of a colloidal suspension in real 
 space using confocal microscopy .
four snapshots of their system containing 
 approximately 4000 particles are shown in figure 6 .
this was a significant step , because nucleation 
 had previously been investigated indirectly , using the structure factor 
 obtained from light - scattering experiments , for example . in their 
 study
, they were able to directly measure the size of a critical nucleus 
 for the first time .
achieving sufficient temporal and spatial resolution 
 at the same time is possible thus far only for colloidal systems ( for 
 more details about experimental techniques , see section 1.2 ) .
they found that the nucleus 
 was rather aspherical with a rough surface ; both of these effects 
 are completely neglected in cnt .
note that aspherical nuclei also 
 appear in lj systems , for example ( see section 2.2.1 ) .
in addition , a random hexagonal close - packed 
 ( rhcp ) structure for the hard spheres was observed , in good agreement 
 with auer and frenkel .
this is interesting , 
 because slightly different systems such as soft spheres and lennard - jones 
 particles seem to favor body - centered - cubic ( bcc ) stacking .
however , 
 gasser et al.s study did not resolve the discrepancy between 
 experimental and simulated nucleation rates , as their results were 
 in agreement with earlier small - angle light - scattering experiments .
red ( large ) and blue ( small ) spheres show crystal- and 
 liquid - like particles , respectively .
the size of the observed volume 
 is 58 m by 55 m by 20 m , containing about 4000 
 particles .
after shear melting of the sample , snapshots were taken 
 after ( a ) 20 , ( b ) 43 , ( c ) 66 , and ( d ) 89 min .
much of the subsequent work focused 
 on trying to resolve this discrepancy 
 between experiments and simulation .
found experimental 
 evidence supporting a two - step crystallization process ( see section 1.1.2 ) in hard - sphere 
 systems .
other systems such as proteins and molecules in solution 
 ( see section 2.5 ) 
 were well - known at that time to crystallize through a more complex 
 mechanism than that assumed by cnt .
even for hard - sphere systems , 
 two - step nucleation processes were reported before 2006 ; the occurrence of this mechanism was attributed to details of the 
 polydispersity of the hard spheres , however .
the new insight provided 
 by schpe et al . in 2006 and 2007 was that the two - step nucleation 
 process is general , and as such , it does not depend on either polydispersity 
 or volume fraction . in 2010 , simulations performed by schilling et 
 al . supported these experimental findings . 
 using unbiased mc simulations , schilling et al .
not 
 even the simplest model system seemed to follow the traditional picture 
 assumed in cnt .
could this two - step mechanism explain why the computational 
 rates disagreed with experiments ?
at 
 first , it seems like a tempting explanation , because auer and frenkel had to introduce order parameters to calculate 
 absolute nucleation rates .
such a conclusion , however , automatically 
 presupposes a reaction pathway , which might not necessarily match 
 the nucleation pathway taken in experiments .
showed in the same year , however , that very 
 different computational approaches [ brute - force md , us , and ffs , which 
 we described earlier ( see section 1.3.2 ) ] led to the same nucleation rates , 
 all in agreement with auer and frenkel .
they therefore concluded that the discrepancy between simulations 
 and experiments did not lie in the computational approach employed 
 by auer and frenkel .
they offered two possible explanations , one being 
 that hydrodynamic effects , completely neglected in the simulations , 
 might play a role and the other being possible difficulties in interpreting 
 the experiments .
schilling et al . tried 
 to address one of the key issues when comparing experiments with simulations : 
 uncertainties and error estimation . whereas the determination of the 
 most characteristic quantity in hard - sphere systems , the volume fractions , 
 is straightforward for simulations ,
experimentalists are confronted 
 with a more difficult task in this case . the typical error in determining 

the volume fraction experimentally is about 0.004 , which translates 
 into an uncertainty in the nucleation rate of about an order of magnitude . 
 upon taking these considerations into account ,
the authors concluded 
 that the discrepancy can be explained by statistical errors and uncertainties .
does this mean that the past 10 years of research tried to explain 
 a discrepancy that is actually not there ?
rightfully pointed out that , whereas the rates 
 between experiments and simulations coincide at high volume fraction , 
 they still clearly disagree in the low - volume - fraction regime .
no 
 simple rescaling justified by statistical uncertainty could possibly 
 resolve that discrepancy . in their article
, they also addressed a 
 different issue . in a computational study in 2010 , kawasaki and tanaka obtained , by means of brownian dynamics , nucleation rates in good agreement with experiments , 
 contrary to the nucleation rates computed by auer and frenkel using 
 brute - force md .
it should be noted that 
 kawasaki and tanaka did not use a pure hard - sphere potential , but 
 used a weeks  chandler  andersen potential instead .
was 
 the approximation of a hard - sphere system , something that can never 
 be fully realized in experiments , the problem all the time ?
what filion 
 et al . showed is that different computational approaches ( brute - force 
 md , us , and fss ) all lead to the same nucleation rates , all of them 
 in disagreement with what kawasaki and tanaka found . through a detailed 
 evaluation of their approach and
that of kawasaki and tanaka , they 
 concluded that their rates are more reliable .
the discrepancy was 
 back on the table , where it still remains and is as large as ever . for a detailed comparison between experimental and computational 
 rates ,
the message we want to convey here is that 
 the disagreement between simulations and experiments in the simplest 
 system still persists today .
it is worth mentioning that this fundamental 
 disagreement between simulations and experiments is not unique to 
 colloids .
other systems such as water ( sections 2.4.1 and 2.4.2 ) and 
 molecules in solution ( section 2.5 ) also show discrepancies of several orders of magnitude 
 in nucleation rates .
this long - standing debate is of great relevance 
 to all investigations dealing with systems modeled using any flavor 
 of hard - sphere potential .
a notable example in this context is the 
 crystallization of proteins , which are usually treated as hard spheres . 
 despite basically neglecting most of the complexity of these systems , 
 this substantial approximation
has allowed for a number of computational 
 studies that , although outside the scope of this review , certainly contributed 
 to furthering our understanding of the self - assembly of biological 
 particles . beyond hard 
 spheres ,
the lennard - jones liquid is 
 a widely studied model system that does just that .
it can be seen 
 as the natural extension of the hard - sphere model , to which it becomes 
 equivalent when the strength of the attractive interactions goes to 
 zero .
lj liquids were first introduced in 1924 , and since then , they have been the subject of countless 
 computational studies .
lj potentials allow for exceedingly fast md 
 simulations , and a wide range of thermodynamic information is available 
 for them , such as the phase diagram and the interfacial free energy .
the stable structure of the lj system up to   is a face - centered - cubic ( fcc ) 
 crystal ; 
 slightly less stable in free energy is a hexagonal - close - packed ( hcp ) 
 structure , which , in turn , is significantly more stable then a third 
 body - centered - cubic ( bcc ) phase . with his 
 study of liquid argon in 1964 , rahman reported what is probably the 
 first lj md simulation .
his findings 
 showed good agreement with experimental data for the pair distribution 
 function and the self - diffusion coefficient , thus demonstrating that 
 lj potentials can properly describe noble elements in their liquid 
 form at ambient pressure .
to 
 the best of our knowledge , nucleation of lj liquids was investigated 
 for the first time in 1969 by de wette et al . and
in 1976 by mandell et al . for two - dimensional and three - dimensional systems , respectively .
early simulations investigating the condensation of lj vapors into a liquid already 
 indicated a substantial discrepancy with cnt rates .
it is worth noticing 
 that the order parameter for crystal - like particles presented by ten 
 wolde et al .
fostered a considerable 
 amount of later work devoted to improving the order parameters customarily 
 used to describe crystal nucleation from the liquid phase ( see , e.g. , 
 ref ( 254 ) ) .
in 2008 , 
 kalikmanov et al . compared cnt and cdft 
 ( see section 1.1 ) 
 simulations with condensation data for argon .
they found that cnt 
 spectacularly failed to reproduce experimental condensation rates , 
 underestimating them by up to 26 orders of magnitude .
this disagreement 
 triggered a number of computational studies aimed at clarifying the 
 assumption of the sphericity of the critical nucleus within the freezing 
 of lj liquids . by embedding pre - existing spherical clusters into supercooled 
 lj liquids , bai and
li found values of the 
 critical nucleus size in excellent agreement with cnt within a broad 
 range of temperatures .
however , these results have been disputed by 
 the umbrella sampling simulations of wang et al .
, for example , as well as the path - sampling investigation 
 of moroni et al . in both cases
in addition , moroni et al . pointed out that the critical nucleus 
 size is determined by a nontrivial interplay between the shape , the 
 size , and the degree of crystallinity of the cluster .
such a scenario 
 is clearly much more complex than the usual cnt picture , as it violates 
 the capillarity approximation ( see 1.1.1 ) . 

nonspherical nuclei were also observed by trudu et al . , who extended the conventional cnt formula to 
 account for ellipsoidal nuclei .
such a tweak gave much better estimations 
 of both the critical nucleus size and the nucleation barrier .
recall 
 that the shape of the critical nuclei can be observed experimentally 
 in very few cases ( see sections 1.2 and 2.1 ) .
however , at 
 very strong supercooling , things fell apart because of the emergence 
 of spinodal effects ( see section 1.1 ) .
note that cnt fails at strong supercooling even 
 without the occurrence of spinodal effects , as the time lag ( transient 
 time ) needed for structural relaxation into the steady - state regime 
 results in a time - dependent nucleation rate .
for instance , huitema et al . showed 
 that incorporating the time dependence into the kinetic prefactor 
 yields an improved estimate of nucleation rates .
in fact , by embedding 
 extensions to the original cnt framework , one can , in some cases , 
 recover a reasonable agreement between simulations and experiments 
 even at strong supercooling .
as an example , peng et al . also showed that including enthalpy - based terms 
 in the formulation of the temperature dependence of   substantially improves 
 the outcomes of 
 cnt .
another aspect that 
 has been thoroughly addressed within the crystal nucleation of lj 
 liquids is the structure of the crystalline clusters involved .
the 
 mean - field theory approach of klein and leyvraz suggests a decrease of the nucleus density as well as an 
 increase of the bcc character when moving toward the spinodal region . 

these findings were confirmed by the umbrella sampling approach of 
 ten wolde et al . , who reported a bcc shell surrounding fcc cores .
furthermore , wang 
 et al . showed that the distinction between 
 the crystalline clusters and the surrounding liquid phase falls off 
 as a function of t .
in fact , the free energy 
 barrier for nucleation , computed by means of umbrella sampling simulations 
 ( see section 1.3.2 ) , was found to be on the order of kbt at t = 52% .
in addition , 
 the nuclei undergo substantial structural changes toward nonsymmetric 
 shapes , a finding validated by the metadynamics simulations of trudu 
 et al .
the same authors investigated 
 the nucleation mechanism close to the critical temperature for spinodal 
 decomposition ,   ( see section 1.1.4 ) , where the free energy basin corresponding 
 to the liquid phase turned out to be ill - defined , that is , already 
 overlapping with the free energy basin of the crystal .
such a finding 
 suggested that , below , there is no free energy barrier for nucleation , 
 indicating that the liquid is unstable rather than metastable and 
 that the crystallization mechanism has changed from nucleation toward 
 the more collective process of spinodal decomposition ( see section 1.1.4 and figure 1 ) .
insights 
 into the interplay between nucleation and polymorphism have been provided 
 by the simulations of ten wolde et al .
, among others , suggesting that , 
 within the early stages of the nucleation process , the crystalline 
 clusters are bcc - like , later turning into fcc crystalline kernels 
 surrounded by bcc shells .
performed 
 a cdft study to determine the difference between the free energy barrier 
 for nucleation required for the creation of an fcc or bcc critical 
 nucleus .
in addition , the difficulty for nucleation of the three different 
 crystal orientations for fcc was ranked ( 100 ) > ( 110 ) > ( 111 ) . 

these 
 studies confirm the presence of a two - step mechanism ( see section 1.1.2 ) and the 
 validity of ostwald s step rule for the lj model . as we will see later ( e.g. , homogeneous ice nucleation , section 2.4.1 ) , nucleation 
 through metastable phases has also been observed for more complicated 
 liquids .
important contributions regarding polymorph control during 
 crystallization were made by desgranges and delhommelle , who investigated nucleation 
 under different thermodynamic conditions . by keeping the temperature 
 constant and altering the pressure , they were able to influence the 
 amount of bcc particles .
this reached up to a point where the nucleus 
 was almost purely bcc - like .
calculation of the bcc  liquid line 
 in the phase diagram showed that these nucleation events occurred 
 in the bcc existence domain .
additionally , the transformation from 
 fcc to hcp during crystal growth , well after the critical nucleus 
 size has been reached , was studied by changing the temperature at 
 constant pressure . as depicted in figure 7 , at t = 10% , a small 
 number of hcp atoms were observed surrounding the fcc core , whereas 
 at t = 22% , much larger hcp domains formed 
 within the crystallite , suggesting that the conversion from hcp to 
 fcc is hindered at higher temperatures . on a final note ,
we emphasize 
 that many findings related to polymorphism are often quite dependent 
 on the choice of the order parameters employed .
this issue is not 
 limited to lj systems , and it is especially important when dealing 
 with similarly dense liquid and crystalline phases ( e.g. , metallic 
 liquids ) , where order parameters usually struggle to properly distinguish 
 the different crystalline phases from the liquid . in particular , it remains to be seen whether the fractional 
 bcc , fcc , and hcp contents of the lj nuclei that we have discussed 
 will stand the test of the last generation of order parameters .
cross section 
 of postcritical crystalline clusters of 5000 lj particles 
 for t = ( a ) 10% and ( b ) 22% .
fcc- , hcp- , and 
 bcc - like particles are depicted in gray , yellow , and red , respectively . 
 at t
= 22% substantial hcp domains form within 
 the crystallite , whereas at t = 10% , hcp particles 
 can be found almost exclusively on the surface of the fcc core . reprinted 
 with permission from ref ( 247 ) .
heterogeneous 
 crystal nucleation has also been investigated for a variety of lj 
 systems .
used 
 umbrella - sampling simulations ( see section 1.3.2 ) to calculate the free energy barrier 
 for heterogeneous nucleation of an lj liquid on top of an ideal impurity , 
 represented by a single fcc ( 111 ) layer of lj particles . by explicitly 
 varying the lattice spacing of the substrate , asub
, they calculated   as a 
 function of asub  aequi , where aequi is the lattice 
 spacing of the equilibrium 
 crystalline phase .
aequi = 0 , whereas 
 for large values of asub  aequi , nucleation proceeds within the bulk of 
 the supercooled liquid phase .
these findings support the early argument 
 of the zero lattice mismatch introduced by turnbull and vonnegut to justify the striking effectiveness of agi 
 crystals in promoting heterogeneous ice nucleation .
in fact , in several 
 situations , one can define a disregistry or lattice mismatch  
 as7values of  close to or even 
 equal to 
 zero have often been celebrated as the main ingredient that makes 
 a crystalline impurity particularly effective in promoting heterogeneous 
 nucleation .
however , the universality of this concept has been severely 
 questioned in the past few decades , as we shall see in section 2.4.2 for heterogeneous 
 ice nucleation .
nonetheless , it seems that the argument regarding 
 zero lattice mismatch can hold for certain simple cases , as demonstrated 
 by mithen and sear , who studied heterogeneous 
 nucleation of lj liquids on the ( 111 ) and ( 100 ) faces of an fcc crystal 
 by means of ffs simulations ( see section 1.3.2 ) .
they reported a maximum in the heterogeneous 
 nucleation rate for a small , albeit nonzero , value of  ( see figure 8) . the difference 
 between their study and that of wang et al . is simply that many more values of  were taken into account 
 by mithen and sear , thus allowing the 
 maximum of   to be determined 
 more precisely . on a different 
 note ,
dellago et al . performed tis simulations 
 ( see section 1.3.2 ) to investigate heterogeneous crystal nucleation of lj supercooled 
 liquids on very small crystalline impurities .
they found that even 
 tiny crystalline clusters of just 10 lj particles can actively 
 promote nucleation and that the morphology of the substrate can play 
 a role as well .
specifically , whereas fcc - like clusters were rather 
 effective in enhancing nucleation rates , no substantial promotion 
 was observed for icosahedrally ordered seeds .
nucleation rates computed 
 with the ffs method for a rigid hexagonal 
 surface of lj atoms in contact with a lj liquid .
potentials 1 and 
 2 describe the interaction between substrate and liquid and differ 
 only slightly by the value of  they use .
these results show that the maximum 
 in the nucleation rate occurs at nonzero values of the lattice mismatch 
 . reprinted with permission from ref ( 271 ) .
mc simulations performed by page 
 and sear have demonstrated that confinement 
 effects can be of great relevance 
 as well .
they computed heterogeneous nucleation rates for a lj liquid 
 walled inbetween two flat crystalline planes characterized by a certain 
 angle sub .
a maximum of   was 
 found for a specific value of sub , boosting the 
 nucleation rate by several orders of magnitude 
 with respect to the promoting effect of a flat crystalline surface . 

in addition , different values of sub led to the 
 formation of different crystalline polymorphs .
recently probed 
 the influence of structured and structureless lj potential walls or 
 nucleation rates .
both types of wall were found to increase the temperature 
 at which nucleation occurs .
we shall see in section 2.4.2 that the interplay between the morphology 
 of the substrate and the strength of the liquid 
md simulations 
 of lj liquids are computationally cheap , making them the perfect candidates 
 to examine how finite - size effects impact crystal nucleation .
the 
 seminal work of honeycutt and andersen took into account up to 1300 lj particles at , which turned out to be too 
 few particles 
 to completely rule out the effects of periodic boundary conditions . 

in fact , the authors suggested that extra care had to be taken because 
 of the diffuseness of the interface between the supercooled liquid 
 phase and the crystalline nucleus , which can induce an artificial 
 long - range order in the system , leading to a nonphysically high nucleation 
 rate .
these findings are particularly relevant , as the critical nucleus 
 size at this t value is on the order of just 
 a few tens of particles , representing a tiny fraction of the whole 
 system .
only a few years later , swope and andersen investigated the same effects by taking into account up 
 to 10 lj particles subjected to the same strong supercooling 
 as probed by honeycutt and andersen . according to their large - scale md simulations ,
. this outcome must be 
 carefully pondered , as currently , the vast majority of simulations 
 dealing with crystallization of realistic systems can not afford to 
 take into account system sizes 3 orders of magnitude larger than the 
 size of the critical nucleus .
examined the nucleation of an lj liquid in a wide range of temperatures 
 ( 70140 k ) .
although nonphysical instantaneous crystallization 
 was observed for systems on the order of 500 particles , simulation 
 boxes containing about 10000 particles seemed to be free from finite - size 
 effects .
recently described a novel class of finite - size 
 effects unrelated to periodic boundary conditions .
in fact , they showed 
 that the equilibrium density of critical nuclei , , can effectively 
 influence the absolute value of nucleation rates .
specifically , at 
 very strong supercooling , the critical nuclei will on average form 
 very shortly after the transient time , whereas at mild t , the stochastic nature of nucleation will lead to a consistent 
 scatter of the nucleation times .
in other words , in the latter scenario , 
 either exceedingly large systems must be taken into account , or a 
 sizable number of independent simulations must be performed to deal 
 with the long tails of the distribution of nucleation times .
examples 
 include the sutton  chen potentials for several metals and the tosi  fumi potential for molten salts such as nacl .
terms accounting 
 for the directionality of covalent bonds have been included , for example , 
 in the stillinger 
weber potential for si ; the bond order potentials of tersoff for si , gaas , and ge ; and the reactive potential of brenner for carbon - based systems .
another class of 
 interatomic potential is based on the concept of local electronic 
 density and includes , for instance , the finnis  sinclair potentials for metallic systems , the whole family of the embedded - atom - method 
 ( eam ) potentials , and the glue potential for au and al .
many of these potentials are still incredibly 
 cheap in terms of computer time , thus allowing for large - scale , unbiased 
 md simulations .
recently , massively parallel md runs succeeded in 
 nucleating supercooled liquid al and 
 fe using an eam potential and a finnis
as up to 10 atoms were taken 
 into account , actual grain boundaries were observed , providing unprecedented 
 insight in to the crystal growth process .
the nucleation of bcc fe 
 crystallites and the evolution of the resulting grain boundaries at 
 different temperatures can be appreciated in figure 9a .
the sizable dimension of the simulation 
 boxes ( 50 nm ) allowed nucleation events to be observed within 
 hundreds of picoseconds , and grain coarsening ( i.e. , the process by 
 which small crystallites end up incorporated into larger ones ) is 
 also clearly visible .
mere visual inspection of the nucleation trajectories 
 depicted in figure 9a suggests different nucleation regimes as a function of temperature . 

in fact , the same authors calculated a temperature profile for the 
 nucleation rate , shown in figure 9b , that demonstrates the emergence of a maximum of   value characteristic 
 of diffusion - limited 
 nucleation ( see section 1.1.1 ) .
crystal nucleation of supercooled fe by means of large - scale 
 md 
 simulations .
yellow circles highlight 
 small crystalline grains doomed to be incorporated into the larger 
 ones later on because of grain coarsening .
( b ) nucleation rate as 
 a function of temperature . reprinted with permission from ref ( 290 ) .
a field that has greatly 
 benefited from md simulations is the crystallization 
 of metal clusters , as nicely reviewed by aguado and jarrold .
for instance , it is possible to probe the interplay 
 between the size of the clusters and the cooling rate upon crystal 
 nucleation and growth . in this context , shibuta reported three different outcomes for supercooled liquid 
 mo nanoparticles modeled by means of a finnis  sinclair potential , 
 namely , the formation of a bcc single crystal , a glassy state , or 
 a polycrystalline phase . in some cases , nucleation rates obtained 
 from simulations were consistent with cnt , as in the case of ni nanodroplets , for which nucleation events were again observed 
 by means of brute - force md simulations using the sutton  chen 
 potential .
the influence of the redox potential on the nucleation 
 process has also been investigated .
milek and zahn employed an enhanced flavor of the eam potential to study 
 the nucleation of ag nanoparticles from solution .
they established 
 that the outcome of nucleation events is strongly influenced by the 
 strength of the redox potential , able to foster either a rather regular 
 fcc phase or a multitwinned polycrystal .
similar to what was done 
 for lj liquids , the effects of confinement were assessed for au nanodomains 
 modeled using the glue potential by pan and shou . according to their findings ,
l and chen instead investigated surface - layering - induced crystallization of 
 ni  si nanodroplets using a modified eam potential .
it seems 
 that , for this particular system , atoms proximal to the free surface 
 of the droplet assume a crystalline - like ordering on very short time 
 scales , thus triggering crystallization in the inner regions of the 
 system .
no such effect has been reported in the case of surface - induced 
 crystallization in supercooled tetrahedral liquids such as si and 
 ge , as investigated by li et al . through 
 ffs simulations employing both tersoff and stillinger 
the presence of the free surface facilitates crystal nucleation 
 for this class of systems as well , but surface layering was not observed . 
 instead
, the authors claimed that the surface reduces the free energy 
 barrier for nucleation as it introduces a pressure - dependent term 
 in the volume free energy change expected for the formation of the 
 crystalline clusters .
the situation is quite different for surface - induced 
 ice nucleation , at least according to the coarse - grained mw model 
 of molinero and moore .
in fact , haji - akbari 
 et al . recently investigated ice nucleation 
 in free - standing films of supercooled mw water using both ffs and 
 us , finding that , in these systems , crystallization is inhibited in 
 the proximity of the vapor
very recently , 
 gianetti et al . extended the investigation 
 of haji - akbari et al . to the crystallization 
 of a whole family of modified stillinger 
weber liquids with 
 different degrees of tetrahedrality  , locating a crossover 
 from surface - enhanced to bulk - dominated crystallization in free - standing 
 films as a function of . another seminal study by li et al .
, again using ffs , focused on homogeneous ice 
 nucleation within supercooled mw water nanodroplets , where nucleation 
 rates turned out to be strongly size dependent and in general consistently 
 smaller ( by several orders of magnitude ) than the bulk case .
ffs was 
 also applied by li et al . to examine 
 homogeneous nucleation of supercooled si .
ffs has also been successful 
 in predicting homogeneous crystal nucleation rates in molten nacl , 
 modeled using a tosi  fumi potential by valeriani et al .
large discrepancies between their results and 
 experimental nucleation rates can be appreciated when cnt is used 
 to extrapolate the calculations to the milder supercooling probed 
 by the actual measurements .
given that the authors obtained consistent 
 results using two different enhanced sampling methods , this study 
 hints again at the many pitfalls of cnt .
a unique 
 example of a class of materials for which nucleation can be effectively 
 addressed by brute - force md simulations is given by so - called phase - change 
 materials .
these systems are of great technological 
 interest as they are widely employed in optical memories ( e.g. , dvd - rw ) 
 and in a promising class of nonvolatile memories known as phase - change 
 memory , based on the fast and reversible 
 transition from the amorphous to the crystalline phase .
although crystal 
 nucleation in amorphous systems , especially metallic and covalent 
 glasses , is beyond the scope of this review , we refer the reader to 
 the excellent work of kelton and greer for a detailed introduction . here
we just note that in phase change 
 memories the amorphous phase is often heated above the glass transition 
 temperature , so that crystal nucleation occurs within the supercooled 
 liquid phase .
phase - change materials used in optical and electronic 
 devices are typically tellurium - based chalcogenide alloys ( see ref ( 305 ) ) . the family of the pseudobinary 
 compounds ( gete)x(sb2te3)y represents a prototypical system . 

although both the structure and dynamics of these systems are far 
 from trivial , nucleation from the melt takes place on the nanosecond 
 time scale for a wide range of supercooling .
thus , with phase - change materials , we have a great opportunity to 
 investigate nucleation in a complex system by means of brute - force 
 md simulations .
we note that the crystallization of these systems 
 has been extensively characterized by different experimental techniques 
 [ particularly tem and afm ( see section 1.2 ) ; the crystallization kinetics has also 
 been recently investigated by means of ultrafast - heating calorimetry and ultrafast x - ray imagining ] , but because of the exceptionally high nucleation rates , 
 it is difficult to extract information about the early stages of the 
 nucleation process .
unfortunately , phase - change materials require ab 
 initio methods or sophisticated interatomic potentials with first - principles 
 accuracy .
in fact , several attempts have been made to study nucleation 
 in phase - change materials by ab initio md in very small systems .
although these studies provided useful insights into the nucleation 
 mechanism , severe finite - size effects prevented the full characterization 
 of the crystallization process .
the limited length and time scales 
 typical of first - principles calculations were recently outstripped 
 in the case of the prototypical phase - change material gete by the 
 capabilities of a neural network interatomic potential .
such potentials allow for a computational speedup 
 of several orders of magnitude compared to conventional ab initio 
 methods while retaining an accuracy close to that of the latter .
although nucleation rates have not yet been 
 calculated , detailed investigations of homogeneous and heterogeneous 
 nucleation have already been reported .
for instance , 
 as shown in figure 10 , a single - crystalline nucleus formed in a 4000-atom model of supercooled 
 liquid gete in the 625675 k temperature regime within a few 
 hundred picoseconds . on the same timescale ,
several nuclei appeared 
 below 600 k , suggesting that the free energy barrier for nucleation 
 is vanishingly small for this class of materials just above the glass 
 transition temperature .
this is because of the fragility of the supercooled liquid , which displays a 
 substantial atomic mobility even at large supercoolings .
thus , in this particular case , the kinetic 
 prefactor   ( see eq 5 ) is not hindered that much by the strong 
 supercooling , 
 whereas the free energy difference between the liquid and the crystal   ( see eq 1 ) skyrockets as expected , leading to the exceedingly 
 high nucleation rates characteristic of these materials . fast crystallization 
 of supercooled gete by means of md simulations 
 with neural - network - derived potentials .
the number of crystalline 
 nuclei larger than 29 atoms at different temperatures in the supercooled 
 liquid phase is reported as a function of time ( notice the exceedingly 
 small time scale at strong supercooling ) .
two snapshots at the highest 
 and lowest temperatures showing only the crystalline atoms are also 
 reported .
at high temperature , a single nucleus is present , whereas 
 several nuclei ( each one depicted in a different color ) appear at 
 low temperature .
copyright 2013 american 
 chemical society . in conclusion , whereas 
 md simulations have by no means exhausted 
 the field of crystal nucleation of atomic liquids , they have certainly 
 provided insight into a number of interesting systems and paved the 
 way for the study of more complex systems , as we shall see in the 
 following sections .
ice 
 nucleation impacts many different areas , ranging from aviation to biological cells and earth s 
 climate .
it is therefore not surprising 
 that a considerable body of work has been carried out to understand 
 this fundamental process .
we can not cover it all here ; instead , we 
 give a general overview of the field , starting with a discussion of 
 nucleation rates .
this allows us to directly compare experiments and 
 simulations and to identify strengths and weaknesses of different 
 approaches .
experimental nucleation 
 rates have been measured over a broad range of temperatures , most 
 often with micrometer - sized water droplets so as to avoid heterogeneous 
 nucleation . in figure 11 , we bring together nucleation rates obtained from various experiments , 
 along with computed nucleation rates . compilation of homogeneous 
 nucleation rates for water , obtained 
 by experiments and simulations .
the x axis shows 
 the supercooling with respect to the melting point of different water 
 models or 273.15 k for experiment .
the y axis shows 
 the logarithm of the nucleation rate in m  s. rates obtained with computational approaches are 
 shown as solid symbols ; experimental rates are shown as crossed symbols . 
 for each computational study , the computational approach and the water 
 force field used
is not included in this graph , because their 
 study was conducted at a small supercooling ( 20 k ) , which resulted 
 in a very low estimated nucleation rate far outside this plot ( it 
 would correspond to 83 on the y axis ) .
taborek performed measurements with different setups , 
 namely , using sorbitan tristearate ( sts ) and sorbitan trioleate ( sto ) 
 as surfactants .
data for the graph were taken from refs ( 124 , 299 , 301 , and 322335 ) . accessing nucleation rates from 
 md simulations became feasible 
 only in the past few years as a result of advances in force fields 
 ( such as the coarse - grained mw potential ) 
 and enhanced sampling techniques described earlier ( see section 1.3.2 ) .
these 
 methods have therefore been widely used for studies of not only homogeneous 
 but also heterogeneous nucleation ( see section 2.4.2 ) . from the comparison of experimental 
 and computational nucleation rates reported in figure 11 ,
first , 
 nucleation rates vary hugely with supercooling , by a factor of more 
 than 10 .
second , nucleation rates differ substantially 
 ( approximately 10 orders of magnitude ) between simulations ( solid 
 symbols ) and experiments ( crossed symbols ) at relatively small supercoolings 
 ( 3050 k ) . at larger supercoolings , the agreement appears 
 to be slightly better , even though very few simulations have been 
 reported at very strong supercooling .
the third striking feature is 
 that , whereas the experimental results agree well with each other 
 ( within 12 orders of magnitude ) , the computational rates differ 
 from each other by a factor of approximately 10 . what is the cause of disagreement between different computational 
 approaches ?
part of the reason is certainly that different water models 
 lead to different rates ; see , for example , espinosa et al . yet , even if the same water model is employed , 
 the rates do not agree with each other very well .
an early 
 study by moore and molinero succeeded 
 in calculating the avrami exponent for the crystallization 
 kinetics of ice from brute - force md simulations at very strong supercooling , 
 obtaining results remarkably similar to experiment . however , mw nucleation rates turned out to be far less encouraging . 

in fact , li et al . and reinhardt and 
 doye both performed simulations using 
 the mw model , obtaining nucleation rates that differed by about 5 
 orders of magnitude .
the only major difference was the enhanced sampling 
 technique employed , ffs by li et al . and us by reinhardt and doye . 

the statistical uncertainties of the two approaches ( 12 orders 
 of magnitude ) are much smaller than the 5-orders - of - magnitude discrepancy 
 between the two studies .
it was also shown that the two methods agree 
 very well with each other for colloids , for example ( see section 2.1 ) .
the use of different computational approaches therefore 
 also seems to be unlikely as the source of the disagreement .
what 
 the cause is remains elusive . because we can not cover all of 
 the work shown in figure 11 in detail here , we now discuss 
 just two studies
first , that of sanz et al . , which agrees best with the experimental rates .
the authors used 
 the tip4p/2005 and tip4p / ice water models in combination with seeded 
 md simulations ( see section 1.3.2 . for more details ,
seeding involves considerably 
 more assumptions than , for example , us or ffs . in particular
, the 
 approach assumes a cnt - like free energy profile , although it does 
 not usually employ the macroscopic interfacial free energy .
furthermore , 
 the temperature dependence of key quantities such as   and   ( see section 1.1.1 ) is approximated .
nevertheless , the 
 agreement between their nucleation rates and experiment seemingly 
 outperformed other approaches . in a more recent work , espinosa et 
 al .
however , it should be noted that the 
 good agreement between the nucleation rates reported in refs ( 51 and 326 ) and the experimental data could originate from error cancellation . 

in fact , whereas the rather conservative definition of crystalline 
 nucleus adopted in these works will lead to small nucleation barriers 
 ( and thus to higher nucleation rates ) , the tip4p family of water models 
 is characterized by small thermodynamic driving forces to nucleation , which , in turn , results in smaller nucleation 
 rates .
the second work we briefly discuss here is the very recent 
 study 
 ( 2015 ) of haji - akbari and debenedetti .
the authors directly calculated the nucleation rate at 230 k of 
 an all - atom model of water ( tip4p / ice ) using a novel ffs sampling 
 approach .
this was a tour de force , 
 but strikingly , their rates differed from experiment by about 11 orders 
 of magnitude .
the authors noted that this might be as close as one 
 can actually get to experiment with current classical water models . 

this is because of the extreme sensitivity of nucleation rates to 
 thermodynamic properties such as   and , which , according to cnt , enter exponentially 
 ( section 1.1.1 ) 
 in the definition of .
for instance , 
 an uncertainty of only 67% 
 for   at 235 k leads to an 
 error of about 9 orders 
 of magnitude in .
experimental 
 estimates for  range from 25 to 35 mn / m ; computational estimates range from about 20 to 35 mn / m . as 
 another example , haji - akbari and debenedetti explicitly quantified the extent to which the tip4p / ice model underestimates 
 the free energy difference   between the crystalline and liquid phases 
 and found that the mismatch between   and   alone 
 leads to an overestimation of the 
 free energy barrier for nucleation of about 60% , which translates 
 into nucleation rates up to 9 orders of magnitude larger .
in fact , 
 taking into account such a discrepancy brings the results of haji - akbari 
 and debenedetti within the confidence interval of the experimental 
 data .
thus , it is clear that we simply do not know some key quantities 
 accurately enough to expect perfect agreement between simulations 
 and experiments .
in addition to issues of modeling water / ice 
 accurately , finite - size 
 effects can be expected to also play a role [ as they do with lennard - jones 
 systems ( section 2.2 ) and molecules in solution ( section 2.5 ) ] . only recently was this issue addressed 
 explicitly for ice nucleation by english and tse in unbiased simulations with the mw model .
they were able 
 to simulate systems containing nearly 10 million water molecules on 
 a microsecond time scale and found that larger systems favor the formation 
 of crystallization precursors compared to smaller ones .
interestingly , 
 lifetimes of the precursors were found to be less sensitive to system 
 size .
a quantitative understanding of finite - size effects on nucleation 
 rates remains elusive nevertheless . in summary , it can be said 
 that , in terms of accurate nucleation 
 rates , experiments are still clearly superior to simulations .
however , 
 the advantage of simulations is that the nucleation mechanism can 
 also be obtained , which , at present , is not possible with experiments , 
 although femtosecond x - ray laser spectroscopy might be able to partially 
 overcome this limitation in the near future . in 2002
were the first to report a nucleation event 
 in an unbiased simulation based on an all - atom model of water ( tip4p ) . 

their landmark work opened the doors to the study of ice nucleation 
 at an atomistic level .
they found that nucleation took place once 
 a sufficient number of long - lived hydrogen bonds were formed with 
 a nucleus of ice .
recent evidence suggests that , most likely , their 
 nucleation trajectory was driven by finite - size effects .
together with the simulations of vrbka and jungwirth , also affected by severe finite - size effects , the work of matsumoto et al . remains , to date , 
 the only seemingly unbiased md simulation observing homogeneous ice 
 nucleation with an all - atom force field .
what really enabled 
 the community to investigate ice formation at a molecular level was 
 the development of the coarse - grained mw potential for water in the early 2010s . using unbiased md simulations 
 based on the mw force field , moore and molinero in 2011
provided evidence that , in the supercooled 
 regime around the homogenous nucleation temperature , th , the fraction of 4-fold - coordinated water molecules 
 increases sharply prior to a nucleation event . in a separate work , 
 the same authors suggested that , at 
 very strong supercooling , the critical nucleus is mostly made of cubic 
 ice , which subsequently evolves into a mixture of stacking - disordered 
 cubic and hexagonal ice layers . in the same year , li et al . identified another structural motif that might 
 play a role in ice nucleation .
they consistently observed a topological 
 defect structure in growing ice nuclei in their ffs simulations based 
 on the mw representation of water .
this defect , depicted in figure 12a , can be described 
 as a twin boundary with 5-fold symmetry , and it has also been observed in nucleation simulations of tetrahedral liquids 
 simulated with the stillinger  weber potential , on which the 
 mw coarse - grained model was built . ( a ) formation of a topological defect 
 with 5-fold symmetry during 
 homogeneous ice nucleation .
the snapshots ( i  iv ) show the time 
 evolution of the defect structure , indicated by black dashed lines . 

reprinted with permission from ref ( 322 ) ( copyright 2011 royal 
 society of chemistry ) , in which li et al
( b ) nucleation 
 of an ice cluster forming homogeneously from i0-rich precritical 
 nuclei .
water molecules belonging to ic , ih , 
 a clathrate - like phase , and i0 are depicted in yellow , 
 green , orange , and magenta , respectively .
( iii ) the crystalline cluster 
 evolves in a postcritical nucleus , formed by an ic - rich 
 core surrounded by an i0-rich shell .
( iv ) the same postcritical 
 nucleus as depicted in iii , but only particles with 12 or more connections 
 ( among ice - like particles ) are shown .
the color map refers to the 
 order parameter q12 specified in ref ( 325 ) , from which this image 
 was reprinted with permission ( copyright 2014 nature publishing group ) . 
 the unbiased md simulations on which the analysis is based feature 
 10000 mw molecules . in 2012 ,
another significant 
 leap in understanding the nucleation 
 mechanism of ice from a structural point of view was made by combining 
 experimental and computational techniques .
specifically , malkin et al . showed that ice forming homogeneously 
 is stacking - disordered ( the corresponding ice structure was called 
 isd ) , meaning that it is made out of cubic and hexagonal 
 ice layers stacked in a random fashion . in 2014 , two studies 
 substantiated the potential relevance of precursor 
 structures prior to ice formation .
liquid transition in supercooled water in a molecular 
 model of water ( st2 ) . in their study
, 
 the authors sampled the energy landscape of supercooled water and 
 found two metastable liquid basins corresponding to low - density ( ldl ) 
 and high - density ( hdl ) water .
the appealing idea behind the transition 
 from hdl to ldl prior to ice nucleation is that ldl is structurally 
 closer to ice than hdl .
note that the existence of two metastable 
 liquid basins was not a general finding : the mw model does not have 
 a basin for ldl , for example .
indeed , 
 the presence of this liquid  liquid phase transition is a highly 
 debated issue .
another conceptually similar 
 idea is ice formation through ice 
 0 ( i0 ) , proposed by russo et al . instead of a liquid 
liquid phase transition that transforms 
 water into another liquid state prior to nucleation , the authors proposed 
 a new ice polymorph ( i0 ) to bridge the gap between supercooled 
 water and ice .
i0 is a metastable ice polymorph and is 
 structurally similar to the supercooled liquid .
it has a low interfacial 
 energy with both liquid water and ice ic / ih . 

russo et al . therefore proposed i0 to bridge liquid water 
 to crystalline ic / ih .
indeed , the authors found 
 i0 at the surface of growing ice nuclei in md simulations ; 
 we show part of a nucleation trajectory in figure 12b .
furthermore , they showed that the shape 
 of the nucleation barrier is much better described by a core  shell - like 
 model ( ic / ih core surrounded by i0 ) compared to the classical nucleation model .
this is important , 
 because it suggests that models that are based solely on cnt assumptions 
 might not be appropriate for describing homogeneous ice nucleation .
however , the emergence of i0 has not yet been reported 
 by any other nucleation study , including the recent work of haji - akbari 
 and debenedetti that we previously mentioned 
 in the context of nucleation rates . in that work
, the authors performed a topological analysis of the nuclei , 
 validated by the substantial statistics provided by the ffs simulations . 

as depicted in figure 13 , the majority of nuclei that reached the critical nucleus size contained 
 a large amount of double - diamond cages ( ddcs , the building blocks 
 of ic ) , whereas nuclei rich in hexagonal cages ( hcs , the 
 building blocks of ih ) had a very low probability to overcome 
 the free energy barrier for nucleation .
in addition , even postcritical 
 nuclei had a high content of ddcs , whereas hcs did not show any preference 
 to appear within the core of the postcritical nuclei .
this evidence 
 is consistent with the findings reported in ref ( 323 ) and in contrast with 
 the widely invoked scenario in which a kernel of thermodynamically 
 stable polymorph ( in this case , ih ) is surrounded by a 
 shell of a less stable crystalline structure ( in this case , ic ) .
( left ) a typical double - diamond cage ( ddc , blue ) and a 
 hexagonal 
 cage ( hc , red ) , the building blocks of ic and ih , respectively .
( right ) temporal evolution of an ice nucleus from 
 ( i , ii ) the early stages of nucleations up to ( v , vi ) postcritical dimensions , 
 as observed in the ffs simulations of haji - akbari and debenedetti .
about 4000 water molecules , modeled with the 
 tip4p / ice potential , were considered 
 in the npt ensemble at t  
 40 k. one can clearly notice the abundance of ddcs throughout the 
 whole temporal evolution .
in contrast , hc - rich nuclei have only a 
 marginal probability to cross the nucleation barrier ( see text ) . reprinted 
 with permission from ref ( 327 ) .
copyright 2015 national academy of sciences . in the past few years , the understanding of homogeneous 
 ice nucleation 
 has improved dramatically .
we now have a good understanding of the 
 structure of ice that forms through homogeneous nucleation , stacking - disordered 
 ice .
furthermore , there is very good agreement ( within 2 orders of 
 magnitude ) between experimental nucleation rates in a certain temperature 
 range .
computational methods face the problem of being very sensitive 
 to some key thermodynamic properties ; the nucleation rates they predict 
 are therefore less accurate . on the other hand ,
they allow us to study 
 conditions that are very challenging to probe experimentally , and 
 they also provide insight into the molecular mechanisms involved in 
 the crystallization process .
as 
 mentioned in the previous section , homogeneous ice nucleation becomes 
 extremely slow at moderate supercooling .
we do not , for example , have to wait for 
 temperatures to reach 30 c before we have to use a deicer 
 on our car windows .
in fact , the formation of ice in nature occurs 
 almost exclusively heterogeneously , thanks to the presence of foreign 
 particles .
these ice - nucleating agents facilitate the formation of 
 ice by lowering the free energy barrier for nucleation ( see figure 1 ) . indeed , the work 
 of sanz et al . , in which homogeneous 
 ice nucleation was studied using seeded md ( see section 1.3.1 ) , found rates so low 
 at temperatures above t = 20 k that they concluded 
 that all ice nucleation above this temperature must occur heterogeneously . 

homogeneous nucleation is still of great importance in atmospheric 
 processes and climate modeling , as under certain conditions , both 
 heterogeneous and homogeneous nucleation are feasible routes toward 
 the formation of ice in clouds , as reported in ref ( 354 ) , for example . in 
 addition to the challenges
( both computational and experimental ) faced 
 when investigating homogeneous ice nucleation , one also has to consider 
 the structure of the water  surface interface and how this impacts 
 the nucleation rate .
generally , the experimental data for the rates 
 and characterization of the interfacial structure come from two different 
 communities : climate scientists have provided much information on 
 how various particles , often dust particles or biological matter such 
 as pollen , affect ice nucleation ( as depicted in figure 14 ) , whereas surface scientists 
 have invested a great deal of effort in trying to understand , at the 
 molecular level , how water interacts with and assembles itself at 
 surfaces ( see , e.g. , ref ( 355 ) ) .
this means that there is a huge gap in our understanding , 
 as the surfaces of the particles used to obtain rates are often not 
 characterized , whereas surface science experiments are generally carried 
 out at pristine , often metallic , surfaces under ultrahigh - vacuum conditions . 

we will see in this section that computational studies have gone some 
 way toward bridging this gap , although there is still much work to 
 be done should we wish to quantitatively predict a material s 
 ice - nucleating efficacy .
potential immersion - mode ice nucleus concentrations , nice , a measure of the efficiency of a given 
 substance 
 to boost heterogeneous ice nucleation , as a function of temperature 
 for a range of atmospheric aerosol species .
note the wide range of 
 nucleating capability for materials as diverse as soot and bacterial 
 fragments over a very broad range of temperatures .
copyright 2012 royal society of chemistry . from a computational 
 perspective , it is the surface science experiments that lend themselves 
 most readily to modeling .
in fact , even relatively expensive computational 
 methods such as density functional theory ( dft ) , which have not featured 
 much in this article , have proven indispensable in furthering our 
 understanding of how water behaves at surfaces , especially when used 
 in conjunction with experiments ( see , e.g. , refs ( 355 and 357 ) for an overview ) . as such , early 
 computational studies focused on understanding how the surface affected 
 the first few layers of water , especially with respect to the concept 
 of lattice mismatch ( see section 2.2 ) , where a surface that has a structure commensurate 
 with ice acts as a template for the crystal .
nutt and co - workers investigated the adsorption structures of water at a model hexagonal 
 surface and at baf2(111 ) using interaction potentials derived 
 from ab initio calculations .
although the surfaces under investigation 
 had structures that matched the basal face of ice well , they found 
 disordered structures of water to be more favorable than ice - like 
 overlayers . using dft , hu and michaelides investigated the adsorption 
 of water on the ( 001 ) face of the clay mineral kaolinite , a known ice - nucleating agent in the atmosphere .
the ( 001 ) surface 
 of kaolinite exposes a pseudohexagonal arrangement of oh groups that 
 were proposed to be the cause of its good ice - nucleating ability .
although they found that a stable ice - like 
 layer could form at the surface , the amphoteric nature of the kaolinite 
 surface , depicted in figure 15 , meant that all of water molecules could participate in four 
 hydrogen bonds , making further growth on top of the ice - like layer 
 unfavorable .
investigated adsorption 
 of water on kaolinite using the clayff and spc / e potentials and grand canonical monte carlo ( gcmc ) .
although some hexagonal 
 patches of water were seen in the contact layer , the overall structure 
 was mostly disordered , and the hexagonal structures that did form 
 were strained relative to those found in ice . also using gcmc , cox 
 et al .
they found that , for atomically flat surfaces , 
 a nominally zero lattice mismatch produced disordered contact layers 
 comprising smaller - sized rings ( i.e. , pentagons and squares ) and observed 
 hexagonal ice - like layers only for surfaces with larger lattice constants .
the 
 amphoteric nature of kaolinite is important to its ice - nucleating 
 ability .
( left ) ice - like contact layers at the kaolinite surface , 
 with the ( a , c ) basal and ( b , d ) prism faces of ice adsorbed on kaolinite , 
 as viewed from the ( a , b ) side and ( c , d ) top .
( right ) adsorption energy 
 of ice on kaolinite when bound through either its basal face ( red 
 data ) or its prism face ( blue data ) for varying numbers of layers 
 of ice .
( open and solid symbols indicate data obtained with a classical 
 force field and with dft , respectively . )
when only the contact layer 
 is present , the basal face structure is more stable than the prism 
 face structure , but as soon as more layers are present , the prism 
 face structure becomes more stable .
this can be understood by the 
 ability of the prism face to donate hydrogen bonds to the surface , 
 and to the water molecules above , through the  dangling  

hydrogen bonds seen in b and d. reprinted with permission from ref ( 371 ) .
2010 , the above types of study were the state 
 of the 
 art for simulations of heterogeneous ice nucleation .
although they 
 provided evidence that properties such as lattice match alone are 
 insufficient to explain a material s ice - nucleating ability , 
 because ice nucleation itself was not directly observed , only inferences 
 could be drawn about how certain properties might actually affect 
 ice nucleation .
yan and patey investigated 
 the effects of electric fields on ice nucleation using brute - force 
 molecular dynamics ( the electric fields were externally applied and 
 were not due to an explicit surface ) .
they found that the electric 
 field needed to act over only a small range ( e.g. , 10  ) and 
 that the ice that formed near the  surface 
was ferroelectric 
 cubic ice , although the rest of the ice that formed above was not . 

cox et al . performed simulations of heterogeneous ice nucleation in which the atomistic natures of both the 
 water and the surface were simulated explicitly , using tip4p/2005 water and clayff to describe kaolinite . despite the fact that 
 the simulations were affected by finite - size effects
, the simulations 
 revealed that the amphoteric nature of the kaolinite is important to ice nucleation . in the liquid
, a strongly bound 
 contact layer was observed , and for ice nucleation to occur , significant 
 rearrangement in the above water layers was required .
it was found 
 that ice nucleated with its prism face , rather than its basal face , 
 bound to the kaolinite , which was unexpected based on the theory that 
 the pseudohexagonal arrangement of oh groups at the surface was responsible 
 for templating the basal face .
cox et al . rationalized the formation 
 of the prism of ice at the kaolinite as being due to its ability to 
 donate hydrogen bonds both to the surface and to the water molecules 
 above ( see figure 15 ) , whereas the basal face maximizes hydrogen bonding to the surface 
 only .
more recent simulation studies , employing 
 rigid and constrained models of kaolinite , have also found the amphoteric 
 nature of the kaolinite surface to be important .
however , the heterogeneous nucleation mechanism of water 
 on clays is yet to be validated by unconstrained simulations unaffected 
 by substantial finite - size effects . as in the case 
 of simulations of homogeneous ice nucleation
, the use of the coarse - grained 
 mw potential has seen the emergence 
 of computational studies that actually quantify the ice - nucleating 
 efficiencies of different surfaces .
recently , lupi et al . investigated ice nucleation at carbonaceous 
 surfaces ( both smooth graphitic and rough amorphous surfaces ) using 
 cooling ramps to measure nonequilibrium freezing temperatures tf
 tf  tfhomo , where tf is the temperature at which 
 ice nucleates in the presence of a surface and tfhomo = 201  
 1 k is the temperature at which homogeneous ice nucleation occurs . 

it was found that the rough amorphous surface did not enhance ice 
 nucleation ( tf = 0 k ) , whereas 
 the smooth graphitic surfaces promoted ice nucleation ( tf = 1113 k ) .
this was attributed to 
 the fact that the smooth graphitic surface induced a layering in the density profile of water above the surface , whereas the rough 
 amorphous surface did not .
lupi and molinero quantified the extent 
 of layering as where (z ) is the density 
 of water at a height z above the surface and 0  (zbulk ) , where zbulk is a height where the density profile is 
 bulk - like . in a subsequent work using the same methodology , lupi and 
 molinero investigated how the hydrophilicity 
 of graphitic surfaces affected ice nucleation .
the hydrophilicity 
 of the surface was modified in two different ways : first , by uniformly 
 modifying the water  surface interaction strength and , second , 
 by introducing hydrophilic species at the surface .
it was found that 
 the two ways produced qualitatively different results : uniformly modifying 
 the interaction potential led to enhanced ice nucleation , whereas 
 increasing the density of hydrophilic species was detrimental to ice 
 nucleation ( although the surfaces still enhanced nucleation relative 
 to homogeneous nucleation ) . it was concluded that hydrophilicity is 
 not a good indicator of the ice - nucleating ability of graphitic surfaces . 
 as for the difference between increasing the hydrophilicity by uniformly 
 modifying of the interaction potential and by introducing hydrophilic 
 species , lupi and molinero again saw that the extent of layering in 
 water s density profile above the surface correlated well with 
 the ice - nucleating efficacy .
the general applicability of the layering 
 mechanism , however , was left as an open question .
cox et al . addressed the question of the general applicability of the layering 
 mechanism by investigating ice nucleation rates over a wider range 
 of hydrophilicities ( by uniformly changing the interaction strength ) 
 on two surfaces with different morphologies : ( i ) the ( 111 ) surface 
 of a face - centered - cubic lj crystal ( fcc-111 ) that provided distinct 
 adsorption sites for the water molecules and ( ii ) a graphitic surface , 
 similar to that of lupi et al .
although 
 it was found that the layering mechanism ( albeit with a slight modification 
 to the definition used by lupi et al . ) could describe the ice - nucleating 
 behavior of the graphitic surface , at the fcc-111 surface , no beneficial 
 effects of layering were observed .
this was attributed to fact that 
 the fcc-111 surface also affected the structure of the water molecules 
 in the second layer above the surface , in a manner detrimental to 
 ice nucleation .
it was concluded that layering of water above the 
 surface can be beneficial to ice nucleation , but only if the surface 
 presents a relatively smooth potential energy surface to the water 
 molecules .
the studies at the carbonaceous and fcc-111 surfaces discussed above hinted that the heterogeneous nucleation mechanism 
 could be very different at different types of surfaces .
although there 
 is experimental evidence that , for example , different carbon nanomaterials 
 are capable of boosting ice nucleation ( see , e.g. , ref ( 378 ) ) , most experiments can 
 only quantify the ice - nucleating ability of the substrates ( see section 1.2 ) .
however , 
 the structure of the water  substrate interface and any insight 
 into the morphology of the nuclei are typically not available , making 
 simulations essential to complement the experimental picture . in this 
 respect , zhang et al . assessed that 
 the ( regular ) pattering of a generic crystalline surface at the nanoscale 
 can strongly affect ice formation .
more generally , the interplay between 
 hydrophobicity and surface morphology was recently elucidated by fitzner 
 et al .
brute - force md simulations of 
 heterogeneous ice nucleation were performed for the mw water model 
 on top of several crystalline faces of a generic fcc crystal , taking 
 into account different values of the water  surface interaction 
 strength , as well as different values of the lattice parameter .
the 
 latter is involved in the rather dated concept of zero lattice mismatch , which we introduced in section 2.2 ( see eq 7 ) and which has been often 
 quoted as the main requirement of an effective ice - nucleating agent . 

however , a surprisingly nontrivial interplay between hydrophobicity 
 and morphology was observed , as depicted in figure 16 .
clearly , neither the layering nor the 
 lattice mismatch alone are sufficient to explain such a diverse scenario . 

in fact , the authors proposed three additional microscopic factors 
 that can effectively aid heterogeneous ice nucleation on crystalline 
 surfaces : ( i ) an in - plane templating of the first water overlayer 
 on top of the crystalline surface ; ( ii ) a first overlayer buckled 
 in an ice - like fashion ; and ( iii ) enhanced nucleation in regions of 
 the liquid beyond the first two overlayers , possibly aided by dynamical 
 effects and/or structural templating effects of the substrate extending 
 past the surface water interface .
in addition , it turned out that 
 different lattice parameters can lead to the nucleation and growth 
 of up to three different faces of ice [ basal , prismatic , and secondary 
 prismatic ( { 1120 } ) ] on top of the very same surface , 
 adding a layer of complexity to the nucleation scenario
. insights 
 into the interplay between hydrophobicity and morphology were also 
 very recently obtained by bi et al . , who investigated heterogeneous ice nucleation on top of graphitic 
 surfaces by means of ffs simulations using the mw model .
among their 
 findings , the authors suggested that the efficiency of ice - nucleating 
 agents can be a function not only of surface chemistry and surface 
 crystallinity but of the elasticity of the substrate as well .
( a ) heat maps representing 
 the values of ice nucleation rates on top of four different fcc surfaces 
 [ ( 111 ) , ( 100 ) , ( 110 ) , ( 211 ) ] , plotted as a function of the adsorption 
 energy , eads , and the lattice parameter , afcc .
the lattice mismatch  with respect 
 to ice on ( 111 ) is indicated below the corresponding graph in panel 
 a. the values of the nucleation rate , , are 
 reported as log10(j / j0 ) , where j0 refers to the 
 homogeneous nucleation rate at the same 
 temperature .
( b ) sketches of the different regions ( white areas ) in 
 ( eads , afcc ) space in which a significant enhancement of the nucleation rate 
 is observed .
each region is labeled according to the face of ih nucleating and growing on top of the surface [ basal , prismatic , 
 or ( 11/200 ) , together with an indication of what 
 it is that enhances the nucleation , where  temp  ,  buck  , 
 and  highe  refer to the in - plane template 
 of the first overlayer , the ice - like buckling of the contact layer , 
 and the nucleation for high adsorption energies on compact surfaces , 
 respectively .
reinhardt and doye used umbrella sampling 
 with the mw model to investigate nucleation at a smooth planar interface 
 and at an ice - like surface .
they found that the flat planar interface 
 did not help nucleate ice and that homogeneous nucleation was the 
 preferred pathway .
one explanation given for this finding was that , 
 as the density of liquid water is higher than that of ice , an attractive 
 surface favors the liquid phase .
it was also noted that the mw potential 
 imposes an energy penalty for nontetrahedral triplets , that removing 
 neighbors at the surface decreases this energetic penalty , and that 
 this reduction in tetrahedrality favors the liquid phase .
cabriolu 
 and li recently studied ice nucleation at graphitic surfaces using 
 forward flux sampling , again with the 
 mw model . under the assumptions that   depends linearly on t and 
 that   does not depend on t ,
cabriolu and li also extracted the values of the contact 
 angle 
 at different temperatures , which , along with the free energy barrier , 
 turned out to be consistent with cnt for heterogeneous nucleation 
 ( see section 1.1.3 ) .
although intriguing , the generality of this finding to surfaces 
 that include strong and localized chemical interactions remains an 
 open question .
we have seen that , for both homogeneous and heterogeneous 
 nucleation , using the coarse - grained mw model has greatly enhanced 
 our ability to perform quantitative , systematic simulation studies 
 of ice nucleation .
we must face the fact , however , that this approach 
 will further our understanding of heterogeneous ice nucleation only 
 so far . as discussed for kaolinite ,
an explicit treatment 
 of the hydrogen bonds is essential in describing heterogeneous ice 
 nucleation .
in addition , the mw model ( as well as the majority of 
 fully atomistic water models ) can not take into account surface - charge 
 effects .
surfaces can polarize water molecules in the proximity of 
 the substrate , alter their protonation state , and even play a role 
 in determining the equilibrium structure of the liquid at the interface . 
 in light of recent studies
, it seems that these 
 effects can heavily affect nucleation rates of many different systems . 

as discussed , enhanced 
 sampling techniques such as umbrella sampling and forward flux sampling have been applied to heterogeneous ice nucleation with the mw model , 
 and we have seen the latter applied successfully to homogeneous nucleation 
 with an all - atom model of water ; the 
 computational cost , however , was huge .
although the presence of an 
 ice - nucleating agent should help reduce this cost , the parameter space 
 that we wish to study is large , and systematically studying how the 
 various properties of a surface affect ice nucleation requires the 
 investigation of many different surfaces .
simulating heterogeneous ice nucleation under realistic conditions 
 does not mean just mild supercooling ; we also need realistic models 
 of the surfaces that we wish to study !
most studies of kaolinite have 
 considered only the planar interface , even though , in nature , kaolinite 
 crystals have many step edges and defects .
ice nucleation at agi was 
 also recently studied , although bulk truncated 
 structures for the exposed crystal faces were used . in the case of 

agi(0001 ) , this is problematic , as the wurtzite structure of the crystal 
 means that this basal face is polar and likely to undergo reconstruction .
furthermore , agi is photosensitive , and it 
 has been shown experimentally that exposure to light enhances its 
 ice - nucleating efficacy , suggesting 
 that structural motifs at the surface very different from those expected 
 from the bulk crystal structure are important .
the development of 
 computational techniques to determine surface structures , along with 
 accurate force fields to describe the interaction with water , will 
 be essential if we are to fully understand heterogeneous ice nucleation .
understanding 
 crystal nucleation from solution is a problem of great practical interest , 
 influencing , for instance , pharmaceutical , chemical , and food processing 
 companies .
being able to obtain a microscopic description of nucleation 
 and growth would allow the selection of specific crystalline polymorphs , 
 which , in turn , can have an enormous impact on the final product .
an ( in)famous case illustrating the importance 
 of this issue is the drug ritonavir , originally 
 marketed as solid capsules to treat hiv .
this compound has at least 
 two polymorphs : the marketed and thoroughly tested polymorph ( pi ) and a second more stable crystalline phase pii that appeared after pi went to market .
pii is basically nonactive as a drug because of a much lower solubility 
 than pi . as such and , most importantly , because of the 
 fact that pii had never been properly tested , ritonavir 
 was withdrawn from the market in favor of a much safer alternative 
 in the form of gel capsules .
many other examples could be listed , as various environmental factors ( such 
 as the temperature , the degree of supersaturation , the type of solvent , 
 and the presence of impurities ) can play a role in determining the 
 final polymorph of many classes of molecular crystals .
thus , it is 
 highly desirable to pinpoint a priori the conditions leading to the 
 formation of a specific polymorph possessing the optimal physical / chemical 
 properties for the application of interest .
the term nucleation from solution encompasses a whole range of systems , 
 from small molecules in aqueous or organic solvents to proteins , peptides , 
 and other macromolecular systems in their natural environment .
these 
 systems are very diverse , and a universal nucleation framework is 
 probably not applicable to all of these cases .
the review by dadey 
 et al . discusses the role of the solvent 
 in determining the final crystal .
many aspects of the nucleation of 
 solute precipitates from solution were recently reviewed by agarwal 
 and peters . in this section ,
a central issue 
 with md simulations of nucleation from solution 
 is the choice of order parameters able to distinguish different polymorphs . 

many of these collective variables have been used in enhanced - sampling 
 simulations ( see section 1.3.2 ) .
md simulations of nucleation from solution are particularly challenging 
 because of finite - size effects due to the nature of the solute / solvent 
 system . in the nvt and npt ensembles , where md simulations of nucleation are usually 
 performed
however , 
 the ratio between the numbers of solute molecules in the crystalline 
 phase and in the solution varies during the nucleation events , leading 
 to a change in the chemical potential of the system .
this occurrence 
 has negligible effects in the thermodynamic limit , but it can substantially affect the outcomes of , for example , 
 free - energy - based enhanced - sampling simulations .
simulations of models 
 containing a large number ( 1010 ) of 
 molecules can alleviate the problem , although this is not always the case .
an analytic correction to the 
 free energy for npt simulations of nucleation of 
 molecules from solution was proposed in refs ( 394 and 404 ) on the basis of a number of previous works ( see , e.g. , refs ( 152 , 400 , and 405 ) ) and applied 
 later in ref ( 406 ) as 
 well .
alternative approaches include seeded md simulations ( see section 1.3.1 ) and simulations mimicking the grand canonical ensemble ( vt ) , where the number of constituents 
 is not a constant and the number of molecules in  in this case 
it is worth noticing that nucleation 
 of molecules in solution is a challenging playground for experiments 
 as well .
for instance , quantitative data about nucleation of ionic 
 solutions are amazingly hard to find within the current literature . 

this is in stark contrast with the vast amount of data covering , for 
 example , ice nucleation ( as illustrated in section 2.4 ) . among the countless 
 organic compounds ,
urea molecules can be regarded as a benchmark for 
 md simulation of nucleation from solution .
this is because urea is 
 a system of great practical importance that ( i ) displays fast nucleation 
 kinetics and ( ii ) has only one experimentally characterized polymorph . 

early studies by piana and co - workers focused on 
 the growth rate of urea crystals , which turned out to be consistent 
 with experimental results .
years later , the inhibition of urea crystal 
 growth by additives was investigated by salvalaglio et al .
the investigation of the early stages of nucleation was tackled 
 only recently by salvalaglio et al . for 
 urea molecules in aqueous and organic ( ethanol , methanol , and acetonitrile ) 
 solvents . in these studies ,
the resulting free energies , modified for finite - size effects related 
 to the solvent , suggested that different 
 solvents lead to different nucleation mechanisms . whereas a single - step 
 nucleation process is favored in methanol and ethanol , a two - step 
 mechanism ( see section 1.1.2 )
emerges for urea molecules in acetonitrile and water , as 
 depicted in figure 17 . in this case , the initial formation of an amorphous  albeit 
 dense 
note that , according to the free energy surface reported in figure 17a , the amorphous 
 clusters ( configurations 2 and 3 in figure 17a , b ) are unstable with respect to the liquid 
 phase ; that is , they are not metastable states having their own free 
 energy bases , but rather , they originate from fluctuations within 
 the liquid phase .
this evidence , together with the fact that the transition 
 state ( configuration 4 in figure 17a , b ) displays a fully crystalline core , prompts the 
 following , long - standing question : if the critical nucleus is mostly 
 crystalline and the amorphous precursors are unstable with respect 
 to the liquid phase , can we truly talk about a two - step mechanism ? 

reference ( 394 ) suggests 
 the terms ripening regime two - step nucleation when 
 dealing with stable amorphous precursors and crystallization - limited 
 two - step nucleation when the amorphous clusters are unstable 
 and the limiting step is the formation of a crystalline core within 
 the clusters .
salvalaglio et al . also 
 observed two polymorphs ( pi and pii ) in the 
 early stages of the nucleation process .
pi corresponds 
 to the experimental crystal structure and is the most stable structure 
 in the limit of an infinite crystal .
pii , however , is more stable for small crystalline clusters . 
 in agreement with the ostwald rule ( see section 2.2 ) , the small crystalline clusters that 
 initially form in solution are of the pii type , and the 
 subsequent conversion from pii to pi seems to 
 be an almost - barrierless process .
( a ) free - energy surface ( fes ) associated 
 with the early stages 
 of nucleation of urea in aqueous solution , as obtained by salvalaglio 
 et al . from a well - tempered metadynamics 
 simulation of 300 urea molecules and 3173 water molecules , within 
 an isothermal  isobaric ensemble at p = 1 bar 
 and t = 300 k ( simulation s2 in ref ( 406 ) with a correction term 
 to
the free - energy included to represent the case of a constant supersaturation 
 of 2.5 ) .
the contour plot of the fes is reported as a function of 
 the number of molecules belonging to the largest connected cluster 
 ( n , along the ordinate ) and the number of molecules 
 in a crystal - like configuration within the largest cluster ( no , along the abscissa ) .
note that n  no by definition and that cnt 
 would prescribe that the evolution of the largest cluster in the simulation 
 box is such that n = no ( i.e. , only the diagonal of the contour plot is populated ) .
the 
 presence of an off - diagonal basin provides evidence of a two - step 
 nucleation of urea crystals from aqueous solutions .
this is further 
 supported by the representative states sampled during the nucleation 
 process , shown in panel b. urea molecules are represented as blue 
 spheres , and red connections are drawn between urea molecules falling 
 within a cutoff distance of 0.6 nm of each other .
an approach similar to that employed in ref ( 413 ) was used to investigate 
 crystal nucleation of 1,3,5-tris(4-bromophenyl)benzene molecules in 
 water and methanol .
these simulations showed the emergence of prenucleation 
 clusters , consistent with recent experimental results based on single - molecule real - time transmission 
 electron microscopy ( smrt - tem ; see section 1.2 ) .
the formation of prenucleation clusters 
 in the early stages of nucleation from solution has been observed 
 in several other cases . this is of great relevance , as cnt is not able to account for two- 
 ( or multi- ) step nucleation .
md simulations have been of help in several 
 cases , validating or supporting a particular mechanism .
for instance , 
 md simulations provided evidence for two - step nucleation in aqueous 
 solutions of -glycine and n - octane ( or n - octanol ) solutions of d-/l - norleucine . sodium chloride 
 ( nacl )
nucleation from supersaturated brines represents an interesting 
 challenge for simulations , as the system is relatively easy to model 
 and experimental nucleation rates are available .
the first simulations 
 of nacl nucleation date back to the early 1990s , when ohtaki and fukushima performed brute - force md simulations using 
 very small systems ( 448 molecules including water molecules and ions ) 
 and exceedingly short simulation times ( 10 ps ) .
thus , the 
 formation of small crystalline clusters that they observed was most 
 likely a consequence of finite - size effects .
more recently , the tps 
 simulations of zahn suggested that the 
 centers of stability for nacl aggregates consist of nonhydrated na ions octahedrally coordinated with cl ions , 
 although the results were related to very small simulations boxes 
 ( containing 310 molecules in total ) .
tentative insight into 
 the structure of the crystalline clusters 
 came with the work of nahtigal et al .
, featuring simulations of 4132 molecules ( 4000 water molecules and 
 132 ions ) in the 6731073 k range for supercritical water at 
 different densities ( 0.170.34 g / cm ) .
they reported 
 a strong dependence of the crystalline cluster size distribution on 
 the system density , with larger clusters formed at lower densities . 

the emergence of 
 amorphous precursors was also reported in the work of chakraborty 
 and patey , who performed large - scale md 
 simulations featuring 56000 water molecules and 4000 ion pairs in 
 the npt ensemble .
the spc / e model was used for water , and the ion parameters were those used 
 in the opls force field .
their findings 
 provided strong evidence for a two - step mechanism of nucleation , where 
 a dense but unstructured nacl nucleus is formed first , followed by 
 a rearrangement into the rock salt structure , as depicted in figure 18a . on a similar 
 note , metadynamics simulations performed by giberti et al . using the gromos force field for the ions and the spc / e model for water suggested the emergence of a wurtzite - like polymorph 
 in the early stages of nucleation .
this precursor could be an intermediate 
 state along the path from brine to the nacl crystal .
however , alejandre 
 and hansen pointed out a strong sensitivity 
 of the nucleation mechanism on the choice of the force field . ( a ) snapshots 
 from an md simulation of crystal nucleation of nacl 
 from aqueous solution . the simulations , carried out by chakraborty 
 and patey , involved 56000 water molecules 
 and 4000 ion pairs ( concentration of 3.97 m ) in the npt ensemble .
all na ( black ) and cl ( yellow ) ions within 2 nm of a reference na ion ( larger 
 and blue ) are shown , together with water molecules ( oxygen and hydrogen 
 atoms in red and white , respectively ) within 0.4 nm from each ion . 
 from the relatively homogeneous solution ( 3 ns ) ,
this fluctuation in the concentration of the 
 ions leads to a subsequent ordering of the disordered cluster ( 10 
 ns ) in a crystalline fashion ( 30 ns ) , consistently with a two - step 
 nucleation mechanism .
( b ) comparison of nacl nucleation rates , , as a function 
 of the driving force for 
 nucleation , reported as 1/(/kbt ) .
red points and blue and gray ( continuous ) 
 lines were estimated by three different approaches in the simulations 
 of zimmermann et al .
experimental data 
 obtained employing an electrodynamic levitator trap ( na et al . ) , an efflorescence chamber ( gao et al . ) , and microcapillaries ( desarnaud et al . )
are also reported , together with a tentative 
 fit ( fitexp , dotted line ) .
note the substantial ( up to about 30 orders of magnitude ) 
 discrepancy between experiments and simulations .
in fact , very recent simulations by zimmermann et al . demonstrated that the gromos force field overestimates 
 the stability of the wurtzite - like polymorph .
the authors employed 
 a seeding approach within an nvt setup for which 
 the absence of depletion effects was explicitly verified .
the force fields used were those developed 
 by joung and cheatham for na and cl and spc / e for water , which provide reliable solubilities and accurate chemical 
 potential driving force . using a methodology 
 introduced in ref ( 193 ) ,
a thorough investigation of 
 the latter demonstrated that the limiting factor for n , which , in turn , strongly affects the kinetics of 
 nucleation ( see section 1.1.1 ) , is not the diffusion of the ions within the solution 
 but is instead the desolvation process needed for the ions to get 
 rid of the solvent and join the crystalline clusters .
moreover , zimmermann 
 et al . evaluated the nucleation free 
 energy barrier as well as the nucleation rate as a function of supersaturation , 
 providing three estimates using different approaches .
the results 
 are compared with experiments in figure 18b , showing a substantial discrepancy as 
 large as 30 orders of magnitude .
interestingly , experimental nucleation 
 rates are much smaller than what is observed in simulations , contrary 
 to what has been observed for colloids , for example ( see section 2.1 ) .
we stress 
 that the work of zimmermann et al . employed state - of - the - art computational 
 techniques and explored nacl nucleation under different conditions 
 using a variety of approaches .
the fact that these tour de force simulations 
 yielded nucleation rates that differed significantly from experiments 
 casts yet another doubt on the possibility of effectively comparing 
 experiments and simulations .
however , it must be noted that zimmermann 
 et al . assumed a value of about 5.0 
 molnacl / kgh2o for the nacl solubility 
 in water , as proposed in ref ( 432 ) .
( 3.64 molnacl / kgh2o ) and more recently by mester 
 and panagiotopoulos ( 3.71 molnacl / kgh2o ) .
, once again demonstrating the severe sensitivity of nucleation rates 
 to any of the ingredients involved in their calculations . on 
 a final note ,
we stress that many other examples of molecular 
 dynamics simulations looking at specific aspects of crystal nucleation 
 from solution exist in the literature . for instance , a recent study 
 by anwar et al .
describes secondary 
 crystal nucleation , where crystalline seeds are already present within 
 the solution .
the authors suggest , for a generic solution represented 
 by lennard - jones particles , a ( secondary ) nucleation mechanism enhanced 
 by the existence of pncs ( see section 1.1.2 ) .
kawska et al . stressed instead the importance of proton transfer within 
 the early stages of nucleation of zinc oxide nanoclusters from an 
 ethanol solution .
the emergence of similar ripening processes , selecting 
 specific crystalline polymorphs , for example , according to the effect 
 of different solvents is still fairly unexplored but bound to be of 
 great relevance in the future .
finally , several computational studies 
 have dealt with the crystallization of calcium carbonate , which was 
 recently reviewed extensively in ref ( 18 ) and thus , together with the broad topic of crystal 
 nucleation of biominerals , is not discussed in here .
natural gas 
 hydrates are crystalline compounds in which small gas molecules are 
 caged ( or enclathrated ) in a host framework of water 
 molecules . as natural gas molecules ( e.g. , methane , ethane , propane ) 

are hydrophobic , gas hydrates are favored by conditions of high pressure 
 and low temperature , and are found to occur naturally in the ocean 
 bed and in permafrost regions . with 
 exceptionally high gas storage capabilities and the fact
that it is 
 believed that gas hydrates exceed conventional gas reserves by at 
 least an order of magnitude , there is 
 interest in trying to exploit gas hydrates as a future energy resource . 

although gas hydrates might potentially play a positive role in the 
 energy industry s future , they are currently considered a hindrance : 
 if mixed phases of water and natural gas are allowed to cool in an 
 oil pipeline , then a hydrate can form and block the line , causing 
 production to stall .
understanding the mechanism(s ) by which gas hydrates 
 nucleate is likely to play an important role in the rational design 
 of more effective hydrate inhibitors .
there are two 
 main types of natural gas hydrates : structure i ( si ) , which has a 
 cubic structure ( space group pm3n ) , and structure ii ( sii ) , which also has a cubic structure ( space 
 group fd3m ) .
( there is also 
 a third , less common type , sh , which has a hexagonal crystal structure , 
 but we do not discuss this structure any further here . ) structurally , 
 the water frameworks of both si and sii hydrates are similar to that 
 of ice ih , with each water molecule finding itself in an 
 approximately tetrahedral environment with its nearest neighbors . 

unlike ice ih , however , the water framework consists of 
 cages , with cavities large enough to accommodate a gas molecule . between the si and sii hydrates , there exist three types of cages , 
 which are denoted 56 depending on the numbers of five- and six - sided faces that 
 make up the cage .
for example , common to both the si and sii hydrates 
 is the 5 cage , where the water molecules sit on the vertices 
 of a pentagonal dodecahedron .
along with 5 cages , the 
 si hydrate also consists of a 56 cages , which 
 have two six - sided faces and 12 five - sided faces : there are two 5 cages and six 56 cages in the unit 
 cell .
the sii hydrate , on the other hand , has a unit cell made up 
 of 16 5 cages and eight 56 cages . 
 because of the larger size of the 56 cage , 
 the sii structure forms in the presence of larger guest molecules 
 such as propane , whereas small guest molecules such as methane favor 
 the si hydrate .
( this is not to say that small guest molecules are 
 not present in sii , just that the presence of larger guest molecules 
 is necessary to stabilize the larger cavities . )
the si , sii , and sh 
 crystals structures are shown in figure 19 , along with the individual cage structures . 

further details regarding the crystal structures of natural gas hydrates 
 can be found in ref ( 439 ) .
crystal structures of the si , sii and sh gas hydrates , along with 
 the corresponding cage structures .
first , sloan and co - workers proposed the labile cluster hypothesis ( lch ) , which essentially describes 
 the nucleation process as the formation of isolated hydrate cages 
 that then agglomerate to form a critical hydrate nucleus .
second , 
 the local structure hypothesis ( lsh ) was proposed 
 after umbrella sampling simulations by radhakrishnan and trout suggested that the guest molecules first arrange 
 themselves in a structure similar to the hydrate phase , which is accompanied 
 by a perturbation ( relative to the bulk mixture ) of the water molecules 
 around the locally ordered guest molecules . for the same reasons as 
 already outlined elsewhere ( see section 1.2 )
, it is experimentally challenging the 
 verify which , if either , of these two nucleation mechanisms is correct . 

what we will see in this section is how computer simulations of gas 
 hydrate nucleation have been used to help shed light on this process .
although not the first computer simulation study of natural gas 
 hydrate formation ( see , e.g. , refs ( 444447 ) ) , one the most influential simulation works on gas hydrate formation 
 is that of walsh et al . , in which methane 
 hydrate formation was directly simulated under conditions of 250 k 
 and 500 bar .
it was found that nucleation proceeded through the cooperative 
 organization of two methane and five water molecules into a stable 
 structure , with the methane molecules adsorbed on opposite sides of 
 a pentagonal ring of water molecules .
this initial structure allowed 
 the growth of more water faces and adsorbed methane , until a 5 cage formed .
after persisting for 30 ns , this 5 cage opened when two new water molecules were inserted into the 
 only face without an adsorbed methane molecule , on the side opposite 
 to that where several new full cages were completed .
this opening 
 of the original 5 cage was then followed by the relatively 
 fast growth of methane hydrate .
after 240 ns , the original 5 cage transformed 
 into a 56 cage , a structure not found in 
 any equilibrium hydrate structure .
walsh et al . also found that 5 cages dominated , in terms of abundance , during the early 
 stages of nucleation .
56 cages ( which along 
 with the 5 cages comprise the si hydrate ) were the second 
 most abundant , although their formation occurred approximately 100 
 ns after that of the initial 5 cages .
a significant amount 
 of the larger 56 cages that are found in 
 the sii hydrate was also observed , which was rationalized by the large 
 number of face - sharing 5 cages providing an appropriate 
 pattern .
the 56 cages were also observed 
 in an abundance close to that of the 56 cages . 

the final structure can be summarized as a mixture of si and sii motifs , 
 linked by 56 cages .
a similar structure had 
 previously been reported as a result of hydrate growth simulations .
( a  c ) a pair of methane molecules is 
 adsorbed on either side of a single pentagonal face of water molecules . 

partial cages form around this pair , near the eventual central violet 
 methane molecule , only to dissociate over several nanoseconds .
( d , e ) 
 a small cage forms around the violet methane , and other methane molecules 
 adsorb at 11 of the 12 pentagonal faces of the cage , creating the 
 bowl - like pattern shown .
( f , g ) the initial central cage opens on the 
 end opposite to the formation of a network of face - sharing cages , 
 and rapid hydrate growth follows .
( h ) a snapshot of the system after 
 hydrate growth shows the fates of those methane molecules that made 
 up the initial bowl - like structure ( other cages not shown ) .
. provided useful insight 
 into the hydrate nucleation mechanism , the 
 conclusions were based on only two independent nucleation trajectories . 
 soon after the publication by walsh et al .
water 
 model based on mw water under conditions 
 of 210 k and 500 atm ( the melting point of the model is approximately 
 300 k ) .
owing to the reduced computational cost of the coarse - grained 
 model , they were also able to study a much larger system size than 
 walsh et al .
( 8000 water and 1153 guest molecules vs 2944 water and 512 guest molecules ) . in agreement with walsh
, the initial stages of 
 the nucleation mechanism were also dominated by 5 cages , 
 and a mixture of si and sii motifs connected by 56 cages was observed .
it was also observed that solvent - separated 
 pairs of guest molecules were stabilized by greater numbers of guest 
 molecules in the cluster . as gas hydrates
are composed of solvent - separated 
 pairs of guest molecules as opposed to contact pairs , this suggests 
 a resemblance to the lsh , where the local ordering of guest molecules 
 drives the nucleation of the hydrate .
, however , also 
 found a likeness to the lch : clusters of guest molecules and their 
 surrounding water molecules formed long - lived blobs that slowly diffused 
 in solution .
these blobs could be considered large analogues of the 
 labile clusters proposed in the lch . through analysis of their simulation 
 data ,
jacobsen et al . concluded that the blob is a guest - rich precursor 
 in the nucleation pathway of gas hydrates with small guest molecules 
 ( such as methane ) .
note that the distinction between blobs and the 
 amorphous clathrate is that the water molecules have yet to be locked 
 into the clathrate hydrate cages in the former .
sketch of the nucleation mechanism of methane hydrates 
 proposed 
 in ref ( 451 ) .
clusters 
 of guest molecules aggregate in blobs , which transform into amorphous 
 clathrates as soon as the water molecules arrange themselves in the 
 cages characteristic of crystalline clathrate , which eventually form 
 upon the reordering of the guest molecules  and thus of the 
 cages  in a crystalline fashion .
note that the difference between 
 the blob and the amorphous clathrate is that the water molecules have 
 yet to be locked into clathrate hydrate cages in the former .
copyright 2010 american chemical society . both the work of walsh et al . and that of jacobson 
 et al . suggest 
 that amorphous hydrate structures are involved in the nucleation mechanism , 
 although both studies were carried out under high driving forces . 
 in ref ( 453 )
, jacobsen 
 and molinero addressed the following two questions raised by the above 
 studies : how could amorphous nuclei grow into a crystalline form ? 

are amorphous nuclei precursors intermediates for clathrate hydrates 
 under less forcing conditions ? by considering the size - dependent melting 
 temperature of spherical particles using the gibbs 
thomson 
 equation , jacobson and molinero found for all temperatures that the 
 size of the crystalline critical nucleus was always smaller than that 
 of the amorphous critical nucleus , with the two becoming virtually 
 indistinguishable in terms of stability for very small nuclei of 15 
 guest molecules ( i.e. , under very forcing conditions ) . from a thermodynamic 
 perspective , this would suggest that nucleation would always proceed 
 through a crystalline nucleus .
the observation of amorphous nuclei , even at temperatures as high 
 as 20% supercooling , hints that their formation might be favored for 
 kinetic reasons . employing the cnt expression for the free energy 
 barrier suggested that the amorphous nuclei could be kinetically favored 
 up to 17% supercooling if a/x = 
 0.5 , where a and x are the surface 
 tensions of the liquid - amorphous and liquid - crystal structures , respectively . 

jacobson and molinero estimated x  36 mj / m and 16 < a < 32 mj / m , 
 so it is certainly plausible that amorphous precursors are intermediates 
 for clathrate hydrates under certain conditions .
the growth of clathrate 
 hydrates from amorphous and crystalline seeds was also studied , where 
 it was found that crystalline clathrate can grow from amorphous nuclei . 
 as the simulation led to fast mass transport
, the growth of postcritical 
 nuclei was relatively quick , and the amorphous seed became encapsulated 
 by a ( poly)crystalline shell . under conditions where an amorphous 
 nucleus forms first because of a smaller free energy barrier but diffusion 
 of the guest species becomes a limiting factor , it is likely that 
 small nuclei would have long enough to anneal to structures of greater 
 crystallinity before growing to the macroscopic crystal phase .
it thus appears that gas hydrates might exhibit a multistep nucleation 
 process involving amorphous precursors for reasonably forcing conditions , 
 but for temperatures close to coexistence , it seems that nucleation 
 should proceed through a single crystalline nucleus . by assuming a 
 cnt expression for the free energy ( as well as the total rate ) , knott 
 et al . used the seeding technique
( see section 1.3.1 ) to compute 
 the nucleation rate for si methane hydrate with relatively mild supersaturation 
 of methane , in a manner similar to that of espinosa et al . for homogeneous ice nucleation as discussed 
 in section 2.4.1 .
they found vanishingly small homogeneous nucleation rates of 10 nuclei cm s , meaning that , even with all of earth s ocean waters , the 
 induction time to form one crystal nucleus homogeneously would be 
 10 years !
knott et al . therefore concluded that , 
 under mild conditions , hydrate nucleation must occur heterogeneously .
compared 
 to homogeneous nucleation , the heterogeneous nucleation of gas hydrates 
 has been little studied .
liang et al . investigated the steady - state growth of a h2s hydrate 
 crystal in the presence of silica surfaces , finding that the crystal 
 preferentially grew in the bulk solution rather than at the interface 
 with the solid .
they also observed that , in one simulation , local 
 gas density fluctuations of the dissolved guest led to the spontaneous 
 formation of a gas bubble from solution , which was located at the 
 silica interface .
this had two effects on the observed growth : ( i ) 
 the bubble depleted most of the gas from solution , leading to an overall 
 decrease of the crystal growth rate , and ( ii ) because of the location 
 of the guest bubble , the silica surface effectively acted like a source 
 of gas , promoting growth of the crystal closer to the interface relative 
 to the bulk .
investigated the heterogeneous nucleation 
 of co2 hydrate in the presence of a fully hydroxylated 
 silica surface , first in a two - phase system where the water and co2 were well - mixed and then in 
 a three - phase system where the co2 and water were initially 
 phase - separated .
for the two - phase system , 
 the authors reported the formation of an ice - like layer at the silica 
 surface , above which a layer composed of semi-5 cage - like 
 structures mediated the structural mismatch between the ice - like contact 
 layer and the si hydrate structure above .
in the three - phase system , 
 nucleation was observed at the three - phase contact line , along which 
 the crystal nucleus also grew .
this was attributed to the stabilizing 
 effect of the silica on the hydrate cages , plus the requirement for 
 the availability of both water and co2 . in a later work , 

bai et al . investigated the effects 
 of surface hydrophilicity ( by decreasing the percentage of surface 
 hydroxyl groups ) and crystallinity on the nucleation of co2 hydrate .
they found that , in the case of decreased hydrophilicity , 
 the ice - like layer at the crystalline surface vanished , replaced instead 
 by a single liquid - like layer upon which the hydrate directly nucleated . 
 whereas shorter induction times to nucleation at the less hydrophilic 
 surfaces were reported ,
although certainly an interesting observation , 
 as only a single trajectory was performed for each system , studies 
 in which multiple trajectories are used to obtain a distribution of 
 induction times would be desirable , and as the hydrate actually appears 
 to form away from the surface in all cases , a full comparison of the 
 heterogeneous and homogeneous rates would also be a worthwhile pursuit .
there have also been a number of studies investigating the potential 
 role of ice in the nucleation of gas hydrates .
pirzadeh and kusalik performed md simulations of methane hydrate 
 nucleation in the presence of ice surfaces and reported that an increased 
 density of methane at the interface induced structural defects ( coupled 
 58 rings ) in the ice that facilitated the formation of hydrate 
 cages .
nguyen et al . used md simulations 
 to directly investigate the interface between a gas hydrate and ice 
 and found the existence of an interfacial transition layer ( itl ) between the two crystal structures .
the water molecules in 
 the itl , which was found to be disordered and two to three layers 
 of water in thickness , had a tetrahedrality and potential energy intermediate 
 between those of either of the crystal structures and liquid water . 

the authors suggested that the itl could assist the heterogeneous 
 nucleation of gas hydrates from ice by providing a lower surface free 
 energy than either of the ice  liquid and hydrate  liquid 
 interfaces .
differential scanning calorimetry experiments by zhang 
 et al . found ice and hydrate formation 
 to occur simultaneously ( on the experimental time scale ) , which was 
 attributed to the heterogeneous nucleation of ice , which , in turn , 
 facilitated hydrate formation .
poon and peters provide a possible explanation for ice acting as a heterogeneous 
 nucleating agent for gas hydrates , aside from the structural considerations 
 of refs ( 460 and 462 ) : at a growing 
 ice front , the local supersaturation of methane can be dramatically 
 increased , to the extent that induction times to nucleation are reduced 
 by as much as a factor 10 .
computer simulations 
 of hydrate nucleation have certainly contributed 
 to our understanding of the underlying mechanisms , especially in the 
 case of homogeneous nucleation .
one fairly consistent observation 
 across many simulation studies ( e.g. , refs ( 444 , 445 , 447 , 448 , 451 , 453 , and 464 ) ) suggests that some type of ordering of dissolved guest molecules 
 precedes the formation of hydrate cages .
another is that amorphous 
 nuclei , consisting of structural elements of both si and sii hydrates 
 form when conditions are forcing enough .
nevertheless , open questions 
 still remain . in particular , the prediction that homogeneous nucleation 
 rates are vanishingly small under mild conditions
emphasizes the need to better understand heterogeneous 
 nucleation . to this end , enhanced sampling techniques such as ffs , 
 which was recently applied to methane hydrate nucleation at 220 k 
 and 500 bar , are likely to be useful , 
 although directly simulating nucleation under mild conditions is still 
 likely to be a daunting task .
another complicating factor is that , 
 aside from the presence of solid particles , the conditions from which 
 natural gas hydrates form are often highly complex ; for example , in 
 an oil or gas line , there is fluid flow , and understanding how this 
 effects the methane distribution in water is likely to be an important 
 factor in determining how fast gas hydrates form . in this respect ,
we have described only a fraction of the many 
 computer simulation 
 studies of crystal nucleation in supercooled liquids and solutions . 

still , we have learned that md simulations have dramatically improved 
 our fundamental understanding of nucleation .
for instance , several 
 studies on colloidal particles ( see section 2.1 ) provided evidence for two - step nucleation 
 mechanisms , and the investigation of lj liquids yielded valuable insights 
 into the effects of confinement ( see section 2.2 ) .
in addition , the investigation of more 
 realistic systems has provided outcomes directly related to problems 
 of great relevance .
for example , the influence of different solvents 
 on the early stages of urea crystallization ( see section 2.5 ) has important consequences 
 in fine chemistry and in the fertilizer industry , and the molecular 
 details of clathrate nucleation ( see section 2.5 ) could help to rationalize and prevent 
 hydrate formation in oil or natural gas pipelines .
thus , it is fair 
 to say that md simulations have been and will remain a powerful complement 
 to experiments .
however , simulations are presently affected 
 by several shortcomings , 
 which hinder a reliable comparison with experimental nucleation rates 
 and limit nucleation studies to systems and/or conditions often far 
 from those investigated experimentally .
these weaknesses can be classified 
 in two main categories : ( i ) limitations related to the accuracy of 
 the computational model used to represent the system and ( ii ) shortcomings 
 due to the computational techniques employed to simulate nucleation 
 events.(i)in 
 an ideal world , ab initio calculations 
 would be the tool of the trade . unfortunately , in all but a handful 
 of cases such as the phase - change materials presented in section 2.3 , the time - scale 
 problem makes ab initio simulations of crystal nucleation unfeasible 
 ( see figure 5 ) .
as 
 this will be the status quo for the next few decades , we are forced 
 to focus our efforts on improving the current classical force fields 
 and on developing novel classical interatomic potentials .
this is 
 a fundamental issue that affects computer simulations of materials 
 as a whole . although this is not really an issue for nucleation of 
 simple systems such as colloids ( section 2.1 ) , things start to fall apart when dealing 
 with more realistic systems ( see , e.g. , sections 2.5 and 2.6 ) and become 
 even worse in the case of heterogeneous nucleation ( see , e.g. , section 2.4.2 ) , as the 
 description of the interface requires extremely transferable and reliable 
 force fields .
machine learning techniques such as neural network potentials ( see section 2.3 and refs ( 467 and 468 ) ) are emerging 
 as possible candidates to allow for classical md simulations with 
 an accuracy closer to first - principles calculations , but the field 
 is constantly looking for other options that are capable of bringing 
 simulations closer to reality.(ii)the limitations of the computational 
 techniques currently employed to study crystal nucleation are those 
 characteristic of rare - events sampling .
brute - force md simulations 
 ( see section 1.3.1 ) allow for an unbiased investigation of nucleation events , but the 
 time - scale problem limits this approach to very few systems , typically 
 very distant from realistic materials ( see , e.g. , sections 2.1 and 2.2)although notable exceptions exist ( see section 2.3 ) .
it is also 
 worth noticing that , whereas brute - force md is not able to provide 
 a full characterization of the nucleation process , useful insight 
 can still be gained , for example , into prenucleation events .
enhanced sampling techniques ( see section 1.3.2 ) are rapidly evolving and have the potential 
 to take the field to the next level .
however , free energy methods 
 as they are do not give access to nucleation kinetics and , in the 
 case of complex systems ( see , e.g. , sections 2.4.1 and 2.5 ) , are 
 strongly dependent on the choice of the order parameter . on the other 
 hand , in light of the body of work reviewed , it seems that path - sampling 
 methods can provide a more comprehensive picture of crystal nucleation . 

however , at the moment , these techniques are computationally expensive , 
 and a general implementation is not available yet , although consistent 
 efforts have recently been put in place .
we believe that the development 
 of efficient enhanced sampling methods specific to crystal nucleation 
 is one of the crucial challenges ahead . in 
 an ideal world , ab initio calculations 
 would be the tool of the trade .
unfortunately , in all but a handful 
 of cases such as the phase - change materials presented in section 2.3 , the time - scale 
 problem makes ab initio simulations of crystal nucleation unfeasible 
 ( see figure 5 ) .
as 
 this will be the status quo for the next few decades , we are forced 
 to focus our efforts on improving the current classical force fields 
 and on developing novel classical interatomic potentials .
this is 
 a fundamental issue that affects computer simulations of materials 
 as a whole .
although this is not really an issue for nucleation of 
 simple systems such as colloids ( section 2.1 ) , things start to fall apart when dealing 
 with more realistic systems ( see , e.g. , sections 2.5 and 2.6 ) and become 
 even worse in the case of heterogeneous nucleation ( see , e.g. , section 2.4.2 ) , as the 
 description of the interface requires extremely transferable and reliable 
 force fields .
machine learning techniques such as neural network potentials ( see section 2.3 and refs ( 467 and 468 ) ) are emerging 
 as possible candidates to allow for classical md simulations with 
 an accuracy closer to first - principles calculations , but the field 
 is constantly looking for other options that are capable of bringing 
 simulations closer to reality .
the limitations of the computational 
 techniques currently employed to study crystal nucleation are those 
 characteristic of rare - events sampling .
brute - force md simulations 
 ( see section 1.3.1 ) allow for an unbiased investigation of nucleation events , but the 
 time - scale problem limits this approach to very few systems , typically 
 very distant from realistic materials ( see , e.g. , sections 2.1 and 2.2)although notable exceptions exist ( see section 2.3 ) .
it is also 
 worth noticing that , whereas brute - force md is not able to provide 
 a full characterization of the nucleation process , useful insight 
 can still be gained , for example , into prenucleation events .
enhanced sampling techniques ( see section 1.3.2 ) are rapidly evolving and have the potential 
 to take the field to the next level .
however , free energy methods 
 as they are do not give access to nucleation kinetics and , in the 
 case of complex systems ( see , e.g. , sections 2.4.1 and 2.5 ) , are 
 strongly dependent on the choice of the order parameter . on the other 
 hand , in light of the body of work reviewed , it seems that path - sampling 
 methods can provide a more comprehensive picture of crystal nucleation . 

however , at the moment , these techniques are computationally expensive , 
 and a general implementation is not available yet , although consistent 
 efforts have recently been put in place .
we believe that the development 
 of efficient enhanced sampling methods specific to crystal nucleation 
 is one of the crucial challenges ahead . at the moment
, simulations of crystal nucleation of complex liquids 
 are restricted to small systems ( 1010 particles ) , most often under idealized conditions .
for instance , 
 it is presently very difficult to take into account impurities or , 
 in the case of heterogeneous nucleation , defects of the substrate . 

indeed , defects seem to be ubiquitous in many different systems , such 
 as ice , hard - sphere crystals , lj crystals , and organic crystals as 
 well .
defects are also often associated with polymorphism , but possibly 
 because of the inherent difficulties in modeling them ( or in characterizing 
 them experimentally ) , they are under - represented in the current literature . 

these are important aspects that almost always impact experimental 
 measurements and that should thus be included in simulations as well
. 
 in general , simulations of nucleation should allow us not only to 
 provide microscopic insight but also to make useful predictions and/or 
 to provide a general understanding to be applied to a variety of systems . 

these two ambitious goals are particularly challenging for simulations 
 of heterogeneous nucleation . in light of the literature we have reviewed 
 in this work ,
we believe that much of the effort in the future has 
 to be devoted to ( i ) enabling atomistic simulations of heterogeneous 
 nucleation dealing with increasingly realistic interfaces and ( ii ) 
 obtaining general , maybe non - material - specific trends able to point 
 the community into the right direction , even at the cost of sacrificing 
 accuracy to a certain extent .
on the other hand , we hope that the 
 body of work reviewed here will inspire future experiments targeting 
 cleaner , well - defined systems by means of novel techniques , possibly 
 characterized by better temporal and spatial resolution .
improving 
 on the current limitations of the computational models and techniques 
 would enable simulations of much larger systems over much longer time 
 scales , with a degree of accuracy that would allow a fruitful comparison 
 with experiments .
we think this should be the long - term objective 
 for the field . up to now
, the only way to connect simulations and 
 experiments has been through the comparison of crystal nucleation 
 rates , which even now still exhibit substantial discrepancies for 
 every single class of systems we have reviewed .
this is true not only 
 for complex liquids such as water ( see section 2.4.1 ) but even for model systems such as 
 colloids ( section 2.1 ) .
this , together with the fact that , in some cases , even experimental 
 data are scattered across several orders of magnitude , suggests that 
 we are dealing with crystal nucleation in liquids within a flawed 
 theoretical framework .
it 
 is thus no wonder that every 
 aspect of this battered theory has been criticized at some point . 

for instance , the emergence of two - step ( or even multistep ) mechanisms 
 for nucleation has been reported for many different systems ( see sections 2.1 , 2.2 , 2.5 , and 2.6 ) and can not be easily embedded in cnt as it is , although 
 several improvements on the original cnt formulation have appeared 
 within the past decade ( see section 1.1.2 ) . nonetheless , cnt is basically the 
 only theory invoked by both experiments and simulations when dealing 
 with crystal nucleation from the liquid phase .
cnt is widely used 
 because it offers a simple and unified picture for nucleation and 
 it is often very useful .
however , as demonstrated by both experiments 
 and simulations , even the basic rules governing the formation of the 
 critical nucleus can change dramatically from one system to another . 

thus , we believe that any sort of theoretical universal approach , 
 a brand new cnt , so to say , will be unlikely to significantly further 
 the field . indeed , we fear that the same reasoning will hold for the 
 computational methods required .
we can not think of a single enhanced 
 sampling technique capable of tackling the complexity of crystal nucleation 
 as a whole .
the interesting but uncomfortable truth is that each class 
 of supercooled liquids often exhibits unique behavior , which , in turn , 
 results in specific features ruling the crystal nucleation process . 

thus , it is very much possible that different systems under different 
 conditions could require different , ad hoc flavors of cnt .
although 
 the latter have been evolving for decades , we believe that a sizable 
 fraction of the new developments in the field should aim at producing 
 particular flavors of cnt , specifically tailored to the problem at 
 hand . in conclusion , it is clear that md simulations have proven 
 themselves 
 to be of the utmost importance in unraveling the microscopic details 
 of crystal nucleation in liquids . we have reviewed important advances 
 that have provided valuable insights into fundamental issues and diverse 
 nucleation scenarios , complementing experiments and furthering our 
 understanding of nucleation as a whole .
we feel that the ultimate 
 goal for simulations should be to get substantially closer to the 
 reality probed by experiments and that , to do so , we have to sharpen 
 our computational and possibly theoretical tools . in particular , we 
 believe that the community should invest in improving the classical 
 interatomic potentials available as well as the enhanced sampling 
 techniques currently used , enabling accurate simulations of crystal 
 nucleation for systems of practical relevance .
 the nucleation of 
 crystals in liquids is one of nature s 
 most ubiquitous phenomena , playing an important role in areas such 
 as climate change and the production of drugs . 
 as the early stages 
 of nucleation involve exceedingly small time and length scales , atomistic 
 computer simulations can provide unique insights into the microscopic 
 aspects of crystallization . in this review , we take stock of the numerous 
 molecular dynamics simulations that , in the past few decades , have 
 unraveled crucial aspects of crystal nucleation in liquids . 
 we put 
 into context the theoretical framework of classical nucleation theory 
 and the state - of - the - art computational methods by reviewing simulations 
 of such processes as ice nucleation and the crystallization of molecules 
 in solutions . 
 we shall see that molecular dynamics simulations have 
 provided key insights into diverse nucleation scenarios , ranging from 
 colloidal particles to natural gas hydrates , and that , as a result , 
 the general applicability of classical nucleation theory has been 
 repeatedly called into question . 
 we have attempted to identify the 
 most pressing open questions in the field . 
 we believe that , by improving 
 ( i ) existing interatomic potentials and ( ii ) currently available enhanced 
 sampling methods , the community can move toward accurate investigations 
 of realistic systems of practical interest , thus bringing simulations 
 a step closer to experiments . 
Introduction
Selected Systems
Future
Perspectives