Das Sh.P. Statistical Physics of Liquids at Freezing and Beyond
Подождите немного. Документ загружается.
2.4 Applications to other systems 101
Fig. 2.7 A schematic diagram showing the spatial variations of B
m
(z) order parameters describing
the interfacial density profile given in eqn. (2.4.33) for an interface of width z. Adapted from
Curtain (1989).
c
American Physical Society.
the interfacial density function n
0
(x) depicted above the order parameters all begin their
decay from the point z = z
0
. The set of parameters for the optimum density distribu-
tion in the interface region is therefore reduced to
0
, and ν. In principle z
0
can also
be treated as a free parameter in characterizing the density profile. The function B
0
is
another variational parameter representing the change of average density in the interfa-
cial region on going from n
s
to n
l
. In the present context B
0
is chosen to be the same
as the order parameter corresponding to the shortest RLV K
1
, i.e., B
0
≡ B
1
. With the
above parametrization of the density function n
0
(x), the weight function w(x − x
) of
the WDA is determined by finding the self-consistent solution of the basic equation of the
WDA (2.2.6) or equivalently of (2.2.7) in the Fourier space. The weighted density ¯n(x) is
obtained using the relation (2.2.2) and the excess free energy F
ex
is computed from (2.2.5).
Note that, instead of using the WDA approach described above, the excess free energy can
be computed with a direct functional expansion (Oxtoby and Haymet, 1982). However, this
requires one to adopt the square-gradient approximation. While this assumption of slow
density variation allows an analytic solution of the problem, in the present context this
can be avoided in the WDA formulation through the direct minimization of free energy as
outlined above. The square-gradient approximation, which involves a different paramet-
ric function n
0
(x) for the interfacial density, will be discussed in the next chapter while
we treat the nucleation of a spherical droplet of crystal in the liquid melt. The DFT dis-
cussed here is strictly at the mean-field level. Hence the role of fluctuations, the existence
of long-wavelength capillary waves, and roughening transitions are ignored in the present
treatment.
The grand potential for the interface is minimized with respect to the parameters {
0
,ν}
of the density function. In Fig. 2.8 the results obtained for a hard-sphere system in terms
102 The freezing transition
Fig. 2.8 The equilibrium plane-averaged interfacial density ρ(z) defined in eqn. (2.4.35) for the
(100) f.c.c.–liquid interface. Solid line, the theoretical result for the hard spheres showing four non-
bulk layers; dashed line, simulation results for soft spheres showing six to seven non-bulk layers.
The centers of the theoretical and simulated profiles have been approximately aligned. From Curtain
(1989).
c
American Physical Society.
of the plane averaged density profiles at the interface are shown. The average density is
defined as
n
0
(z) =
1
A
dx
dy n
0
(x), (2.4.35)
where A is the area of the interface.
The interface considered here is a (100) plane of the crystal. The value of z
0
appearing in
the definition of n
0
(x) is set by choosing it to be at the middle of the two layers of the crys-
tal interfaces. The Percus–Yevick approximation is used throughout for the structure factor
and equation of state. This also allows use of the coexistence conditions stated in Table 2.2.
The surface energy is given by γ
srf
= /A, where A is the area of the interface. The
numerical computation indicates that the situation of an interface with four atomic lay-
ers each of width a/2 (where a is the cubic-lattice spacing) and the exponent ν = 0.25
corresponds to the state of lowest free energy (Curtain, 1989). The surface free energy is
2.4 Applications to other systems 103
Fig. 2.9 Layer-by-layer contributions to the surface free energy βγ σ
2
of the equilibrium hard-sphere
(111) f.c.c.–liquid interface. The interface spreads between layer 16 and layer 19. From Curtain
(1989).
c
American Physical Society.
computed to be γ
srf
= 0.66 ± 0.02 in units of βσ
2
, where σ is the hard-sphere diameter.
The (111)-plane-averaged interface is also found to have the lowest surface free energy
corresponding to the same conditions, namely a width of four layers, each having width
a/
√
3, and an exponent ν = 0.25. The surface free energy γ
srf
= 0.63±0.02 is nearly equal
to that of the (100) plane. An analysis of the layer-by-layer contribution to the surface free
energy of the (111) interface is presented in Fig. 2.9. The largest contribution to the free
energy occurs in the middle of the interface. There is also an asymmetry in the free-energy
contributions, with the near-crystal region costing more free energy than that near the liq-
uid side. The results of the density-functional model for the interfaces show very little
anisotropy, with the surface free energies differing by less than 5%. However, it is useful to
note that the two-parameter form of the density variation depicted above does not include
anisotropy. The above results on the structure of the interfaces are comparable to those
obtained from computer simulations of similar systems. For soft-sphere systems with an
r
−12
interaction potential the (100) (Cape and Woodcock, 1980) and (111) (Tallon, 1986)
planes and for the Lennard-Jones system the (100) and (111) planes all obtain interfaces
that consist of six or seven atomic layers of density variation. The surface free energy γ
srf
for the hard-sphere system computed above can be shown to be in reasonable agreement
with simulation results for the soft-sphere systems by expressing the simulation results in
terms of an equivalent hard-sphere system (Curtain, 1989). The above method of studying
the interface is extended to the Lennard-Jones systems by computing the contribution due
to the attractive part of the interaction in the mean-field approximation. Broughton and
Gilmer (1986) obtained the surface tension from simulation of the Lennard-Jones system,
which showed that the (111), (100), and (110) faces have the same surface free energies.
The work mentioned is discussed in Chapter 3, in which we treat nucleation.
104 The freezing transition
The basic step in the weighted-density-functional theories involves a mapping of the
inhomogeneous system (crystal) into a coarse-grained system that has inhomogeneity over
a different scale (WDA) or has constant density (MWDA). This procedure represents a
reduction in number of the associated degrees of freedom at the shortest length scales. The
basic interaction in the original and coarse-grained system, however, remains the same
here. The present treatment of the first-order freezing transition is in some way anal-
ogous to the renormalization-group (RG) approach (Wilson and Kogut, 1974; Wilson,
1975) which is applied to the theory of critical phenomena of second-order phase tran-
sitions. For first-order transitions application of the RG has not been common (Nienhuis
and Nauenberg, 1975).
Appendix to Chapter 2
A2.1 Correlation functions for the inhomogeneous solid
We demonstrate below that the two-point function c
(2)
(x
1
, x
2
) for the inhomogeneous state
is related to the pair correlation in the fluid. The one-body potential c
(1)
(x
1
, n
0
(x)) defines
the inhomogeneous density n(x
1
) through the relation (2.1.33) or equivalently through
(2.1.31). By taking a functional derivative of (2.1.31) with respect to n
0
(x
2
) we obtain the
following result for the two-point function c
(2)
(x
1
, x
2
;n
0
(x)):
c
(2)
(x
1
, x
2
;n
0
(x)) ≡
δ
δn
0
(x
2
)
c
(1)
(x
1
;n
0
(x)) = β
δu(x
1
)
δn
0
(x
2
)
−
δ(x
1
− x
2
)
n
0
(x
1
)
. (A2.1.1)
However, using the definition (2.1.9) of the equilibrium density n
0
(x
1
) and taking one more
functional derivative of it with respect to the field u(x
2
), we obtain the following result for
the average of the product of the microscopic densities ˆn at two points,
β
−1
δn
0
(x
1
)
δu(x
2
)
=ˆn(x
1
) ˆn(x
2
)
0
− n
0
(x
1
)n
0
(x
2
), (A2.1.2)
where the microscopic density ˆn(x) is defined as in (2.1.6) and the subscript 0 on the angu-
lar bracket refers to the equilibrium average with the grand-canonical ensemble. The first
term on the RHS of (A2.1.2) can be expressed in terms of one- and two-point contributions
as follows:
ˆn(x
1
) ˆn(x
2
)
0
=
/
α,β
δ(x
1
− x
α
)δ(x
1
− x
β
)
0
0
=
/
α=β
δ(x
1
− x
α
)δ(x
1
− x
β
)
0
0
+ δ(x
1
− x
2
)ˆn(x
1
)
0
= n
0
(x
1
)n
0
(x
2
)g(x
1
, x
2
) + n
0
(x
1
)δ(x
1
− x
2
), (A2.1.3)
105
106 Appendix to Chapter 2
where we have used the definition (1.2.58) for the pair correlation function introduced in
Chapter 1.Using(A2.1.3), the above relation (A2.1.2) reduces to the form
β
−1
δn
0
(x
1
)
δu(x
2
)
= h(x
1
, x
2
)n
0
(x
1
)n
0
(x
2
) + n
0
(x
1
)δ(x
1
− x
2
). (A2.1.4)
We now use the two functional derivatives in eqns. (A2.1.1) and (A2.1.4) in the functional
identity (u(x) and n
0
(x) being uniquely related to each other)
dx
3
δu(x
1
)
δn
0
(x
3
)
δn
0
(x
3
)
δu(x
2
)
=
δu(x
1
)
δu(x
2
)
= δ(x
1
− x
2
), (A2.1.5)
to obtain the following integral relation of the two-point kernel c
(2)
(x
1
, x
2
) with h(x
1
, x
2
) =
g(x
1
, x
2
) − 1:
h(x
1
, x
2
) = c
(2)
(x
1
, x
2
) +
dx
3
{c
(2)
(x
1
, x
3
)}h(x
3
, x
2
)n
0
(x
3
). (A2.1.6)
For the isotropic homogeneous liquid with constant density the relation (A2.1.6) reduces to
the standard Ornstein–Zernike relation (1.2.89) presented in Chapter 1. The kernel function
c
(2)
appearing in eqn. (2.1.37) is now evaluated in terms of the pair correlation function
g(r) or the static structure factor S(k) for the isotropic liquid. These structure functions are
known from various integral-equation theories of the liquid state (discussed in Chapter 1),
from experiments or from computer simulations, and form the starting point of the DFT for
the freezing transition. Similarly, any one of the higher-order generalizations of the direct
correlation function (as defined in eqn. (2.1.36))(Barrat et al., 1987, 1988) is also related
to the corresponding multi-particle correlations in the fluid.
A2.2 The Ramakrishnan–Yussouff model
We consider here the deduction and solution scheme of the density-functional equations in
the Ramakrishnan–Yussouff model. We begin with the form of the inhomogeneous density
function given by
n
0
(x) = n
l
1 + η +
∞
m=1
A
m
e
iK
m
· x
(A2.2.1)
defining the crystalline state in terms of the order parameters η and {A
m
}. By taking a
Fourier transform of (2.1.37) the following set of equations for the order parameters η and
A
m
is obtained:
1 + η = e
˜c
0
η
dx
V
exp
∞
m=1
˜c
m
A
m
e
iK
m
· x
, (A2.2.2)
A
j
= e
˜c
0
η
dx
V
e
−iK
j
· x
exp
∞
m=1
˜c
m
A
m
e
iK
m
· x
. (A2.2.3)
A2.2 The Ramakrishnan–Yussouff model 107
To simplify the discussion of the above equations, we define the integrals κ
m
in terms of
the rescaled amplitudes ξ
m
=˜c
m
A
m
as
κ
j
=
dx
V
e
−iK
j
· x
exp
⎡
⎣
s=0
ξ
s
e
iK
s
· x
⎤
⎦
. (A2.2.4)
Equations (A2.2.2) and (A2.2.3) are now obtained in the simplified forms
1 + η = e
˜c
0
η
κ
0
, (A2.2.5)
ξ
j
˜c
j
=
[
1 + η
]
κ
j
κ
0
. (A2.2.6)
Equations (A2.2.5) and (A2.2.6) together with the definition (A2.2.4) form a set of coupled
equations for the order parameters η and A
m
of the freezing transition.
At the phase-transition point the grand potentials for the two phases are equal. By
substituting the ansatz (A2.2.2) into the RHS of (2.1.41), the difference of the thermo-
dynamic potential between the two phases is obtained in terms of the order parameters η
and A
m
. Using for the integral the result
dx
1
dx
2
c
(2)
(x
1
, x
2
;n
l
)δn
0
(x
1
)δn
0
(x
2
) = n
l
V
⎡
⎣
˜c
0
η
2
+
m=0
ξ
2
˜c
m
⎤
⎦
(A2.2.7)
and some trivial algebra, we obtain from (2.1.41) the following expression for :
= n
l
V
⎡
⎣
( ˜c
0
− 1)η +
1
2
˜c
0
η
2
+
1
2
m=0
ξ
2
m
˜c
m
⎤
⎦
. (A2.2.8)
It is useful to note here that the RLV |K
j
| is expressed in units of the inverse of the lattice
constant a and the direct correlation function is expressed as a function of |K
j
|σ , hence
the ratio σ/a enters the numerical computation. Since
σ
a
=
(1 + η)n
l
σ
3
3
1/3
(A2.2.9)
the input structural entity c
j
on the RHS of (A2.2.8) is dependent on η.
The next step in the calculation is to use a set of RLVs for the crystalline structure.
The choice of a proper set of RLVs is related to the symmetry of the chosen crystal lattice
whose stability is being tested. For an f.c.c. lattice, for example, the first set of eight nearest-
neighbor RLVs is given by
K
1
=
2π
a
{±1, ±1, ±1}, (A2.2.10)
which is a set of eight vectors, each of length |K
1
|=(2π/a)
√
3. In general the RLV K
m
is expressed in terms of its components,
K
m
=
2π
a
{±m
x
, ±m
y
, ±m
z
}, (A2.2.11)
108 Appendix to Chapter 2
where m
x
, m
y
, and m
z
are integers. The length of this vector is |K
m
|=(2π/a)
:
m
2
x
+ m
2
y
+ m
2
z
. We adopt a notation in which the RLV K
m
denotes a set of vectors
K
mα
, α = 1, 2,...,2p (p is an integer), forming a shell of radius K
m
. For all the RLVs
belonging to a particular shell s (say) the common amplitude is denoted by A
s
. Hence all
the ξ
i
belonging to this shell are assumed to be the same and denoted by ξ
s
. Using this
symmetry of the different RLVs belonging to a particular shell s, we simplify the sum
5
s=0
ξ
s
e
iK
s
· x
(appearing on the RHS of eqn. (A2.2.4)) and obtain the following result for
κ
0
corresponding to the case K
0
= 0:
κ
0
=
dx
V
exp
s
ω
s
(x)ξ
s
, (A2.2.12)
where
ω
m
(x) =
α
exp[iK
mα
· x]≡
α
exp[−iK
mα
· x]. (A2.2.13)
Thus, for example, with {K
1α
} being the first set of eight lattice vectors denoted by eqn.
(A2.2.10), we obtain that
ω
1
(x) = 8 cos(π x)cos(π y)cos(πz) (A2.2.14)
and so on. We obtain κ
0
using the result (A2.2.13), giving that
κ
0
=
1
0
dx
1
0
dy
1
0
dz exp
s
ω
s
(x)ξ
s
, (A2.2.15)
where the volume integral is performed over one primitive cell of the lattice and the index
s in the sum denotes the different concentric shells in the RLV space. For computing κ
j
for
j = 0, we make use of the relation
κ
j
=
∂
∂ξ
j
κ
0
, (A2.2.16)
which follows from using the symmetry of the RHS of (A2.2.4) for all the RLVs {K
j
}
belonging to a particular shell. With a suitably chosen set of RLVs extending to a certain
number of shells, the corresponding amplitudes κ
j
and η are obtained from numerical
solution of eqns. (A2.2.5) and (A2.2.6) together with the defining relations (A2.2.4) and
(A2.2.16). The average solid-state density n
s
and hence the lattice constant (assuming the
solid to be free of defects) is obtained in terms of the fractional density change,
n
s
= n
l
(1 + η) =
4
a
3
(A2.2.17)
(for the f.c.c. structure). The lattice constant thus follows self-consistently from the mini-
mization of the free-energy functional.
The results are listed in Table A2.4.
A2.3 The weighted-density-functional approximation 109
Table A2.4 Freezing parameters for hard spheres in the
Ramakrishnan–Yussouff theory, using the Percus–Yevick structure
factor (with Verlet–Weiss correction) for the uniform liquid state.
From Haymet and Oxtoby (1986).
Number of Lattice vectors Order Lindemann
RLVs N
n
(±i, ±j, ±k) parameter A
s
ratio L
0
1 0 0.0907
8 (1, 1, 1) 1.0436 0.0668
6 (2, 0, 0) 1.0299 0.0660
12 (2, 2, 0) 0.9780 0.0644
24 (3, 1, 1) 0.9427 0.0635
8 (2, 2, 2) 0.9309 0.0633
6 (4, 0, 0) 0.8857 0.0629
24 (3, 3, 1) 0.8540 0.0626
24 (4, 2, 0) 0.8444 0.0624
24 (4, 2, 2) 0.8053 0.0620
8 (3, 3, 3) 0.7767 0.0618
24 (5, 1, 1) 0.7775 0.0617
12 (4, 4, 0) 0.7324 0.0615
48 (5, 3, 1) 0.7081 0.0612
24 (4, 4, 2) 0.6995 0.0612
6 (6, 0, 0) 0.7018 0.0610
24 (6, 2, 0) 0.6701 0.0608
24 (5, 3, 3) 0.6464 0.0608
24 (6, 2, 2) 0.6398 0.0607
8 (4, 4, 4) 0.6103 0.0606
A2.3 The weighted-density-functional approximation
Equations (2.2.1) and (2.2.2), which form the basis of the WDA, can in fact be justified
from general considerations (Ashcroft, 1996). In order to demonstrate this we start from
the following definition of functional derivatives of the excess free energy:
F
ex
[n
λ+λ
(x)]−F
ex
[n
λ
(x)]=
dx
δF
ex
[n]
δn
λ
(x)
n
λ
(x), (A2.3.1)
wherewehaveusedforthedeviationn
λ
(x) ≡ n
λ+λ
(x) − n
λ
(x). The parameter λ
associated with the density function n
λ
(x) specifies the nature or the “path” along which the
density is being varied. It follows from the definition (2.1.36) that the functional derivative
on the RHS of (A2.3.1) is the one-particle distribution function c
(1)
[x;n
λ
(x)]. A suitable
expression for the excess free energy F
ex
is obtained by formally integrating the density
along the “path” n
λ
(x) = λn
0
(x) from λ = 0 to 1. Along this “path”
n
λ
(x) =
∂n
λ
∂λ
dλ = n
0
(x)dλ. (A2.3.2)
110 Appendix to Chapter 2
The expression (A2.3.1) for the excess free energy reduces to the form
F
ex
=−
dx
1
0
dλ c
(1)
[x;λn
0
(x)]n
0
(x). (A2.3.3)
In obtaining the result (A2.3.3) we also used the fact that the excess free energy F
ex
is zero
in the limit of zero density. In an exactly similar way, by applying the definition (2.1.36)
for the case n = 2 at the appropriate density, the following formal expression for the
one-particle function c
(1)
[x, n
0
(x)] is obtained:
c
(1)
[x;n
0
(x)]=c
(1)
[x;¯n
0
(x)]+
1
0
dλ
dx
c
(2)
[x
1
, x
2
;n
λ
(x)]n
λ
(x), (A2.3.4)
where the functional integration in (A2.3.4) is performed along the “path”
n
λ
(x) =¯n
0
(x) + λ
{n
0
(x) −¯n
0
(x)} (A2.3.5)
from λ
= 0to1,¯n
0
(x) being some suitable density function to be further specified below.
Here c
(1)
[x;¯n
0
(x)] is the one-particle direct correlation function of a homogeneous liquid
of density ¯n
0
(x). On substituting the exact expression (A2.3.4) for the one-particle function
c
(1)
into (A2.3.3) we obtain for the excess free energy
F
ex
=
dx n
0
(x) f
0
[¯n
0
(x)]
− β
−1
dx
dx
n
0
(x)
dλ
dλ
c
(2)
[x
1
, x
2
;λn
λ
(x)]n
λ
(x), (A2.3.6)
where the function f
0
( ¯n
0
) in the first term represents the excess free energy of the homo-
geneous liquid at density ¯n
0
. Using the definition for n
λ
along this chosen path (A2.3.5),
F
ex
is obtained as
F
ex
=
dx n
0
(x) f
0
[¯n
0
(x)]
− β
−1
dx
dx
n
0
(x)
!
n
0
(x
) −¯n
0
(x)
"
W[x, x
;n
0
(x), ¯n
0
(x)], (A2.3.7)
where we have adopted the following definition for the function W:
W[x, x
;n
0
(x), ¯n
0
(x)]≡
1
0
dλ{λ}
1
0
dλ
c
(2)
(x, x
;λn
λ
). (A2.3.8)
The density in the argument of c
(2)
is λn
λ
(x) = λ{¯n
0
(x) + λ
[n
0
(x) −¯n
0
(x)]}.Ifwe
demand that the excess free energy is given by only the first term on the RHS of (A2.3.7)
in a form similar to (2.2.1), then the second term must go to zero. Using this condition,