A fingerprint of surface-tension anisotropy in the free-energy cost of nucleation

# A fingerprint of surface-tension anisotropy in the free-energy cost of nucleation

Santi Prestipino111Corresponding author. E-mail: sprestipino@unime.it, Alessandro Laio222E-mail: laio@sissa.it, and Erio Tosatti333E-mail: tosatti@sissa.it Università degli Studi di Messina, Dipartimento di Fisica e di Scienze della Terra, Contrada Papardo, I-98166 Messina, Italy
International School for Advanced Studies (SISSA) and UOS Democritos, CNR-IOM, Via Bonomea 265, I-34136 Trieste, Italy
The Abdus Salam International Centre for Theoretical Physics (ICTP), P.O. Box 586, I-34151 Trieste, Italy
July 20, 2019
###### Abstract

We focus on the Gibbs free energy for nucleating a droplet of the stable phase (e.g. solid) inside the metastable parent phase (e.g. liquid), close to the first-order transition temperature. This quantity is central to the theory of homogeneous nucleation, since it superintends the nucleation rate. We recently introduced a field theory describing the dependence of on the droplet volume , taking into account besides the microscopic fuzziness of the droplet-parent interface, also small fluctuations around the spherical shape whose effect, assuming isotropy, was found to be a characteristic logarithmic term. Here we extend this theory, introducing the effect of anisotropy in the surface tension, and show that in the limit of strong anisotropy once more develops a term logarithmic on , now with a prefactor of opposite sign with respect to the isotropic case. Based on this result, we argue that the geometrical shape that large solid nuclei mostly prefer could be inferred from the prefactor of the logarithmic term in the droplet free energy, as determined from the optimization of its near-coexistence profile.

###### pacs:
64.60.qe, 68.03.Cd, 68.35.Md

## I Introduction

When a homogeneous, defect-free bulk system is brought across a first-order phase boundary, it may survive in its metastable state even for a long time, until the stable phase spontaneously nucleates Kelton (); Kashchiev (). The nucleation process has attracted much attention over the years, both from a fundamental point of view as well as for its great practical interest. To mention but one example, a better control of crystal nucleation in protein solutions could help hinder protein condensation which is at the heart of several human pathologies Vekilov (). Thermal fluctuations continuously sprout droplets of the stable phase inside the metastable mother phase. Small droplets dissolve, for the gain in volume free energy fails to compensate the loss in surface free energy. Occasionally, a droplet is sufficiently large that it is favorable for it to grow. Once this happens, the solid nucleus expands until the whole liquid crystallizes. Quenching the system deeper and deeper lowers the nucleation barrier until the point where the barrier vanishes (kinetic spinodal limit). Beyond this threshold, nucleation ceases and the phase transition occurs through spinodal decomposition and coarsening (i.e., uniformly throughout the material). Classical nucleation theory (CNT) Volmer (); Farkas (); Becker () provides the simplest theoretical framework in which the initial stage of the phase transformation can be described. In this theory, an isolated droplet is schematized, regardless of its size, as a sphere of bulk solid, separated from the liquid by a sharp interface with a constant free-energy cost per unit area (“capillarity approximation”). This gives rise to a (Gibbs) free-energy difference between the supercooled liquid system with and without a solid cluster, that is

 ΔG(V)=−ρs|Δμ|V+(36π)1/3σV2/3, (1.1)

where is the cluster volume, is the difference in chemical potential between solid and liquid, and is the bulk-solid number density. The droplet grows if it exceeds a critical size corresponding to the maximum (), which thus provides the activation barrier to nucleation Bagdassarian ().

The cluster free energy can be accessed numerically via the statistics of cluster size, through which the validity of Eq. (1.1) for specific model interactions can be directly tested. We recently showed Prestipino () that the accuracy of CNT is less than satisfactory in estimating the size probability distribution of clusters, especially the smaller ones, implying that interface-tension estimates based on the use of CNT are systematically in error. We then proposed a more detailed field theory of the nucleation barrier, based on the assumption that clusters are soft and not sharp, and can deviate mildly from the spherical shape (“quasispherical” approximation). If the solid-liquid interface tension is taken to be isotropic, the volume dependence of the Gibbs free energy of a cluster is of the Dillmann-Meier form Dillmann (),

 ΔG(V)=−ρs|Δμ|V+AV2/3+BV1/3+C−79kBTlnVa3, (1.2)

where , and can all be expressed as explicit functions of the “microscopic” parameters entering a Landau free energy, and is a microscopic length. It turned out that the numerical profiles of in a few test cases and at various supersaturations are better reproduced by this theory.

Here we critically reconsider the most severe assumption made in that derivation, namely the isotropy of the solid-liquid interface tension. We show that the theory introduced in Ref. Prestipino () can be extended relaxing this important approximation, and that the results change. Starting once again from a Landau-like theory, we derive an interface Hamiltonian, that allows estimating the probability of observing a cluster of any shape and size. The angular dependence of the interface tension is taken into account by terms that depend on the local orientation of the cluster surface. Within this framework, we calculate in the limit of strong surface anisotropy and compare it with the isotropic case. For large anisotropy, the cluster free energy still retains at large size a logarithmic term, however with a prefactor of opposite sign to the isotropic one. On account of this, we suggest that the nominal shape of large solid nuclei could be guessed from the optimization of the actual close to coexistence. Looking for a numerical exemplification, we conducted 3D Monte Carlo simulations of the Ising model extracting for clusters of variable size , at various distances from coexistence. Although we could not really attain sizes where the anisotropic shape effects are heavy, we do detect evidence that the sign is as expected for large anisotropy.

The paper is organized as follows. We start in Section II by relaxing the approximation of an infinitely sharp cluster interface, with the introduction of a Landau free energy. From that, an effective sharp-interface Hamiltonian is derived in Section III, as an intermediate step to building up a field theory for isotropic surfaces where small shape fluctuations are allowed (Section IV.A). Eventually, this leads to a modified-CNT expression of . In Section IV.B, the issue of interface anisotropies is addressed, and we show by examples how the dependence of the interface free energy on the local surface normal affects the formation energy of a large cluster. Next, in Section V, we check our theory against old and fresh Monte Carlo simulation data for the nucleation barrier to magnetization reversal in the 3D Ising model above the roughening temperature. While confirming that CNT is not generally adequate to fit the numerical data, this analysis also gives a quantitative measure of the errors made with CNT and demonstrates their cancellation in the more general theory. Finally, our conclusions are presented in Section VI.

## Ii Diffuse interface: Landau theory

The main assumption behind CNT is that of a sharp and spherical cluster surface. A way to relax this approximation is through the introduction of a scalar, non-conserved order-parameter (OP) field (“crystallinity”) which varies smoothly from one phase to the other. Hence, the solid-liquid interface becomes diffuse in space, even though only on a microscopic scale. In practice, may be thought of as the local value of the main Fourier coefficients of the crystal-periodic one-body density , i.e., those relative to the reciprocal-lattice vectors which are closest in modulus to the point where the liquid structure factor reaches its maximum Shih (). Otherwise, may be identified with the parameter discriminating between solid and liquid in an ansatz like

 n(x)=(ϕπ)3/2∑Re−ϕ(x−R)2=ρs∑Ge−G2/(4ϕ)eiG⋅x, (2.1)

assuming a specific crystal symmetry and an overall number density .

Across the solid-liquid interface, is no longer constant and, for a system with short-range forces, the thermodynamic cost of the interface may be described through the free-energy functional Cahn (); Cahn2 (); Fisher (); Harrowell (); Shen ()

 G[ϕ;^n]=∫d3x{c(^n)2(∇ϕ)2+κ(^n)2(∇2ϕ)2+g(ϕ(x))}, (2.2)

where are stiffness parameters dependent on the interface orientation as defined by the unit normal and is the specific Landau free energy of the homogeneous system, taken the bulk liquid as a reference. In Eq. (2.2), besides the customary square-gradient term, also a square-laplacian term appears. This is the next-to-leading isotropic term in the gradient expansion of the Landau free-energy density Kardar (). Even though being a fourth-order gradient term, it is however only second-order in the order parameter, and this places it on the same footing as the square-gradient term (hence, potentially relevant). We shall see below that, without such a term, the bending rigidity (i.e., the coefficient of in Eq. (3.15) below) would simply be zero. Below the melting temperature , shows, besides the liquid minimum, also a second and deeper solid minimum. Exactly at coexistence, the two minima are equal, falling at in the bulk solid and at in the bulk liquid, which means that while otherwise.

When boundary conditions are applied such that for , a planar interface orthogonal to is forced to appear in the system. The corresponding OP profile is the stationary solution of (2.2) that satisfies the boundary conditions:

 c(^n)ϕ′′0−κ(^n)ϕ′′′′0=dgdϕ(ϕ0;T=Tm),withϕ0(−∞)=ϕs0andϕ0(+∞)=0. (2.3)

From now on, we simplify the notation by dropping any reference to in , and . Equation (2.3) can be simplified by multiplying both sides by and integrating by parts. We thus arrive at a new boundary value problem:

 κϕ′0ϕ′′′0=c2ϕ′20+κ2ϕ′′20−g(ϕ0),withϕ0(−∞)=ϕs0andϕ0(+∞)=0. (2.4)

Obviously, represents the free-energy cost of the interface at .

At temperature below coexistence, the absolute minimum of falls at for . This can be described by

 g(ϕ)=c2ϕ2+c3ϕ3+c4ϕ4+… (2.5)

with (), all other coefficients being constant.

For the remaining part of this Section, we will assume that and do not depend on . Under this condition, a large solid cluster can be assumed to be spherical, with a OP profile described by  Fisher (), provided the center of is at . From this ansatz, in Prestipino () we derived an expression for the cluster free energy,

 ΔG(R)=4πR2σL(1−2δLR+ϵLR2)−43πR3ρs|Δμ|, (2.6)

in terms of quantities () which depend linearly on the supersaturation . Equation (2.6) resembles the CNT expression, Eq. (1.1), with the crucial difference that the interface free energy is now a function of both and :

 σ(R;T)=σL(1−2δLR+ϵLR2). (2.7)

Exactly of this form is the tension of the equilibrium interface between a liquid droplet and the vapour background in the Lennard-Jones model, as being extracted from the particle-number histogram in grand-canonical simulations of samples of increasing size Block (). At coexistence, the solid-liquid interface tension and the Tolman length Tolman () are given by:

 σm≡σL(Tm)=∫+∞−∞dz[cϕ′20(z)+2κϕ′′20(z)]; δm≡δL(Tm)=−∫+∞−∞dzz[cϕ′20(z)+2κϕ′′20(z)]∫+∞−∞dz[cϕ′20(z)+2κϕ′′20(z)] (2.8)

A nonzero occurs if and when is asymmetric around zero, as is generally the case for the interface between phases of a different nature (see Appendix A). Summing up, Eq. (2.6) describes the corrections to CNT which arise by replacing the assumption of a sharp solid-liquid interface with a more realistic finite width, in the case of isotropic surface tension and Tolman length.

## Iii Shape fluctuations: the interface Hamiltonian

A real cluster may be spherical only on average. Far from being static, clusters fluctuate widely away from their mean shape Filion (); Zykova (). To describe fluctuations, we switch from a description in terms of the crystallinity OP to another in which the cluster shape itself rises to the role of fundamental variable. We begin by deriving a coarse-grained, purely geometrical Hamiltonian for the cluster surface directly from the microscopic free-energy functional (2.2), under the assumption of small deviations of the interface from planarity. The outcome is a Canham-Helfrich (CH) Hamiltonian Canham (); Helfrich (), containing spontaneous-curvature and bending penalty terms in addition to interface tension.

For the present derivation, we build on Refs. Kogon (); Kassner (). Other attempts to derive an effective interface Hamiltonian from a mean-field density functional are described in Napiorkowski (); Segovia-Lopez (). Let the cluster “surface” be depicted as a closed mathematical surface embedded in three-dimensional space and let be the parametrization (coordinate patch) of an infinitesimal piece of . We switch from 3D cartesian coordinates, , to new coordinates (tangential and normal to ) by the transformation

 r=R(u,v)+ζˆn(u,v), (3.1)

where

 ˆn(u,v)=Ru∧Rv|Ru∧Rv| (3.2)

is the unit normal to . For a patch that deviates only slightly from planarity, we may adopt a free energy , thus arriving at the surface Hamiltonian

 Hs[Σ]=∫dudvdζJ{c2(∇ϕ0(ζ))2+κ2(∇2ϕ0(ζ))2+g(ϕ0(ζ))} (3.3)

with . In order to make Eq. (3.3) simpler, it is convenient to view the patch as parametrized in terms of orthonormal, arc-length coordinates, i.e., and all over the patch. Although this construction is rigorously possible only for surfaces having zero Gaussian curvature (Abate (), we can reasonably expect that only small errors of order are made for quasiplanar interfaces. With this caution in mind, we go on to get (see Appendix B):

 ∂r∂u = (1−ζκ(1)n)Ru−ζτgRv; ∂r∂v = −ζτgRu+(1−ζκ(2)n)Rv; ∂r∂ζ = ˆn, (3.4)

where and are the normal curvatures of the - and -lines respectively, and is the geodetic torsion. From Eqs. (3.4), we readily derive the metric tensor ,

 (3.5)

and the Jacobian,

 J=(1−ζκ(1)n)(1−ζκ(2)n)−ζ2τ2g=√g, (3.6)

being the determinant of (3.5). Considering that covariant and contravariant components of a vector are built by projecting it on the bases and , respectively, we can calculate the gradient of a scalar field and the divergence of a vector field A in local coordinates as follows:

 ∇ϕ=∂ϕ∂qαgαβ∂r∂qβand∇⋅A=1√g∂∂qα(√gAα), (3.7)

being the inverse of (3.5). In particular,

 ∇ϕ(ζ)=ϕ′(ζ)ˆnand∇2ϕ(ζ)=ϕ′′(ζ)+ϕ′(ζ)∇⋅ˆn, (3.8)

where

 ∇⋅ˆn=1√g(−κ(1)n−κ(2)n−2ζτ2g). (3.9)

Finally, the mean and Gaussian curvatures of the patch are given by

 H=12∇⋅ˆn∣∣∣ζ=0=−12(κ(1)n+κ(2)n) (3.10)

and

 K=ˆn⋅(∂ˆn∂u∧∂ˆn∂v)=κ(1)nκ(2)n−τ2g. (3.11)

Hence, 1) the mean curvature, which is defined only up to a sign depending on our convention on the orientation of , is half the sum of the two normal curvatures relative to any orthogonal parametrization, i.e., not necessarily the two principal curvatures; 2) since is the product of the two principal curvatures, the geodetic torsion must vanish when the coordinate lines are also lines of curvature.

We are now in a position to simplify Eq. (3.3). Upon using Eq. (2.4) to eliminate in favor of , and inserting Eqs. (3.6), (3.8), (3.10), and (3.11), we eventually get

 Hs = ∫dudvdζ(1+2ζH+ζ2K){cϕ′20(ζ)+32κϕ′′20(ζ) (3.12) + κ2(ϕ′′0(ζ)+ϕ′0(ζ)2H−2ζτ2g1+2ζH+ζ2K)2−κ(ϕ′0(ζ)ϕ′′0(ζ))′⎫⎬⎭.

We now argue that, to a first approximation, any term of order higher than and can be discarded. Moreover, since . Lastly, the geodetic torsion vanishes if we perform a change of integration variables (that is, a change of parametrization) such that the coordinate lines are also lines of curvature Abate2 (). In the end, we are left with the classic Canham-Helfrich Hamiltonian for fluid membranes:

 Hs=∫ΣdS(¯¯¯a+¯¯bH+¯cH2+¯¯¯dK), (3.13)

with the following explicit expressions for the coefficients:

 ¯¯¯a = ∫+∞−∞dζ[cϕ′20(ζ)+2κϕ′′20(ζ)]; ¯¯b = 2∫+∞−∞dζζ[cϕ′20(ζ)+2κϕ′′20(ζ)]; ¯¯c = 2κ∫+∞−∞dζϕ′20(ζ); ¯¯¯d = ∫+∞−∞dζ{ζ2[cϕ′20(ζ)+2κϕ′′20(ζ)]−κϕ′20(ζ)}. (3.14)

A few remarks are now in order: 1) and are reparametrization invariants, hence no ambiguity arises from the arbitrariness of the parametrization used. 2) The above derivation actually applies for just one patch. However, upon viewing as the union of many disjoint patches, the Hamiltonian (3.13) holds for the whole as well. 3) As anticipated, the coefficient of the term in (3.13) could be different from the quoted one since a parametrization in terms of orthonormal coordinates does not generally exist. However, as far as we only allow for clusters with the topology of a sphere, takes the constant value of by the Gauss-Bonnet theorem and the term in can be dropped. Upon comparing the definition of and in Eqs. (3.14) with Eqs. (2.8), we can rewrite Eq. (3.13) in the form (restoring everywhere the dependence upon interface orientation):

 Hs=∫ΣdS(σm(^n)−2σm(^n)δm(^n)H+2λ(^n)H2), (3.15)

where (we note that under the same hypotheses for which Eq. (A.15) holds). 4) The term linear in is related to the spontaneous curvature of , , which is proportional to the Tolman length . A nonzero value of yields a difference in energy between inward and outward interface protrusions, thus entailing a non-zero . The additional fact that in systems, such as the Ising model, where the symmetry is perfect between the two phases then , has long been known Fisher ().

We point out that Eq. (3.15) retains the same form as in the isotropic case Prestipino (). In the general anisotropic case, the dependence of the Hamiltonian parameters on the interface normal is through the constants and , and the function .

## Iv The cluster free energy in two extreme cases: isotropic and strongly anisotropic interface tension

Considering that every single realization of the profile of the cluster surface should be sampled in equilibrium with a weight proportional to , it is natural to define a volume-dependent cost of cluster formation through

 ΔG(V)=−ρs|Δμ|V+Fs(V) (4.1)

with

 Fs(V)=−kBTlnZs(V)=−kBTln{a3∫DΣe−βHsδ(V[Σ]−V)}. (4.2)

In the above expression of the constrained partition function , is a microscopic length of the system, is the volume enclosed by the closed surface , and a yet-to-be-specified integral measure.

While the calculation of for a realistic form of -dependent parameters in (3.15) is certainly possible numerically once the admissible surfaces have been parametrized in terms of a basis of eigenfunctions, some restrictions are to be made in practice if we want to make analytical progress. In the following, we examine two limiting cases for , according to whether it is constant or strongly anisotropic. In general, a strongly anisotropic is typical of e.g. systems where melting is very strongly first order, implying very sharp and thus direction dependent solid-liquid interfaces, such as for example in the case of alkali halides Zykova2 (). That brings about a non-spherical cluster shape through the prescription that the surface free energy be the minimum possible for the given cluster volume . The same condition is responsible for a spherical shape when the interface free energy is isotropic.

### iv.1 Isotropic interfaces

If , and in Eq. (3.15) do not depend on , the shape of a cluster is on average spherical. We here compute the free energy (4.2) assuming small deviations from this shape.

Neglecting overhangs and liquid inclusions, let be the equation of in spherical coordinates. We assume only small deviations from a sphere, i.e., , with  Milner (). Then, we expand in real spherical harmonics,

 ϵ(θ,ϕ)=∞∑l=1l∑m=−lxl,mYl,m(θ,ϕ), (4.3)

and we agree to ignore, from now on, all terms beyond second-order in the coefficients . With these specifications, we obtain approximate expressions for the area of and its enclosed volume, as well as for the mean curvature . Upon inserting this form of in terms of the into Eq. (4.2), we are left with the evaluation of a Gaussian integral. While we refer the reader to Appendix C for all the technicalities, we here quote the result of the calculation. The free energy cost of cluster formation for large is

 ΔG(V) = −ρs|Δμ|V+(36π)1/3σQSV2/3−(384π2)1/3σQSδQSV1/3 (4.4) + 4πσQSϵQS−76kBTln((36π)1/3(Va3)2/3),

where , and can be read in Eq. (C.21). The above formula is strictly valid only near coexistence, where the various assumptions beneath its derivation are expected to hold true. We have thus found that the surface free energy has a form consistent with the Dillmann-Meier ansatz, with -dependent parameters , and that are different (even at !) from the corresponding ones in Landau-theory , and , and with a universal logarithmic correction to the mean-field form of . This term is responsible for the well known exponential prefactor to the nucleation rate Guenther ().

### iv.2 Anisotropic interfaces

We now consider an interface tension of the form:

 σ(^n)=σ100[1+M(^n4x+^n4y+^n4z−1)2] (4.5)

with , written in terms of the cartesian components of the outer normal to the cluster surface. In the infinite- limit, the equilibrium crystal shape is a cube, though rectangular cuboids are also admissible, though not optimal, shapes (they arise at non-zero temperatures). The terms in Eq. (3.15) beyond the first are singular in the limit; however, they would contribute to the surface free energy if were large but not infinite, see more in Appendix D. In the same Appendix we show that the asymptotic, large- free-energy cost of cluster formation is given by:

 ΔG(V)=−ρs|Δμ|V+6σPPV2/3+12νPPV1/3+kBTln(6(Va3)2/3)+const. (4.6)

with . Similarly to the isotropic case, in the cluster free energy (4.6) both a logarithmic term and an offset are added to the classical CNT expression of for a cubic cluster of side . The Tolman term in Eq. (4.6) only appears if we envisage an energy penalty, that is per unit length, also for the edges.

More generally, in all the anisotropic-nucleation models examined in Appendix D, the consideration of clusters of same type but unequal edges/semiaxes provides for “breathing” fluctuations of the surface that determine the appearance of a logarithmic term in . In fact, for all such models, the analytically computed is asymptotically given, as in Eq. (4.6), by the CNT expression – as written for the respective symmetric shape – plus subleading terms in the form of a Tolman term, a universal logarithm ( in dimensions), and a negative offset. The value of is for rectangles and 1 for both cuboids and ellipsoids. This is to be contrasted with the quasispherical-cluster case, where by Eq. (4.4). Apparently, the value of is sensitive to both the space dimensionality and the number of independent parameters that are needed to describe the cluster shape, in turn crucial to determine the entropy contents of the surface degrees of freedom (for a quasispherical cluster, this number of parameters goes to infinity with ). In short, a large anisotropy in the interface tension has the overall effect of drastically reducing the spectrum of thermal fluctuations of cluster shape. The reduction cancels the entropy gain which these fluctuations produced in the isotropic case.

This attractive prediction is a difficult one to fully validate numerically at present. A logarithmic correction to CNT can only be detected if we push the numerical investigation of so close to coexistence as to make the Dillmann-Meier form exact for all but the smallest clusters, and that is still a difficult task (see more in the next Section). In the near future, with faster computers becoming available, we can imagine that it will be possible to directly probe the cluster geometry through the optimization of the logarithmic prefactor in an ansatz of the kind (4.4) or (4.6), and thus choose among the many cluster models on the market the one which is most appropriate to the problem at hand.

## V Numerical assessment of the theory

We now critically consider if there are signatures of the degree of anisotropy of the interface free energy in the free-energy cost of cluster formation for a specific instance of microscopic interaction.

We first recall how the work of formation of a -particle cluster is calculated from simulationsi tenWolde (); Reiss (); Bowles (). Given a criterion to identify solid-like clusters within a predominantly liquid system of particles, the average number of -clusters is given, for , by , where is the chemical potential of the liquid and is the Gibbs free energy of the -cluster, including also the contribution associated with the wiggling of the cluster center of mass within a cavity of volume (observe that CNT estimates as , where is the chemical potential of the solid and a geometrical factor). For rare clusters, it thus follows that . This equation is then taken to represent the work of cluster formation for all . Maibaum Maibaum () has shown that the same formula applies for the Ising model.

However, for quenches that are not too deep, the spontaneous occurrence of a large solid cluster in the metastable liquid is a rare event. This poses a problem of poor statistics in the Monte Carlo (MC) estimation of , which is overcome through e.g. the use of a biasing potential that couples with the size of the largest cluster. In practice, this keeps the system in the metastable state for all the ’s of interest. By properly reweighting the sampled microstates one eventually recovers the ordinary ensemble averages. This umbrella-sampling (US) method was used in Refs. tenWolde (); Pan () to compute for the Lennard-Jones fluid and the 3D Ising model, respectively. The main obstacle to the calculation of by US is the necessity of performing the identification of the largest cluster in the system after every MC move. This problem can be somewhat mitigated by the use of a hybrid MC algorithm Gelb (), which in our case reduced the simulation time by a factor of about 20.

A low-temperature Ising magnet where the majority of spins point against the applied field probably yields the simplest possible setup for the study of nucleation. Along the first-order transition line of the model, where two (“up” and “down”) ferromagnetic phases coexist, the interface (say, (100)) between the two phases undergoes a roughening transition at a certain . The up-down interface tension at coexistence is strongly anisotropic close to zero temperature; moreover, it is either singular or smooth according to whether is below or above . Strictly speaking, the interface tension is anisotropic also above , though less and less so when approaching the critical temperature from below Mon (); Hasenbusch1 (). Exactly at the interface tension critically vanishes Rottman (). When a sample originally prepared in the “down” phase is slightly pushed away from coexistence by a small positive field and thus made metastable, the critical droplet of the “up” phase is expected to be less and less spherical as decreases.

With the 3D Ising model as a test system, we carried out a series of extensive US simulations, computing the cluster free energy relative to the nucleation process of magnetization reversal for a fixed , slightly above the roughening temperature of the (100) facet ( Hasenbusch2 ()), and for a number of values of the external field (, in units). Two up spins are said to belong to the same cluster if there is a sequence of neighboring up spins between them; the counting of clusters was done with the Hoshen-Kopelman algorithm Hoshen (). The absolute value of was determined through a standard MC simulation of the system with all spins down, with no bias imposed on the sampling of the equilibrium distribution. We point out that, at the chosen temperature, the Ising surface tension is barely anisotropic Mon (), which would exclude a net preference for either the spherical or the cubic shape. Furthermore, we are sufficiently far away from not to worry about the percolation transition of geometric clusters which was first described in MuellerKrumbhaar (). This event, which would invalidate the assumption (at the heart of the conventional picture of nucleation) of a dilute gas of clusters, is still far away here.

Coherently with the physical picture at the basis of our theory, we verified for all the considered that clusters close to critical indeed contain the vast majority of up spins in the system. A sample of the critical cluster for is shown in Fig. 1. Looking at this picture, it is hard to say whether this particular realization of the critical cluster resembles more a sphere or a cube. When moving to , a spherical shape is eventually preferred over the cube far above , whereas the opposite occurs much below .

In Fig. 2, the ratio of the surface free energy to the area of the cluster surface is reported as a function of , and the data are fit using the functions (4.4) and (4.6) (we stress that different expressions apply for on the left and right panels of Fig. 2, i.e., and respectively; accordingly, the spherical ’s would typically turn out a factor larger than the cubic ’s). Both fits are based on three parameters, namely , , and , which enter in a different way in Eqs. (4.4) and (4.6). However, the dependence on is similar for the two fitting functions, except for the numerical factor in front of the (parameter-free) term. Looking at Fig. 2, it appears that the quality of the “cubic” fit is slightly better than that of the “spherical” fit, in line with the fact that, for , the Ising surface tension is moderately anisotropic. Clearly, at the nucleus is neither spherical nor cubic, and one may object that neither of the fits would actually be meaningful. We nonetheless argue that, within the uncertainty associated with the finite value in the simulations, the better one of the fits will correspond to the regular shape which is closest to that of the real nucleus, thus giving a qualitative indication of the prevailing isotropic or anisotropic character of the solid-liquid interface tension. When going to smaller and smaller , and provided is sufficiently above , we expect that the “spherical” fit would eventually become better than the “cubic” fit.

## Vi Conclusions

In order to estimate from nucleation the solid-liquid interface free energy of a substance, two indirect routes are available: one is through the measurement of the solid nucleation rate as a function of temperature (see e.g. Li (); Li2 (); Franke ()), the other is via the free energy of solid-cluster formation in a supercooled-liquid host, as determined for example in a numerical simulation experiment for a system model. In both cases, the theoretical framework of classical nucleation theory (CNT) has routinely been employed to extract . This is far from satisfactory, as discussed at length in Ref. Prestipino () and in many other papers, due to the neglected cluster interface-tension dependence on both the droplet volume and the supersaturation .

Concentrating on the expression of the cluster formation energy as a function of and , we gave here an extension of the modified CNT theory first introduced in Prestipino (), now including anisotropy, which is important when only a few interface orientations survive in the equilibrium average cluster shape. We showed that, also in this case, a universal non-CNT term is found in the asymptotic expression of the surface free energy versus volume, so long as an infinity of regular shapes is allowed to occur. However that term has now a different prefactor with respect to the quasispherical case. In particular, the sign is positive for large anisotropy and negative for vanishing anisotropy. The sign of that prefactor, which we surmise is related to the amount of surface entropy developed by cluster shape fluctuations, is proposed as the imprinted signature of the geometrical shapes most preferred by the nucleation cluster – negative for spherical or very isotropic shapes, positive for nearly polyhedric or anyway very anisotropic shapes. For the 3D Ising model slightly above the (100) roughening transition temperature, the detected sign suggests cubic rather than spherical cluster symmetry for moderate supersaturation/external field. Much more work and larger simulation sizes should be needed in the future in order to verify the expected change of sign of the term as spherical shapes will be approached closer and closer to the coexistence line when the temperature is quite larger than (though still far from the critical region).

Acknowledgements

This project was co-sponsored by CNR through ESF Eurocore Project FANAS AFRI, by the Italian Ministry of Education and Research through PRIN COFIN Contract 2010LLKJBX004, by SNF Sinergia Project CRSII2_136287/1, and by EU ERC Advanced Grant 320796.

## Appendix A Calculation of σm and δm

In this Appendix, we provide approximate expressions for the quantities and in Eqs. (2.8) for a specific model of homogeneous-system free energy in the functional (2.2).

Once the exact OP profile of the planar interface has been determined for the given , the explicit values of and , and of can be computed. While and are strictly positive quantities, the sign of is not a priori definite. A special but sufficiently general case of function is the following:

 g(ϕ;T=Tm)=c20ϕ2(1−ϕϕs0)2[1+(γ5+2γ6)ϕϕs0+γ6ϕ2ϕ2s0] (A.1)

with for and for . Equation (A.1) is the most general sixth-degree polynomial which admits two non-equivalent minimum valleys at 0 and , and no further negative minimum between them. For this , the differential equation (2.4) is still too difficult to solve in closed form for generic , even when . Hence, we decided to work perturbatively in , and .

At zeroth order, i.e., , corresponding to theory, the solution to (2.4) is

 ¯¯¯ϕ0(z)=ϕs02{1−tanh(z−Cℓ)} (A.2)

with and arbitrary . We fix by requiring that the interface is centered at (hence ). Then, by still keeping , we switch on and search for a second-order solution to Eq. (2.4) in the form

 ϕ0(z)=¯¯¯ϕ0(z)+κcℓ2χ1(z)+(κcℓ2)2χ2(z). (A.3)

We thus arrive at the two equations:

 c¯¯¯ϕ′0χ′1−g′0(¯¯¯ϕ0)χ1=cℓ2(¯¯¯ϕ′0¯¯¯ϕ′′′0−12¯¯¯ϕ′′20) (A.4)

and

 c¯¯¯ϕ′0χ′2−g′0(¯¯¯ϕ0)χ2=cℓ2(¯¯¯ϕ′0χ′′′1+χ′1¯¯¯ϕ′′′0−¯¯¯ϕ′′0χ′′1)−c2χ′21+g′′0(¯¯¯ϕ0)2χ21, (A.5)

where

 g0(ϕ)=c20ϕ2(1−ϕϕs0)2. (A.6)

By requiring that is centered at we obtain

 χ1(z)=ϕs0cosh2(z/ℓ)(2tanhzℓ−zℓ) (A.7)

and

 χ2(z)=ϕs0cosh2(z/ℓ)(32tanh3zℓ−12zℓtanh2zℓ−8tanhzℓ+2(zℓ)2tanhzℓ−3zℓ). (A.8)

Hence, we find since the function is even. Actually, the result is valid at any order in when (see below). Up to second order in , the values of and are given by:

 σm = [1+25κcℓ2−3835(κcℓ2)2]cϕ2s03ℓ; ϵm = [π2−612+(265−π23)κcℓ2+(1566175−4π23)(κcℓ2)2]ℓ2. (A.9)

Next, we take , and all non-zero and of the same order of magnitude, and search for a first-order solution to (2.4) in the form

 ϕ0(z)=¯¯¯ϕ0(z)+γ5ψ1(z)+γ6ξ1(z)+κcℓ2χ1(z). (A.10)

Upon inserting (A.10) into Eq. (2.4), we obtain two independent equations for and , namely

 c¯¯¯ϕ′0ψ′1−g′0(¯¯¯ϕ0)ψ1=¯¯¯ϕ0g0(¯¯¯ϕ0)ϕs0 (A.11)

and

 c¯¯¯ϕ′0ξ′1−g′0(¯¯¯ϕ0)ξ1=⎛⎜⎝2¯¯¯ϕ0ϕs0+¯¯¯ϕ20ϕ2s0⎞⎟⎠g0(¯¯¯ϕ0), (A.12)

while is still given by Eq. (A.7). The solutions to Eqs. (A.11) and (A.12) such that each term of (A.10) separately meets the requirement of being centered at zero are the following:

 ψ1(z)=−ϕs08cosh2(z/ℓ)(1−ln2+zℓ−lncoshzℓ) (A.13)

and

 ξ1(z)=−ϕs08cosh2(z/ℓ)[3(1−ln2)+3zℓ−3lncoshzℓ−12tanhzℓ]. (A.14)

Upon plugging the by now specified in the integrals defining , and , we eventually obtain the formulae:

 σm = (1+14γ5+1320γ6+25κcℓ2)cϕ2s03ℓ,δm=548(γ5+3γ6)ℓ,and ϵm = [π2−612−π2−648γ5−(17π2240−12)γ6+(265−π23)κcℓ2]ℓ2. (A.15)

We thus see that is generically non-zero and may be of both signs.

In conclusion, we give a proof that vanishes identically for

 g(ϕ)=c20ϕ2(1−ϕϕs0)2, (A.16)

whatever is (a different argument can be found in Fisher ()). Let be a solution to Eq. (2.4) obeying the boundary conditions

 ϕ(−∞)=ϕs0,ϕ(+∞)=0,ϕ′(±∞)=ϕ′′(±∞)=…=0. (A.17)

There is an infinite number of such solutions, differing from each other by a simple translation. Let us first prove that is also a solution to (2.4). We have:

 g(˜ϕ(z))=g(ϕ(−z));˜ϕ′(z)=ϕ′(−z);˜ϕ′′(z)=−ϕ′′(−z);˜ϕ′′′(z)=ϕ′′′(−z). (A.18)

We thus see that

 κ˜ϕ′(z)˜ϕ′′′(z)−c2˜ϕ′2(z)−κ2˜ϕ′′2(z)+g(˜ϕ(z))= κϕ′(−z)ϕ′′′(−z)−c2ϕ′2(−z)−κ2ϕ′′2(−z)+g(ϕ(−z))=0, (A.19)

since Eq. (2.4) is satisfied by for any . Hence, obeys the differential equation (2.4). Moreover, like , also satisfies the conditions (A.17). This is not yet sufficient to conclude that and are the same function since they could differ by a translation along . However, if among the infinite possibilities the one is selected such that , then and the two functions coincide: , implying

 ϕ(z)+ϕ(−z)=ϕs0foranyz. (A.20)

Upon differentiating (A.20) with respect to we find that and the function is even. This is enough to conclude that (the interface is centered in 0). Differentiating (A.20) once more, we obtain and is an odd function of (while is even). As a result, and the proof is complete.

## Appendix B Derivation of Eq. (3.4)

Let be a curve in parametrized by the arc length and denote the Frenet trihedron in . Note that we are using nearly the same symbol for the normal to () and for the normal vector to in (n), though the two vectors are generally distinct. Now consider the Darboux frame with (the unit normal to in ), and . Clearly, by a convenient rotation around , n and b are carried to N and B, respectively. Calling the rotation angle,

 ⎛⎜⎝TNB⎞⎟⎠=⎛⎜⎝1000cosαsinα0−sinαcosα⎞⎟⎠⎛⎜⎝tnb⎞⎟⎠. (B.1)

Using the Frenet-Serret formulae, namely