Heterogeneous Nucleation of a Droplet Pinned at a Chemically Inhomogeneous Substrate: A Simulation Study of the Two-dimensional Ising Case
Marta L. Trobo Ezequiel V. Albano Kurt Binder












  Heterogeneous nucleation is studied by Monte Carlo simulations and
  phenomenological theory, using the two-dimensional lattice gas model
  with suitable boundary fields. A chemical inhomogeneity of length
  b at one boundary favors the liquid phase, while elsewhere the
  vapor is favored.  Switching on the bulk field H_b favoring the
  liquid, nucleation and growth of the liquid phase starting from
  the region of the chemical inhomogeneity is analyzed. Three regimes
  occur: for small fields, H_b < H_b^crit, the critical droplet
  radius is so large that a critical droplet having the
  contact angle θ_c required  by  Young's equation in the
  region of the chemical inhomogeneity
  does not  yet  “fit” there, since the baseline length of the
  circle-cut sphere droplet would exceed b.
  For H_b^crit < H_b <H_b^∗, such droplets fit inside
  the inhomogeneity, and are indeed
  found in simulations with large enough observation times, but
  these droplets remain pinned to the chemical inhomogeneity when
  their baseline has grown to the length b. Assuming that these
  pinned droplets have a circle cut shape and effective contact
  angles θ_eff in the regime θ_c < θ_eff < π/2,
  the density excess due to these droplets can be predicted, and is
  found to be in reasonable agreement with the simulation results.
  On general grounds one can predict that the effective contact
  angle θ_eff as well as the excess density of the droplets,
  scaled by b^2, are functions of the product b H_b, but do
  not depend on both variables separately. Since the free energy
  barrier for the “depinning” of the droplet (i.e., growth of
  θ_eff to π - θ_c ) vanishes when θ_eff
  approaches  π/2, in practice  only angles θ_eff up
  to about θ_eff^max≃ 70^∘  were observed.
  For larger fields (H_b > H_b^*) the droplets nucleated at the chemical inhomogeneity
  grow to the full system size. While the relaxation time for the
  growth scales as τ_G  ∝ H_b^-1, the nucleation
  time τ_N scales as lnτ_N ∝ H_b^-1. However,
  the prefactor in the latter relation, as evaluated for our
  simulations results, is not in accord with an extension of the
  Volmer-Turnbull theory to two-dimensions, when the theoretical
  contact angle θ_c is used.


Heterogeneous Nucleation of a Droplet Pinned at a Chemically Inhomogeneous Substrate:

A Simulation Study of the Two-dimensional Ising Case
    
        
        
         
               Marta L. Trobo 

               Instituto de Física de Líquidos y Sistemas Biológicos (IFLYSIB). CCT-CONICET La Plata,

                UNLP. Calle 59 Nro. 789, (1900) La Plata, Argentina and

                 Departamento de Ciencias Básicas, Facultad de Ingeniería,
                   Universidad Nacional de La Plata, Argentina.

                                      Ezequiel V. Albano 

                               Instituto de Física de Líquidos y Sistemas Biológicos (IFLYSIB). CCT-CONICET La Plata,

                UNLP. Calle 59 Nro. 789, (1900) La Plata, Argentina and

                 Departamento de Física, Facultad de Ciencias  Exactas,
                                  Universidad Nacional de La Plata, Argentina.

                       Kurt Binder

                                 Institut für Physik, Johannes Gutenberg-Universität Mainz

                 Staudinger Weg 7, D-55099 Mainz, Germany,
    December 30, 2023
==========================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================



§ INTRODUCTION


When thermodynamic variables such as temperature T, pressure p, or
external fields (e.g. a magnetic field H) are varied, discontinuous 
changes in the state of matter can occur. Examples of such so called 
first-order phase transitions <cit.> range from the condensation 
of water, melting of ice, crystal formation in solidifying melts, etc., 
to the magnetization reversal of ferromagnetic devices. These phenomena 
are of great importance not only for condensed matter physics, but also 
for the atmospheric sciences, geosciences and material science, as well 
as for numerous technical applications. However, a common feature of 
all these phase changes is that they are triggered by nucleation events, 
i.e. on the background of the old (metastable) phase a nanoscopically small 
nucleus of the new phase needs to be formed, and such nucleation phenomena 
are rare events since a high free energy barrier needs to be 
crossed <cit.>. Actually, for most conditions of practical 
interest the spontaneous formation of nuclei by statistical fluctuations, 
i.e the so called "homogeneous nucleation", involves too high barriers and 
cannot happen. In contrast, "heterogeneous nucleation" at defects, 
e.g. condensation of water droplets on dust particles in the atmosphere, 
or surface-induced crystallization starting at the walls of a container, etc., 
occurs much more frequently.
Also, processes such as the formation of dew droplets on 
car windows or plants are familiar from everyday life <cit.>.
However, the nanoscopic size of the nucleus, which 
typically contains only a few hundred of atoms or molecules, is a stumbling
block already for the theoretical description of homogeneous 
nucleation <cit.>. The large variety of defects that can 
cause heterogeneous nucleation makes a comprehensive description even more 
difficult, see e.g. <cit.>. Thus theoretical work on heterogeneous
nucleation is rather 
scarce <cit.>, 
while theoretical work on homogeneous nucleation is abundant, see 
e.g.<cit.> for an overview of work done for Ising/lattice gas models.

In this paper, we reconsider the problem of heterogeneous nucleation on
flat substrates by focusing on a chemically inhomogeneous surface, where 
nucleation preferentially occurs in a region of finite (nanoscopic) linear dimension
b, c.f. Figure  <ref>. There are several motivations for such a choice: 
i) chemically structured surfaces are useful for many applications in
nanotechnology such as the fabrication of nanodevices, for the processing 
of nanoscopic amounts of matter ("lab in a chip"), etc. <cit.>. 
ii) Since b is comparable to the size of the nucleus formed in a single 
heterogeneous nucleation event, it is straightforward to study the
characteristics of such isolated nucleation events by computer simulation in 
great detail. In contrast, when one studies nucleation in a system with a 
homogeneous macroscopic surface, surface-attached nuclei can occur anywhere on 
this substrate, and one easily reaches conditions where several nuclei form 
and compete during their growth. This case is familiar from studies of
homogeneous nucleation (e.g. <cit.>) and observing 
the lifetime of metastable states in this limit of multi-nuclei formation 
and growth allows only rather indirect conclusions on the nucleation rates 
and nucleation barriers. 
iii) Due to the fact that dust or soot particles at which nucleation happens 
in the atmosphere are often of μ m size only and have irregular shapes and 
need not be chemically homogeneous, it is of practical interest to study cases 
where conditions favorable for heterogeneous nucleation are limited 
to regions of nanoscopic extent. Of course, also other types of localized
defects are suitable to study isolated heterogeneous nucleation events, 
e.g. in square lattices with free boundary conditions nucleation starts 
at the corners of the square <cit.>.
 
Of course, for this problem there occurs a challenging interplay between 
surface effects due to the substrate and interfacial effects of the material 
forming the nucleus, also statistical fluctuations and finite-size effects 
play a role. Therefore, the development of an analytical theory for the 
treatment of such problems is very difficult <cit.>. 
In the present work, we hence restrict our attention to an approach by 
Monte Carlo simulation <cit.> of a simple model, namely the Ising/lattice 
gas model. As has been discussed elsewhere<cit.>, even this simple model 
presents severe difficulties due to the incomplete knowledge of the 
anisotropy of the interfacial tension between bulk coexisting phases. 
Also understanding of wetting phenomena, e.g. 
contact angles of macroscopic sessile droplets at walls <cit.>, 
is a problem in the three-dimensional case (see e.g. <cit.>). 
Thus, we focus here on the lattice gas/Ising model in d = 2 dimensions, 
for which both bulk and interfacial phenomena including wetting 
behavior <cit.> are well understood. 
Recently, we have already studied droplets pinned at chemically heterogeneous 
substrates under conditions of bulk phase coexistence <cit.>. 
This knowledge also is useful for the present work where we consider such
droplets under out-of-equilibrium conditions, and the dynamics of the
resulting growth process during nucleation events.

The outline of our paper is as follows:
in Section II, we precisely characterize the model and the simulation method, 
in Section III, we present our results for wall-attached precursors droplets, 
which are in metastable equilibrium during the "observation time" of the
simulation. Furthermore, we discuss a scaling description in terms of the 
variables b (spatial extent of the inhomogeneity) and bulk magnetic 
field H_b (which characterizes the "distance" of the metastable state from 
phase coexistence that occurs in the bulk at H_b = 0, of course). 
In Section IV, we present a phenomenological theory of pinned metastable droplets,
which have a baseline b and a non-equilibrium contact angle controlled by the
bulk magnetic field. We discuss the stability limit of these droplets, where they depin
from the chemical inhomogeneity and grow beyond it, causing a fast phase transformation.
In Section V, we analyze the dynamics of nucleation events, characterizing 
both the growth process of a single nucleus from nanoscopic to macroscopic
sizes, and the distribution of nucleation times. Also, we compare these
results to previous findings for single-droplet nucleation in the 
bulk <cit.>. 
Finally, Section VI summarizes our conclusions. 
The extension of the classical theory of heterogeneous nucleation 
of Volmer and Turnbull <cit.> to the two-dimensional case is given in Appendix A.
In Appendix B the depinning of droplets from a chemical heterogeneity in d = 3 dimensions is
briefly discussed.



§ MODEL AND SIMULATION DETAILS


The chosen model is similar to our previous work <cit.> where a two-dimensional
Ising/lattice gas model on the square lattice in L × M geometry was
studied at phase coexistence (bulk field H_b =0). We apply periodic boundary
conditions in x-direction only, while free boundaries are used in
y-direction, and the Ising spins in the first (i = 1) and last
(i = L) rows experience boundary fields, c.f. Figure  <ref>.

At the upper boundary, we choose a 
homogeneous boundary field H_w1 = -0.225 throughout, which hence favors
the minus spins, S(i,j) = -1; in lattice gas language, the local density
variable ρ(i,j) = (1 + S(i,j))/2 at a lattice site with coordinates 
i (i=1,...,L) in y-direction and j (j=1,...,M) in x-direction then is
ρ(i,j) = 0. At the lower boundary, i=1, we choose the boundary field 
inhomogeneous: from the site j = (M-b)/2 + 1 to the site 
j = (M+b)/2 we choose a positive boundary field H_w3 = +0.90, to favor
along a line of length b (b is a large odd integer) the liquid phase.
Outside this region, the boundary field in the first row is H_w2 = -0.90
throughout, to ensure for the chosen total linear dimensions that in the
absence of a bulk field H_b the vapor phase (or phase with negative
magnetization, respectively) is the stable phase, for all choices of b 
that were considered.


The temperature T is measured in units of the critical temperature T_cb of the bulk; 
i.e. <cit.> k_BT_cb/J=2/ln (1+ √(2)) ≈ 2.27, where k_B is Boltzmann's constant, 
and J the exchange constant, respectively. Boundary and bulk fields H_w1, H_w2, H_w3,
H_b are measured in units of J. Thus, the used Hamiltonian is


    ℋ =-J ∑_i=1^L ∑_j=1^M S (i,j) [S(i+1,j) + S(i-1,j) + S(i,j +1 ) + S(i, j-1 ) ]/2 
          -H_w1∑_j=1^M S (L,j ) - ∑^M_j=1 H_w (j ) S(1,j ) -
    H_b ∑_i=1^L ∑_j=1^M S (i,j)    ;    S (i,j ) = ±1   ,


where H_w(j) = H_w2 for 1 ≤ j ≤ (M-b)/2 and (M+b)/2 + 1 ≤ j ≤ M, while H_w (j)=H_w3 
for (M-b)/2 + 1 ≤ j ≤ (M +b)/2. Also, S(i,j)=0 is taken for missing neighbors.

The chosen value of H_w3=0.90 leads to a wetting critical temperature 
t_w=T_w/T_cb≃ 0.4866 <cit.>. Furthermore heterogeneous nucleation 
is studied for T < T_w, so that the order within the bulk domains 
is almost perfect, and the correlation length in the bulk is of the order of the lattice
spacing. For all temperatures T < T_w we have a nonzero contact angle θ_c < 90^∘
in the region where H_w3 acts, while due to the antisymmetric choice H_w2 = -H_w3
the contact angle is π -θ_c in the region where H_w2 acts.

For all simulations we choose M = 453 and L = 300 to make sure that there are no finite size effects 
associated with interfacial fluctuations. Also, in order to make it easier for
the reader to establish the connection to fluid droplets, we will describe all 
our results in terms of local densities defined via ρ(i,j)=(⟨ S(i,j) ⟩ +1 )/2.
Monte Carlo simulations were carried out with standard single spin-flip
algorithms <cit.>, which as is well known <cit.>, can be
interpreted as a simulation of a stochastic process, where (in a lattice gas
terminology) particles are randomly adsorbed at, or desorbed from, the sites
of the lattice; the rates of these processes satisfy the detailed balance 
condition with the Hamiltonian, Equation (<ref>). Note that the density in the
considered lattice gas is not a conserved variable, of course; hence the
physical situation that is simulated is a two-dimensional substrate in
equilibrium with an ideal gas reservoir, at the specified temperature 
and chemical potential, corresponding to the chosen bulk magnetic field.
Of course, the Monte Carlo process has no intrinsic time units for the rates 
of Monte Carlo moves: so time is measured in units of "attempted Monte Carlo
steps per spin [MCS]" rather than any physical time units.
Simulations are started with an initial condition where all spins are taken 
as S(i,j)=-1, i.e.  ρ(i,j) = 0 in accord with the nonwet ground state 
of the system. Then, the system is equilibrated up to  5× 10^6 
MCS by taking H_b = 0, and subsequently is “quenched” to 
H_b > 0 in order to observe nucleation. Runs are performed for additional 2× 10^7 
MCS, and averages are taken after disregarding 1×
10^7 MCS. 



§ WALL ATTACHED PRECURSOR DROPLETS IN METASTABLE EQUILIBRIUM.


In our system, at least in the limit M →∞, L →∞,
the vapor phase is the true equilibrium phase only for H_b ≤ 0. But when we consider 
a situation where equilibrium in the bulk has been established for H_b = 0,
and at the time of the Monte Carlo sampling that we take as the origin of time (τ = 0) the field
H_b is instantaneously switched to a small positive value, the vapor phase in the bulk 
may reach a state of metastable equilibrium, with a "lifetime" larger than the
observation time τ_obs of the simulation. This metastability can be understood 
qualitatively already by the classical theory of heterogeneous nucleation due to 
Volmer and Turnbull <cit.>, see the Appendix A.

In this region where the time τ_N to observe nucleation satisfies the condition 
τ_N ≫τ_obs, one may observe the formation of wall-attached precursor droplets
in the part of the sample where the surface field H_w3 acts, Figure <ref>.
Note that we restrict attention to temperatures T distinctly lower than
T_w, i.e. when partial wetting of the wall occurs. To avoid confusion, we stress that
we denote the boundaries at i = 1 and i = L in our system as "walls" although they are
one-dimensional lines only, Figure 1. If we would choose T > T_w, i.e. the case of 
complete wetting of the wall for b →∞, there would no longer   
occur any nucleation barrier, the liquid wetting layer that occurs then at the wall 
already for H_b = 0 would immediately grow by increasing its thickness as soon as  H_b > 0.   
For finite b the interface of the liquid droplet still is pinned to the walls
near the sites where the boundary field changes from H_w3 to H_w2.
In our previous work <cit.> we have verified a prediction of 
Jakubczyk et al. <cit.> for the case H_b = 0 and T > T_w,
based on the solid-on-solid (SOS) model in terms
of the interface Hamiltonian approach for the excess density Δρ
due to the droplet, given by 


    Δρ=(ρ^ coex_ℓ - ρ^ coex_v ) 1/4 b^3/2√(π/2 Σ (T))   , 
       b →∞   .


In Equation (<ref>), ρ^ coex_ℓ and ρ^ coex_v are the coexisting
liquid and vapor densities in the bulk, respectively. Note that
ρ^ coex_ℓ - ρ^ coex_v = m^ coex, the spontaneous magnetization of the
Ising model <cit.>.
Also, Σ (T) is the interfacial 
stiffness of the Ising model <cit.>.  Equation (<ref>) implies that the average droplet shape
is a semi-ellipse, with small axis proportional to b^1/2 for b →∞. For large b  even
this small axis can easily exceed the critical droplet radius R^* for homogeneous nucleation in 
the bulk (see Appendix A), and thus it is plausible that for
large but finite b the nucleation barrier will be very 
small.  Although a study of the phase transformation for  T > T_w and not so large b may be 
interesting in its own right, we here fixed the attention to  T < T_w, where the linear dimension
in the y-direction of
the wall-attached droplet for H_b = 0 remains finite for b →∞. We have also 
shown <cit.> that the density in the mid-point j_half = (M+1)/2 of the inhomogeneity of the 
boundary is compatible with an exponential decay with the distance y from the wall (y is only 
defined at the discrete lattice indices i, of course),   


    ρ (i, j_half) = A_0(t) exp (-y/ξ_ (b,t)) ,    H_b = 0    .


Note that in the context of various theoretical concepts on interfaces
  in this constrained geometry the use of a continuum description (in terms of coordinates
  x, y) is mandatory; but for a precise characterization of the simulation
  setup on the lattice, discrete indices of the lattice points (j in x-direction,
  i in y-direction) have also to be used.
Here, A_0(t) is an amplitude factor, t = T/T_cb must be less than
t_w = T_w/T_cb, and the decay length  ξ_ (b,t) converges to the standard transverse
correlation length ξ_ (t) familiar from the theory of critical wetting <cit.> 
in d = 2 dimensions, with ξ_ (t) ∝ (t_w - t)^-1. In order to avoid 
critical fluctuations associated with the second-order wetting transition, we consider
here only temperatures distinctly lower than t_w, namely t = 0.30, 0.35, and 0.40   
(recall that t_w = 0.4866 for our choice of h_w3 <cit.>).



When we now consider metastable phases with small but nonzero H_b > 0, we find that
the actual droplet configurations are strongly fluctuating (Figure <ref>) even though the 
average profile ρ(i, j = j_half) is a smooth function (Figure <ref>(a)). We note that
for small enough H_b (such as H_b ≤ 0.02) the decay of the profile with the distance 
y from the boundary still is compatible with Equation (<ref>), and both the amplitude 
A_0(t) and the decay length  ξ_ (b,t) increases with H_b only very slowly.
However, when we reach the apparent limit of metastability H_b^max, where for the 
chosen observation time τ_obs = 2 × 10^7 MCS nucleation becomes observable,
which in the case of Figure <ref>(a) is H_b^max≈ 0.033, the behavior changes:
In the regime 0.02 ≤ H_b≤ 0.032 the profiles change from a simple
exponential decay to sigmoidal, and A_0(t)  moves towards ρ^ coex_ℓ, which
on the scale of Figure <ref>(a) is indistinguishable from unity, while ρ^ coex_v is not 
distinguishable from zero here. Thus we have defined an effective droplet height h_eff(t,H_b)
by measuring the distance y of the inflection point from the coordinate origin (Figure <ref>(b)). 
The increase of h_eff(t,H_b) with H_b is clearly faster than a straight line through 
the origin. An intriguing feature is the fact that the ratio h_eff(t,H_b) / b
seems to be a function of the scaling combination b H_b only (Figure <ref>(c)).



In order to characterize the droplet shape more precisely, also density profiles 
ρ(i, j) in the x-direction parallel to the boundary were taken (Figure <ref>). It is
clear that in the shown example the droplet has considerable extent in the x-direction,
comparable to b, as long as the distance i from the boundary is clearly less than 
h_eff. Thus we have introduced a characteristic length of the wall-attached droplet in the 
x-direction, defining its "baseline" b_eff as the area below the density profile for i = 1,
namely 

    b_eff = ∑_j = (M-b)/2 + 1^(M +b)/2ρ(i = 1, j)   .




Figure <ref>(a) shows a plot of b_eff/b for three choices of b as a function of
H_b. One can see that for  H_b ≪ H_b^max  b_eff/b is small, but saturates
at unity when H_b approaches H_b^max. Again the data for the three choices of 
b almost superimpose on a master curve when one plots b_eff/b  as a function 
of the product bH_b (Figure <ref>(b)); we shall discuss this scaling
behavior of the characteristic
lengths h_eff and b_eff below.



When the ratio  b_eff/b is no longer small, one can still use density profiles such as
those shown in Figure <ref> to evaluate the effective length of the droplet parallel
to the boundary but evaluated at a distance i > 1 from the wall (e.g. by performing 
the corresponding summation of Equation (<ref>)) and characterize the average 
shape of the droplet (note that this procedure is equivalent to the construction of contours 
of constant density), see Figure <ref>.  The slope of these contours can be used to extract estimates
for the effective contact angle θ_eff(H_b) from the part of the contours at 
small i values, where these contours representing the coarse grained interface positions
reach the wall. Figure <ref> exploits this idea for the case t = 0.40, b = 51, and several
choices of H_b. By fitting straight lines to these contours in the region close to the wall
(i.e. for i = 1, 2, 3), estimates of θ_eff(H_b) can be extracted (Figure <ref>(a)). It is 
seen that for choices of H_b for which b_eff (H_b) is distinctly smaller than b,     
θ_eff(H_b) is essentially independent of H_b, and of the order of  
θ_eff≈ 10 ± 2^∘ in this case. However, when b_eff (H_b) starts to saturate
at b, the ratio h_eff / b as well as the effective contact angle  θ_eff both start
to increase rather distinctly. As Figure <ref>(b) demonstrates, this increase of θ_eff 
starts at b H_b ≈ 1. Metastable precursor droplets are found up to angles of 
θ_eff^max≈ 70^∘ when H_b reaches H_b^max. 





A related conclusion can be drawn with even less ambiguity, since it does not require 
an analysis of the shape of the precursor droplet, when we simply record the excess density
Δρ in the system due to the droplet (Figure <ref>). In fact, Figure <ref>(a)
shows plots of Δρ versus H_b at three choices of t for b = 51,
while Figure <ref>(b) shows plots of Δρ versus H_b for t = 0.40 and
various choices of b, as indicated. Note that in both cases the states where the 
transition to the liquid in the simulated system has occurred are included, namely 
all the data points with  Δρ = 1 implies that the whole simulation box 
is filled uniformly by a liquid. The excess density of the precursor droplet in these
plots is normalized by dividing the "excess mass" Δ m contained in the precursor 
droplet by the total number of lattice sites, N = L × M = 300 × 453 = 135900.
Recall that the excess density is the difference between the
density measured in the presence of the inhomogeneity in the wall and that obtained for an homogeneous wall.
Figure <ref>(a) implies that the value of Δ m  where the precursor droplet is large
enough to trigger nucleation of the liquid phase in the system depends only weakly on temperature;
it is Δ m ≈ 400 for b = 51, and the choice H_w2 = 0.90 that was made here.  
But again there is a pronounced dependence of Δρ (or Δ m, respectively) 
at the transition on the choice of b (Figure <ref>(b)). 



At this point it is interesting to make contact with the classical theory of heterogeneous
nucleation as formulated in d = 3 dimensions by Volmer and Turnbull <cit.>.
A simple adaptation of this theory to d = 2 (see Appendix A) implies that a critical droplet 
on a homogeneous substrate with contact angle  θ_c involves an excess density mass
given by

    Δ m = (ρ^ coex_ℓ - ρ^ coex_v ) (R^*)^2π f_VT(θ_c)     .

Here, R^* is the critical radius of homogeneous nucleation, which in the classical
theory simply is 

    R^* = f_int/2 ( ρ^ coex_ℓ - ρ^ coex_v ) H_b   ,
 
with f_int being the interfacial tension between coexisting bulk phases
separated by a flat interface. As discussed in the Appendix A, Equations (<ref>, <ref>)
ignore the anisotropy of the interfacial tension, assuming a circular droplet in the 
bulk, and a circle cut shape of the droplet at the boundary, where the coarse-grained interface 
makes an angle θ_c. The Volmer-Turnbull function f_VT(θ_c) measures the reduction
of the droplet area of the circle cut relative to the full circle, and in this approximation (see Appendix A)
is given by 

    f_VT(θ) = 1/π( θ - sin(2θ)/2)    .

    
The upper curve in Figure <ref>(b) shows Δρ = Δ m /LM as
obtained according to 
Equations (<ref>-<ref>) by using the "observed" contact angle
θ_c = θ_eff^max≈ 70^∘ just at the
transition to the liquid phase. However, by using a theoretical estimate of
the contact angle in thermal equilibrium obtained by means of a SOS  approximation θ_c 
due to Abraham et al. <cit.> would imply θ_c≈ 34^∘, predicting hence
distinctly smaller critical droplets (Figure <ref>(b)). However, for both choices of  
θ_c our results (Figures <ref>, <ref>) imply that b_eff is of the same order as b
when nucleation occurs. Since the baseline length b_drop^* of the circle cut critical droplet 
satisfies the geometrical relationship

    b_drop^* = 2 R^* sin(θ_c)    ,
 
we can eliminate R^* in Equations (<ref>, <ref>) in
favor of b_drop^* = b, and in this way the two theoretical curves in Figure <ref>(b) were
obtained. Figure <ref>(a) shows then the scaling plot of  Δρ /b^2 versus
b H_b, validating the idea that  Δ m scales like b^2 and is a function of 
the product b H_b, as b_eff and h_eff (Figures <ref>(c) and <ref>(b)).  
Both figures  <ref> and <ref> show that two regimes exist: for bH_b ≲ 1,
θ_eff is very small, l_eff < b (Figure <ref>), Δρ increases
slowly with bH_b, and the droplet density profile (Figure <ref>(a) ) decays with the
distance from the inhomogeneity like an exponential. All these properties change for
bH_b ≳ 1. Note that  <ref>(a) includes both the regime
bH_b ≲ 1, where scaling is not expected hold, as well as the true
scaling regime bH_b ≳ 1. Of course, for H_b = 0 in the nonwet regime the
excess density Δρ only scales like Δρ∝ b; thus in the
regime of small bH_b a crossover from Δρ∝ b to Δρ∝ b^2
with increasing bH_b must occur. 




In order to give a physical interpretation of the simple scaling of all quantities (b_eff,
h_eff, and Δρ) with the product b H_b, we recall the description of droplets at
chemical inhomogeneous substrates in d = 2 dimensions in terms of the interface Hamiltonian
proposed by Jakubczyk et al.<cit.>. In this description in the spirit of a Solid on
Solid (SOS) model, the problem is described by a one-dimensional degree of freedom, namely 
the distance  y = ℓ (x) of the (locally sharp) interface from the boundary at 
y = 0 (involving a continuum approximation). So the effective (coarse-grained) Hamiltonian
is, absorbing a factor 1/k_BT here, 

    ℋ[ℓ(x) ] =∫^+M/2_-M/2  dx [Σ (T)/2(d ℓ/dx)^2 + V(x, ℓ) ] ,
 
where both fluctuations in the bulk and overhangs of the interface are neglected,
V(x, ℓ) being the effective potential acting on the interface. Recall
  that in the SOS treatment the interfacial stiffness Σ (T) <cit.> of the
  one dimensional interface is considered, instead of the actual interfacial
  tension f_int of the Ising model <cit.>. For the considered situation,
for |x| > b/2 the boundary at y = 0 strongly favors the vapor phase, so we have
essentially  ℓ (x) = 0 there, as one can verify from Figure <ref>. So, Equation (<ref>) 
can be reduced to   

    ℋ[ℓ(x) ] =- ∫^+b/2_-b/2  dx [Σ (T)/2(d ℓ/dx)^2 + V(x, ℓ) ] .

Note that Equations (<ref>, <ref>) also assume |d ℓ/dx| ≪ 1
everywhere, an assumption that is somewhat questionable in view of the actual snapshots of the
interfacial configurations (Figure <ref>), at least near x = ±b/2;
however since no actual calculations on the basis of Equations (<ref>, <ref>) are done here,
this problem does not matter.

Now the effective potential  V(x, ℓ) can be written as 

    V(x, ℓ) = [V_0(x, ℓ) - ℓ(ρ^ coex_ℓ - ρ^ coex_v ) H_b] /k_BT    ,

where V_0(x, ℓ) is the potential binding the interface to the wall for H_b = 0.
Only this latter case has been considered by Jakubczyk et al.<cit.>. Applying a field 
H_b > 0 favors the liquid phase, and thus the potential decreases proportional to ℓ H_b.

Now the key observation is that the dependence on b is elucidated when we rescale all 
distances by b, namely 

    x = b x '   ,   ℓ = b ℓ '    ,
    
which yields ℋ = bℋ' with

    ℋ' [ℓ'(x')  ] =- ∫^+1_-1  dx' [Σ (T)/2(d ℓ'/dx')^2 + V(x', ℓ') ] ,

with 

    V(x', ℓ') = V_0(x', ℓ') - ℓ' (ρ^ coex_ℓ - ρ^ coex_v ) b H_b    .

Now the partition function needs to be evaluated as a path integral,

    𝒵 =  ∫𝒟ℓ' exp(-bℋ' )      ,

and from Equations (<ref>)-(<ref>) we conclude that the boundary excess
free energy Δ F = -k_B T ln (𝒵) due to the droplet depends on the variables b, t,
and H_b in the following scaled form, f(t, bH_b) being the free energy density per length unit along 
the boundary 

    Δ F =  b f(t, bH_b)     .


The excess density due to the droplet is obtained from Equation (<ref>) via a 
derivative with respect to H_b, i.e.

    Δρ =  b^2M̃(t, b H_b)     ,

where M̃ is the resulting scaling function of the excess mass. Equation (<ref>)
hence justifies the choice of scaling variables for Figure <ref>(a).

This scaling property is subtle, of course, due to the requirement of metastable equilibrium:
it is implied also by Equation(<ref>), that there cannot be for H_b > 0 a true equilibrium 
at any finite value of ℓ, so Equation (<ref>) makes sense only for a suitably constrained 
partition function.

So the droplets studied so far can only be found in a suitable "window" of observation times
τ_obs. In fact, τ_obs must be large enough to allow that the wall-attached droplet
reaches local equilibrium in spite of the slow and sluggish  fluctuations of the interface 
configuration (Figure <ref>). But at the same time, τ_obs must be small enough that
nucleation events (where the droplet grows fast to the full size of the system, see Section V)
still are negligible.

This consideration is exemplified in Figure <ref>(b): here a log-log
plot of the excess density in the system versus H_b is shown, for four
choices of  τ_obs. For very small H_b, such as H_b = 0.0025 the
dependence on  τ_obs is negligible since the time  τ_N needed to
nucleate is astronomically large, and the wall attached droplet is very tightly
bound to the wall (cf. Figure <ref>(a)), so it is rather easily
equilibrated.  However, for  H_b = 0.01 we see that data for
 τ_obs = 2 × 10^5 and τ_obs = 2 × 10^6 perfectly
agree, nucleation is not yet possible; but the result for 
τ_obs = 2 × 10^4 is clearly smaller, this observation time was 
insufficient to sample fluctuations such as those seen in Figure <ref> 
exhaustively. For  H_b = 0.018 however, there is also a systematic
difference between  τ_obs = 2 × 10^5 and τ_obs = 2 ×
10^6: for the latter time, nucleation typically has occurred, while for the
former time, the metastable boundary-attached droplet still is visible.

Since Figure <ref> suggests that we can (for the choices of b used
here) observe metastable boundary attached droplets up to b_eff≈ b,
we have also tested as a possible hypothesis that these boundary-attached
droplets with  b_eff = b can be described by the Volmer-Turnbull theory of
heterogeneous nucleation (see the Appendix A). This theory assumes that the
critical droplet causing nucleation is a cut from a sphere (circle in our d =
2 case) with radius R^* [Equation (<ref>)], the angle of the
sphere cut with the straight line representing the boundary being 
the contact angle θ_c. Geometry then implies b_eff = b as quoted in 
Equation (<ref>), and combining Equations (<ref>), (<ref>)
yields a relationship between b and H_b

    b_drop^* = b = f_int sin(θ_c) / (ρ^ coex_ℓ - ρ^ coex_v )  H_b^crit   ,

where  H_b^crit is the prediction of standard theory of
  heterogeneous nucleation for the
  critical field at the onset of nucleation of the liquid phase when
  the length of the baseline of the droplet is b_drop^* = b.
This result is plotted in Figure <ref>, using for f_int
the Onsager result <cit.> for the interface tension of a straight
interface oriented perpendicular to the lattice axis. For the contact angle    
θ_c we use results derived for Abraham et al. <cit.> in the SOS
approximation, namely 

    tan [θ_c(T_w,T)] = sinh [2 (K - K_w ] / [ cosh(K) - cosh[2(K -
    K_w)] ]     ,

where K = J/k_B T and  K_w = J/k_B T_w, respectively.



For our choice of H_w3 implying T_w = 0.4866 T_cb, Equation
(<ref>) yields θ_c ≈ 67,11^∘, 50,05^∘ and
34,11^∘, for t = T/T_cb = 0.30, 0.35 and 0.40, respectively. 
The resulting straight lines in the log-log plot 
for b versus H_b [Equation (<ref>)] are compared 
to the estimates for the actual critical field H_b^max where the 
onset of nucleation of the liquid phase has been observed
in the simulations [Figure <ref>]. 
For each temperature, these lines correspond to the theoretical 
conditions where nucleation on the inhomogeneity of length b becomes possible.
The actual data points included in Figure <ref> separate the region of parameters where 
metastable precursor droplets are found, and the liquid phase that 
already takes the whole simulation box.   
It is seen that the actual critical fields, measured for t = 0.40 
always are larger than the predictions
based on Equations (<ref>), (<ref>). So, if we would fit the numerical data
to equation (<ref>) we would obtain
the  "observed" contact angle θ_c^max = θ_eff^max≈ 67.5^∘, in
excellent agreement with our previous estimations, e.g. showing that 
metastable precursor droplets are found up to angles of 
θ_eff^max≈ 70^∘ (Figure <ref>(b)), as well
from the scaling plot of Figure <ref> that also 
yields θ_eff^max≈ 70^∘. The main reason for the difference between the
field H_b^crit, predicted by the standard theory for heterogeneous nucleation and
  defined from equation (<ref>) and calculated by using the contact angle θ_c
obtained by means of the SOS approximation (equation (<ref>)), and the actual critical
field found in the simulations H_b^max is that in
the regime H_b^crit < H_b < H_b^max the droplets nucleated
at the inhomogeneity are pinned, see Section IV, i.e., their baseline cannot grow
beyond b. However, it should be stressed that within this regime the
area of the droplets actually grows by simultaneously increasing the contact angle and
decreasing their radius. Only for H_b > H_b^max droplets “depin” and further growth is possible,
with b_drop > b and θ = π - θ_c, see also below.
These pinned droplets should not be mistaken for the droplets described by the
standard Volmer-Turnbull theory of heterogeneous nucleation, as discussed in 
the Appendix A.
In Section IV, we shall attempt a theoretical estimation of the field H_b^max. We also
note that for H_b < H_b^crit only subcritical nuclei (R < R^*) can form
on the inhomogeneity, i.e. transient fluctuations occur whose average effect shows
up in the exponentially decaying density profiles for H_b ≤ 0.022
in figure <ref>(a).

Furthermore, it is worth discussing that Equation
(<ref>) is not expected to be quantitatively accurate for several reasons: (i) The
interface tension for a straight interface f_int is used here, neglecting
possible corrections due to the curvature of the droplet interface. (ii) Due
to the anisotropy of the interfacial tension of the Ising lattice model, the
actual shape of a large droplet of the liquid coexisting with surrounding
vapor is not a circle, it rather resembles a square with
rounded corners at low temperatures <cit.> (see also the
largest droplet in Figure <ref>, left-hand side panel). For heterogeneous
nucleation, the droplet shape resulting from the appropriate Winterbottom
construction <cit.> is then nontrivial to find, and the
Volmer-Turnbull theory as presented in the Appendix A needs to be extended to
account for this anisotropy. For not very large droplets, also 
the “point” where the droplet interface meets the boundary can play a role,
modifying Equation(<ref>) further, in analogy with the effect of the
line tension of the sphere-cap shaped droplet
on the contact angle in d = 3 dimensions <cit.>.

In view of all these shortcomings of the existing theories, a 
more quantitative analysis of our numerical
data for the boundary-attached droplets (Figures <ref> - <ref>)
suffers from the incomplete knowledge of both the droplet shape and the equilibrium
contact angle. Nevertheless, we attempt a phenomenological analysis of pinned
droplets and their depinning in the next section.





§ DROPLETS PINNED AT CHEMICAL INHOMOGENEITIES AND THEIR “DEPINNING TRANSITION”.


When one deals with heterogeneous nucleation at homogeneous substrates, the nucleation
barrier Δ F_het^* (Equation (<ref>)) corresponding to the droplet having the
critical radius R^* (equation (<ref>) or (<ref>), respectively) is all what matters:
when such a droplet (of circle cut shape, with contact angle θ_c) corresponding
to the top of the free energy Δ F_drop(R) (equation (<ref>)) has been formed by a
(rare) statistical fluctuation, with 50 % probability this drop will grow with time τ
after the nucleation event. For small fields H_b the growth velocity is small, and then
we have “local equilibrium” of the growing droplet at the contact line; this means, at
any instant of time growing droplets with R > R^* are still described by equation (<ref>),
and, in particular, their contact angle has the equilibrium value θ_c.

However, this description cannot apply when we have a substrate with a chemical inhomogeneity
of extent b (c.f. figure <ref>). We assume here conditions (corresponding to our actual
choice of the boundary fields H_w1, H_w2, and H_w3) where nucleation rates in the region
where wall fields H_w1, H_w2 act are negligibly small; so only nucleation within the region
of the chemical inhomogeneity needs to be considered, i.e. circle-cut shaped droplets with baseline
b_drop = 2 R sin(θ_c) (equation (<ref>)) smaller than b. Such droplets can grow
at constant contact angle with time only until b_drop = b_drop^* = b, and then
either get pinned and grow
in area and angle up to some nontrivial values, which we shall study 
in this section or they “depin” and grow with baseline b_drop > b, b_drop = 2 R sin(π - θ_c).

Thus we turn to an analysis of the regime where b_drop≃ b.
In this regime we have to use Equation (<ref>)
for the area of the droplet, and hence write the free energy of the droplet as

    Δ F_drop = constant +  f_int 2 R θ - 2 m_coex H_b R^2[θ - 1/2 sin(2 θ)]    ,
 
where the constant is fixed by the requirement that for θ = θ_c and b_drop = b
the previous expression for Δ F_drop (equation (<ref>)) results,
i.e. constant = -b f_int cos(θ_c). Note
that now R is not R^* but rather R = b/(2 sin(θ)) from geometry (see
Figure <ref> in Appendix A).
Thus we obtain, eliminating R in favor of b/(2 sin(θ)),

    Δ F_drop /b f_int = -cos(θ_c) + θ/sin(θ) - m_coex H_b b/2 f_int [θ/sin^2(θ) - cos(θ)/sin(θ)]    .
 

Now the angle θ is found from the condition

    ∂/∂θ (Δ F_drop/bf_int) = 0     ,
 
which after simple algebra yields the minimum of the free energy for  

    sin(θ_min) = m_coex H_b b /f_int   ,   θ_c < θ_min < π/2     ,
 
while the angle θ_max = π - θ_min also is a solution of equation (<ref>),
but corresponds to the maximum of the free energy. Note that the condition θ_c < θ_min
has been added since equation (<ref>) makes sense only for θ > θ_c. 



Using this result it is straightforward to evaluate the area A of the droplet and hence 
the excess density Δρ due to the droplet.  Figure <ref> a) presents 
plots of θ versus Δρ /b^2, comparing Equation (<ref>) to the simulation
data in the regime 20^∘  < θ < 70^∘. The agreement
is reasonable, in particular since no adjustable
parameter whatsoever is present. Note that the knowledge of the contact angle θ_c is not needed
here (apart from defining the range on which this relationship should be used).

Figure <ref> b) presents a plot of θ versus b H_b. The dashed 
horizontal line shows the estimate of θ_c from the SOS approximation Equation (<ref>). The 
actual variation of θ with b H_b (Figure <ref>b ) should be compared to this 
figure only for b H_b > 1, since for small fields, where the excess mass due to the droplet is small,
the assumptions of the above quasi-macroscopic analysis clearly are inapplicable. At least, for   
b H_b ≥ 1.5 the prediction is close to the observations from the simulations (Figure <ref> b),
furthermore the values of θ_eff^max obtained for the larger inhomogeneities,
i.e. 51 ≤ b ≤ 101, are in full agreement with the theoretical result given
by equation (<ref>) taken
f_int = 1.7987 for t =  0.40 <cit.> (f_int is taken in units of J).    



It really is illuminating to plot the free energy given by equation (<ref>) as function
of θ for various values of b H_b, as shown in figure <ref>. In fact, both the
minimum and the maximum of the free energy can clearly be observed. Furthermore, the
free energy difference between the maximum and the minimum (Δ F_2 / b f_int)
monotonically decreases when b H_b increases, and finally vanish
where θ_min and θ_max merge at θ_min = θ_max = π/2.
A more quantitative evaluation of the free energy barrier can be performed  
by reinserting Equation (<ref>) in the free energy  Δ F_drop (equation (<ref>)),
obtaining the following relationships for the minimum

    Δ F_drop^min /b f_int = -cos(θ_c) + 1/2[θ_min/sin(θ_min) + cos(θ_min)]    ,

and the maximum

    Δ F_drop^max /b f_int = -cos(θ_c) + 1/2[π - θ_min/sin(θ_min) - cos(θ_min)]    ,

respectively. Then, the difference between the maximum and the minimum is given by

    Δ F_2 / b f_int = (π/2 - θ_min)/sin(θ_min) - cos(θ_min)      ,
   
which provides the height of the barrier preventing that the system can move
from the angle θ_c to the angle π - θ_c, which is needed for the droplet to
subsequently grow
increasing its baseline beyond b at fixed angle π - θ_c. Expanding equation (<ref>)
in terms of the angle π/2 - θ_min = α, one gets 
Δ F_2 / b f_int≃2/3 (α^3).
For θ_min = 70^∘ (α = 20^∘) this leads to a barrier (in units of k_B T) of about 2.6.

Actually, when m_coex b H_b approaches f_int, then the angles where
the minimum and maximum of the free energy occur, merge at θ = 90^∘
(c.f. figure <ref> for b H_b = 1.8). However, already at a smaller field
(keeping b constant) the barrier caused by the free energy maximum, Δ F_2
given by equation (<ref>), will be small enough so that the second
nucleation event by which the angle grows from θ_c to π - θ_c
can take place. Note that the analytical formula for the barrier can also be extracted
from equations (<ref>), (<ref>), and figure <ref> shows that long
before θ reaches 90^∘  it will be of order of a few k_B T only. This argument  also
explains why the temperature dependence of the apparent angle (close to
70^∘) where the depinning transition occurs is rather weak (see Table I):
the scale for the barrier is simply set by b f_int, and this quantity
does not vary strongly with T for the choices we have made.
Thus when this barrier is small enough the instability that would occur
for θ = 90^∘ (where θ_min and θ_max merge) is preempted by the
jump of the angle θ from θ_c to π - θ_c .



On the other hand, one can change the height of the barrier by around
one order of magnitude just by taking a fixed temperature (t = 0.4),
such as f_int = 1.79873, but varying the length of the inhomogeneity
13 ≤ b ≤ 101. In this way an increment of the
apparent angle of about 20^∘ is observed, as shown in Table II.
It is also obvious from equation (<ref>) that solutions for θ_min,
corresponding to pinned droplets, exist only for m_coex H_b b /f_int < 1 :
for larger fields droplet growth with time is not hindered by any barriers,
after they have been nucleated.



The smooth variation of θ from θ_c to π - θ_c with increasing droplet 
area as predicted by Lipowsky et al. <cit.> is a special consequence of the canonical 
ensemble, where the droplet volume (in d = 3) or the droplet area (in d = 2) is taken as a
fixed independent variable. In contrast,  
only part of the variation is realizable as a metastability effect in the grandcanonical 
ensemble, where H_b is given. According to the theory outlined above, metastable
pinned droplets should exist only up to a "spinodal" where 
θ_min = π/2. The spinodal field is then H_b^spin = f_int/b m_coex,
such that  H_b^spin / H_b^crit = 1/sin(θ_c)
(see also equations  (<ref>) and (<ref>)).
However, it should be kept in mind that for systems with short-range
interactions "spinodals" are a somewhat ill defined concept <cit.>
and cannot be reached in practice.
In the context of nucleation phenomena, the present case
of a grandcanonical ensemble is the physically meaningful choice, of course.

For the understanding of the results observed in the simulations, it is hence
crucial to consider the combined effects of the primary nucleation event of the
wall-attached droplet and a further growth of this droplet.

We have made the hypothesis, that for this growth a "local equilibrium" assumption
holds, in particular near the point where the droplet-vapor  interface
meets the substrate.
This implies, for the case where the length b_drop of the growing droplet still is
less than b, that we have θ = θ_c for the contact line of the growing droplet
(see Figure <ref> in Appendix A, top panel). However, for droplets that have b_drop > b, we have
θ = π - θ_c (see Figure <ref> in Appendix A, lower panel).
For fields H_b < H_b^spin = f_int/(b m_coex) metastable pinned droplets are predicted,
and a barrier Δ F_2 for the “depinning” of these droplets could be estimated (Equation (<ref>)).
The actual limit of stability of metastable pinned droplets, as seen in Figures <ref>, <ref>, is
somewhat smaller than H_b^spin: this happens because when the barrier Δ F_2
is small, it can be overcome in a second nucleation event.

Of course, a perfect quantitative agreement of the predictions 
based on our phenomenological theory for pinned droplets with 
the corresponding simulation results should not be expected: (i) 
the mean-field like treatment of equations (<ref>)-(<ref>)
disregards the huge statistical fluctuations that are 
present (Figure <ref>), (ii) the anisotropy of the interfacial
free energy should lead to some deviations of the shape of 
the pinned droplets from the circle cut, which should cause
some systematic deviations from the free energy plotted in 
figure <ref>. Also, the curvature of the interface 
may modify the effective surface tension.



§ NUCLEATION KINETICS AND DROPLET GROWTH.


Already in earlier work on studies of homogeneous nucleation
in bulk Ising models (see e.g.  <cit.> for a recent review) it has
been shown for conditions where the phase transformation is caused
by nucleation and growth of a single droplet one needs to distinguish
two very different time scales: The typical lifetime of the metastable
state, which is then simply inversely proportional to the nucleation rate; and
the time needed for the nucleated droplet to grow and essentially occupy the
total volume of the simulation box. However, often these processes are somewhat
confused by the crossover to the regime where during the phase transformation
many droplets are nucleated in different parts of the system and the lifetime of the
metastable state then is limited by this competitive growth of many droplets. This
latter regime is dominant when the simulation volume is relatively large
and H_b is not so small, so nucleation becomes relatively easy <cit.>.

In the present work, conditions were chosen such that homogeneous nucleation
is not observable at all, and heterogeneous nucleation is restricted to the boundary
region of length b.  It then is rather straightforward to follow the growth of the single
droplet (Figure <ref>), and it turns out that the time intervals between
the snapshots of the growing droplet are indeed very small in comparison with the
nucleation times. In order to give further insight on the involved times
as well as on the growing and nucleation process of the droplets, Figure 
<ref> shows plots of the time evolution of both the total excess density Δρ
due to the droplet (upper panel)
and the linear density excess  Δρ_ measured in the direction perpendicular to the wall
just at the center of the droplet. Each curve is the average over several hundred
individual time evolutions of the system.

In Figure <ref>(a) one can roughly
estimate the average time required by the system to achieve the metastable
state (τ_MS), which for the case shown (i.e. t = 0.40 and b = 51) 
is τ_MS≃ 10^4 MCS. Choosing H_b = 0.030 no nucleation events
are detected during the observation time (τ_obs = 10^6), and the curve
remains flat after achieving the metastable state. For H_b = 0.032 few
nucleation events are detected, and each of them shows up as an upward step
in the corresponding plot. The height of each individual step simply is the
inverse of the number of runs, since in each run when nucleation has occurred the
droplet grows fast (on the scale of τ_N) to fill the available area.
Here, one can estimate the typical growth time (τ_G)
required for each already nucleated droplet to expand over the whole sample,
namely τ_G ≃ 10^4 MCS. However, for  H_b = 0.034 nucleation is
dominant and one can estimate τ_N ≃ 4× 10^4 (also by discounting
τ_MS one can get τ_N ≃ 3× 10^4 MCS). On the other hand,
the time evolution of the density per unit length as measured in the
direction perpendicular to the sample (Figure <ref>(b)), which shows
the development of the droplet in that direction, is fully
consistent with the above discussed scenario. Furthermore, here one can also
estimate the average time elapsed between the onset of nucleation and the 
achievement of a full liquid phase covering the whole sample, 
τ_LP≃ 8× 10^5 MCS. Note that this averaged time results from 
the contribution of many growing events of already nucleated droplets, 
occurring at different times over a wide time interval, 
(actually 318 events for H_b = 0.034) each of them having a short lifetime 
of the order of τ_G ≃ 10^4 MCS, as already discussed.

The individual nucleation events seen in figure <ref> actually all relate to
overcoming the barrier Δ F_2 discussed in Figure <ref>, since the fields H_b
all are slightly below the stability limit H_b^spin.  

For a more quantitative analysis, we have also recorded both the nucleation time
distribution (P(τ_N)) and the growth time distribution (P(τ_G)) for
the case t = 0.40,  b = 51, and for different choices of H_b, as shown 
in  Figures <ref> (a) and (b), respectively.


These choices all refer to H_b > H_b^spin, and hence for them the barrier
Δ F_2 does no longer occur.
A simple comparison of
both figures indicates that the characteristic times, as estimated from the
location of the peaks of the distributions, are roughly of the same order
for larger fields (H_b ≥ 0.08), while τ_N > τ_G in the opposite
limit. In fact, figure <ref> (a) shows the monotonic 
increase of the ratio τ_G / τ_N when it is plotted versus H_b,
spanning the range 0.25 ≤τ_G / τ_N ≤ 1.



Also, Figure <ref> (b) shows a plot of τ_N versus
1/H_b. Nucleation theory predicts  

    ln(τ_N) ∝Δ F^*/k_B T = π/2 f_int^21/(ρ^ coex_ℓ - ρ^ coex_v )  H_bf_VT/k_B T   .

However, the curvature of the log-linear plot indicates that only part of the
chosen region of fields is in the regime where the barrier due to the
heterogeneous nucleation on the inhomogeneity controls the kinetic exclusively;
in fact, when H_b approaches H_b^spin a slowing down related to the
barrier Δ F_2 that occurs for H_b < H_b^spin may be present. 
The best fit of the data of Figure <ref> (b), within the linear regime,
yields Δ F^* H_b /k_B T = 0.111.
This number is smaller than the theoretical expectation given by
π/2 f_int^21/(ρ^ coex_ℓ - ρ^ coex_v )f_VT/k_B T
(see equation (<ref>)),
obtained by taking f_int = 1.7987 <cit.>, which yields
Δ F^* H_b /k_B T = 0.231 for θ_c = 34^∘
in the Volmer-Turnbull factor (Equation (<ref>)).
It is a subtle issue to understand where this discrepancy 
of about a factor two in the effective barrier height comes from.
On inmediate thought concerns the curvature dependence of the 
interfacial free energy f_int(R). In d = 2 there is 
evidence from field theoretical calculations <cit.>, Monte
Carlo simulations of cluster-size distributions <cit.>,
and analysis of the two-phase coexistence <cit.> that

    f_int(R)/f_int(∞) = 1 + 5/4π f_int(∞)ln(R)/R +
            const/R   ,

where the constant in the last term on the right-hand side
is non-universal, while the prefactor 5/4π of the logarithmic 
term is universal. If only this correction would be taken into account,
the interfacial tension be enhanced by a factor 1 + 0.221  ln(R)/R,
which for typical values of R (e.g. R = 16) is an enhancement of about 4%.
Neither the magnitude nor the sign of this effect can account for the observed 
discrepancy. Actually a more plausible assumption is that our estimate  of the 
contact angle θ_c and hence the factor f_VT(θ_c) is an overestimate.
Since  f_VT(θ_c) ≈4/3 πθ^3 a decrease of θ_c
by a few degrees already suffies to reduce  f_VT(θ_c) by a factor of two,
e.g.  f_VT(27^∘)/f_VT(34^∘) ≃ 0.51.
Moreover the equation for f_VT(θ_c) holds only for circle-cut shaped 
droplets, and the effect of anisotropy cuasing somewhat non-circular
shapes (see figure <ref> in the Appendix A) of the droplet
on f_VT(θ_c) still needs to be clarified. Also, the curvature
of the plot shown in figure <ref> b) may be taken as an indication
that it is questionable whether the asymptotic region where the theory 
holds has been reached.  Thus clearly the conclusion emerges that 
in spite of the simplicity of the Ising model still more work is needed 
to understand there heterogeneous nucleation quantitatively.  
  
Figure <ref> c) shows a log-log plot of τ_G versus
H_b to show that the growth time τ_G scales inversely with H_b, as expected.
On the other hand, figure <ref>
shows log-linear plots of  P(τ_N) versus  τ_N, for different
choices of the length of the inhomogeneity b and the bulk field H_b,  
demonstrating an exponential distribution for the long times, as theoretically
expected <cit.>.



Of course, it would be interesting to explore the kinetics of heterogeneous nucleation
systematically for a wide range of b, but due to excessive needs for computer
time this has not been attempted. So, we have only studied one other choice, b = 25.
For b = 25 one has H_b^spin/J ≃ 0.060, so the
measurements performed slightly above H_b^spin indicate that τ_N 
depends strongly on H_b, and decrease for larger values of H_b, as
expected. The fact that τ_N for b= 51 and H_b/J = 0.034 is smaller than 
for the case b = 25 and  H_b/J = 0.064 means that one has smaller
nucleation times for larger inhomogeneities, since the larger the inhomogeneity
is, the less tendency is found that the nucleated droplet gets pinned.




§ CONCLUSIONS.


In this work, we have considered the effect of a chemical inhomogeneity on heterogeneous
nucleation on a flat substrate. As a generic case, we have considered the simple Ising
lattice gas system in d=2 dimensions, where the flat substrate is just a straight line,
and the chemical inhomogeneity is represented by a positive boundary field H_w3
(favoring the liquid phase of the lattice gas) acting on a part of this line of length b,
while in the remaining part of this lower boundary of the system a boundary field
H_w2 = - H_w3  is applied  favoring the vapor phase of the lattice gas.
For the Monte Carlo simulations of our model, we choose a boundary with finite
length M and periodic boundary condition in the direction parallel to this boundary,
while in the direction perpendicular to this boundary a finite linear dimension L is used,
and at the upper boundary a negative field H_w1 = H_w2/4 acts, to stabilize the vapor
as a bulk phase of the system in the absence of a bulk field, H_b = 0.
For conditions of partial wetting, the density inhomogeneity in the lattice gas caused
by the chemical inhomogeneity then extends only over a distance of the order of one lattice
spacing, for H_b = 0 (Fig.<ref> a)). However, when  a small field H_b > 0 is applied,
the vapor phase chosen as the initial state of the system becomes metastable, and the
structure of the density inhomogeneity caused by the chemical inhomogeneity, Figs. <ref>-<ref>,
as well as the decay rate of the metastable state due to nucleation of a boundary-attached
droplet and its growth Figs. <ref>-<ref>, are the subjects of investigation.
Conditions are chosen such that neither homogeneous nucleation in the bulk nor nucleation
starting in the boundary regions favoring the vapor phase can ever be observed.

We perform for each choice of temperature, width b of the chemical inhomogeneity, and fields H_b,
many hundred equivalent Monte Carlo runs, differing by the pseudorandom numbers used to realize
the time evolution of the Monte Carlo sampling process. By using over 20 × 10^6 Monte
Carlo steps (MCS) per site, we are able to reliably  estimate various relaxations times and
their distributions (Figs. <ref>-<ref>) over 6 decades of time. The initial stages of the
relaxation process are characterized by the equilibration of the metastable state,
after the field H_b has been switched on at time τ = 0, taking a time τ_MS.
If H_b is small enough, e.g. H_b  ≤ 0.030 for the choice t=0.40, b=51, no decay of the
metastable state is observed, which implies that the nucleation time τ_N exceeds the observation
time τ_obs.

If we would study heterogeneous nucleation on a chemically homogeneous boundary
of linear dimension M, the nucleation time (for the regime of fields where single-droplet
nucleation matters) would be related to the nucleation rate J^het by τ_N = (M J^het)^-1).
When a nucleation event has occurred, it takes a time τ_G for the critical droplet to grow
until the whole (finite) system has transformed; only when τ_N  ≫τ_G is the
phase transformation triggered by single nuclei the dominant process. In the regime where
τ_N and τ_G are comparable the simultaneous growth of multiple nucleated droplets
needs to be considered, making separate estimations of τ_N and τ_G difficult.
By choosing our geometry with a chemical inhomogeneity, we extend the regime where the
transformation triggered by single nuclei is the dominant process: note that the baseline of
the critical droplet is b_drop^∗ = 2 R^∗sin(θ_c), where R^∗ is the critical droplet
radius and θ_c the contact angle, assuming droplets of circle-cut shape;
only when b_drop≪ b, phase transformations affected by nucleation of multiple droplets and
their competitive growth could matter. The detailed analysis of our observed phase transformation
events has allowed us a separate analysis of the distributions of τ_N and τ_G;
as theoretically expected, the growth rate of supercritical droplets is proportional to H_b,
and hence τ_G ∝ 1/H_b (Fig. <ref>c) ), while τ_N varies exponentially
with 1/H_b, lnτ_N  ∝ 1/H_b (Fig.<ref>b)), as expected from nucleation
theory in d=2 dimensions. Unfortunately, only a very small range of H_b, much less than
a decade (Fig.<ref> b)), is available when the time scales for nucleation and growth
are well separated.

Very interesting behavior was found for the metastable regime, where during observation
time τ_obs no phase transformation occurs. On general grounds one can predict that
then the chemical inhomogeneity causes an excess density Δρ in the system,
which exhibits a scaling behavior Δρ = b^2 M̃(t, bH_b), Eq.(<ref>),
the effective droplet height h_eff/b similarly is a function of the product b H_b
only (Fig. <ref> c)), as well as the effective contact angle θ_eff (Fig. <ref> b)).
We hence identified two regimes: for very small values of H_b such that b_drop > b,
nucleation of droplets with the contact angle θ_c  “ preferred “ by  the chemical
inhomogeneity still is geometrically impossible, it does not matter how large observation
times are chosen. In fact, in this regime critical droplets would have the shape as
shown in Fig. <ref> (lower part), their baseline b_drop extending beyond b and
the contact angle being π - θ_c, but the corresponding nucleation barriers correspond
to astronomically large nucleation times, and hence are of no interest here.
In this regime, thermal fluctuations allow only the occasional formation of subcritical nuclei
with R < R^∗, of circle cut shape with contact angle θ_c.  Thus the average effect
of such fluctuations is measured by the scaling function M̃(t, b H_b) for H_b < H_b^crit,
where H_b^crit can be estimated as H_b^crit = sin(θ_c) f_int/(b m_coex),
when we ignore anisotropy effects on the interfacial tension in the lattice gas model.
In any cases H_b^crit is the smallest field where a droplet with the correct contact angle
θ_c fits to the chemical inhomogeneity. For H_b^crit < H_b < H_b^spin we may
encounter pinned droplets, having a baseline of length b, with contact angles exceeding
the equilibrium value, θ_c < θ < π/2 (cf. Fig. <ref>). Assuming that
these droplets still have circle cut shape,  we have predicted that these droplets become
unstable for H_b^spin = f_int/(b m_coex), and we have obtained an approximation
for M̃(t, bH_b) in this regime (Figure <ref> a)).  Note that further growth of
the droplets with baseline b_drop > b requires that the contact angle grows up to
π - θ_c, and in the regime H_b^crit < H_b < H_b^spin this is hindered by
a free energy barrier Δ F_2, see Figure <ref> and Equation (<ref>).
Since this barrier is only of the order of a few k_BT when θ  has reached about
70^∘, pinned droplets with 70^∘ < θ < π/2 actually were not observed;
droplets with a shape as sketched in the lower part of Figure <ref> then appear in a second nucleation
event, and grow to complete the phase transformation. So the actual limit of metastability
H_b^max of pinned droplets is somewhat smaller than H_b^spin, e.g. (for b = 51)
H_b^max b ≃ 1.71 while H_b^spin b ≃ 1.7987.


Figure <ref> Caption. Sketches summarizing qualitatively the different 
regions of wall-attached droplets for the case of an inhomogeneity of linear
extension b identified in this work.
Case a) shows the range where the bulk field H_b is very small,
namely H_b < H_b^crit, ( H_b^crit = f_intsin(θ_c)/(b m_coex)).
Then the critical droplet having a contact angle θ_c would require 
a baseline larger than b, since then R^* exceeds b/(2 sin(θ_c)),
see the case H_b = 0.030 in part d). No nucleation then is possible,
and due to the average effect of subcritical droplets a density excess 
Δρ∝ b ξ_ occurs on the inhomogeneity.
Case b) shows the regime H_b^crit < H_b < H_b^spin (H_b^spin = f_int/(b m_coex)),
where critical droplets of radius R^* and contact angle θ_c with baseline
b_drop^* < b are nucleated, and their radius grows until their baseline is equal to b.
Then these droplets can lower their free energy further by increasing  their
contact angle from  θ_c to  θ_min. These metastable 
pinned droplets are characterized by sin(θ_min) = m_coex H_b b/f_int.
The case H_b = 0.075 in part d) illustrate the corresponding Δ F(R) in the regime
where θ = θ_c.
Case c) shows the behavior for H_b^spin < H_b, e.g. the case  H_b = 0.20 in part d),
where the critical droplet nucleates with contact angle θ_c has such a
small radius  R^* and corresponding baseline, that after growth to the baseline b     
the increase of the contact angle is no longer pinned, and when the contact angle 
π - θ_c has been reached, the droplet can grow further with this contact angle 
and increasing thereby its contact line b_drop(τ) beyond b with increasing time.
Thus, the two critical fields  H_b^crit and H_b^spin simply correspond to the cases
R^* = b/(2 sin(θ_c)) and R^* = b/2, respectively. Note that in panel c)
the length of the inhomogeneity in the upper sketch was taken a factor two larger 
than in the medium and lower sketches for the sake of clarity. 
Panel d) shows plots of the excess free energy relative to the wall without droplet
Δ F(R) vs R (see equation (<ref>), as obtained for different fields
corresponding to the regimes shown in panels a)-c). Theoretical curves are obtained by taken
θ_c = 35^∘ and f_int = 1.7987 in order to illustrate the expected 
behavior for t = 0.40. The chosen fields are then suitable to describe the
case of an inhomogeneity of extension b = 21, such that  b/(2 sin(θ_c)) = 18.31,
see the vertical dashed line. The values of R^* corresponding to the selected
fields are shown along the horizontal axis.  
More details in the text.



Thus, an unexpectedly rich behavior concerning nucleation at a chemically
  inhomogeneous substrate has been found within the context of the grandcanonical ensemble
  used in our calculations (i.e. when the pressure of the fluid or equivalently the magnetic field
  of the Ising model is given as a control parameter). We have shown that this scenario differs from that
  corresponding to the canonical ensemble (i.e. when the volume in d = 3 or the area
  in d = 2 of the droplets is taken as a control parameter). In order to acquaint the reader
  with a clear description of the relevant findings reported in this paper, we have summarized
  and discussed our results in figure <ref>. In this way we addressed the relevant regimes
  encuntered in our study performed in the framework of the grandcanonical ensemble:
  a) The regime H_b < H_b^crit where no nucleation is possible. b) The regime 
  H_b^crit < H_b < H_b^spin, where the droplet grows with contact angle θ_c until
  its baseline matches the length of the inhomogeneity, and then subsequently grows by
  keeping its baseline constant but increasing the contact angle. Finally, the regime c)
  corresponds to larger fields H_b > H_b^spin that lie beyond the stability limit,
  so that the droplets can grow with baseline larger than the length of the inhomogeneity
  and contact angle π - θ_c. Furthermore, all these three regimes are
  properly identified with the corresponding free energy functions F(R) shown in panel
  d) of figure <ref>.

















§ HETEROGENEOUS NUCLEATION ON HOMOGENEOUS AND CHEMICALLY INHOMOGENEOUS
SUBSTRATES: PHENOMENOLOGICAL THEORY.


We consider here a fluid in d = 2 dimensions exposed to a 
(one-dimensional) flat wall, under conditions of partial wetting.
For simplicity, we shall consider explicitly only the nearest-neighbor 
lattice gas (Ising) model and the wall is then oriented along the x-direction 
of the lattice; the Ising spins adjacent to this wall then have a row of 
missing neighbors but experience a surface field H_w3 chosen such that the liquid 
phase of the lattice gas model is favored; unlike Figure <ref> we consider
first the case where this boundary field is homogeneous independent of the coordinate x.
In the semi-infinite system at zero bulk field H_b = 0, the vapor phase can be stable
(due to a suitable boundary field at y = L →∞), and near the boundary
where H_w3 acts only a liquid film of the average thickness of order ξ_
is stabilized; ξ_is the perpendicular correlation length of the wetting 
transition, and in the regime of partial wetting this is a microscopic length, 
i.e. of the order of the lattice spacing in the lattice gas model.

Note that in the grandcanonical ensemble (temperature T and  bulk field
H_b chosen as control variables) liquid droplets in full thermal equilibrium
can exist neither in the bulk nor attached to the wall, even if a nonzero
bulk field favoring the liquid phase is switched on.  However, metastable droplets
can exist only for a finite "lifetime", and we shall address this case using concepts
of the theory of heterogeneous nucleation.  Thus, it is useful to recall that stable
liquid droplets of macroscopic size do exist when we consider the canonical ensemble,
choosing the density of the vapor ρ in between the density of coexisting vapor
(ρ_v) and liquid (ρ_ l) phases.  The area fraction of the liquid then 
is given by the lever rule, X = (ρ - ρ_v)/(ρ_ l - ρ_v), and 
the shape of the liquid domain does depend on the choice of boundary conditions.
In our case, for small X we obtain a wall attached droplet having the shape of
a circle cut with baseline b_drop given in terms of the droplet radius R and contact 
angle θ_c by 

    b_drop = 2 R sin (θ_c)

and the contact angle for an isotropic interface tension f_int between the 
coexisting vapor and liquid phases is given by Young’s equation, namely

    f_int cos(θ_c) = f_wv - f_w l   ,

where f_wv, f_w l are the surface excess free energy   
densities of the vapor phase (f_wv) and liquid phase (f_w l)
due to the wall. 
Actually, Equation (<ref>) is valid for a fluid in continuous space,
but not for the lattice gas model, where the interface tension f_int(θ)
depends on the angle θ between the interface normal and the x axis
of the lattice. Then Equation (<ref>) needs to be replaced by   

    f_int(θ_c) cos(θ_c) - sin(θ_c) df_int(θ_c)/dθ_c|_θ = θ_c  = f_wv - f_w l   .
 

While in the isotropic case straightforward geometric considerations yield the area (A) of the 
circle cut as 

    A = R^2 (θ_c -1/2 sin(2 θ_c)) = 1/4 b_drop^2 (θ_c/sin^2(θ_c)  -cos(θ_c)/sin(θ_c))        ,

and the length of the vapor-liquid interface line is

    l_ lv = 2 R θ_c = b_dropθ_c/sin(θ_c)   ,

finding the droplet shape for the anisotropic case is less straightforward.

In the bulk this problem is solved in terms of the well known Wulff construction, which for
the d = 2 Ising model can be worked out explicitly, and the shape of the droplet is given by
the equation <cit.>

    cosh(x̃) + cosh(ỹ) = cosh(2J/k_B T)/tanh(2J/k_B T)      ,

where x̃ and ỹ are the x, y coordinates of the curve describing
the droplet shape. Equation (<ref>) interpolates smoothly between a square shape
(for T → 0) and a circle (for T → T_cb). When we inscribe a circle
that touches the actual shape at x̃ = 0 and at ỹ = 0, it has a radius R_in
given by

    R_in = arcosh [cosh(2J/k_B T)/tanh(2J/k_B T) -1 ]      .

As an example we hence plotted Equation (<ref>), 
in figure <ref>, together with the inscribed circle of radius R_in  
evaluated at t = T/T_cb = 0.40, since for this choice of the reduced temperature
most of our simulations were made.

We found that the deviations from the spherical
shape are already rather minor, and this justifies our neglect of these anisotropy
effects, at least as a first approximation. 
The solution of Equation (<ref>), reduces to the equation
of a circle near T_cb, where x̃→ 0 and ỹ→ 0 and hence

    x̃^2 + ỹ^2 = 2 cosh(2J/k_B T)/tanh(2J/k_B T)  - 4    ,

recalling that cosh(2J/k_B T_cb) = 2^1/2, sinh(2J/k_B T_cb) = 1, and hence 
R_in→ 0 as well. The solution for the wall-attached droplet then is given by the
Winterbotton construction <cit.>, i.e. we have to cut the droplet shown in figure <ref>
by an horizontal straight line such that the angle of the tangent is θ_c as given by
Equation (<ref>). The linear dimensions ỹ_cut and R_in then follow from the 
condition that the area above the cut yields the desired area fraction X.

However, for the sake of simplicity we shall ignore these anisotropy effects in the following,
working with droplets of circular shape only. But even then there is one fundamental problem: there
is no physical reason for the x-coordinate of the center of mass of the droplet to coincide
with the origin of the coordinate system. In fact, this center of mass coordinate can be anywhere
on the x-axis when the boundary field is homogeneous, independent of x. Even in the
inhomogeneous case the droplets are only on average centered in the middle of the inhomogeneity,
as e.g. can be qualitatively observed in the snapshots of Figure <ref>. This
fact creates a translational entropy contribution k_B T ln(M) for the droplet, where
1 ≤ x ≤ M in our finite lattice of length M in the homogeneous case.
Similar translational entropy contributions 
are known to hamper the numerical study of interfacial free energies <cit.>. Thus, in 
a straightforward simulation study of the present problem the droplet would diffuse along the 
x-axis and its density profile ρ(x,y) would be completely smeared out, until only the 
average translationally invariant density profile ρ_av(y) is left, containing little
information on the droplet. Thus a "demon" would be needed to constrain the sampling of 
configurations such, that in each microstate of the system that is sampled the droplet center of mass
has its x-coordinate in the origin. Practical implementation of such a constraint 
is not completely trivial, since the size and the shape of the droplet due to their nanoscale 
dimensions is strongly fluctuating (c.f. figure <ref>).

We now consider the main subject of interest of the present paper, namely a boundary condition 
of the type shown in Figure <ref>, where the surface field H_w3 acts only over a distance 
b along the x-axis, while in the remaining boundary a field H_w2 = -|H_w3| acts,
and hence the contact angle θ^'_c = π - θ_c applies.

In d = 3 dimensions in the canonical ensemble, this situation has already been 
considered by Lipowsky et al. <cit.>.

They pointed out that three regimes
need to be distinguished, namely: (I) b_drop < b; (II) b_drop = b; and (III) b_drop > b, 
see Figure <ref> adapted to our d = 2 dimensional case. Here the area A of the
wall attached droplet is the control parameter that is varied: For sufficiently small area taken 
by the liquid baseline b_drop that will result from A and θ_c via Equation (<ref>) will
be in the regime I, and the x-coordinate of the center of mass of the droplet can be
anywhere in the interval from x = -(b-b_drop)/2  to   x = +(b-b_drop)/2. Unlike Lipowsky   
et al. <cit.> we do not assume that the droplet is exactly centered at x = 0,
the center of the inhomogeneity of the wall, which is our origin. This center certainly is
the most probable position, but there will be a broad probability distribution for this
center of mass coordinate, and when we consider the average density profile  ρ(x,y) 
obtained by convoluting the density profile of the droplet with baseline b_drop and contact angle
θ_c with this probability distribution, a density distribution ρ_ave(x,y) must result
that is considerably flattened in comparison with ρ(x,y).  From ρ_ave(x,y) one would obtain 
an effective contact angle θ_c^eff that clearly will be much smaller than the correct one,
if b_drop≪ b. This entropic effect was disregarded  by Lipowsky et al. <cit.>, but clearly 
must be present in our simulations and thus hampers their interpretation.
It is tempting to associate the small values of θ_eff in figures <ref>  and <ref> observed
for l_eff < b with this flattened profiles due to the fluctuations in the
center of mass position of small wall-attached droplets.

The most interesting situations of course, are found when b_drop as given by Equations (<ref>),
(<ref>) has reached the value  b_drop = b: then the prediction is that further increase of A 
does not cause a further growth of b_drop. Rather, what happens is a growth of the contact angle 
θ of the droplet from the value θ_c given by the Young's equation to a larger value,
satisfying an equation analogous to Equation (<ref>), namely

    A = 1/4 b^2 (θ/sin^2(θ)  -cos(θ)/sin(θ) )    ,
    θ_c < θ < θ_c^'   .

Thus in a sense the interface between liquid and vapor is pinned at the points 
x = ± b/2 when A has increased up to the value where Equation (<ref>) yields
θ =  θ_c^' (= π - θ_c, in our case),
depinning of the interface from the inhomogeneities of the boundary occurs, and 
θ stays at θ_c^', while b_drop > b. Again Lipowsky et al. <cit.>
have assumed that the x-coordinate of the center of mass of the droplet is at x = 0, but
we maintain that again fluctuations will occur. However, the 
region b_drop in between the two contact points of the interface will always encompass the 
region of the inhomogeneity, from x = -b/2 to x = +b/2, and the average position of the 
center of mass of the droplet will hence have x-coordinate x = 0.

When we now turn to the description in the grandcanonical ensemble, we note that a correspondence
to the droplet configurations discussed for the canonical ensemble can exist only when the droplet
configurations in the grandcanonical ensemble are still metastable.

For the problem without boundary inhomogeneity we have the standard problem of heterogeneous
nucleation at the wall. The free energy cost of the forming droplet is written as the 
excess free energy relative to the wall without droplet, namely

    Δ F_drop = -2 m_coex H_b A + Δ F_int

where m_coex is the spontaneous magnetization of the Ising model,
and A is given by Equation (<ref>) and Δ F_int becomes,

    Δ F_int = 2 R f_intθ_c + (f_wl - f_wv) 2 R sin(θ_c) = 
    2 f_int R[θ_c - 1/2 sin(2 θ_c)]      ,
 
where Equation (<ref>) was used. Equations (<ref>) and (<ref>) yield

    Δ F_drop = [θ_c - 1/2 sin(2 θ_c)][- 2 m_coex H_b R^2 + 2 f_int R].
 

Minimizing Δ F_drop with respect to R yields

    R^* = f_int / (2 m_coex H_b)     , Δ F^*_het = Δ F^*_homo f_VT (θ_c)     ,
 
where Δ F^*_homo is the standard result for the free energy barrier 
against homogeneous nucleation in d = 2 dimensions  

    Δ F^*_homo = π/2 f_int^2 /(m_coex H_b)      ,
 
and f_VT (θ_c) is the analog of the well-known Volmer-Turnbull
function in d = 2 dimensions, given by

    f_VT (θ_c) = 1/π  [θ_c - 1/2 sin(2 θ_c)]    .


Note that f_VT (θ_c) ≈2/3 πθ_c^3 for θ_c → 0,
when complete wetting begins.  It turns out, of course, that use of R^* = f_int / (2 m_coex H_b)
in Equation (<ref>) yields b_drop≪ b only for rather large fields. All the data where the metastable
droplets are encountered do not fall in this regime, as expected.




§ PINNED DROPLETS: COMPARING THE TWO- AND THE THREE-DIMENSIONAL CASES.


While the numerical simulation work exclusively has addressed the case of a
two-dimensional system with a one-dimensional boundary where the positive surface
field (favoring the liquid phase of the lattice gas) acts on a length b,
it is also instructive to consider the three-dimensional case, where the positive
surface field acts on a circular heterogeneity with radius r.
For the sake of clarity, the geometry of the pinned droplet is sketched in Figure <ref>.


The radius of curvature of the sphere-cap shaped droplet is R. Then 


    r = R sin (θ),       h = R(1 - cos(θ))     ,

  
where it is convenient to express all quantities in terms of the
height h of the droplet above the substrate. The angle θ that the droplet
makes with the substrate can be in the range 


    θ_c ≤θ≤π - θ_c      ,

  
where θ_c is the contact angle given by Young's equation. Notice that
only for angles in the quoted range droplets with basal radius r exist; however,
only for θ≤π/2 such droplets are metastable, while for
π /2 ≤θ≤π - θ_c they are  unstable.

Now the volume of the sphere cap is


    V = π h/6 (3 r^2 + h^2)    ,
   
also the basis surface is π r^2, while the upper surface is 


    A_u = π (r^2 + h^2)    .
 

So, the free energy difference of the droplet of height h
relative to a disk-shaped droplet of radius r and height h = 0 is
(the choice of this reference state is arbitrary, of course)


    Δ F =  f_int A_u  - 2 m_coex H_b V     ,

  
where H_b is the bulk field. Then, by using Equations (<ref>)
and (<ref>) one obtains 


    Δ F =  f_intπ (r^2 + h^2)  -  m_coex H_b π/3 h (3 r^2 + h^2)    .


It is convenient to find the extrema of Δ F simply as a function of h; then


    (∂ (Δ F) / ∂ h)_H_b = 0     ,
 

yields


    h^2 - 2h f_int/m_coex H_b  + r^2 = 0        ,
  

such that in terms of H̃_b = H_b r m_coex / f_int
one finds for H̃_b < 1 two solutions, namely


    h/r = H̃_b^-1±√((H̃_b^-2 - 1))   .


The minus sign yields the free energy minimum, corresponding to
the pinned droplet, while the plus sign corresponds to a surface
free energy maximum, and the corresponding angle θ can be read off from

    tan(θ/2) = 1 - cos(θ)/sin(θ) = h/r   .

  
The limiting case H̃_b = 1 means h/r = 1, θ = π/2, i.e.
a semispherical droplet. The free energy function can be written as

    Δ F /  f_intπ  r^2 = 2/3 (H̃_b^-2±H̃_b^-1√((H̃_b^-2 - 1))∓H̃_b √((H̃_b^-2 - 1)) )     ,

and hence the barrier for the depinning transition of the pinned droplet becomes

    Δ F_2 /  f_intπ  r^2 = 4/3H̃_b (H̃_b^-2 -1)^3/2   .
   
From this calculation it is obvious that the mathematics in d = 3
is even simpler than in d = 2, since the use of h instead of the angle
θ makes the description of  Δ F  (c.f. Equation (<ref>))
very simple. In d = 2, equation (<ref>) also holds, but equation (<ref>)
shows that both θ and sin(2θ) enter in the free energy expression,
so no simple formula for Δ F(h) in d = 2 can be written down. When one works
out h/r and Δ F_2 in both d = 2 and d = 3, one notes a very similar behavior:
near the point  H̃_b = 1 the barrier vanishes like
(1 - H̃_b^-1)^3/2, i.e. with a vanishing slope, and h/r reaches the semicircle or semisphere
configuration with a square-root cusp.

Figure <ref> shows a comparison of the free energy barrier of the
depinning transition of pinned droplets for the two- and the three-dimensional cases. 


Acknowledgments: E.V.A. is grateful to the Alexander von Humboldt foundation and 
to the Deutsche Forschungsgemeinschaft (DFG, SFB TRR 146) for partial support of his research 
stays at the Institut für Physik of the Johannes Gutenberg Universität Mainz.
Also, M.L.T. and E.V.A are grateful to the CONICET and UNLP (Argentina) for financial
support.



99

1h L. D. Landau and E. M. Lifshitz, Statistical Physics, (Pergamon Press, Oxford, 1958).

2h K.Binder, Rep. Progr. Phys. 50, 783 (1987)

3h A. C. Zettlemoyer (ed.), Nucleation, (M. Debber, New York, 1969). 

4h P. G. Debenedetti, Metastable Liquids, (Princenton University Press, Princenton, 1996).

5h D. Kashchiev, Nucleation: Basic Theory with Applications, (Butterworth-Heinemann, Oxford, 2000). 

6h S. Balibar and J. Villain (eds.) Nucleation, C. R. Phys. 7,
fascicule 9-10, (special issue, 2006).

6hb D. Beysens, C. R. Phys. 7, 1082 (2006); 
D. Beysens, V. Pruvost, and B. Pruvost, J. Arid Env. 135, 90 (2016). 

7h H. Biloni, in Physical Metallurgy. (R. W. Cahen and P. Haasen, eds.) p. 477 (North-Holland, Amsterdam, 1983).

7hb M. Volmer, Kinetik der  Phasenbildung,
(Th.Steinkopff, Dresden and Leipzig, 1939).

8h D. Turnbull, J. Appl. Phys.21, 1022 (1950).

9h D. Turnbull, J. Chem. Phys.18, 198 (1950).

10h D. Turnbull and B. Vonnegut, Ind. Eng. Chem. 44, 1292 (1952).

11h N. H. Fletcher, J. Chem. Phys.29, 572 (1958).

12h N. H. Fletcher, J. Chem. Phys. 31, 1136 (1959).

13h N. H. Fletcher, J. Chem. Phys. 38, 237 (1963).

14h V. Talanquer and D. W. Oxtoby,  J. Chem. Phys. 104, 1483 (1996).

15h V. Talanquer and D. W. Oxtoby,  J. Chem. Phys. 114, 2793 (2001).

16h T. V. Bykov and X. C. Zheng,  J. Chem. Phys. 117, 1851 (2002).

17h R. P. Sear, J. Phys. Condens. Matter 14, 3693 (2002).

18h S. Auer and D. Frenkel, Phys. Rev. Lett. 91, 015703 (2003).

19h A. Cacciuto, S. Auer and D. Frenkel, Nature 428, 404 (2004).

20h A. Cacciuto, and D. Frenkel, Phys. Rev. E. 72, 041604 (2005).

21h E. A. Ustinov, and D. D. Do, J. Phys. Chem. 109, 11653 (2005).

22h R. P. Sear, J. Phys. Chem. B. 110, 4985 (2006).

23h T. V. Bykov and X. C. Zheng,  J. Chem. Phys. 125, 144515 (2006).

24h A. J. Page and R. P. Sear, Phys. Rev. Lett. 97, 065701 (2006).

25h R. P. Sear, J. Phys.: Condens. Matter. 19, 033101 (2007).

26h D. Winter, P. Virnau, and K. Binder, J. Phys.: Condens. Matter. 21, 464118 (2009).

27h D. Winter, P. Virnau, and K. Binder, Phys. Rev. Lett. 103, 225703 (2009).

28h D. Deb, A. Winkler, P. Virnau, and K. Binder, J. Chem. Phys. 136, 134710 (2012).
 
29h F. Schmitz, P. Virnau, and K. Binder, Phys. Rev. E. 87, 053302 (2013).

30h K. S. Singha, P. K. Das, and B. Maiti. J. Chem. Phys. 142, 104706 (2015).

31h K. Binder, and P. Virnau, J. Chem. Phys. 145, 211701 (2016).

32h K.D. Sattler (ed.) Handbook of Nanophysics. Principles and Methods. CRC Press, Boca Raton, 2011.
(Wiley, New York 1995)

33h Th.  Ondarçuhu and J.-P. Aimé (eds.)
Nanoscale Liquid Interfaces: Wetting, Patterning, and Force Microscopy at the Molecular Scale.
Pan Stanford Publ. Pte. Ltd., Stanford 2013.

34h K. Binder, and D. Stauffer, Adv. Phys. 25, 343 (1976).

36h P. A. Rikvold, H. Tomita, S. Miyashita, and S. W.Sides, Phys. Rev. E. 49, 5080 (1994).

37h B. A. Berg, and S. Dubey, Phys. Rev. Lett. 100, 165702 (2008).
  
37hb E. N. M. Cirillo, and J. L. Lebowitz, 
J. Stat. Phys. 90, 211 (1998).

38h M. Rauscher and S. Dietrich, Ann. Rev. Mater. Sci. 38, 143 (2008)

39h M. Rauscher and S. Dietrich, in Ref. <cit.>, Chapter 11.

40h S. Dietrich, M. Rauscher and M. Napiorkowski, in Ref. <cit.>, Chapter 3.

41h K. Binder and D.W. Heermann, 
Monte Carlo Simulation in Statistical Physics. An Introduction. 5th Edition (Springer, Berlin, 2010)

42h D. Bonn, J.Eggers, J. Indekeu, J. Meunier, and E. Rolley, Rev. Mod. Phys. 81, 739 (2009).

43h K. Binder, D. P. Landau, and S. Wansleben. Phys. Rev. B 40, 6971 (1989).

44h E. V.Albano, and K. Binder. Phys. Rev. Lett. 109, 036101 (2012).

45h P. Bryk, and K. Binder. Phys. Rev. E. 88, 030401 (2013).

46h L. Onsager, Phys. Rev. 65, 117 (1944).

47h C.N. Yang, Phys. Rev. 85, 808 (1952).

48h D.B. Abraham, Phys. Rev. Lett. 44, 1165 (1980).

49h D.B. Abraham, in Phase Transitions and Critical Phenomena, Vol. 10
(C. Domb and J.L. Lebowitz, eds.) Academic Press, London, 1986, Ch. 1.

50h D.B. Abraham, J. De Coninck and F. Dunlop, Phys. Rev. B. 39, 4708 (1989).

51h M.E. Fisher, J. Stat. Phys. 34, 667 (1984).

52h X.-T. Wu, D.B. Abraham and J. O. Indekeu,  Phys. Rev. Lett. 116, 046101 (2016). 

53h M. L. Trobo, E. V. Albano, and K. Binder, Phys. Rev. E. 93, 052805 (2016).

54h K. Brendel, G. T. Barkema, and H. van Beijeren, Phys. Rev. E. 71, 031601 (2005). 

35h H. Tomita, and S. Miyashita, Phys. Rev. B. 46, 8886 (1992).  

55h P. Jakubczyk, M. Napiórkowski, and A.O. Parry, Phys. Rev. E 74, 031608 (2006).

56h P. Jakubczyk and M. Napiórkowski, J. Phys. A: Math. Theor. 40, 2363 (2007).

57h M. P. A. Fisher, D. S. Fisher, and J. D. Weeks, Phys. Rev. Lett. 48, 368 (1982).

60h R. P. K. Zia, and J. E. Avron, Phys. Rev. B 25, 2042 (1982).

61h C. Rottman, and M. Wortis, Phys. Rep. 103, 59 (1984).
  
62h W. L. Winterbottom, Acta Metall. 15, 303 (1967).

90 L. Schimmele, M. Napiórkowski and S. Dietrich, J. Chem. Phys.
  127, 164715 (2007). 

66h P. Lenz and R. Lipowsky. Phys. Rev. Lett. 80, 1920 (1998).

67h R. Lipowsky, P. Lenz and P. S. Swain. Coll. & Surf. 161, 3 (2000).

68h R. Lipowsky, M. Brinkmann, R. Dimova, T. Franke, J. Kierfeld and X. Zhang,
  J. Phys. Cond. Matt. 17, S537 (2005).

x1 N. J. Günther, D. A. Nicole and D. J. Wallance,
J. Phys. A: Math. Gen. 13, 1755 (1980).  

x2 S. Ryu and W. Cai. Phys. Rev. E. 81, 030601(R) (2010).

x3 A. Tröster, F. Schmitz, P. Virnau and K. Binder.
J. Phys. Chem. B. (Submitted). 

64h F. Schmitz, P. Virnau, and K. Binder, Phys. Rev. Lett. 112, 126701 (2014).
  
65h F. Schmitz, P. Virnau, and K. Binder, Phys. Rev. E. 90, 012128 (2014).
  


































































































