A Curie-Weiss Theory of the Continuum Widom-Rowlinson Model
Yuri Kozitsky Mikhailo Kozlovskii
math-ph math.MP





Heusler compounds – how to tune the magnetocrystalline anisotropy
    H. C. Herper
    January 19, 2018
=================================================================




A version of the continuum Widom-Rowlinson model is introduced and
studied. It is a two-component gas of point particles placed in
𝐑^d in which like particles do not interact and unlike
particles contained in a given vessel of volume V repel each other
with intensity a/V. This model is thermodynamically equivalent to
a one-component gas with multi-particle interaction. For both
models, a rigorous theory of a phase transition is presented and the
ways of its construction in the framework of the grand canonical
formalism are outlined.





                             





§ INTRODUCTION



The rigorous theory of phase transitions in continuum particle
systems has got much more modest results as compared to its
counterpart dealing with lattices, graphs, etc. In fact, there exist
only few models with local interactions in which the existence of a
liquid-vapor phase transition was rigorously proved. One of them is
the model introduced in <cit.> by B. Widom and J. S. Rowlinson in
which the potential energy of n point particles located at x_1,
…, x_n ∈𝐑^d is set to be θ[ W(x_1 , … ,
x_n) - n], where θ>0 is a parameter and W is the volume of
the area ∪_i=1^n B(x_i) covered by the balls of unite volume
centered at the corresponding particles. This model is
thermodynamically equivalent to a two-component system with binary
interactions in which the interaction between unlike particles is a
hard-core repulsion and is zero otherwise. In <cit.>, D. Ruelle
proved that this model undergoes a phase transition of first order
whenever d≥ 2. Later on, this result was extended in
<cit.>, see also <cit.> for a review. However, these
theories give not too much for understanding the details of the
phenomenon. No rigorous results are available on the behavior at the
phase-transition threshold. The very existence of such a threshold
remains unknown for this model. At the same time, for a number of
lattice models the mean field results essentially improve
understanding phase transitions in the corresponding models with
local interactions <cit.>. It is then quite natural to
develop the mean field theory of phase transitions also in continuum
systems. For the Widom-Rowlinson model, the first attempt to do this
was undertaken already in <cit.>. Assuming that the
particles are distributed in a given vessel “at random" the authors
heuristically deduced an equation of state (eq. (7.4) in <cit.>)
which manifests a first order phase transition. Unfortunately, such
and similar heuristic mean-field theories (sometimes called
naive, cf. <cit.>) are not free from
mathematical inconsistencies and other drawbacks of mathematical
nature that diminish the value of their results. The first steps in
developing a rigorous mean-field theory of a liquid-vapor phase
transition were made by J. L. Lebowitz and O. Penrose in <cit.>
where Kac-type interaction potentials were employed to obtain (in
the canonical formalism) an explicit form of the pressure-density
dependence that describes a first order phase transition. Another
way of developing a rigorous mean field theory is to use Curie-Weiss
interaction potentials, see <cit.>. For lattice systems,
this way was formulated as a coherent mathematical theory in the
framework of which thermodynamic phases are constructed as
probability measures on the spin-configuration spaces, see, e.g.,
<cit.>. The aim of this letter is to report the
results of a rigorous study of an analog of the Widom-Rowlinson
model with Curie-Weiss interaction potentials made in the grand
canonical formalism.





§  THE MODEL



We deal with infinite systems of point particles placed in the space
𝐑^d. If the particles do not interact, their state of
thermal equilibrium (phase) is a Poisson probability measure P_z.
The only parameter characterizing this measure is activity
z=e^μ, where μ = (physical chemical potential)/ k_B T,
k_B and T>0 being Boltzmann's constant and absolute temperature,
respectively. For a vessel Λ⊂𝐑^d of finite
volume V and a nonnegative integer n, the measure P_z assigns
the probability

    P_z(Γ_Λ,n) = (z V)^n/ n!exp( - z
    V)

to the event Γ_Λ,n: “Λ contains n
particles". Let now point particles of two types, 0 and 1, be placed
in the same space 𝐑^d. If they do not interact, their
state of thermal equilibrium is the Poisson measure P_z_0 ,z_1 =
P_z_0⊗ P_z_1 such that  the event Γ_Λ,
n_0×Γ_Λ, n_1 has probability

    P_z_0,z_1(Γ_Λ,n_0×Γ_Λ, n_1) =
    P_z_0(Γ_Λ,n_0) · P_z_1(Γ_Λ,n_1),

where P_z_i(Γ_Λ,n_i), i=0,1, are as in (<ref>)
and the event is “Λ contains n_0 particles of type 0
and n_1 particles of type 1". In this case, the particle
densities are ϱ_i = z_i = e^μ_i, i=0,1, and the
equation of state reads

    p = ϱ_0 + ϱ_1 = e^μ_0 + e^μ_1,

where p is the pressure in the system.

For interacting particles, phases are constructed as limits Λ→𝐑^d of local Gibbs measures P_z^Λ,Φ
(P_z_0,z_1^Λ,Φ for two-component systems) describing
the portion of the particles contained in the vessel Λ and
interacting with the energy Φ, see, e.g., <cit.>. In this
work,  we introduce two models that – like the Widom-Rowlinson
model – can be considered as two versions of the same model. The
first one is a two-component gas of point particles in
𝐑^d. For a vessel Λ⊂𝐑^d of
volume V, unlike particles contained in Λ repel each other
with intensity a/V>0, whereas like particles do not interact.
Thus, the potential energy of the collection of n_0 particles of
type 0 located at x_1^0, …, x_n_0^0∈Λ and of
n_1 particles of type 1 located at x_1^1, …, x_n_1^1∈Λ is

    Φ_V (x_1^0, …, x_n_0^0; x_1^1, …, x_n_1^1) =
     ∑_i=1^n_0∑_j=1^n_1a/V = a/Vn_0 n_1.

The grand canonical partition function of this collection then is

    Ξ_V (a, μ_0 , μ_1)
         =    ∑_n_0, n_1=0^∞1/n_0! n_1 !∫_Λ^n_0∫_Λ^n_1exp( μ_0 n_0 + μ_1 n_1 -
    a/Vn_0 n_1)d x_1^0⋯ d
    x_n_0^0 d x_1^1⋯ d x_n_1^1
       =    ∑_n_0, n_1=0^∞V^n_0 + n_1/n_0! n_1 !exp(
    μ_0 n_0 + μ_1 n_1 - a/Vn_0 n_1).

Here the interaction parameter a>0 and the chemical potentials
μ_i∈𝐑, i=1,2, include the reciprocal temperature
β and thus are dimensionless. The second our model is a
one-component system of point particles interacting as follows. For
a vessel Λ of volume V, the potential energy of the
collection of n particles located at x_0, … , x_n∈Λ
is set to be

    Φ_V (x_1 , … , x_n) = V θ[ 1 -
    exp( - a/V
     n)].

Here θ>0 is a parameter, similar to that in <cit.>
mentioned above. Then the corresponding grand canonical partition
function is

    Ξ_V (a, μ , θ)      =    ∑_n=0^∞1/n!∫_Λ^nexp( μ n - V θ[ 1 -
    exp( - a/V
     n)]
    )d x_1 ⋯ d
    x_n
       =    ∑_n=0^∞V^n/n!exp(
    μ n  - V θ[ 1 - exp( - a/V
     n)] ) 
       =    exp(
    - Vθ) Ξ_V(a,μ,lnθ).

The latter equality can readily be derived by summing out in
(<ref>) over n_1. The dependence of the pressure p in the
two-component system (resp. p in the one-component
system) on a and μ_i, i=0,1 (resp. on a, θ and
μ) is then obtained in the thermodynamic limit

    p = p(a,μ_0, μ_1)    =    lim_V→ +∞1/VlnΞ_V(a,μ_0,μ_1),
    p =  p(a,  θ, μ)    =    lim_V→ +∞1/VlnΞ_V(a,θ,μ),

which by (<ref>) yields p= p - θ. Thus, the
particle density ϱ in the one-component system and the
density ϱ_0 of the particles of type 0 in the two-component
system are related to each other by

    ϱ = ∂p/∂μ = ∂
    p/∂μ_0|_μ_0=μ,  μ_1 = lnθ =
    ϱ_0.

Therefore, the study of the one-component system amounts to studying
its two-component counterpart by employing (<ref>) and (<ref>).





§  THE RESULTS


In view of (<ref>), we will deal with three thermodynamic
variables a, μ_0, μ_1, and hence with the phase space

    ℱ= { (a, μ_0, μ_1): a ≥ 0, μ_0 , μ_1 ∈𝐑}.
 Our main result is the statement that ℱ
can be divided into three disjoint subsets, i.e., presented as

    ℱ= ℛ∪𝒞∪ℳ.

A point (a,μ_0, μ_1) belongs to one of these subsets according
to the number of global maxima of the function

    E(y) = f(a,μ_0 +y) + f(a,μ_1 - y) - y^2/2a,

where

    f(a,x)  =  a/2[ u(a,x)]^2 + u(a,x),    x∈𝐑.

Here u is a special function, which can be expressed through
Lambert's W-function<cit.>  as follows

    u(a,x) = 1/a W ( a e^x),      x ∈𝐑.

For a fixed a>0, the function 𝐑∋ x ↦ u(a, x)
can be obtained as the inverse to (0, +∞)∋ u ↦ x(u) =
a u + ln u, by which one gets that

    u(a,x) exp[a u(a,x)  ] = e^x,
          ∂/∂ x u(a, x)  =   u(a, x)/1 + a u(a,
      x).

By (<ref>) and (<ref>) it follows that the points of global
maximum of E are also its local maximum points and hence can be
obtained from the equation

    y = a u(a, μ_0 + y ) - a u(a, μ_1-y),    y∈𝐑.

The single-phase domain ℛ= ℱ∖
(𝒞∪ℳ), cf. (<ref>), consists of all
those (a, μ_0, μ_1) for each of which there exists a unique
global maximum of E (at some y_*). Then the unique phase
existing at this point (a, μ_0, μ_1)  is the Poisson state
P_z̃_0 , z̃_1, see (<ref>) and (<ref>), where

    z̃_0 = u(a,μ_0 + y_*) ,   z̃_1 = u(a,μ_1 -
    y_*).

The set

    ℳ := { (a, μ, μ) : a>0, μ > 1 - ln a}

consists of phase coexistence points, and

    𝒞:={ ( a , 1- ln a , 1 - ln a): a >0}

is the line of critical points. For (a,μ_0, μ_1) ∈ℳ (i.e., for μ_0 = μ_1 = μ > 1- ln a), the
function in (<ref>) has two equal maxima located at ±y̅(a,μ), where the order parameter y̅(a,μ)>0 is the
unique solution of the equation ψ (y ) = μ - (1-ln a) with

    ψ (y ) = y + y/e^y -1 - 1 + lny/e^y -1,   
    y>0 .

For small y>0, we have that ψ(y) = y^2/24 + o(y^2). Thus,

    y̅(a, μ) = √(24(μ-(1-ln a))) + o(μ-(1-ln a))
 for
small positive μ - (1-ln a). For y̅(a, μ)>0, there
exist two phases P_z̃^+ , z̃^- and
P_z̃^- , z̃^+, where

    z̃^± = u(a, μ±y̅(a,μ)).

The equation of state of the two-component system has the following
form, cf. (<ref>),

    p = a ϱ_0 ϱ_1 + ϱ_0 + ϱ_1,

where the densities ϱ_i are, cf.  (<ref>),

    ϱ_0 = u(a, μ_0 + y_*),   ϱ_1 = u(a, μ_1 -
      y_*).

Note that each ϱ_i depends on both μ_0, μ_1 and
ϱ_i = ∂ p/∂μ_i, i=0,1, cf.
(<ref>). Note also that the densities satisfy

    ϱ_0 = exp(μ_0 - a ϱ_1 ),   ϱ_1 =
    exp(μ_1 - a ϱ_0 ).

That is, due to the repulsion both densities are smaller than in the
free case a=0.

Let us turn now to the ground states of the two-component model
which one obtains by passing to the limit a→ +∞. To this
end, we consider ϱ_i, i=0,1, as differentiable functions
of a defined in (<ref>). Let ϱ̇_i, i=0,1, stand
for the corresponding a-derivatives. Differentiating both sides of
each equality in (<ref>) after some calculations we get

    ϱ̇_0 - ϱ̇_1 = a ϱ_0 ϱ_1/1 -
    a^2 ϱ_0 ϱ_1( ϱ_0 - ϱ_1 ).

The denominator here is positive by the fact that y_* used in
(<ref>) is the point of local maximum of E given in (<ref>).
For μ_0
> μ_1, we have that ϱ_0
> ϱ_1 for all a>0. Indeed, assuming
ϱ_0= ϱ_1 for some a>0, we then get by (<ref>)
that e^μ_0=e^μ_1, which contradicts the assumed inequality
μ_0
> μ_1.
Thus, by (<ref>) ϱ_0 - ϱ_1 is an increasing
function of a, which yields ϱ_0 - ϱ_1 ≥ϱ_0
- ϱ_1|_a=0 = e^μ_0 - e^μ_1. By (<ref>) and the
latter estimate we then get

    ϱ_1     =    ϱ_0 exp(- (μ_0 -μ_1) - a (ϱ_0 - ϱ_1)
     ) 
       ≤   ϱ_0 exp(- (μ_0 -μ_1) - a (e^μ_0 -
    e^μ_1)
     ).

Since ϱ_0 ≤ e^μ_0, see (<ref>), by (<ref>) we
get that aϱ_1 → 0, and hence ϱ_1 → 0, as a→
+∞. At the  same time, ϱ_0 ≥ϱ_1 +
(e^μ_0-e^μ_1), which by (<ref>) yields that ϱ_0
→ e^μ_0 as a→ +∞. By (<ref>) and (<ref>) we
thus conclude that the two-component system has two ground states:
P_z_0,0 and P_0,z_1. In each of them, there is only one free
component.

Turn now to the one-component system. Its equation of state reads

    p = a θϱ e^- a ϱ + ϱ - θ(1 - e^- a ϱ),

that can be obtained by (<ref>) and the formula p= p-
θ. Here ϱ is a function of μ∈𝐑
obtained from (<ref>), i.e., ϱ= u(a, μ+y_*). It is
increasing and continuous whenever θ≤ e/a. For θ
> e/a, ϱ makes a jump at μ=lnθ with one-sided
limits lim_μ→lnθ± 0ϱ = u(a, μ±y̅(a,μ)). That is, the system undergoes a first-order phase
transition with the increment of the density ϱ =
y̅(a,μ)/a. The expression in (<ref>) can also be used to
define p as a function of ϱ. Namely,
p is as in (<ref>) for ϱ≤z̃^-_0
and ϱ≥z̃^+_0, and p≡p_* := a z̃^+_0 z̃^-_0 +
z̃^+_0+ z̃^-_0 - θ for
ϱ∈[z̃^-_0 , z̃^+_0].   Note that the
equation of state in (<ref>) with a=1 formally coincides with
that found heuristically in <cit.>, in which, however, the
horizontal part p≡p_* should be found
from the Maxwell rule, see <cit.> for more detail.

The part of the phase diagram of the two-component system in the
plane in ℱ with constant a>0 is presented in Fig
<ref>.

Points from the grayed area correspond to the existence of three
solutions of (<ref>), one of which is y_*. Note that y_* >0
for μ_0 > μ_1. At the boundaries of this area (symmetric under
 μ_0 ↔μ_1), (<ref>) has only two
solutions. The upper branch of the boundary is described by the
equation 
    η = √(ξ^2 -1) + ln(ξ - √(ξ^2 -1)),
 where ξ = (μ_0 + μ_1 )/2 + ln a and η =
(μ_1 - μ_0)/2. For all points from the complement to the grayed
area, (<ref>) has only one solution. Note that y_* =0 for
μ_0 = μ_1 < 1 - ln a.








§ DERIVING THE RESULTS


Here we outline the way of deriving the results presented in the
preceding section. In <cit.>, these results are formulated and
proved as mathematical statements.

By means of the identity 
    - a/Vn_0 n_1 = -
    a/2Vn^2_0 - a/2Vn^2_1 + a/2V(n_0 - n_1)^2,

and by a standard Gaussian formula we rewrite (<ref>) as

    Ξ_V (a, μ_0 , μ_1) = √(V/2π a)∫_-∞^+∞exp(V E_V (y) ) d y,

with

    E_V(y) = f_V (a, μ_0 +y) + f_V(a, μ_1-y) - y^2/2a.

Here f_V is defined by the following relation

    exp( V f_V (a, x)) = ∑_n=0^∞V^n/n!exp( x n - a/2V n^2 ),

and thus is an infinitely differentiable function of x∈𝐑 for each fixed a>0 and V>0. Then so is E_V as a
function of y∈𝐑. In view of (<ref>) we have to find
the large V asymptotic in (<ref>). To this end we employ a more
advanced version of Laplace's method as E_V depends on V.
According to <cit.>, this amounts to proving that
the third x-derivative of f_V is bounded as V→ +∞,
which is the most challenging aspect of the theory. Taking the
x-derivative of both sides of (<ref>) we obtain that u_V (a,
x) := ∂ f_V(a,x) /
∂ x satisfies, cf. (<ref>),

    u_V (a,x) exp( ∫_0^1 u_V ( a, x- a/V t)
    dt)= exp( x -  a/2V).

On the other hand, also by (<ref>) we get that

    u_V (a,x)  = ⟨ n ⟩ /V = 1/V∑_n=1^∞ n
    π_n,

with 
    π_n =   V^n/n!exp( n x - an^2/2V) /∑_n=0^∞V^n/n!exp( n x - an^2/2V).

Then the consecutive x-derivatives of u_V are 
    u'_V (a, x) =
    1/V⟨ (n-⟨ n⟩)^2 ⟩,   u”_V (a, x)
    = 1/V⟨ (n-⟨ n⟩)^3 ⟩.
 These
formulas allow for obtaining uniform in V upper bounds for u'_V
(a, x) and |u”_V (a, x)|. Thereafter, we apply the mentioned
version of Laplace's method by which the problem of calculating p
in (<ref>) is reduced to finding the global maxima of E given in
(<ref>). This yields, in particular, that p is given by
(<ref>). We also prove that f_V → f, u_V → u and u'_V
→ u' as V→∞, where f, u and u' are as in
(<ref>), (<ref>) and (<ref>), respectively. The proof that
the limiting states are P_z̃_0, z̃_1 with
z̃_i given in (<ref>) is done by showing that the
correlation functions of the local Gibbs measures, see
<cit.>, converge as V→ +∞ to those of
P_z̃_0, z̃_1. Here we also employ the result of
<cit.> mentioned above, as well as the convergence of the
derivatives of f_V just mentioned. Then the representation of
ℱ in (<ref>) with ℳ and 𝒞 as
in (<ref>) and (<ref>), respectively, are obtained by studying
the global maxima of E.







§  CONCLUDING REMARKS


We recall that the parameters a, V, μ_0 and μ_1 in the
expression for Ξ_V in (<ref>) are dimensionless. By chosen
them we, in fact, fix some metric of the habitat space
𝐑^d. The change of this metric (scale) can be done by
passing to V^α := α V for some scale parameter
α>0. This, of course, leads to the change ϱ_i →ϱ_i^α = ϱ_i/α, i=0,1, and also to a→
a^α = α a and p→ p^α = p/α, cf. (<ref>)
and (<ref>), respectively. Then by (<ref>) we obtain that the
rescaled densities and the interaction parameter a^α satisfy

    ϱ_0^α = exp( μ_0^α - a^αϱ_1^α),     ϱ_1^α = exp(
    μ_1^α - a^αϱ_0^α),

with μ_i^α = μ_i - lnα, i=0,1. By this we
conclude that the description of the phases corresponding to the
points (a, μ_0, μ_1)∈ℱ is scale-invariant. Thus,
one can choose α = 1/a and hence consider the plane in
ℱ with a=1. Then the results corresponding to a
general point (a, μ_0, μ_1)∈ℱ can be obtained
from those obtained for a=1 by the rescaling as just described.
Note, however, that in view of taking the thermodynamic limit V→
+∞, the role of the scale of V may get be lost. Then
considering a as an interaction parameter helps to reveal it.
Moreover, from our analysis it follows that the phase transition
occurs for an arbitrarily small interaction a. Also for the
following reasons it might be worth to deal with general values of
a: (i) to be able to pass to the free case a=0; (ii) to obtain
the ground states in the limit a→ +∞ as described above;
(iii) to get clues on the phase transition in the original
Widom-Rowlinson model in which the metric is rather fixed by the
choosing the radius of the hard-core repulsion. In this model the
interaction energy of the collection of n_0 particles of type zero
located at x_1^0, … , x_n_0^0 with n_1 particles of type
one located at x_1^1, … , x_n_1^1 is written in the form,
cf. <cit.>,

    Φ(x_1^0, … , x_n_0^0; x_1^1, … , x_n_1^1) =
    ∑_k=1^n_0∑_l=1^n_1∞1_|x_k^0 -
    x_l^1|≤ r_d(x_k^0, x_l^1),

where 1_|x_k^0 - x_l^1|≤ r_d is the corresponding
characteristic function and r_d is the radius of the ball in
𝐑^d of unit volume. In contrast to our Φ_V given in
(<ref>) Φ takes values either zero or infinity, and hence no
parameter like our a can be associated with this model. Regarding
its thermodynamic phases the following is known. There exist μ_*,
μ^* ∈𝐑 such that μ_*<μ^* and: (a) for μ_0
=μ_1 =μ<μ_*, there exists only one thermodynamic phase; (b)
for μ_0 =μ_1 =μ>μ^*, there exist at least two
thermodynamic phases. In contrast to our description, nothing is
known about the threshold that might exist on the interval [μ_*,
μ^*]. Having in mind that for spin systems on 𝐙^d with
local (properly normed) interactions the mean-field theory of a
phase transition becomes exact in the limit d→ +∞, see
<cit.>, one might speculate that
also for the model in (<ref>) there exists a threshold value
μ_c = 1 - ln a_ WR (d)∈ [μ_*, μ^*] with a_
WR(d) satisfying a_ WR(d) → a_ WR as d→ +∞
for some a_ WR>0. If this is true, then our model with a
particular value a=a_ WR can be considered as a mean-field
limit of the original Widom-Rowlinson model.




ll
WR Widom B and Rowlinson J S, 1970 New model for the study of liquid-vapor phase transition, J. Chem.
Phys. 52 1670–1684

R Ruelle D, 1971 Existence of phase transition in a continuous classical system, Phys. Rev. Lett. 27 1040–1

Kot Chayes J T, Chayes L and Kotecky R, 1995 The analysis of the Widom-Rowlinson model by stochastic geometric
methods, Commun. Math. Phys. 172 551–569

GHM Georgii H-O, Häggström O and Maes C, 2000 The random geometry of equilibrium phases, in: Phase Transitions
and Critical Phenomena vol 18, ed C Domb and J L Lebowitz (New
York: Academic) pp 1–142

Biskup Biskup M, Chayes L and Crawford N, 2006 Mean-field driven first-order phase transitions in systems with
long-range interactions, J. Stat. Phys. 122 1139–1193

Simon Simon B, 1993 The Statistical Mechanics of Lattice
Gases. I, Princeton University Press,

LP Lebowitz J L and Penrose O, 1966 Rigorous treatment of the Van Der Waals-Maxwell theory of the liquid-vapor
transition, J. Math. Phys. 7 98–113


EN Ellis R S and Newman C N, 1978 The statistics of
Curie-Weiss models, J. Stat. Phys. 19 149–161


E Ellis R S, 2006 Entropy, Large Deviations, and
Statistical Mechanics, 2nd edition, Grundlehren der mathematischen
Wissenschaften, Springer


Zag Georgii H-O, Miracle-Solé S, Ruiz J and Zagrebnov
V A, 2006 Mean-field theory of the Potts gas, J. Phys A: Math.
Gen. 39 9045–9053

Ku  Külske Ch and Opoku A, 2008 Continuous spin mean-field models:
limiting kernels and Gibbs properties of local transforms, J.
Math. Phys. 49 125215–31


W Coreless R M, Gonnet G H, Hare D E G, Jeffrey D J and
Knuth D E, 1996 On the Lambert W function, Advanced in
Computational Mathematics 5 329–359

KK Kozitsky Y and Kozlovskii M, 2018 A phase transition in a Widom-Rowlinson model with Curie-Weiss
interaction, ArXiv:1801.01428

Fed Fedoryuk M V 1989 Asymtotic Methods in Analysis in Analysis I: Integral Representations and Asymptotic
Methods eds.
 Evgrafov M A and Gamkrelidze R V Encyclopaedia of Mathematical Sciences vol 13 (Springer-Verlag Berlin Heidelberg) 83-191

Ruelle  Ruelle D, 1999 Statistical Mechanics. Rigorous
Results World Scientific





