Do the contact angle and line tension of surface-attached droplets depend on the radius of curvature?
Subir K. Das, Sergei A. Egorov, Peter Virnau, David Winter, Kurt Binder
aa r X i v : . [ c ond - m a t . s o f t ] J a n Do the contact angle and line tension of surface-attached droplets depend on theradius of curvature?
Subir K. Das , Sergei A. Egorov , , Peter Virnau , David Winter , and Kurt Binder Theoretical Sciences Unit, Jawaharlal Nehru Centre for Advanced Scientific Research, Jakkur, Bangalore, 56004 India Department of Chemistry, University of Virginia, Charlottesville, USA Leibniz Institut f¨ur Polymerforschung Dresden,Hohe Strasse 6, D-01069 Dresden, Germany and Institut f¨ur Physik, Johannes Gutenberg-Universit¨at Mainz, Staudinger Weg 9, 55099 Mainz, Germany
Results from Monte Carlo simulations of wall-attached droplets in the three-dimensional Isinglattice gas model and in a symmetric binary Lennard-Jones fluid, confined by antisymmetric walls,are analyzed, with the aim to estimate the dependence of the contact angle (Θ) on the dropletradius ( R ) of curvature. Sphere-cap shape of the wall-attached droplets is assumed throughout. Anapproach, based purely on “thermodynamic” observables, e.g., chemical potential, excess densitydue to the droplet, etc., is used, to avoid ambiguities in the decision which particles belong (ordo not belong, respectively) to the droplet. It is found that the results are compatible with avariation [Θ( R ) − Θ ∞ ] ∝ /R , Θ ∞ being the contact angle in the thermodynamic limit ( R = ∞ ).The possibility to use such results to estimate the excess free energy related to the contact line ofthe droplet, namely the line tension, at the wall, is discussed. Various problems that hamper thisapproach and were not fully recognized in previous attempts to extract the line tension are identified.It is also found that the dependence of wall tensions on the difference of chemical potential of thedroplet from that at the bulk coexistence provides effectively a change of the contact angle ofsimilar magnitude. The simulation approach yields precise estimates for the excess density dueto wall-attached droplets and the corresponding free energy excess, relative to a system withouta droplet at the same chemical potential. It is shown that this information suffices to estimatenucleation barriers, not affected by ambiguities on droplet shape, contact angle and line tension. PACS numbers: e-mail address: [email protected]
I. INTRODUCTION
Fluid nano-droplets on planar substrates are important in the context of heterogeneous nucleation [1–5] and ofmuch recent interest [6–8], given that these have diverse applications in connection with microfluidics, nanofluidics,lithography, etc. [9, 10]. Their properties in (metastable) equilibrium are controlled by a subtle interplay of varioussurface tensions and the excess free energy attributed to the three phase contact line, the so-called line tension ( τ )[11–19] (see Fig. 1 for schematic depiction). However, theoretical and conceptual aspects of line tension are stillcontroversially discussed (e.g. [14–17]) and attempts to estimate it experimentally have often led to results whosevalidity is still debated [13].In this work we shall study various aspects related to wall-attached droplets via computer simulations, as previouslydone, see e.g. Refs. [20–28]. Unlike our previous works [22–25], we pay attention to the fact that for nano-droplets,in heterogeneous context, the contact angle (Θ) may depend on the droplet radius, e.g. because of the line tension[19, 29, 30]. This dependence needs to be taken into account self-consistently when one tries to estimate the line tension.Even if the line tension is disregarded, the contact angle of a nano-droplet may differ from its macroscopic value (Θ ∞ )for other reasons [15]. This problem can be avoided, at least in the framework of models possessing the particle-holesymmetry, like the Ising lattice gas system, between bulk coexisting phases, by considering a simulation geometry(shown in Fig. 2) [28] involving inclined planar interfaces between “antisymmetric” walls. We have demonstrated theviability of this approach for estimating contact angles also by lattice-based self-consistent field theory (SCFT), seee.g. [31], for a symmetric binary mixture between two antisymmetric walls. However, as stated above, the objectiveof the present work is to explore the curvature dependence.The sketch drawn in Fig. 1(a) implicitly uses a description appropriate for macroscopic droplets, such as waterdroplets on a car window visible by the naked eye: the droplet surface is drawn infinitely thin, and the contactline where it hits the wall is well-defined. However, it is doubtful that such a description can still be applied atthe nanoscale without ambiguities resulting from the fact that on molecular scales interfaces are diffuse, there is nowell-defined contact line then, and the nano-droplet shape may differ from the sphere-cap shape appropriate in themacroscopic limit (in the absence of gravity). In this work, we shall present a simulation approach that providesthe possibility to extract properties such as contact angles, line tensions, etc., for liquid nanodroplets. At the sametime, we shall demonstrate that it is possible to obtain estimates for the free energy barrier for the (heterogeneous)nucleation of such droplets and these estimates do not suffer from the ambiguities resulting from situation when themacroscopic description invoked in Fig. 1 (a) is taken too literally.In Sec. II we shall concisely summarize the theoretical background of our work. In Sec. III we briefly describe ourapproach for the Ising lattice gas model, before showing that the resulting estimates for the contact angles indeeddiffer significantly from previous work [22, 23], exhibiting a reasonably pronounced dependence on the droplet radius R . This dependence is qualitatively consistent with standard predictions [19, 29] of a correction of order τ /R , whenline tension estimates for planar interfaces [28] (see Fig. 2) are inserted. Due to anisotropy effects from the lattice,quantitative conclusions from our work are delicate, however. In addition, the lack in accuracy of the data [22, 23]does not allow us to obtain significantly improved estimates for the line tension. Sec. IV presents correspondingresults for a symmetric binary ( A, B ) Lennard-Jones mixture, along with a description of the model and methods.Again quite noticeable differences with respect to the approximate approach of the previous work [25] are found. Inthis case also the results can be interpreted as an effective radius dependence of the contact angle, as mentioned above.Sec. V then discusses estimates for the δµ -dependence of contact angle, obtained via a mean-field type approach [31],coming via the wall tensions, for a lattice model of a binary mixture. Finally we summarize our findings in Sec. VI. II. THEORETICAL BACKGROUND
We now discuss the situation sketched in Fig. 1(a) in more detail. Note that in the thermodynamic limit, whenthe droplet radius of curvature R → ∞ , stable phase-coexistence between a liquid droplet (at density ρ coex ℓ ) andsurrounding vapor (at density ρ coex v ) is possible in the bulk. Then, of course, the chemical potential difference( δµ = µ − µ coex ) between the (liquid) droplet and vapor (at bulk coexistence) is identically zero. As opposed to thespherical shape in the bulk, the droplet has sphere-cap shape at the surface (or wall) and the related contact angle isgiven by Young’s equation [12, 32] γ ℓv cos Θ ∞ = γ wv ( δµ ) − γ wℓ ( δµ ) , δµ = 0 . (1)Here γ ℓv is the tension related to a planar liquid-vapor interface, while γ wv ( δµ ) and γ wℓ ( δµ ) are the surface tensionsof vapor and liquid against the wall [28].In the situation of interest, however, where R is finite, we need to distinguish Θ from Θ ∞ . If we consider the dropletof Fig. 1 (a) to be nanoscopic, the chemical potential difference δµ is positive [see the qualitative demonstration inFig. 1 (b)]. The wall tensions γ wv ( δµ ) and γ wℓ ( δµ ), in that case, will, in general, differ from their counterparts for δµ = 0 – see Ref. [15]. In the framework of Monte Carlo (MC) simulations, it is easily possible to compute thesewall tensions with the variation of δµ (or of densities ρ v and ρ ℓ , respectively) [33, 34]. At least for the model studiedin Ref. [33, 34] the dependence of γ wv and γ wℓ on density is rather pronounced. Since γ wv decreases with density,while γ wℓ increases, this effect does not cancel out in the difference γ wv − γ wℓ , which matters in Young’s equation.In addition, the excess free energy due to the line tension [11–19] plays crucial role. Some authors [29, 35], hence,argued that Eq. (1) should be replaced by γ ℓv cos Θ = γ wv ( δµ ) − γ wℓ ( δµ ) − τR sin Θ . (2)Eq. (2) seems to be almost obvious as a consequence of the mechanical force balance at the contact line [19].However, the real situation is not as simple, given that all interfaces are diffuse objects. In that situation, the notionof a contact line (particularly its length) may be influenced by the use of somewhat arbitrarily chosen “dividingsurfaces” for the interfaces [14]. Even for a spherical droplet in the bulk, the “equimolar dividing surface” needs tobe distinguished from the “surface of tension” [12]. Naturally, the length of the contact line, for a sphere cap, willdepend on the choice of a dividing surface [14, 15]. In any case, when the dependence of the wall tensions on δµ canbe neglected, Eq. (2) can be reduced to the simpler result [19] τR = γ ℓv sin Θ(cos Θ ∞ − cos Θ) . (3)Now we remind the reader of the basic concepts of the theory of heterogeneous nucleation [1–4]. There one doesnot consider a small finite size system with linear dimensions comparable to the droplet radius as drawn in Fig. 1.Systems of interest, in fact, are macroscopic ones, being in a metastable vapor phase ( δµ > δµ , rather than the density ρ in thesystem, is chosen as an independent variable. We emphasize at the outset that this will not be done in the simulationsof the later sections. We will use the density, rather than the chemical potential, as the given independent variable.Then one asks the question, assuming that a wall-attached sphere-cap shaped droplet with radius of curvature R and contact angle Θ forms, what would be its Gibbs free energy cost ∆ F ( R ). Note that this question refers toa somewhat hypothetical situation, since such droplets are intrinsically unstable, out of equilibrium, and hence theuse of equilibrium-type descriptions is somewhat questionable. However, in the end one is mostly interested in theposition R ∗ of the free energy maximum and its value ∆ F ( R ∗ ) only.Hence, the Gibbs excess free energy of the droplet in Fig. 1 (a), relative to a state of the system at density ρ v [seeFig. 1 (b)], without a droplet, becomes (for small δµ )∆ F ( R ) = − V d δµ ( ρ coex ℓ − ρ coex v ) + Aγ ℓv + A b ( γ wℓ ( δµ ) − γ wv ( δµ )) + ℓτ. (4)where ℓ = 2 πR sin Θ is the length of the contact line. Furthermore, simple geometric considerations yield the volume V d , upper surface area A and basal area A b of the droplet as V d = 4 π R f V T (Θ); f V T (Θ) ≡
14 (2 + cos Θ)(1 − cos Θ) , (5) A = 2 πR (1 − cos Θ) , (6)and A b = R π sin Θ , (7)where f V T is the Volmer-Turnbull (VT) function [1, 36]. The first term on the right hand side of Eq. (4) can bemotivated by considering that the volume contribution to ∆ F is typically written as − V d ∆ p , where ∆ p is the pressuredifference between the coexisting phases. The linear expansion of ∆ p , by recognizing that density is a derivativeof pressure with respect to the chemical potential, provides this form. A further assumption here is that the liquidand vapor densities, for the overall box density ρ box , are very close to the corresponding coexistence densities in thethermodynamic limit.In this paper, we shall ignore the possibility that the vapor-liquid interfacial tension γ ℓυ contains a curvaturecorrection [12, 37–45] of order 1 /R , which would make the contribution of the term related to the upper surface areain Eq. (4) of the same order ( ∝ R ) as the line tension term [14] (when R dependence in τ is ignored). In fact, for themodels that we shall explicitly study, viz., the Ising lattice gas and a symmetrical binary Lennard-Jones fluid, no suchcurvature correction is expected – the leading correction there is proportional [41, 42, 46, 47] to R − . Here we shallassume throughout that the radii of interest are large enough such that this curvature dependence of γ ℓυ [12, 37–45]can be neglected. This fact we will justify later. For the lattice gas, however, it is a severe approximation to ignorethe change in interfacial free energies due to lattice anisotropy (see e.g. Ref. [28] for a discussion). To avoid thisproblem, one should, in fact, choose for the latter model the vicinity of the critical point. Critical-like fluctuationsthen will get rid of the anisotropy, allowing one to have spherical droplets in the bulk and sphere-cap shaped dropletson planar walls. As a passing remark, for γ wl ( δµ ) and γ wv ( δµ ) we do not discard here the possibility of the leadingcorrection being proportional to 1 /R . In fact, we will see that δµ ∼ /R which brings a correction linear in 1 /R tothe wall tensions.When we consider a large spherical droplet in the bulk, we have the simpler result [cf. Eq. (4)]∆ F sphere = − π R δµ ( ρ coex ℓ − ρ coex v ) + 4 πR γ ℓv . (8)Extremizing this expression, viz., by equating ∂ ∆ F sphere /∂R | R ∗ to zero, one obtains the standard theoretical ex-pressions for the critical radius R ∗ and associated nucleation free energy barrier ∆ F ∗ hom related to the homogeneousnucleation [1–4]. These are R ∗ = 2 γ ℓv δµ ( ρ coex ℓ − ρ coex v ) , ∆ F ∗ hom = 4 π R ∗ γ ℓv . (9)When the contribution of the line tension to the droplet free energy is negligible, the standard result for the barrier,∆ F ∗ het , against heterogeneous nucleation would result [1, 36]. If we assume that the contact angle takes its macroscopicvalue Θ ∞ , then ∆ F ∗ het = ∆ F ∗ hom f V T (Θ ∞ ) = 4 π R ∗ γ ℓv f V T (Θ ∞ ) , (10)the same VT function f V T (Θ ∞ ) [see Eq. (5)], as in the case of volume, describing the reduction of the barrier incomparison with the homogeneous case.In our previous simulation studies, attempting to extract the line tension [22, 23, 25], the δµ -dependence of thewall excess free energies as well as the change of the contact angle due to line tension [see Eqs. (2) and (3)] or otherreasons, were neglected. Thus, Eq. (10) was simply replaced by an expression∆ F ∗∞ het = 4 π R ∗ γ ℓυ f V T (Θ ∞ ) + 2 πR ∗ τ sin Θ ∞ , (11)that is believed to be correct only in the limit R ∗ → ∞ . However, taking Eq. (3) into account, one rather finds [19],minimizing the free energy at constant droplet volume, that∆ F ∗ het = 4 π R ∗ γ ℓυ f V T (Θ) + πR ∗ τ sin Θ . (12)Note that part of the correction due to the line tension is accounted for by replacing Θ ∞ by Θ in f V T . This is seenby expanding Eq. (3) in leading order in Θ − Θ ∞ and considering the first term on the right hand side of Eq. (12): τRγ ℓυ ≃ sin Θ ∞ (Θ − Θ ∞ ) , (13)and 4 π R ∗ γ ℓυ f V T (Θ) ≃ π R ∗ γ ℓυ f V T (Θ ∞ ) + πR ∗ τ sin Θ ∞ . (14)Accepting Eq. (12), one can show that for large R ∗ the error made by Eq. (11) scales like τ :∆ F ∗ het − ∆ F ∗∞ het ≃ π cos Θ ∞ γ ℓυ sin Θ ∞ τ . (15)Clearly, it is desirable to avoid Eq. (11), which for small R ∗ , as used in previous simulations, is unjustified.One purpose of the present work hence is to reconsider the simulations presented in Refs. [22, 23, 25] and attemptan analysis where both τ and Θ are extracted from the simulation data by avoiding the use of Eq. (11) completely.Ref. [19] contains a detailed derivation of Eq. (12), promising “a rigorous thermodynamic formulation”. However,this derivation also ignores, to some extent like that of Eq. (11), both a possible dependence of τ on R and Θ, andthe possible δµ -dependence of the wall tensions, by the step from Eq. (2) to Eq. (3). In view of these approximations,in our proposed approach we keep in mind that the estimates should be done in such a way that it does not matter towhat extent the difference Θ − Θ ∞ is due to the line tension or due to the possible δµ -dependence of the wall tensions.A second purpose of the present work is to demonstrate that from a simulation of wall-attached droplets in equilib-rium in finite systems in the canonical ensemble, as indicated in Fig. 1, one can obtain ∆ F ∗ het directly, as a functionof δµ , without the need to use the above equations. The idea is to study the density excess ρ box − ρ v (see Fig. 1b)and the associated free energy difference between the two systems (with and without droplet) at the same δµ . III. ISING MODEL SIMULATIONS AND THEIR ANALYSIS
We study the nearest-neighbor Ising model, on a simple cubic lattice, at a temperature k B T /J = 3 . J beingthe exchange constant and k B the Boltzmann constant. At this temperature, effects due to critical fluctuations arestill negligible, given that the critical temperature is [48] k B T c /J ≃ .
51. At the same time, the temperature is highenough so that the anisotropy effects on the interface are sufficiently small and thus, the droplet shape is almostspherical [49].As depicted in Fig. 1 (a), we choose a geometry restricted in the z -direction with linear dimension L z . At thebottom layer ( n = 1), perpendicular to the latter direction, a short-range ( δ -function) positive surface field H acts.An antisymmetric boundary condition is created by putting a surface field of equal magnitude, but of opposite sign,at n = L z . On the other hand, periodic boundary conditions are applied in x and y directions, for which the lineardimensions L x and L y equal L . We emphasize that the present approach, as described below, does not rely on a“microscopic” identification of which spins do belong to a droplet or not. While such an identification would bepossible [49], the strong fluctuations in the interfacial region [see Fig. 3 (a) for a picture in a bulk system] make suchmicroscopic approach inconvenient. Before getting into the situation with a wall-attached droplet we provide a briefdiscussion with respect to the bulk.In this work we employ a lattice variant of the Widom particle insertion method [50], which is used to determinethe chemical potential as a function of density [Fig. 3(c)]. This way nuclei of various shapes and sizes in a system canbe probed. Integration of this curve with respect to density yields the free energy. Alternatively, one could employa successive umbrella sampling technique [51, 52] by varying the overall composition in the system in a quasistaticmanner [22, 23]. The method allows one to obtain the probability distribution for the density in the full range of thelatter. From there one could obtain a free energy profile that contains information on the bulk coexistence values ofdensities as well as excess free energies for interfaces and lines. The chemical potential difference, δµ , for a particularstructure and size with respect to the bulk coexistence, could then be obtained from this free energy profile, by takingderivative with respect to density.The basic principle one needs to apply then is that a liquid droplet of density ρ ℓ coexists with a vapor of density ρ v at the same value of δµ , as depicted in Fig. 1 (b). This, for a particular overall particle density, ρ box , inside thebox, identifies ρ ℓ and ρ v . The latter allows one to estimate the volume contribution of free energy, subtraction ofwhich from the total value provides the excess free energy at ρ box . In addition, the corresponding radius R , in a bulksystem of volume V box , can be extracted from the lever rule as ρ box = ρ v + ( ρ ℓ − ρ v ) 4 πR V box , (16)for a spherical structure of the liquid droplet. This, of course, implicitly implies the use of an equimolar dividingsurface between vapor and liquid so that there are no excess particles attributed to the interface. Since, for any choiceof ρ box we know the appropriate value of δµ and R is straightforwardly calculated, the relation δµR vs. R can beconstructed. The corresponding result in Fig. 3 (b), from bulk, implies an inverse relation between δµ and R for largeenough R . We make use of this relation in the situation with antisymmetric walls as well. In the confined case, withsurface fields ± H , we have ρ box = ρ v + ( ρ ℓ − ρ v ) 4 πR f V T (Θ)3 V box . (17)Of course, a crucial assumption of this approach is that also for radii R that are finite, but much larger than thelattice spacing, the deviation of the shape of a droplet from a sphere-cap can be neglected. Only when this assumptionis accurate, one can write the volume of the wall-attached droplet as [1, 36] 4 πR f V T (Θ) /
3. Assuming that Eq. (17)is accurate, for any given choice of ρ box we know δµ [cf. Fig. 3 (c), for simulation results from confined systems], aswell as the radius R from bulk simulation [cf. Fig. 3 (b)]. Since at the obtained value of δµ , ρ v and ρ ℓ are also known,we can use this information to obtain f V T (Θ) from Eq. (17). This allows us to calculate hence the contact angle Θfrom Eq. (5). E.g., for the case H = 0 and R = 10, this analysis yields Θ ≈ . o , while Θ ∞ = 90 o . Althoughthis difference Θ ∞ − Θ is relatively small, it must not be neglected, contrary to the previous studies [22, 23, 25]. Weemphasize that for this method of estimation of the contact angle of wall-attached sphere-cap shaped droplets for achosen value of δµ (or R ) from Eq. (17), at a chosen value of H /J , neither the knowledge of Θ ∞ nor the knowledgeof the line tension τ is required.At this point we mention that the assumptions implied by Eq. (17) can be avoided if we only like to know theexcess number of particles N exc = L L z ( ρ box − ρ v ) , (18)due to the wall-attached droplet in the simulation box. This is of interest, since the “volume term” of nucleationtheory can also be written as N exc δµ . We will make use of this later.Fig. 4 now shows some central results obtained in the present paper, viz., estimates for Θ as a function of 1 /R , forvarious choices of H /J , as indicated in the figure. Note that the values Θ ∞ , shown for 1 /R = 0, are from the studyof Block et al. [28], for planar inclined interfaces in the Ising model. For H > ∞ are included:the smaller values are obtained from Young’s equation, Eq. (1), while the upper ones are from a “first principle’s”approach [28]. The latter method takes the effects of lattice anisotropy on the interfacial free energy, γ ℓυ , of the latticegas model, into account. This lattice effect is expected to become negligible when the temperature approaches thecritical value, but turns out to yield a difference of a few degrees for the contact angles at the considered temperature k B T /J = 3 .
0. In any case, it is remarkable to observe that in the accessible range of radii (8 ≤ R ≤
20) the contactangles of the droplets are always smaller than both the estimates for Θ ∞ . Roughly, the variation of Θ with 1 /R iscompatible with a linear behavior, as expected from Eq. (3), when R → ∞ [see Eq. (13)].As a first qualitative test of Eq. (3), in Fig. 4 we plot the result of this equation, using for γ ℓυ the estimate γ ℓυ ( R → ∞ ) /k B T = 0 . τ (Θ) the estimates due to Block et al. [28] for planar inclined interfaces that arehitting a flat wall. It is seen that the resulting curves that use Θ ∞ from Eq. (1) always fall slightly below the actualdata for Θ, while the other predictions on Θ ∞ , taking the anisotropic feature of the interfacial tension into account,overestimate it. Clearly, the effects due to this anisotropy, which also is expected to cause a slight but systematicdeviation from sphere shape [49], preclude a more precise test of Eq. (3).Next we consider the numerical computation of the excess free energy of the droplet, making use of a thermodynamicintegration procedure ∆ f s ( ρ box ) = ρ box Z ρ v δµ ( ρ ) dρ, (19)using data such as shown in Fig. 3 (c). Eq. (19) is the excess free energy due to the droplet per lattice site,in the canonical ensemble. Note that Eqs. (4) and (8) refer to the grand-canonical ensemble. In the canonicalensemble the volume term proportional to δµ is not present. So we pick up from ∆ f s ( ρ coex ) just the surface and linetension contributions. Since there is a one-to-one correspondence between ρ box and δµ and hence R , the quantity L L z ∆ f s ( ρ box ) can be re-interpreted in terms of the surface excess free energy F s ( R ). Recently, it has been recognizedthat the translational entropy of the droplet in the simulation box should not be ignored [45], if one compares dataobtained from simulations [Eq. (19)] with theoretical estimates. For a wall attached droplet, we hence should add k B T ln( L ) to ∆ f s ( ρ box ) to obtain F s ( R ). For L = 40 this amounts to about 7 . k B T .In previous work it was assumed that F s ( R ) can be interpreted as [22, 23] F s ( R ) = 4 πR γ ℓυ f V T (Θ ∞ ) + 2 πRτ sin Θ ∞ , (20)which leads to the barrier for heterogeneous nucleation – see Eq. (11). For an understanding on the effects ofcurvature dependence of γ ℓv , by neglecting the correction due to the line tension, the surface part of the free energyof the wall-attached droplet simply could be written as F s ( R ) = 4 πR f V T (Θ) γ ℓυ ( R ) . (21)At k B T /J = 3 .
0, the data of Winter et al. [22, 23] are compatible with γ ℓυ ( R ) ≃ γ ℓυ ( ∞ )1 + ( a/R ) , (22)with γ ℓυ ( ∞ ) /k B T = 0 .
444 and a =1.85 lattice units. Note that for R ≥ F sphere ( R ) = − πR ( ρ coex ℓ − ρ coex υ ) δµ + 4 πR γ ℓυ ( ∞ ) − πγ ℓυ ( ∞ ) a . (23)Hence the result for R ∗ [see Eq. (9)], resulting from ∂ [∆ F sphere ( R )] /∂R ] | R ∗ = 0, is not truly affected by this correctiondue to the R -dependence of γ ℓυ ( R ).Naively one might correct Eq. (20), by rewriting it after taking the actual contact angles into account, as F s ( R ) = 4 πR f V T (Θ) γ ℓυ ( R ) + 2 πrτ ( r, Θ) , (24)where r = R sin Θ is the radius of the circular contact line. The resulting estimates for τ ( r, Θ) are plotted in Fig. 5,versus 1 /r . Unfortunately, the scatter of the resulting data is large and not systematic, indicating that the statisticalaccuracy with which F s ( R ) is obtained from Eq. (19), as well as the accuracy of Θ, does not suffice for a meaningfulanalysis of the dependence of the line tension on the circular radius r and contact angle Θ. It is seen that the resultsfor a planar interface (from Block et al. [28]) increase slightly from τ ( ∞ , o ) ≈ − .
249 to τ ( ∞ , o ) ≈ − . ∞ is varied. The data for τ ( r, Θ) extracted from the droplets seem to be only very roughly consistent with thistrend, but probably lack the necessary accuracy to allow for a clear conclusion. But it is evident that the dependenceof τ on R and Θ is not very strong. However, if Eq. (12) would be used to analyse the data for F s ( R ), the absolutemagnitude of τ ( r, Θ) would be twice as large, making it completely incompatible with the limiting behavior for r → ∞ (planar interface).At first sight, all these facts are surprising. However, while the droplet surface is a two-dimensional object here,the contact line is one-dimensional, and hence much larger finite-size corrections due to fluctuations can be expected.When, for the sake of analogy, we consider the droplets in the two-dimensional Ising model, where the droplet “surface”also is a one-dimensional line, a much larger size effect than in d = 3 [see Eq. (22)] occurs [39, 45] γ ℓυ ( R ) = γ ℓυ ( ∞ ) + 54 π ln RR + const R . (25)Returning now to the starting point of our discussion, introducing the line tension from the mechanical equilibriumof the contact line, Eq. (2), we note that the wall tensions, to first order in δµ , can be written as γ wv ( δµ ) = γ wv (0) + δµ Γ v , (26)and γ wℓ ( δµ ) = γ wℓ (0) + δµ Γ ℓ , (27)which yield [see Eqs. (1) and (2)] γ ℓυ (cos Θ − cos Θ ∞ ) = δµ ∆Γ − τR sin Θ , (28)where ∆Γ is the difference in the adsorption from the vapor (Γ v ) and the liquid (Γ ℓ ). Using Eq. (9), to replace δµ by R , we note that the first term on the right hand side of Eq. (28) is of the same order as the second term. Still anotherexpression has been derived in Ref. [14] ( r = R sin Θ): γ ℓυ (cid:16) cos Θ − cos Θ ∞ (cid:17) = − τR sin Θ − dτdr − R cos Θ ∞ dτd Θ , (29)where a further correction involving the Tolman length [37] has been omitted. Eq. (2) [or Eq. (28)] is compatiblewith Eq. (29) only when τ ( r, Θ) depends neither on r nor on Θ, which clearly is not true.While our work hence cannot fully clarify the problems relating to the line tension, we do get from our study validand reasonably accurate results for the barrier against heterogeneous nucleation as∆ F ∗ het k B T = − δµk B T N exc + F s , (30)where F s = ∆ f s ( ρ box ) + k B T ln( L ), as noted above. Fig. 6 shows a log-log plot of ∆ F ∗ het /k B T versus δµ/k B T . For H /J = 0 and 0 . F ∗ het , class k B T = 4 π R ∗ γ ℓv f V T (Θ ∞ ) k B T = 16 π (cid:16) γ ℓv k B T (cid:17) ( ρ coex ℓ − ρ coex v ) − (cid:16) δµk B T (cid:17) − f V T (Θ ∞ ) . (31)One sees that the actual barriers, say for H /J = 0, are about a factor 1.6 lower, though the general trend is similar.If there is no line tension effect, one would also have the relation ∆ F ∗ het = F s /
3. For H /J = 0 this estimate also hasbeen included, with the above choice of F s , and one sees that this still amounts to an overestimation. Similar, butsomewhat larger, discrepancies apply to the other choices of H /J as well. In fact the simulation results for ∆ F ∗ het are expected to be smaller than the predictions. This is because the negative values of line tension (see Fig. 5) reducethe barriers.Unfortunately, for each choice of H /J only a rather restricted range of δµ/k B T is accessible: this happens becausewe cannot use data for ρ = ρ box [in Figs. 1(b) and 3(c)] that are close to the peaks in the δµ vs. ρ curve; as is wellknown, these peaks correspond to the droplet evaporation/condensation transition [22, 23, 41]. Likewise, data for toolarge ρ (where in the curves of Fig. 3(c) an inflection point is seen) cannot be used either; there the droplet shapechanges from sphere cap to cylinder cap shape (stabilized by the periodic boundary conditions). But the range ofbarriers that is accessible here (about 10 < ∆ F ∗ het /k B T <
60) is in a range that would be physically significant. Notethat the regime of much larger droplets, for which the theory outlined in Sec. II is presumably more accurate, wouldrelate to much larger barriers which are physically irrelevant.
IV. ANALYSIS OF SIMULATION RESULTS FOR THE SYMMETRIC BINARY FLUID
The inter-particle interaction in the binary (
A, B ) fluid mixture, in this section, is defined in terms of the Lennard-Jones (LJ) potential ( r = | ~r i − ~r j | , ~r i and ~r j being the positions of particles i and j , respectively) U LJ = 4 ε αβ [ (cid:0) σ αβ r (cid:1) − (cid:0) σ αβ r (cid:1) ] , α, β ∈ A, B, (32)with σ AA = σ BB = σ AB = σ, ε AA = ε BB = ε. (33)For the interaction strength, ε AB is chosen differently [25, 46, 53], viz., ε AB = ε/
2, to facilitate liquid-liquid phaseseparation, as opposed to the choice of equal diameter ( σ ) for all combinations of particles. To improve the speed of thesimulations, it is preferable to cut the potential at some distance r = r c , for which we choose 2 . σ . But for moleculardynamics simulations (that were performed [53, 54] for understanding of various equilibrium and nonequilibriumdynamical properties of this model) it became essential to make the potential and also the force continuous at r = r c .This can be achieved by modifying the LJ potential as [25, 46, 53] u ( r ) = U LJ ( r ) − U LJ ( r c ) − ( r − r c ) dU LJ ( r ) dr | r = rc , (34)which we used for the present study.Temperature in this model has the unit k B /ε , all lengths will be measured in unit of σ and the dimensionlessdensity is calculated as ρ box = N σ /V box . Here N is the total number of particles (= N A + N B , N A and N B being thenumbers for A and B particles) and V box = L x L y L z σ , where L α ( α = x, y, z ) is the box length along the Cartesiandirection α . In the following, for convenience, we set k B , ε and σ to unity. We may, however, explicitly use themwhen verification of an expression with respect to dimension is required. We set ρ box = 1 and perform all simulationsat T = 1, which is far below the critical temperature ( T c ≃ . γ AB = 0 . ± . A and B particles. The model has the further advantage that the interfacial tension is fully isotropic,due to the off-lattice character.In order to realize a situation with antisymmetric walls, similar to Fig. 1(a), wall potentials u A ( z ) and u B ( z ), actingon the two types of particles, are chosen as [24, 25, 46] u A ( z ) = 2 πρ h ε r n(cid:16) σz + δ (cid:17) + (cid:16) σL z + δ − z (cid:17) o − ε a (cid:16) σz + δ (cid:17) i , (35) u B ( z ) = 2 πρ h ε r n(cid:16) σz + δ (cid:17) + (cid:16) σL z + δ − z (cid:17) o − ε a (cid:16) σL z + δ − z (cid:17) i , (36)where z (0 ≤ z ≤ L z ), as in the previous section, again is the coordinate perpendicular to the walls. An offset δ = σ/ ε r = ε/
15) on both types of particles, so that no particles can leave the system bypenetrating the walls. Note that periodic boundary conditions are used in the x and y directions, with L x = L y = L .The attractive potential, with prefactor ε a [see the last terms on the right hand sides of Eq. (35) and (36)], acts on A particles only from the wall at z = 0, while for B particles it comes only from the wall at z = L z . This is the analogof the antisymmetric surface field ±| H | in the Ising model (Fig. 1). However, while for the lattice model the wallpotential was chosen to be of strictly short range, here we do not use any cut-off for the wall potentials.Our Monte Carlo (MC) simulations [48] with this model, of course, consider local displacement moves of the particles.For such trial moves we randomly choose particles and change their Cartesian coordinates, again randomly, such thatthe displacement lies in the range [ − σ/ , + σ/ N , V box and T are kept fixed, this provides realization ofthe canonical ensemble [48]. For the purpose of efficiently obtaining the effective free-energy density f ( x A , T ) /k B T , x A (= N A /N ), the variable analogous to ρ for the vapor-liquid case in the previous section, being the concentrationof A particles, in the two-phase coexistence region, we use (successive) umbrella sampling method [25, 51]. Thus, inaddition to using displacement as trial moves, we have used identity switches ( A → B → A ) as well. For such moveswe have chosen a particle randomly, identified its type and changed it. This calls for an introduction of differencebetween chemical potential between A and B particles in the Boltzmann factor, for the execution of the Metropolisalgorithm [48]. However, along coexistence such difference is identically zero. Umbrella sampling MC simulations,in addition to providing information on the coexistence curve, helps obtaining the accurate probability distributionfor the fluctuation of x A over the whole spectrum ( ∈ [0 , x A , can be straightforwardly obtained, for bulk as well as for the confined systems. Atthe chosen conditions, the system in the bulk exhibits two-phase coexistence for x A lying between 0 .
03 and 0 .
97. Forthe confined system, we find x coex A = 0 .
035 [see Fig. 7 (a)].Fig. 7 shows examples for both f ( x A , T ) /k B T [see part (a)] and δµ ( x A , T ) /k B T [see part (b)]. The latter can beobtained from the derivative of f with respect to x A . While in part (a) we include results from three choices of ε a ,including ε a = 0, part (b) contains results only for ε a = 0 .
1, since chemical potential difference, vs. x A , is shown onlyto demonstrate the method. For both the quantities we have included results only for a part of the overall x A range,that is relevant for the analysis. It is clear from Fig. 7 (a) that, as expected, with the increase of ε a the energy barrierdecreases, due to enhancement of favorable wetting condition. The regions between the first two knees of the freeenergy curves correspond to sphere-cap structure of droplets rich in A particles, radius of which increases with theincrease of x A . By moving towards right one will encounter transitions to cylinder-like and slab-like structures, belowthe wetting transition. These are not of interest in this work. In this figure we have outlined the estimation of theconcentration difference ∆ x , from which, via a relation analogous to Eq. (17), the actual contact angle is extracted.Specifically, ∆ x is related to the (wall attached) droplet volume ( V d ) as∆ x = (1 − x coex A ) V d V box . (37)Thus, one can straightforwardly calculate f V T , thus, Θ, as a function of x A , by obtaining information on the volumereduction compared to a droplet in the bulk with same radius. This, of course, requires the knowledge of R . Estimationof the latter we discuss below, though already described in the previous section. In fact, missing information in Sec.III, if any, can be obtained from the details below.Due to the pronounced statistical fluctuations of the chemical potential difference, the estimation of the latter, aswell as of ∆ x , for chosen states x box A , has to be done with care. To facilitate working with smooth curves, we havetaken help of fitting of the simulation data to nonlinear functions. E.g., the ascending branch of δµ vs. x A plot wasfitted to polynomial forms of degree four, while the descending branch was fitted to power laws. For the correspondingradius R needed here, we use the relation between R and δµ/k B T found in the bulk. Here we make use of the factthat droplets of same R exist at same chemical potential difference, irrespective of whether they are attached to thewall or not. Related results are presented in Fig. 8 (a). To avoid the finite-size effects we have considered onlyoverlapping data from different system sizes, as shown. For this purpose also we have used the fitting method – seecaption for details. There, essentially an inverse relation between δµ and R emerges. This fact is consistent with ourdiscussion above [see Eqs. (8) and (9), along with related text] with respect to vapor-liquid transition, concerning thehomogeneous nucleation theory for critical radius and corresponding energy barrier. Even the prefactor ( ≃ .
51) isalmost theoretically perfect, given that for the binary mixture R = 2 γ AB (1 − x coex A ) δµ , (38)where γ AB is A − B interfacial tension.Here we mention that in an earlier work [25] we had obtained R by using Θ( R = ∞ ), calculated from an independentstudy that used a thermodynamic integration method, in the expression V d = (4 / πR f V T (Θ). This and subsequentmethod did not naturally allow us to obtain any curvature dependence of Θ and τ .From Fig. 7 (a) we can extract the actual (surface) excess free energy due to the droplet, from the difference (∆ f )between f ( x A , T ) /k B T for the states at x νA and x box A , as a function of x box A , and thus, of R , given that F s = ∆ f V box .This can be compared with the result one would obtain if the line tension effects were negligible [see lines of varioustypes in Fig. 8 (b)], i.e., with F s /k B T = 4 πR γ AB f V T (Θ), using the actual contact angles found from the knowledgeof ∆ x . The simulation results for F s are presented as symbols in Fig. 8 (b). One sees that the actual surface excessfree energies are clearly smaller, indicating the presence of a negative line tension.The line tension can be calculated from the discrepancies between the continuous lines and the symbols in Fig. 8(b), by using Eq. (11). Our best estimates for both contact angle and line tension are shown in Fig. 9, vs. 1 /R . InFig. 9 (a) we show contact angles from three choices of ε a . Note that for ε a = 0 one expects Θ( R = ∞ ) = 90 ◦ . In thisfigure we have also included the results that one would obtain from Eq. (3). While Eq. (3) gives qualitatively the righttrend, a good quantitative agreement is not obtained. Unfortunately, only a rather small range of values for 1 /R isaccessible for our analysis. Hence, as in the case of the Ising model, rather tentative conclusions on the dependence ofthe line tension on contact angle and radius can be inferred. However, we can conclude rather safely that the contactangle for small wall-attached droplet is always significantly smaller than its asymptotic value. Hence, as mentionedabove, the results for the line tension obtained previously in Ref. [25], for this model, suffer from systematic errorsdue to this change of contact angle for small droplets.We add here that for this model we have an estimate of τ from planar interface available only [24, 25] for ε a = 0.For other values of ε a we have obtained thermodynamic limit value of τ from extrapolations of R -dependent data of τ to R = ∞ . (Note that these plots are slightly different from the Ising lattice gas case for which we have plotted τ as a function of 1 /r , not 1 /R .) This procedure is shown in Fig. 9 (b) – see the right frames. In the left frames of thisfigure we have shown magnified pictures for the radius dependence of contact angle for the nonzero values of ε a . Thevalues of Θ( R = ∞ ) obtained from here are used for the theoretical lines of Fig. 9 (a).Before closing this section, in Fig. 10 we present results for the barriers for heterogeneous nucleation, i.e., weplot ∆ F ∗ het /k B T versus δµ/k B T , for different values of ε a . These results are analogous to Fig. 6 (that correspondsto the lattice gas model). Here also we compare the simulation results with the classical theory (marked as CNT).The discrepancies between theory and simulation again imply negative line tension. Here we have not included theentropic contribution to F s .0 V. MEAN-FIELD TYPE CALCULATION FOR THE DEPENDENCE OF CONTACT ANGLE ON THECHEMICAL POTENTIAL DIFFERENCE
In this section, we return to the Ising model between antisymmetric walls, and reinterpret the two spin values asparticles of type A and B. Being interested in temperatures far below the bulk critical temperature, it is temptingto study contact angles within a mean-field type approximation. This is to obtain clarification to what extent thevariation of the line tension with droplet radius, implied by the results of the previous sections, is a fluctuation effect,or simply an effect due to the dependence of wall free energies on the chemical potential difference between the dropletand its environment.A convenient construction of mean-field theory for inhomogeneous lattice problems is the Scheutjens-Fleer formu-lation [55] of self-consistent field theory (SCFT). While SCFT was originally intended for polymeric systems, it canbe applied to a lattice model for binary mixtures as well. Since SCFT directly yields the free energy of the sys-tem, one can obtain wall excess free energies straightforwardly [31, 55]. Denoting the normalized interaction energy(Flory-Huggins parameter) as χ AB = q [ ε AB − . ε AA + ε BB )] k B T , (39)where q = 6 is the coordination number of the simple cubic lattice, we note that criticality in the bulk would occurfor χ crit AB = 2. So the regime of interest is χ AB >
2. Similarly, one can define the normalized interaction energy of A -particles with S -particles forming the wall as χ AS = q [ ε AS − . ε AA + ε SS )] k B T , (40)and likewise, the normalized interaction energy of B -particles with S -particles ( χ BS ) can be defined by replacing A with B everywhere in the above equation.Using walls, which attract either A -particles with attraction strength χ AS or B -particles with attraction strength χ BS , one can both obtain wall excess free energies at coexistence conditions, and study also antisymmetric systemswith an inclined interface, similar to M¨uller et al. [31]. Using Young’s equation, Eq. (1), one readily obtains thecontact angle Θ ∞ as a function of these energy parameters (see Fig. 11 (a) as an example). One sees that Θ ∞ vanishes when the wetting transition occurs but for χ AB = 2 . ∞ occur. Hence weshall focus on the R -dependence of Θ for this case in the following.For the calculation of wall tensions γ wv ( δµ ) and γ wℓ ( δµ ) from this mean-field theory, one uses a L × L × L z geometry,with periodic boundary conditions invoked in x and y directions, and on a mean-field level only inhomogeneity in z -direction occurs. Thus the implementation of SCFT for this calculation is rather straightforward.In order to carry out a meaningful comparison with the result of the previous sections, we assume chemical potentialdifferences related to a droplet radius according to Eq. (38). For the chosen example (Fig. 11) with χ AB = 2 . γ AB = 0 . − x coex A = 0 .
37. Being only interested in the leading behavior for large R , it suffices toexpand the wall tensions at the coexistence curve according to Eqs. (26) and (27) by identifying the liquid phase with A -rich phase and the vapor phase with the B -rich phase. From Eqs. (38), (26) and (27) we immediately concludecos Θ = γ wB − γ wA γ AB = cos Θ ∞ + 2∆Γ(1 − x coex A ) 1 R , (41)where ∆Γ = Γ A − Γ B , computed at coexistence. Note that Eq. (41) holds beyond mean field, for R → ∞ , if thedependence on wall excess free energies on δµ is the only reason for a difference between Θ and Θ ∞ .In Eq. (41), cos Θ ∞ simply is the result according to Young’s equation, applying for R → ∞ . Clearly, mean-fieldtheory predicts that the dependence of the contact angle on the chemical potential difference δµ , and hence viaEq. (38) indirectly on R , is a rather pronounced effect! This finding confirms the concerns raised by Schimmele andDietrich [15] that one has to be very careful when one wishes to associate the radius dependence of the contact angleof wall attached sphere-cap shaped droplets with the “true” line tension.It is tempting to use a higher-dimensional version of this mean-field theory for computing the line tension [56].However, due to lattice effects (interfaces in mean-field theory stay non-rough up to bulk criticality) this is ratherdifficult to implement, and hence not done here.1 VI. DISCUSSION
In this paper, we have presented a re-analysis of simulation data, that were used in two earlier attempts [22, 23, 25],to extract estimates for the contact angle Θ of wall-attached sphere-cap shaped droplets and the line tension relatedto the corresponding three-phase contact, by paying attention to a dependence on the droplet radius. In the previousworks, a heuristic assumption was made that the contact angle of these droplets is the same as the contact angle Θ ∞ of macroscopic droplets, which can be estimated independently from Young’s equation. As a matter of fact there isno theoretical basis for this assumption, and hence it is clearly desirable to avoid it. We show here that this task canindeed be achieved by the simulation strategy of the previous work [22, 23, 25], by using the explicit knowledge ofthe relation between the droplet radius of curvature R and the chemical potential difference δµ that characterizes theequilibrium between the system containing the droplet and bulk phase coexistence, and the excess density due to thedroplet.However, a crucial approximation which remains is that for the excess volume of the wall-attached droplet a sphere-cap shape is accurate, i.e., corrections to Eq. (5) due to deviations of the droplet shape near the contact line aresufficiently small so that they can be neglected. This is plausible since only a volume of order 2 πRa sin Θ, where a is a length of the order of molecular distances, should be affected. Hence for the volume of the sphere cap this is acorrection of order ( a/R ) . So we expect a correction for our contact angle of order ( a/R ) , while the corrections thatare found and discussed here are of order 1 /R , much stronger than above errors can bring in. In the case of the latticegas model, a further caveat is that the anisotropy of the lattice also causes systematic deviations of droplet shapesfrom being spherical (or sphere-cap, respectively). But at the chosen temperature these anisotropy effects should bevery small.If these assumptions are accepted, direct estimation of Θ as a function of R is possible. We have interpreted thedeviations between Θ and Θ ∞ in terms of line tension effects, notwithstanding the knowledge that in view of thecriticisms raised by Schimmele et al. [14, 15] this is doubtful. For the Ising model, an estimation of the dependenceof the “true” line tension (referring to macroscopic planar interfaces) on contact angle is available [28], and thisdependence is incompatible with previous estimates resulting from droplets [22, 23]. This led to a motivation for thepresent work. However, the present results are hardly compatible with the implicit assumption of [19], viz., the linetension is a true constant, i.e., τ does not depend upon R and Θ. We have, however, provided evidence, already onthe mean-field level, that the dependence of contact angles (even for R → ∞ ) on the chemical potential δµ in thesystem (for finite R , δµ must be nonzero in equilibrium) is of a similar magnitude as the dependence attributed toline tension effects.The above findings agree with concerns raised by Schimmele and Dietrich [14, 15]. Of course, computer simulationssuffer from well-known problems such as finite size effects, statistical errors, incomplete equilibration, etc. [48]. Henceit is not straightforward to clarify all the issues that our work raises. However, other numerical approaches such asdensity functional theory [18] seem to provide a rather strong variation ( τ ∝ R ) for small R , and also raise questions.In conclusion, although the concept of the line tension was introduced by Gibbs more than a century ago [11], itsnumerical estimation by experiment, simulation and theory still is difficult. An interesting point is that for both themodels, lattice gas and binary Lennard-Jones, the variation of the contact angle of droplets seems to be compatiblewith a linear behavior in 1 /R .Finally, we again like to emphasize that our simulation approach yields directly estimates for the barrier ∆ F ∗ het that needs to be overcome in the heterogeneous nucleation events, which are not at all affected by uncertainties in ourknowledge of contact angles and line tensions. These estimates (e.g. Fig. 6 and Fig. 10) do not make any assumptionwhat the actual droplet shapes are, but rely fully on our estimates of the excess density due to the droplet andassociated free energy excess. Acknowledgment
SKD acknowledges the Marie Curie Actions Plan of European Commission (FP7-PEOPLE-2013-IRSES grant No.612707, DIONICOS), International Centre for Theoretical Physics, Trieste, and Johannes-Gutenberg University ofMainz for partial supports.2 L L p.b.c z x(y) wall (favoring vapor) wall (favoring liquid)vapor (density ϱ , δμ ) l center of mass v γ γ ( δμ ) γ ( δμ ) vl wv Rsin Θ z contactline Θ x liquid (density ϱ , δμ )sphere centerR (a) (b) FIG. 1: (a) Simulation geometry used for the study of wall-attached sphere-cap shaped droplets in the Ising lattice gas model,showing schematically a cross section through the droplet center of mass (marked by a cross) in the zx -plane. The lattice haslinear dimensions L in x and y direction and L z in z -direction, with periodic boundary conditions (p.b.c.) in x and y directionsonly. At the lower wall [first plane ( n = 1) of the simple cubic lattice in z -direction] a positive surface field H acts and atthe upper wall ( n = L z ) a negative surface field H L z = − H acts, so that an antisymmetric boundary condition is created.Otherwise, free boundary conditions in z -direction are used (i.e, missing spins in planes n = 0 and n = L z + 1). The radius ofcurvature R and the contact angle Θ of the droplet are indicated. The contact line meeting with the shown plane is indicatedby a full dot (only on the left side). Note that the shown state is in stable equilibrium for a chemical potential µ = µ coex + δµ ,where µ coex is the chemical potential at bulk coexistence, with vapor and liquid densities ρ coex v and ρ coex ℓ , respectively. (b)Chemical potential difference δµ plotted vs. density ρ , for the geometry of part (a) which refers to a choice of the density ρ = ρ box to the right of the peak (related to the droplet evaporation/condensation transition, and is not of interest here). Forfinite size of the radius R in (a), δµ > ρ v and ρ ℓ are enhanced in comparison with theircorresponding values at two-phase coexistence in the bulk. magnetic (cid:1) eld: FIG. 2: Sketch of the geometry used to obtain the “macroscopic” line tension for the Ising model, considering a simple cubiclattice with linear dimensions L x , L y and L z in the x, y , z directions, respectively, in the limit where all these linear dimensionsbecome macroscopically large. While a periodic boundary condition is used in the z -direction, in the x -direction two freesurfaces of linear dimensions L y × L z are used, at which boundary fields H = −| H | at the bottom and H n = + | H | at thetop act. When one uses a generalized anti-periodic boundary condition (GAPBC), see [28], phase coexistence between twodomains with opposite magnetization + m coex and − m coex in the bulk (symbolized by the two double arrows) occurs. Thesedomains are separated by a planar interface that is inclined by an angle Θ with respect to the yz -plane. (a) ( µ − µ c o e x ) R / ( k B T ) (b) ρ -0.100.10.20.3 ( µ − µ c o e x ) / ( k T ) period. H /J=0.0H /J=0.1H /J=0.2H /J=0.3H /J=0.4H /J=0.5 (c) FIG. 3: (a) Snapshot of an Ising model in the bulk at k B T /J = 3 . L × L × L simulation box with L = 40 lattice units,periodic boundary conditions throughout, and ρ box = 0 .
1. Occupied lattice sites are marked by dots. (b) δµR/k B T is plotted vs. R , for k B T /J = 3 .
0. Note that for R ≥ R -dependence of this product any more. The observed averagevalue (0.937(1)), highlighted by a broken horizontal line, slightly exceeds the theoretical value 0.889(1), Eq. (9), predicted fromthe interface tension γ ℓυ /k B T = 0 . ρ coex ℓ − ρ coex v = 0 .
948 was also used. (c) Chemical potentialdifference δµ/k B T is plotted vs. ρ at k B T /J = 3 .
0, for a L × L × L system with periodic conditions and L × L × L z systems[of Fig. 1 (a)] with different choices of the boundary field H /J , as indicated. Here L = L z = 40 throughout. Θ ( d e g ) H /J FIG. 4: Plot of the contact angle Θ of wall-attached droplets, for several choices of H /J (symbols), against 1 /R , obtained byassuming a sphere-cap shape of the droplet and using the data of Fig. 3(c) and Eq. (17). The curves were constructed by usingthe line tension estimates for planar inclined interfaces (extracted in Ref. [28]) in Eq. (3). The broken curves correspond toΘ ∞ obtained from the Young’s equation and the full curves used the values of Θ ∞ that were estimated by accounting for theanisotropy effects on the interfacial free energy (see Ref. [28]). -0.35-0.3-0.25-0.2-0.15-0.1 τ (r , Θ ∞ ) H /J FIG. 5: Plots of the line tension τ , in units of k B T , vs. the inverse circular radius 1 /r , for various choices of H /J , as indicated,for k B T /J = 3 .
0, in the simple cubic Ising model. Data for r − = 0 were obtained for planar interfaces by Block et al. [28]. δµ /k B T10100 ∆ F * h e t / k B T H /J F s /3, H /J=0 classical theory, H /J=0classical theory, H /J=0.5 FIG. 6: Nucleation barrier for k B T /J = 3 . L = 40) in the simple cubic Ising model as a function of chemical potential δµ/k B T . The simulation results for various values of H /J are represented by symbols. For H /J = 0 and 0 . H /J = 0 we have also shown a plot of F s / f [ k B T ] - ε a =0.0 0.1 0.20 0.2 0.4 x A δ µ [ k B T ] - L L z : L=24, L z =12 ε a =0.1 (b)L L z : L=24, L z =12(a) x A XXX X x A x A ∆ x coex ν box FIG. 7: (a) Plots of the effective free energy f ( x A , T ) for three values of ε a , as quoted in the figure, versus x A , the concentrationof A particles during umbrella sampling Monte Carlo simulations of the symmetric binary fluid. The chosen linear dimensionsare L = 24 and L z = 12. Compositions and values of f for two states (one with and the other without a droplet – see below)at the same chemical potential difference, identified by a dotted horizontal straight line [see part (b)], for ε a = 0 .
1, are markedby crosses. (b) Chemical potential difference δµ , in units of k B T , is plotted vs. x A , for the case ε a = 0 . L = 24, and L z = 12.The state with ∆ µ = 0 corresponds to bulk coexistence, x coex A . State x νA [marked with a cross on the left] is the analog of thestate denoted as ρ v in Fig. 1(b) for the vapor-liquid coexistence, and the state x box A is the chosen average concentration wherean A-rich droplet coexists with the surrounding B-rich fluid at the same chemical potential difference. δµ/ k B T46810 R L=18 20 24 26R=1.5109 x -0.998 ;Bulk (L )(a) x = δµ/ k B T F s [ k B T ] - L z =12 L=24 L=30L=30 ε a =0 ε a =0.05 ε a =0.1(b) FIG. 8: (a) Variation of droplet radius R with chemical potential difference in the bulk symmetric binary Lennard-Jones fluid(obtained by using cubic boxes of linear dimension L with periodic boundary conditions throughout; the chosen values of L areindicated in the figure). All data are compatible with the relation R = 1 . µ − η where the effective exponent is η = 0 . σ , and ∆ µ in units of k B T . (b) Effective surface free energy F s (in units of k B T ) is plotted versusthe droplet radius R for three choices of ε a , as indicated. For ε a = 0, L = 24 was used, and L = 30 for the other choices. Thecurves are the corresponding predictions F s ( R, Θ ∞ ) = F sim s ( R ) f V T (Θ ∞ ), with F sim s ( R ) = 4 πR γ AB . Θ ( R ) ( d e g . ) ε a =0 ε a =0.05 ε a =0.1(a) Θ ( R ) ( d e g . ) Θ ( R ) ( d e g . ) τ ( R ) σ / k B T τ ( R ) σ / k B T ε a =0.05 ε a =0.1 ε a =0.1 ε a =0.05(b) FIG. 9: (a) The lines are plots of the contact angle Θ, versus 1 /R , resulting from Eq. (3). For the case ε a = 0, the line tensionestimate (for planar interface) [25] τ σ/k B T = − .
52 is used, while in the other two cases the extrapolated results from part(b) of the figure are used to draw the shown straight lines. Data shown by symbols are the values obtained from the analysisof the MC simulation results. (b) Simulation estimates for contact angles Θ( R ) (left) and line tension τ ( R ) (right) are plottedvs. 1 /R , for ε a = 0 .
05 (upper frames) and ε a = 0 . τ and Θ at R = ∞ from (tentative)possible linear fits. All results are from symmetric binary Lennard-Jones mixture. δµ /k B T10100 ∆ F * [ k B T ] - ε a =0, L=24 ε a =0.05, L=30 ε a =0.1, L=30 L z =120.2CNT ( ε a =0.1) CNT ( ε a =0)F s /3( ε a =0) h e t FIG. 10: Plots of nucleation barrier as a function of chemical potential, for the symmetric binary Lennard-Jones fluid confinedbetween antisymmetric walls. Simulation results (symbols) for a few different values of ε a are included. For ε a = 0 and 0 . ε a = 0 we have included a plot for F s / χ AB Θ χ AS =-0.10 χ AS =-0.15 χ AS =-0.19 Θ χ AS =-0.10 χ AS =-0.15 χ AS =-0.19 FIG. 11: a) Contact angle as a function of the (Flory Huggins) parameter χ AB for the mean-field model, for three values ofthe wall attraction strength χ AS . b) Contact angle, computed from Eq. (41), is plotted vs. 1 /R , for χ AB = 2 . χ AS as shown in part (a). [1] Volmer M 1939 Kinetic der Phasenbildung (Dresden: Th. Steinkopff)[2] Zettlemoyer AC (ed.) 1969
Nucleation (New York: Dekker)[3] Kashchiev D 2000
Nucleation: Basic Theory with Applications (Oxford Butterworth-Heinemann)[4] Kelton KF and Greer AI, 2009
Nucleation (Oxford: Pergamon)[5] Dietrich S 1988
Phase Transitions and Critical Phenomena
Vol XII, ed C. Domb and J.L. Lebowitz (New York: Academic)p1[6] Jamali V, Biggers EG, van der Schoot P and Pasquali M 2017 Langmuir Introduction to Microfluidics (Oxford: Oxford University)[10] Natelson D 2015
Nanostrucutres and Nanotechnology (Singapore: World Scientific)[11] Gibbs JW 1961
The Scientific Papers, Vol 1 (New York: Dover)[12] Rowlinson R and Widom B 1982
Molecular Theory of Capillarity (Oxford: Oxford University)[13] Amirfazli A and Neumann AW 2004 Adv. Colloid Interface Sci. Wetting Behavior of the Liquid-Vapor Phase Coexistence of a Colloid-Polymer Mixture (Mainz: JohannesGutenberg Universit¨at, Diplomarbeit unpublished)[35] Borovka L and Neumann AW 1977 J. Chem. Phys. A Guide to Monte Carlo Simulation in Statistical Physics, 4 th ed. (Cambridge University)[49] Schmitz F, Virnau P and Binder K 2013 Phys. Rev. E Polymers at Interfaces (London Chapmannand Hall)[56] Ibagon I, Bier M, and Dietrich S 2016 J. Phys.: Condens. Matter28