11institutetext: Yasutaka Yamaguchi 22institutetext: Department of Mechanical
Engineering, 2-1 Yamadaoka, Suita 565-0871, Japan 22email:
<EMAIL_ADDRESS>
# Introduction to Molecular-Scale Understanding of Surface Tension

Yasutaka Yamaguchi  Hiroshi Kawamura

###### Abstract

In this article, microscopic understanding of the surface tension are
provided, which needs basic knowledge of thermodynamics, statistical mechanics
as well as continuum mechanics. By introducing the intermolecular interaction
potential and temperature definition, and by showing conceptual pictures
including some results obtained by molecular dynamics simulations, the author
hopes that the target readers of undergraduate level students can find
fascinating aspects of surface tension as the boundary of macroscopic and
microscopic physics.

##  1 Introduction

[scale=.80]./fig01-interface.eps

Figure 1: Schematic for the molecular scale explanation of the surface
tension.

Liquids are supposed to be imcompressible due to their small compressibility,
e.g., its value of water is about $5\times 10^{-10}$ Pa-1. Thus, when liquids
are in an open space or in a closed space with a volume larger than theirs,
they dispose solid-liquid and/or solid-gas interfaces, hence, liquids we see
in our daily life almost always have interface. Interfacial tensions are
termed as the force exerted on such interfaces.

At present, the existence of molecules is commonly accepted, and the mechanism
of the liquid-vapor (LV) or liquid-gas (LG) surface tension, known as a
tensile force at the interface, is usually explained by using a schematic as
exemplified in Fig. 1. As in this figure, the molecules around the liquid-
vapor interface have less partners to interact with than those in liquid bulk,
and the liquid-vapor interface is generally disadvantageous to liquid bulk.
This seems to be simple and intuitively understandable; however, as Ref.
Marchand2011  pointed out, Fig. 1 with the arrow directing outward the
molecules gives an impression that molecules in the liquid bulk and at the
liquid-vapor interface ‘pull’ each other, and the molecule at the interface is
subject to a net force along the direction perpendicular to the interface, not
a force parallel to the interface. This kind of confusing questions are often
posed especially for interfacial phenomena, and indeed, I personally think
this is one of the fascinating aspects of interface physics which includes
macroscopic and microscopic features.

In this article, basic pictures for understanding of the surface tension are
provided with some simulation results, This concept is similar to Refs.
Marchand2011  and deGennes2003 , basically targeting undergraduate level
students, contrary to many other textbooks on this subject; however, I believe
that basic but wide and throughout outlook of thermodynamics, statistical
mechanics and continuum mechanics is needed to understand the physics of
interface, and more importantly, to avoid misleading which also have trapped
professional scientists Gao2009 . Partly related to this point, molecular
dynamics (MD) simulations have become common especially with the development
of high-performance computers and simulation packages including LAMMPS lammps
and GROMACS gromacs  from the beginning of the 21st century, and we can easily
have access to the movies of various (colorful) molecules. The great
scientists in the history could not have such movies, and they had to imagine
basically from macroscopic experimental results; however, this does not mean
that we will not be trapped by the misunderstanding. To advance science, we
have to make use of both the robust theories constructed by the great
scientists in the history and computer simulations to extend our
understanding.

##  2 Brief history

[scale=.80]./fig02-young_balance.eps

Figure 2: Schematic of Young’s equation proposed as the surface lateral force
balance of interfacial tensions exerted on the contact line.

Although there are several opinions about the origin of the study of
capillarity or surface tension deGennes2003 ; Rowlinson2002 ; Ono1980 ,
scientific and mathematical approach toward its understanding through the
modeling of wetting started from the beginning of the 19th century with T.
Young (1773–1829) and P.-S. Laplace (1749–1827). In 1805, Young proposed the
following equation Young1805

$\gamma_{\mathrm{SL}}-\gamma_{\mathrm{SV}}+\gamma_{\mathrm{LV}}\cos\theta=0,$
(1)

where $\gamma_{\mathrm{SL}}$, $\gamma_{\mathrm{SV}}$ and
$\gamma_{\mathrm{LV}}$ are the solid-liquid, solid-vapor, and liquid-vapor
interfacial tensions, repectively, and $\theta$ is the contact angle as the
angle between the SL and LV interfaces. Equation (1) is called Young’s
equation. Although his original article Young1805 , indeed in “essay” style,
has only one table and gives no equation nor figure, Eq. (1) was originally
suggested to express the force balance parallel to the solid surface exerted
on the contact line as shown in Fig. 2. Note that the proposal of Young’s
equation (1) is sometimes referred to as in 1804, because he made a talk about
this topic in this year. Young’s equation (1) was extended to evaluate
wettability of a liquid on a solid surface, known as the Young-Dupré equation:

$S=\gamma_{\mathrm{LV}}(\cos\theta-1),$ (2)

where the spreading coefficient $S$ given by

$S=\gamma_{\mathrm{SV}}-(\gamma_{\mathrm{SL}}+\gamma_{\mathrm{LV}})$ (3)

expresses the wettability: $S=0$ gives $\cos\theta=1$ in Eq. (2) and positive
$S$ means that the liquid completely covers the solid surface.

The year of 1804 or 1805 was before the establishment of thermodynamics: N. L.
S. Carnot (1796–1832) was under ten, and J. P. Joule (1818–1881), W. Thomson
(Lord Kelvin, 1824–1907) R. Clausius (1822–1888), and J. W. Gibbs (1839–1903)
were not born, and this clearly means that Young did not and could not bring
the concept of thermodynamics nor molecular interaction into his modelling of
wetting and surface tension as shown in Fig. 1 Gao2009 . Note that the primary
hydrodynamic description summarized in “Hydrodynamica” by D. Bernoulli
(1700–1782) DBernoulli1738 , which was based on Newtonian mechanics by I.
Newton (1642–1727), was already published at that time, and the Euler
Equations about ideal fluids by L. Euler (1707–1783) as the mathematical
framework of fluid mechanics as well as the Lagrangian mechanics by J.-L.
Lagrange (1736–1813) were already available Fujiwara1995 . Probably, continuum
mechanics including the concept of stress tensor was not available considering
the ages of A.-L. Cauchy (1789–1857), C. L. M. H. Navier (1785–1836) and S. D.
Poisson (1781–1840).

Anyway, by the frontiers including the above-mentioned scientists,
thermodynamics was established in the 19th century, and based on it, the van
der Waals equation of state expressing the phase coexistence was proposed by
J. D. van der Waals (1837–1923). He also wrote an article about surface
tension in 1893 vdW1893 ; Rowlinson1979 , and introduced the concept of
dividing surface and formulated a thermodynamic framework with a liquid-vapor
interface. Just after the tragedy of L. Boltzmann (1844–1906), J. Perrin
(1870–1942) experimentally proved the existence of molecules inspired by the
idea of A. Einstein (1879–1955) Fujiwara1995 , and this enabled the above-
mentioned explanation of the surface tension through the interaction potential
between molecules interacting with a potential well depth, i.e., with a
moderate attraction force for long intermolecular distance range and a strong
repulsion for short distance range described in detail below. Based on this
interaction potential model, G. Bakker (1856–1938), R. C. Tolman (1881–1948),
J. G. Kirkwood (1907–1959), F. P. Buff (1924–2009), etc. constructed the
theoretical framework of the surface tension in thermodynamic equilibrium.
Note that Japanese scientists including S. Ono (1918–1995), S. Kondo
(1922–2014), A. Harashima (1908–1986), etc. largely contributed to the
development. Especially, Ono and Kondo wrote a great review Ono1960 , and they
categorized the approaches of surface tension into thermodynamic and quasi-
thermodynamic ones, statistical mechanical ones and hydrodynamic (also
mentioned as “mechanical” in Ono1980 ) ones. Ono also left an instructive
textbook in Japanese Ono1980 , and he wrote in the introduction of the book:
“surface is very thin but still has thickness, and it is a pity that its
physics is difficult to understand because of this fact.” This sentence
clearly points out the difficulty of surface, and also shows the fascinating
feature of surface physics. Another great textbook by J. S. Rowlinson and B.
Widom Rowlinson1982  is also available; indeed, the history in the present
article mentioned above is based on these books Ono1980 ; Rowlinson1982 .

As described above, it is basically necessary to include the microscopic
concept to fully explain or understand the mechanism of surface tension and
wetting. On the other hand, as we see water droplets on solid surfaces almost
everyday, such a liquid motion affected by the surface tension is very common
and macroscopic, i.e., visible by bare eyes. As Young modelled, it is possible
to understand and predict the macroscopic liquid behavior by simply
considering the surface tension as the force to reduce the surface area
without the above-mentioned microscopic knowledge. Indeed in the middle of
18th century, an experimental approach to measure the surface tension was
proposed by L. F. Wilhelmy (1812–1864). In this method called the Wilhelmy
plate method Wilhelmy1863 , which is commonly used also at present, a solid
plate or a cylinder is vertically immersed into a liquid pool, and the surface
tension is evaluated by measuring the contact angle and the vertical force
exerted on the solid. This standpoint is equivalent to the above-mentioned
hydrodynamics approach that Ono and Kondo categorized, in which the
microscopic physics of interface with non-zero thickness is integrated into
the surface tension on a surface with zero thickness. As the name indicates,
this hydrodynamic implementation matches well with the governing equations of
fluid dynamics, and is especially familiar with the numerical simulation
technique called the computational fluid dynamics (CFD). In practice,
thermodynamic and statistical mechanical approaches cannot be applied for
dynamic systems in principle because these approaches assume static
thermodynamic equilibrium. Anyway, CFD simulations of liquid flow with moving
interface are common at present owing to the development of high-speed
computers from 1980s.

Parallel to the CFD, the (classical) molecular dynamics (MD) method, which
solves the motion of constituent molecules of fluids governed by inter-
molecular potential functions based on the Newtonian mechanics, has also been
developed, and this method enabled the simulations of liquid behavior with
interfaces as a microscopic approach. Indeed, one can obtain an equilibrium
liquid-vapor coexistence system with an interface by locating a certain number
of molecules in a simulation cell, whose intermolecular force is given by a
potential function exemplified in Fig. 3, and by controlling the temperature
of the system. As pointed out in Ref. Rowlinson1982 , this solves the main
difficulty in the statistical mechanical approach of how to analytically
obtain the equilibrium density distribution in such liquid-vapor coexistence
system from the microscopic intermolecular potential function.

##  3 Why phase separation happens and why interface is formed – from
thermodynamics –

[width=7.5cm]./fig03-potential.eps

Figure 3:  Intermolecular interaction between ideal gas molecules and between
real molecules. No intermolecular force acts between ideal gas molecules and
they have no volume. On the other hand, real molecules are generally attracted
by each other for an intermediate interatomic distance and repelled for a
smaller distance around the molecular diameter, and they have almost no
interaction for a longer distance. The minimum point of the potential energy
is the equilibrium distance.

###  3.1 Ideal gas model and temperature

In order to understand the properties of the liquid-vapor interface, a
primitive question of why the constituent molecules try to gather to form
liquid and vapor phases having different densities as illustrated in Fig. 1
instead of forming a single homogeneous phase with a uniform density. Indeed,
a fluid which always stays single homogeneous phase exists as a conceptual
one: an ideal gas. The principal (macroscipic) definition of ideal gases is
that they always obeys the equation of state:

$pV=nRT,$ (4)

where $p$, $V$, $n$, $T$ and $R$ denote the pressure, the volume, the amount
of substance of the gas in moles, the absolute temperature, and the universal
gas constant, respectively. Indeed, Eq. (4) holds under high temperature and
low pressure conditions for real gases, and an ideal gas model is an extended
imaginary substance which obeys Eq. (4) irrespective of the temperature and
pressure. D. Bernoulli was the first to propose the basis for the kinetic
theory of gases in his book of Hydrodynamica DBernoulli1738  in the 18th
century. In this book, he argued that a gas consists of a huge number of
molecules moving in all directions, and their impact on a surface causes the
pressure of the gas, and thus their average kinetic energy determines the
temperature of the gas under the following assumptions:

  * •

The constituent molecules of the gas are infinitesimally small hard spheres,
besides they are subject to elastic collisions among each other and with the
surroundings (container wall).

  * •

There are no attractive or repulsive forces between the molecules apart from
those that determine their point-like collisions, i.e., the only forces
between the gas molecules and the surroundings are impulsive force upon the
point-like collisions.

  * •

The molecules are constantly moving, and as a result of multiple random
collisions, their velocity distribution reaches a certain isotropic random
directions as an equilibrium state.

Almost after a century, the idea was revisited in the 19th century by Clausius
and J. C. Maxwell (1831–1879), etc. From the assumption above, it is possible
to relate the mean square velocity of the constituent molecules
$\left<|\mbox{\boldmath$v$}|^{2}\right>$ with the pressure $p$ with the
density $\rho$ of the gas in three-dimension as

$p=\frac{\rho\left<|\mbox{\boldmath$v$}|^{2}\right>}{3}$ (5)

by considering the sum of impulse exerted from the molecules on the container
wall of unit area per unit time, where
$\left<|\mbox{\boldmath$v$}|^{2}\right>$ is defined as the space, temporal and
molecular average by

$\left<|\mbox{\boldmath$v$}|^{2}\right>\equiv\lim_{t_{\infty}\rightarrow\infty}\frac{1}{t_{\infty}}\int_{0}^{t_{\infty}}\left[\frac{1}{N}\sum_{i=1}^{N}|\mbox{\boldmath$v$}_{i}|^{2}\right]dt$
(6)

for $N$ molecules with the velocity of $i$-th molecule being
$\mbox{\boldmath$v$}_{i}$. Note that $\mbox{\boldmath$v$}_{i}$ satisfies

$\frac{1}{N}\sum_{i=1}^{N}\mbox{\boldmath$v$}_{i}=\mbox{\boldmath$0$}$ (7)

for present equilibrium gases, e.g., confined in a closed container. In case
the container is moving at a constant velocity,
$\left<|\mbox{\boldmath$v$}|^{2}\right>$ must be defined by using the
molecular velocities relative to the group motion. By the way, the relation in
Eq. (5) derived by Clausius is amazing because the speed of invisible
molecules can be estimated only by the two measurable macroscopic values of
the pressure $p$ and the density $\rho$ Fujiwara1995 . By inserting Eq. (5)
into Eq. (4), it follows

$M\left<|\mbox{\boldmath$v$}|^{2}\right>=3RT,$ (8)

where $M$ denotes the mass of gas per mole. Then, let $m$ and $k_{\mathrm{B}}$
be the mass of single molecule and the Boltzmann constant, respectively
defined by

$m\equiv\frac{M}{N_{\mathrm{A}}},\quad
k_{\mathrm{B}}\equiv\frac{R}{N_{\mathrm{A}}}=1.38\times 10^{-23}\ {\rm J/K}$
(9)

using the Avogadro number $N_{\mathrm{A}}$, we obtain a fundamental relation
between the kinetic energy of the constituent molecule and the absolute
temperature:

$\frac{1}{2}m\left<|\mbox{\boldmath$v$}|^{2}\right>=\frac{3}{2}k_{\mathrm{B}}T.$
(10)

Equation (10) can be extended to real molecules as the microscopic definition
of temperature $T$.

###  3.2 Intermolecular potential

The intermolecular interaction between ideal gas molecules can be
schematically illustrated in Fig. 3 with the red line. The horizontal line for
overall distance except at zero distance shown with the vertical line
indicates that the molecules have no interaction force except the hard point-
like collision. In contrast to this ideal gas feature, the interaction
potential between real molecules can be modeled as the blue line in Fig. 3.
Two real molecules within an intermediate intermolecular distance attract each
other whereas a strong repulsive force acts between the two within a smaller
distance around the molecular diameter, and they have almost no interaction
for a longer distance. The distance giving the minimum potential energy of
$-\varepsilon$ is the equilibrium distance.

[width=100mm]./fig04-phasechange-schematic.eps

Figure 4: Schematic of the intuitive understanding of the temperature-
dependent phase change from intermolecular interaction between real molecules.
Left: molecules are trapped in the potential well because the average kinetic
energy corresponding to $k_{\mathrm{B}}T$ is not sufficiently large. Right:
molecules can escape from the potential well with the kinetic energy and
freely moving in a ballistic manner. The left corresponds to the solid or
liquid phase whereas the right is the gas phase.

Now we try to intuitively understand phase change from this interaction
potential of real molecules with a potential well depth and the definition of
absolute temperature $T$ in Eq. (10). Figure 4 shows the schematic for the
understanding of the temperature-dependent phase change from the
intermolecular interaction between real molecules. Suppose two molecules
vibrating around the potential well at a low temperature without total
translational motion as in the left panel, i.e., satisfying Eq. (7). Then, as
indicated by Eq. (10), the two molecules have certain relative average kinetic
energy with a scale of $k_{\mathrm{B}}T$.

At a low temperature with the corresponding kinetic energy sufficiently
smaller than the potential well depth $\varepsilon$, the molecules are trapped
in the potential well. In that case, molecules rarely change their interaction
pairs and if there more than two molecules, they try to maximize the number of
neighboring pairs to minimize the potential energy, and they eventually form a
solid crystal typically with a closed packing structure. Note that thermal
expansion of solid can also be understood from the asymmetric vibration
feature around the equilibrium point; the intermolecular distance at $T=0$ is
apparently at the minimum of the potential well while the time-averaged mean
distance at $T>0$ becomes longer as the temperature increases. This expansion
is not expected if the intermolecular interaction is modeled by a simple
harmonic spring connecting the two molecules, which gives a symmetric
vibration around the point of minimum potential energy.

Under a moderate temperature, molecules can frequently escape from the
potential well and change the pairs without having certain fixed structures.
This can be understood as liquid. At high temperature as shown in the right
panel of Fig. 4, molecules are not trapped in the potential well due to their
sufficient kinetic energy and freely moving in space. This corresponds to the
gas phase.

Hence, the key features of real molecules can be summarized as follows:

  * •

There are short-range repulsive and long range attractive forces between the
molecules to form a potential well.

  * •

The molecules are constantly moving, and as a result of multiple random
collisions, their velocity distribution reaches a certain isotropic random
directions as an equilibrium state.

The second feature is the same as ideal gases and the kinetic energy of the
random motion is equally distributed to each degree of freedom: the
equipartition under thermal equilibrium as a primitive basis of statistical
mechanics. It is easy to imagine why Eq. (4) holds under high temperature and
low pressure conditions for real gases from the intermolecular interaction
potential of real molecules; however, the history of thermodynamics and
statistical mechanics tells that constructing it from ideal gas potential is
not that easy. Indeed, the outstanding idea of D. Bernoulli that gases as
fluid consist of huge number of infinitesimally small molecules, proposed in
the first half 18th century was not further investigated for about a century
until Clausius and Maxwell Fujiwara1995 .

[scale=.80]./fig05-drop-homogen.eps

Figure 5: MD simulation results of simple Lennard-Jones molecules confined in
a calculation cell.

Anyway, now we go back to the interpretation of the surface tension in Fig. 1.
This figure does not mean that the molecular pairs in the liquid bulk pull the
partner molecule each other but they are at a certain mean inter-molecular
distance as indicated in the right panel of Fig. 3. Molecules around the
interface have less partners to interact with than those in the liquid bulk,
and those in the vapor bulk have almost no partners. As mentioned above, with
the present computers, MD simulations can be easily run even with a laptop
computers. Figure 5 shows the snapshot of equilibrium systems consisting of
molecules with their intermolecular interaction described by the Lennard-Jones
potential popularly used as a simple model expressing the above-mentioned
feature:

$\Phi_{\mathrm{LJ}}(r_{ij})=4\varepsilon\left[\left(\frac{\sigma}{r_{ij}}\right)^{12}-\left(\frac{\sigma}{r_{ij}}\right)^{6}\right],$
(11)

where $\varepsilon$ and $\sigma$ denote the potential well depth and diameter
shown in Fig. 3. One thousand molecules are confined in a cubic box sized
$30\sigma\times 30\sigma\times 30\sigma$ with the time-averaged system
temperatures, i.e., kinetic energy given by Eq. (10), set at
$k_{\mathrm{B}}T/\varepsilon$ equal to 0.75 and 1.5 in the left and right
systems, respectively.

Note that the molecules in the systems are constantly moving as described
above, and even if one starts a simulation from arbitrary initial
configuration, for instance, with locating the molecules at grid points, the
molecules spontaneously form corresponding equilibrium states after a certain
time as in the snapshots. This term ‘equilibrium’ in the microscopic scale
roughly means that the molecules are moving with the positions and momenta of
the constituent molecules changing with time, but their apparent feature as a
group is unchanged. For instance, most of the constituent molecules are
gathered to form a spherical liquid droplet with remaining molecules flying as
vapor around the droplet in the left system as an equilibrium state; the
constituent molecules of liquid and vapor change but the volume as well as the
spherical structure of the droplet are unchanged. Similarly, the positions of
momenta of the molecules change, but they always keep vapor phase in the right
system after reaching the equilibrium state spontaneously achieved after a
certain time irrespective of the initial condition. As described later, this
spontaneous change of the system is due to the minimization of the Helmholtz
free energy of the system.

As clearly observed from the figure, liquid-vapor coexistence kept at the
lower temperature of $k_{\mathrm{B}}T/\varepsilon=0.75$ in the left panel is
lost with the temperature rise, and the system is filled with vapor (indeed
supercritical phase) at the higher temperature of
$k_{\mathrm{B}}T/\varepsilon=1.5$. This phase change can be intuitively
understood from the schematic in Fig. 4.

We further think about why the droplet in the left panel takes the spherical
structure to minimize the surface area from an energetic point of view. Let
the potential energy per molecule in the vapor bulk, at the LV interface, and
in the liquid bulk be $e_{\mathrm{V}}$, $e_{\mathrm{LV}}$ and
$e_{\mathrm{L}}$, respectively. Then, they satisfy

$e_{\mathrm{L}}<e_{\mathrm{LV}}\ll e_{\mathrm{V}}\approx 0$ (12)

considering the schematic in Fig. 1 and intermolecular potential in Fig. 3,
which indicates that $e_{\mathrm{V}}$, $e_{\mathrm{LV}}$ are basically
negative whilst $e_{\mathrm{L}}$ is almost zero for molecules in the vapor
phase without partners to interact with. The first inequality implies that the
number of neighbouring molecules is smaller at the interface than in the
liquid bulk. Consequently

$\quad e_{\mathrm{V}}-e_{\mathrm{LV}}>0$ (13)

holds. In addition, let the number of corresponding molecules and volumes
being $N_{\mathrm{L}}$, $N_{\mathrm{LV}}$ and $N_{\mathrm{V}}$, and
$V_{\mathrm{L}}$, $V_{\mathrm{LV}}$ and $V_{\mathrm{V}}$, respectively. Then,
the number of molecules $N$ and volume $V$ of the system write

$N=N_{\mathrm{L}}+N_{\mathrm{LV}}+N_{\mathrm{V}}$ (14)

and

$\displaystyle V=V_{\mathrm{L}}+V_{\mathrm{LV}}+V_{\mathrm{V}},$ (15)

respectively. In the macroscopic scale, we assume $N_{\mathrm{LV}}=0$ and
$V_{\mathrm{LV}}=0$ while the LV interface is a region in the molecular scale
and molecules indeed exist there. We define the corresponding volumes per
molecule as

$v_{\mathrm{L}}=\frac{V_{\mathrm{L}}}{N_{\mathrm{L}}},\
v_{\mathrm{LV}}=\frac{V_{\mathrm{LV}}}{N_{\mathrm{LV}}},\
v_{\mathrm{V}}=\frac{V_{\mathrm{V}}}{N_{\mathrm{V}}}.$ (16)

Assuming that $v_{\mathrm{LV}}$ is approximately the same as $v_{\mathrm{L}}$,
it follows for the volume per molecule that

$0<v_{\mathrm{L}}\approx v_{\mathrm{LV}}\ll v_{\mathrm{V}}.$ (17)

Then, Eq. (15) rewrites

$V=N_{\mathrm{V}}v_{\mathrm{V}}++N_{\mathrm{LV}}v_{\mathrm{LV}}+N_{\mathrm{V}}v_{\mathrm{L}}\approx(N_{\mathrm{L}}+N_{\mathrm{LV}})v_{\mathrm{L}}+N_{\mathrm{V}}v_{\mathrm{V}}.$
(18)

On the other hand, the internal energy $U$ of a system is separated into the
kinetic and potential contributions $U^{\mathrm{kin}}$ and $U^{\mathrm{pot}}$,
respectively as

$U=U^{\mathrm{kin}}+U^{\mathrm{pot}}.$ (19)

Using Eq.(10), the former is given by

$\displaystyle
U^{\mathrm{kin}}\equiv\sum_{i=1}^{N}\frac{1}{2}m|\mbox{\boldmath$v$}_{i}|^{2}=\frac{3}{2}Nk_{\mathrm{B}}T,$
(20)

which is constant under the constant temperature condition. The latter of
potential term is formally given by

$\displaystyle\displaystyle
U^{\mathrm{pot}}\equiv\sum_{i=1}\sum_{j(>i)}\Phi(r_{ij}).$ (21)

This can be approximated using the mean potential energy per molecule as

$\displaystyle\displaystyle U^{\mathrm{pot}}$ $\displaystyle\approx$
$\displaystyle
N_{\mathrm{V}}e_{\mathrm{V}}+N_{\mathrm{LV}}e_{\mathrm{LV}}+N_{\mathrm{L}}e_{\mathrm{L}}$
(22) $\displaystyle=$
$\displaystyle(N_{\mathrm{L}}+N_{\mathrm{LV}})e_{\mathrm{L}}+N_{\mathrm{LV}}(e_{\mathrm{LV}}-e_{\mathrm{L}}).$

From Eqs. (18)-(22) and Inequality (13), it is shown that with keeping
$N_{\mathrm{L}}+N_{\mathrm{LV}}$ constant, reducing the internal energy $U$ in
Eq. (22) is possible by decreasing $N_{\mathrm{LV}}$ without changing the
volume $V$ in Eq. (18). This means that the system internal energy can be
reduced by decreasing the surface area $A_{\mathrm{LV}}$ because
$N_{\mathrm{LV}}$ is apparently proportional to $A_{\mathrm{LV}}$. This leads
to the formation of the spherical droplet in the left panel of Fig. 5 as an
equilibrium shape with the smallest surface area in the three-dimensional
space.

###  3.3 Free energy and entropy

The explanation above looks reasonable; however, another question arises why
the system prefers to take the ‘all vapor’ state with less internal energy $U$
at higher temperatures. The trap model in Fig. (4) merely indicates that the
molecular pairs easily dissociate from each other, but they can have a lower
internal energy even at a higher temperature if they form more intermolecular
pairs. To answer this question, the free energy instead of the internal energy
must be considered. The second law of thermodynamics tells us that the
Helmholtz free energy $F$ of a system given by

$F\equiv U-TS$ (23)

decreases until the system reaches the equilibrium state, i.e.,

$dF\leq 0\ (\mathrm{nonequilibrium}),\quad F=F_{\mathrm{min}}\
(\mathrm{equilibrium})$ (24)

holds for $F$ for a system under constant volume $V$ and constant temperature
$T$ condition. This corresponds to a closed container in contact with a a
constant-temperature heat bath. The key point is that the Helmholtz free
energy in Eq. (23) includes the entropy of the system $S$.

In classical statistical mechanics, the entropy $S$ is related to the “number
of possible equivalent microstates corresponding to a certain macroscopic
state” $\Omega$ as follows:111 Boltzmann’s tombstone bears the inscription of
‘$S=k\log W$’.

$S=k_{\mathrm{B}}\ln\Omega,$ (25)

where a microstate in 3-dimensional system of $N$-molecules is determined by
giving 3$N$-positions and $3N$-momenta of constituent molecules. From Eq.
(25), it is possible to show that the position contribution to entropy $S$ is
the largest for a system with homogeneous density. We will see that in the
following example:

(Example) Suppose a closed equilibrium system with ideal gas molecules in a
box at a constant temperature, and let $N_{\mathrm{left}}$ and
$N_{\mathrm{right}}$ be the numbers of gas molecules in the left-half and
right-half of a box, respectively, and denote the case by
$[N_{\mathrm{left}},N_{\mathrm{right}}]$
($N=N_{\mathrm{left}}+N_{\mathrm{right}}$). Evaluate the number of possible
cases and compare the two for the following three cases, respectively:

  1. 1.

Case of small total number $N_{\mathrm{total}}=200$: $[99,101]$ vs.
$[100,100]$

  2. 2.

Case of a larger total number $N_{\mathrm{total}}=2000$: $[990,1010]$ vs.
$[1000,1000]$

  3. 3.

Case of a huge (Avogadro scale) total number of $N_{\mathrm{total}}=2\times
10^{23}$: $[0.99N_{\mathrm{A}},1.01N_{\mathrm{A}}]$ vs.
$[N_{\mathrm{left}}=N_{\mathrm{A}},N_{\mathrm{right}}=N_{\mathrm{A}}]$

(Answer) The numbers of cases $\Omega$ to separate
$N_{\mathrm{left}}+N_{\mathrm{right}}$ into left and right are given by

  1. 1.

$N_{\mathrm{total}}=200$:

$\displaystyle\displaystyle\Omega_{[99,101]}=\frac{200!}{99!101!},\quad\Omega_{[100,100]}=\frac{200!}{100!100!}$
$\displaystyle\displaystyle\frac{\Omega_{[99,101]}}{\Omega_{[100,100]}}=\frac{100!100!}{99!101!}=\frac{100}{101}$

  2. 2.

$N_{\mathrm{total}}=2000$:

$\displaystyle\displaystyle\Omega_{[990,1010]}=\frac{2000!}{990!1010!},\quad\Omega_{[1000,1000]}=\frac{2000!}{1000!1000!}$
$\displaystyle\displaystyle\frac{\Omega_{[990,1010]}}{\Omega_{[1000,1000]}}=\frac{1000!1000!}{990!1010!}=\frac{991}{1001}\cdot\frac{992}{1002}\cdot\ldots\frac{1000}{1010},$

The ratio above satisfies

$\left(\frac{99}{100}\right)^{10}<\frac{\Omega_{[990,1010]}}{\Omega_{[1000,1000]}}<\left(\frac{100}{101}\right)^{10}\quad\therefore\frac{\Omega_{[990,1010]}}{\Omega_{[1000,1000]}}\approx\left(\frac{99}{100}\right)^{10}$

  3. 3.

$N_{\mathrm{total}}=2\times 10^{23}$:

$\displaystyle\displaystyle\Omega_{[0.99N_{\mathrm{A}},1.01N_{\mathrm{A}}]}=\frac{2N_{\mathrm{A}}!}{(0.99N_{\mathrm{A}})!(1.01N_{\mathrm{A}})!},\quad\Omega_{[N_{\mathrm{A}},N_{\mathrm{A}}]}=\frac{2N_{\mathrm{A}}!}{N_{\mathrm{A}}!N_{\mathrm{A}}!}$
$\displaystyle\displaystyle\frac{\Omega_{[0.99N_{\mathrm{A}},1.01N_{\mathrm{A}}]}}{\Omega_{[N_{\mathrm{A}},N_{\mathrm{A}}]}}=\frac{N_{\mathrm{A}}!N_{\mathrm{A}}!}{(0.99N_{\mathrm{A}})!(1.01N_{\mathrm{A}})!}\approx\left(\frac{99}{100}\right)^{0.01N_{\mathrm{A}}}$
(26)

Equation (26) indicates that the probability to see 2 % density difference
($0.99\rho_{0}$ and $1.01\rho_{0}$) between the left and right for ideal gas
systems approaches to zero for $N_{\mathrm{A}}\gg 1$. In addition, by
inserting Eq. (25) into Eq. (26), it follows for the entropy difference
$\Delta S$ between $S_{[N_{\mathrm{A}},N_{\mathrm{A}}]}$ and
$S_{[0.99N_{\mathrm{A}},1.01N_{\mathrm{A}}]}$ that

$\Delta S\equiv
S_{[N_{\mathrm{A}},N_{\mathrm{A}}]}-S_{[0.99N_{\mathrm{A}},1.01N_{\mathrm{A}}]}\approx
0.01N_{\mathrm{A}}k_{\mathrm{B}}\ln\left(\frac{100}{99}\right)>0.$ (27)

This indicates that the entropy $S$ is the largest for uniform density, and
also that $S$ is proportional to the number of molecules.

Now we go back to the comparison of the Helmholtz free energy in Fig. 5. The
entropy $S$ is larger for a system with homogeneous density than for an
inhomogeneous system, meaning that the right system is advantageous from the
entropy aspect. At low temperature, $U$ contribution in Eq. (23) is dominant
and the system tries to decrease $U$ to minimize $F$, whereas at high
temperature, the entropy contribution $TS$ overcomes the potential
minimization effect. This balance between the potential energy and the entropy
governs the basic mechanism of phase separation as well as the surface
tension. Note that for ideal gas systems, the potential contribution is
constant (zero for the potential in Fig. 3), and the gas molecules tries to
fill the system with a homogeneous density to maximize entropy $S$. Note as
well that if a solid surface exists near the droplet at a temperature for the
liquid-vapor phase separation, and the droplet is thermodynamically more
stable on the surface, a hemispherical droplet is formed on the solid surface,
and the equilibrium shape, i.e., the contact angle $\theta$ in Fig. 2, is
determined so that the total free energy of the system may become minimum
depending on the solid-fluid interaction strength.

##  4 Surface tension

[scale=.50]./fig06-thought_exp-LV.eps

Figure 6: Thought experiment of Bakker’s equation for a flat liquid-vapor
interface. The red arrows denote the pressure, i.e., normal fluid stress with
its sign inverted, acting on the piston. Reprinted with permission from
Yamaguchi et al., J. Chem. Phys. 150, 044701 (2019). Copyright 2019 Author(s),
licensed under a Creative Commons Attribution (CC BY) license 4.0 License.

###  4.1 Bakker’s equation

As seen in the above example, we have to inevitably introduce the relation
between the free energy of a fluid with interface and the local force to
quantitatively evaluate the surface tension. For that purpose, we set a
thought experiment with a flat liquid-vapor interface shown in Fig. 6 for the
basic connection between the thermodynamic work as energy and the surface
tension. Indeed, such flat interface can easily be achieved by MD simulations
by using the periodic boundary condition in the surface lateral directions. In
this thought experiment, one piston is set normal to a flat liquid-vapor
interface, and it covers from $z=z_{\mathrm{L}}$ to $z=z_{\mathrm{V}}$ at the
liquid and vapor bulk regions, respectively, across the plane of the liquid-
vapor interface. Another piston parallel to the interface is set in the vapor
bulk far from the interface. Through simultaneous virtual infinitesimal
displacements of the pistons, only the interface area can be changed without
changing the liquid and vapor volumes, $V_{\mathrm{L}}$ and $V_{\mathrm{V}}$,
respectively. Note that the change of bottom area in the figure is not
considered. Let $l$ be the depth normal to the $xz$-plane, and $\delta x$ be
the corresponding displacement of the side piston, the change of interface
area $A_{\mathrm{LV}}$ writes

$\delta A_{\mathrm{LV}}=l\delta x.$ (28)

In order to let the internal energies before and after the displacement
unchanged both for the liquid and vapor parts, this displacement must be done
under constant temperature. Then, the net minimum mechanical work $\delta W$
exerted from the top and side pistons required for this change with quasi-
static process can be associated with the change in the Helmholtz energy $F$.
Let $\gamma_{\mathrm{LV}}$ be the LV interfacial energy per area, it follows
for the quasi-static change that

$\delta F=\gamma_{\mathrm{LV}}\delta A_{\mathrm{LV}}.$ (29)

Thus, $\gamma_{\mathrm{LV}}$ writes

$\gamma_{\mathrm{LV}}=\left(\frac{\partial F}{\partial
A_{\mathrm{LV}}}\right)_{N,V_{\mathrm{L}},V_{\mathrm{V}},T}.$ (30)

This means that surface tension is intrinsically a thermodynamic force.

We further relate the mechanical stress changing as a continuous function with
the LV interfacial tension $\gamma_{\mathrm{LV}}$. The fluid stress tensor is
not isotropic at phase interfaces even at equilibrium as shown with the red
arrows. Due to the one-dimensional feature of the system along the
$z$-direction, the symmetric stress tensor $\tau$ writes

$\mbox{\boldmath$\tau$}=\left(\begin{matrix}\tau_{xx}&0&0\\\ 0&\tau_{yy}&0\\\
0&0&\tau_{zz}\end{matrix}\right),$ (31)

where the surface-normal components $\tau_{zz}$ is constant in the entire
region because of the force balance in the $z$-direction to be satisfy in the
static equilibrium system for the present $z$-normal flat LV interface system.
This constant value is equal to the saturated vapor pressure
$p^{\mathrm{bulk}}$ with its sign inverted, i.e.,

$\tau_{zz}=-p^{\mathrm{bulk}}.$ (32)

On the other hand, the surface-lateral diagonal components $\tau_{xx}$ and
$\tau_{yy}$ satisfy

$\tau_{xx}(z)=\tau_{yy}(z)\equiv-p_{\mathrm{T}}(z),$ (33)

where the surface tangential pressure is denoted by $p_{\mathrm{T}}(z)$, which
is a unique function of position $z$ and is equal to the isotropic vapor
pressure $p^{\mathrm{bulk}}$ in the liquid and vapor bulks. We set the liquid
and vapor bulk positions at $z=z_{\mathrm{L}}$ and $z=z_{\mathrm{V}}$,
respectively. Using these values, the changes of the volume and the Helmholtz
free energy as the work upon the quasi-static displacement $\delta x$ of the
side piston are expressed by

$\delta V=l\delta x\int_{z_{\mathrm{L}}}^{z_{\mathrm{V}}}dz,$ (34)

and

$dF\equiv\delta W=p^{\mathrm{bulk}}\delta V+l\delta
x\int_{z_{\mathrm{L}}}^{z_{\mathrm{V}}}p_{\mathrm{T}}\left(z\right)dz,$ (35)

respectively. By inserting Eqs. (28) and (32)-(35) into Eq. (30),

$\displaystyle\gamma_{\mathrm{LV}}$ $\displaystyle=$
$\displaystyle\displaystyle\lim_{\delta x\rightarrow
0}\left.\frac{dF}{dA_{\mathrm{LV}}}\right|_{N,V_{\mathrm{L}},V_{\mathrm{V}},T}$
(36) $\displaystyle=$ $\displaystyle
p^{\mathrm{bulk}}\int_{z_{\mathrm{L}}}^{z_{\mathrm{V}}}dz-\int_{z_{\mathrm{L}}}^{z_{\mathrm{V}}}p_{\mathrm{T}}\left(z\right)dz=\int_{z_{\mathrm{L}}}^{z_{\mathrm{V}}}\left[p^{\mathrm{bulk}}-p_{T}\left(z\right)\right]dz$
$\displaystyle=$
$\displaystyle\int_{z_{\mathrm{L}}}^{z_{\mathrm{V}}}\left[\frac{\tau_{xx}(z)+\tau_{yy}(z)}{2}-\tau_{zz}\right]dz$

is derived. Equation (36) is called Bakker’s equation, which serves as the
basic connection between the microscopic anisotropic stress distribution and
the surface tension. Note that the integration range in Eq. (36) is sometimes
denoted by $\int_{-\infty}^{\infty}$; however, the necessary condition for
this range is that $z_{\mathrm{L}}$ and $z_{\mathrm{V}}$ cover the entire
range of the interface with anisotropic fluid stress, i.e., between the liquid
and vapor bulks.

###  4.2 Molecular dynamics simulations

Figure 7:  (Left) snapshot of an equilibrium molecular dynamics simulation
system of Lennard-Jones molecules (argon) with flat solid-liquid and liquid-
vapor interfaces. (Middle) corresponding density distribution. (Right)
distributions of stress diagonal components $\tau_{xx}$ and $\tau_{zz}$.

As mentioned in Section 2, molecular dynamics analysis is a powerful choice to
overcome the theoretical difficulty of obtaining the equilibrium density
distribution in liquid-vapor coexistence systems and resulting stress
distribution to calculate the surface tension by Bakker’s equation. The left
panel of Fig. (7) shows a MD simulation system of Lennard-Jones (LJ) molecules
(argon: $\sigma=0.34$ nm, $\varepsilon=1.65\times 10^{-21}$ J in Eq. (11) and
$m=39.948$ g/mol) with flat solid-liquid (SL) and liquid-vapor (LV)
interfaces, where 2000 argon molecules were confined in a rectangular
simulation cell of $4\times 4\times 20$ nm3. The periodic boundary conditions
were imposed in the surface lateral $x$\- and $y$-directions and the mirror
boundary condition was set on the top boundary. In addition, a solid wall
modelled through an integrated LJ potential was set on the bottom of the
system to let the LJ-liquid be adhered on the solid. The system was
equilibrated at a temperature of $T=100$ K
($k_{\mathrm{B}}T/\varepsilon=0.835$) until the apparent feature of the system
did not change. With such setting, one can easily realize a quasi-one-
dimensional equilibrium MD system with flat LV and SL interfaces normal to the
$z$-direction.

Figure 8:  Schematics of the stress tensor calculation based on the method of
plane (MoP).

It is easy to calculate the density of this equilibrium system by taking the
time-average of the mass in each flat bin volume set normal to the
$z$-direction with a small thickness of $\Delta z$. On the other hand, to
calculate the stress, a method called the Method of Plane (MoP) is used, where
one sets control surfaces in the system instead of the bin volumes
Thompson1984 ; Yaguchi2010 . The local fluid stress tensor
$\mbox{\boldmath$\tau$}(x,z)$ is calculated by setting $z$-normal and $x$\- or
$y$-normal flat bin faces in the system in the Cartesian coordinate. Figure 8
shows the schematics of the stress tensor calculation by the MoP. The fluid
stress tensor component $\tau_{\alpha\beta}$, which expresses the stress in
$\beta$-direction exerted on a surface element with an outward normal in
$\alpha$-direction, is given by kinetic term
$\tau_{\alpha\beta}^{\mathrm{kin}}$ and inter-molecular interaction term
$\tau_{\alpha\beta}^{\mathrm{int}}$ as

$\tau_{\alpha\beta}=\tau_{\alpha\beta}^{\mathrm{kin}}+\tau_{\alpha\beta}^{\mathrm{int}}.$
(37)

In the MoP, the kinetic term on an $\alpha$-normal bin face with an area
$A_{\alpha}$ is calculated by

$\tau_{\alpha\beta}^{\mathrm{kin}}=-\frac{1}{A_{\alpha}}\left<\sum_{i\in\mathrm{fluid},\delta
t}^{\mathrm{across}A_{\alpha}}\frac{\left(2\Theta(\mbox{\boldmath$v$}_{i}\cdot\mbox{\boldmath$e$}_{\alpha})-1\right)m_{i}\mbox{\boldmath$v$}_{i}\cdot\mbox{\boldmath$e$}_{\beta}}{\delta
t}\right>,$ (38)

where $m_{i}$ and $\mbox{\boldmath$v$}_{i}$ denotes the mass and velocity
vector of $i$-th fluid molecule, and $\mbox{\boldmath$e$}_{\alpha}$ and
$\mbox{\boldmath$e$}_{\beta}$ are the unit normal vectors in $\alpha$\- and
$\beta$-directions, respectively. The angle brackets means the time average,
and the summation $\sum_{i,\delta t}^{\mathrm{across}A_{\alpha}}$ is taken for
every fluid molecule $i$ passing through the bin face within a time interval
of $\delta t$, which is equal to the time increment for the numerical
integration. A switching function
$2\Theta(\mbox{\boldmath$v$}_{i}\cdot\mbox{\boldmath$e$}_{\alpha})-1$, which
gives $\pm 1$ depending on the sign of
$\mbox{\boldmath$v$}_{i}\cdot\mbox{\boldmath$e$}_{\alpha}$ implemented through
the Heaviside step function $\Theta$, is included in the RHS of Eq. (38). The
meaning of the kinetic stress tensor in Eq. (38) is shown with an example case
with $\alpha=x$ and $\beta=y$. Suppose that molecule $i$ with its
$x$-directional velocity $v^{x}_{i}$ passes through a $x$-normal plane from
left to right within a time interval $\delta t$ as in the upper molecule in
Fig. 8(a). This is only achieved for positive $v^{x(=\alpha)}_{i}$ satisfying
$\mbox{\boldmath$v$}_{i}\cdot\mbox{\boldmath$e$}_{\alpha}>0$ with
$\mbox{\boldmath$e$}_{\alpha}=(1,0,0)$, which corresponds to
$2\Theta(\mbox{\boldmath$v$}_{i}\cdot\mbox{\boldmath$e$}_{\alpha})-1=1$. Then,
the $y$-directional momentum of the fluid in the left is reduced by
$m_{i}v^{y}_{i}$, i.e., is increased by
$-m_{i}\mbox{\boldmath$v$}_{i}\cdot\mbox{\boldmath$e$}_{\beta}$ with
$\mbox{\boldmath$e$}_{\beta}=(0,1,0)$. This change of momentum is equivalent
to the impulse $\tau_{xy}A\delta t$ with the stress tensor. If the passage is
in the opposite direction from right to left as in the lower molecule in Fig.
8(a), then,
$2\Theta(\mbox{\boldmath$v$}_{i}\cdot\mbox{\boldmath$e$}_{\alpha})-1=-1$ and a
momentum of $m_{i}\mbox{\boldmath$v$}_{i}\cdot\mbox{\boldmath$e$}_{\beta}$ is
given to the fluid in the left. From this formulation, it is known that the
kinetic contribution to surface normal stress $\tau_{\alpha\alpha}$ is always
negative, i.e., results in positive pressure, and this indeed corresponds to
the ideal gas pressure $p$ given in Eq. (5), which is proportional to the
density $\rho$ on the face, and is also proportional to the temperature $T$
from Eq. (10).

On the other hand, the intermolecular interaction term in Eq. (37) is given by

$\tau_{\alpha\beta}^{\mathrm{int}}=-\frac{1}{A_{\alpha}}\left<\sum_{(i,j)\in\mathrm{fluid}}^{\mathrm{across}A_{\alpha}}\left(2\Theta(\mbox{\boldmath$r$}_{ij}\cdot\mbox{\boldmath$e$}_{\alpha})-1\right)\mbox{\boldmath$f$}_{ij}\cdot\mbox{\boldmath$e$}_{\beta}\right>,$
(39)

where $\mbox{\boldmath$r$}_{ij}$ and $\mbox{\boldmath$f$}_{ij}$ denote the
relative position vector $\mbox{\boldmath$r$}_{j}-\mbox{\boldmath$r$}_{i}$ and
force vector exerted on molecule $j$ at position $\mbox{\boldmath$r$}_{j}$
from molecule $i$ at $\mbox{\boldmath$r$}_{i}$, respectively, and the
summation $\sum_{(i,j)\in\mathrm{fluid}}^{\mathrm{across}A_{\alpha}}$ was
taken for all line segments between $\mbox{\boldmath$r$}_{i}$ and
$\mbox{\boldmath$r$}_{j}$ which crossed the bin face. In contrast to the
kinetic contribution, this interaction contribution may give positive surface
normal stress $\tau_{\alpha\alpha}$ on the bin face with an attractive
interaction $\mbox{\boldmath$f$}_{ij}$.

The middle and right panels of Fig. 7 shows the density and surface-normal
stress distributions, respectively. As expected, the LV interface with a large
change in the density has a certain thickness, and in this region with non-
constant density, the normal stress $\tau_{xx}$ tangential to the interface is
different from the bulk stress $\tau_{zz}$, which is constant in the whole
system except near the solid at which the fluid is subject to an external
force from the solid. Note that a remarkable wiggling structures in the
density and stress at bottom is specific to the solid-liquid interface, and
this indeed is an interesting topic especially related to the wetting;
however, the effect vanishes as leaving away from the solid surface, and we
mainly focus on the liquid-vapor interface. It is really amazing that the
frontiers could expect the stress distribution at least around the LV
interface without MD simulations. Also note that even by such great
scientists, quantitatively evaluating the density and stress distributions
only from the intermolecular potential was still difficult, and that we have a
powerful tool of MD simulations to complement the missing piece. Related to
this point, calculating the SL (and SV) interfacial tension from the stress
distribution obtained through the MD simulations toward the understanding of
wetting is a hot topic at present Surblys2014 ; Yamaguchi2019 ; Kusudo2021 ;
however, I will not describe about this in detail considering the scope of
this article.

Figure 9:  The simplest model to explain the stress anisotropy at the
interface. The interaction terms $\tau^{\mathrm{int}}_{zz}$ and
$\tau^{\mathrm{int}}_{xx}$ on the point of interest (black crosses) are
represented by the interaction between two points having densities
$(\rho_{\mathrm{L}},\rho_{\mathrm{V}})$, and
$(\frac{\rho_{\mathrm{L}}+\rho_{\mathrm{V}}}{2},\frac{\rho_{\mathrm{L}}+\rho_{\mathrm{V}}}{2})$,
respectively connected by the double-headed red arrows.

A very simple and intuitive model to explain the stress anisotropy at the
liquid-vapor interface, i.e., $\tau_{zz}<\tau_{xx}$, is given here. Suppose
$z$-normal and $x$-normal planes on a $z$-normal flat LV interface for the
calculation of the normal stress components $\tau_{zz}$ and $\tau_{xx}$,
respectively based on the MoP as shown in Fig. 9. In this model, we assume
average densities $\rho_{\mathrm{L}}$ and
$\rho_{\mathrm{V}}(<\rho_{\mathrm{L}})$ and
$\frac{\rho_{\mathrm{L}}+\rho_{\mathrm{V}}}{2}$ in the liquid and vapor bulks
and on the interface, respectively for simplicity. In addition, we introduce a
mean-field like approach to the MoP here; instead of counting the molecules
passing through the MoP plane for the kinetic stress term
$\tau^{\mathrm{kin}}$ in Eq. (38) or summing up the intermolecular interaction
force for the line segment between two molecules crossing the MoP plane for
the interaction stress term $\tau_{\mathrm{int}}$ in Eq. (39), we evaluate
$\tau^{\mathrm{kin}}$ with the density on the point of interest (black cross)
on the MoP plane, and $\tau^{\mathrm{int}}$ with the two densities on the
points which sandwich the point of interest (double-headed red arrow) in Fig.
9. As easily imagined, to correctly calculate the latter interaction term of
$\tau^{\mathrm{int}}$, complicated space integration for the two points using
a position- and direction-dependent radial distribution function is needed,
and this is indeed the difficulties that the great scientists could not solve
analytically, but here we try to represent the interaction term by the most
probable single distance in the present simplest model. The kinetic terms
$\mbox{\boldmath$\tau$}^{\mathrm{kin}}$ is isotropic in equilibrium systems
and is given by

$\mbox{\boldmath$\tau$}^{\mathrm{kin}}=-\frac{Nk_{\mathrm{B}}T}{V}\mbox{\boldmath$I$}=-\rho\frac{k_{\mathrm{B}}T}{m}\mbox{\boldmath$I$},$
(40)

where $I$ denotes the identity tensor. Equation (40) corresponds to the ideal-
gas pressure, which can be interpreted with Eqs. (10) and (38) as follows:
fast molecules giving a momentum of $m\mbox{\boldmath$v$}_{i}$ in Eq. (38)
pass through the plane more with proportional to $|\mbox{\boldmath$v$}_{i}|$.
From Eq. (40), $\tau^{\mathrm{kin}}_{zz}$ and $\tau^{\mathrm{kin}}_{xx}$ are
identical and those at the LV interface at $z=z_{\mathrm{LV}}$ are expressed
by

$\tau^{\mathrm{kin}}_{zz}(z_{\mathrm{LV}})=\tau^{\mathrm{kin}}_{xx}(z_{\mathrm{LV}})=-\frac{\rho_{\mathrm{L}}+\rho_{\mathrm{V}}}{2}\frac{k_{\mathrm{B}}T}{m}.$
(41)

On the other hand, we assume that the interaction term is represented by a
single line of interaction force of the most probable distance (double-headed
red arrows) which gives an intermolecular force $f_{\mathrm{rep}}(z)$, then
the interaction term $\tau^{\mathrm{int}}_{\alpha\alpha}$ writes

$\tau^{\mathrm{int}}_{\alpha\alpha}(z)\equiv-\rho_{1}\rho_{2}c(z)f_{\mathrm{rep}}(z),$
(42)

where $\rho_{1}$ and $\rho_{2}$ are the density of the two points sandwiching
the point of interest, and $c(>0)$ is a certain positive coefficient for the
model due to this representation. Note that negative $f_{\mathrm{rep}}(z)$
corresponds to attractive force across the MoP plane to give positive normal
stress. Since $\tau_{zz}$ must satisfy Eq. (32), $\tau_{zz}(z_{\mathrm{LV}})$
at the interface $z=z_{\mathrm{LV}}$ given by

$\tau_{zz}(z_{\mathrm{LV}})=\tau_{zz}^{\mathrm{kin}}(z_{\mathrm{LV}})+\tau_{zz}^{\mathrm{int}}(z_{\mathrm{LV}})=-\frac{\rho_{\mathrm{L}}+\rho_{\mathrm{V}}}{2}\frac{k_{\mathrm{B}}T}{m}-\rho_{\mathrm{V}}\rho_{\mathrm{L}}c(z_{\mathrm{LV}})f_{\mathrm{rep}}(z_{\mathrm{LV}})$
(43)

with $(\rho_{1},\rho_{2})=(\rho_{\mathrm{L}},\rho_{\mathrm{V}})$ in Eq. (42),
is equal to $\tau_{zz}(z_{\mathrm{V}})$ in the vapor bulk at
$z=z_{\mathrm{V}}$, corresponding to the saturated vapor pressure. This is
given by

$\tau_{zz}(z_{\mathrm{V}})=\tau_{zz}^{\mathrm{kin}}(z_{\mathrm{V}})+\tau_{zz}^{\mathrm{int}}(z_{\mathrm{V}})=-{\rho_{\mathrm{V}}}\frac{k_{\mathrm{B}}T}{m}-\rho_{\mathrm{V}}^{2}c(z_{\mathrm{V}})f_{\mathrm{rep}}(z_{\mathrm{V}})\approx-{\rho_{\mathrm{V}}}\frac{k_{\mathrm{B}}T}{m},$
(44)

where the interaction term is assumed to be zero for the last approximation
considering that the intermolecular interaction in the vapor is negligible.
Thus, it follows

$c(z_{\mathrm{LV}})f_{\mathrm{rep}}(z_{\mathrm{LV}})=-\frac{\rho_{\mathrm{L}}-\rho_{\mathrm{V}}}{2}\frac{k_{\mathrm{B}}T}{m}<0,$
(45)

indicating that $f_{\mathrm{rep}}(z_{\mathrm{LV}})$ is negative, i.e., the
representative interaction across the MoP plane at the LV interface is
attractive. Assuming as well that
$c(z_{\mathrm{LV}})f_{\mathrm{rep}}(z_{\mathrm{LV}})$ is independent of the
direction at $z_{\mathrm{LV}}$, then the normal stress lateral to the
interface writes

$\left.\tau_{xx}\right|_{z=z_{\mathrm{LV}}}=-\frac{\rho_{\mathrm{L}}+\rho_{\mathrm{V}}}{2}\frac{k_{\mathrm{B}}T}{m}-\left(\frac{\rho_{\mathrm{V}}+\rho_{\mathrm{L}}}{2}\right)^{2}c(z_{\mathrm{LV}})f_{\mathrm{rep}}(z_{\mathrm{LV}}).$
(46)

Using the relation between the arithmetic and geometric means

$\rho_{\mathrm{L}}\rho_{\mathrm{V}}<\left(\frac{\rho_{\mathrm{L}}+\rho_{\mathrm{V}}}{2}\right)^{2},$
(47)

it follows for $\tau_{zz}(z_{\mathrm{LV}})$ and $\tau_{xx}(z_{\mathrm{LV}})$
in Eqs. (43) and (46), respectively that

$\tau_{zz}(z_{\mathrm{LV}})<\tau_{xx}(z_{\mathrm{LV}}).$ (48)

This inequality is a simple model to explain the stress anisotropy shown in
Fig. 7, due to the density anisotropy, and this is the fundamental basis of
Bakker’s equation (36).

Figure 10:  (Left) distributions of the normal stress $\tau_{xx}$ tangential
to the interface and density $\rho$ at three different temperature $T=90,95$
and 100 K enlarged around the LV interface. (Right) liquid-vapor interfacial
tension $\gamma_{\mathrm{LV}}$ and $p^{\mathrm{bulk}}$ at different
temperature $T$.

The left panel of Fig. (10) shows the distributions of the normal stress
$\tau_{xx}$ tangential to the interface and density $\rho$ at three different
temperature $T=90,95$ and 100 K enlarged around the LV interface. As the
temperature increases, the stress difference between the values of
$\tau_{xx}(z)$ at the interface and away from the interface becomes smaller.
The latter away from the interface is equal to and $\tau_{zz}$, which is equal
to the constant stress value away from the LV interface as shown in Fig. 7.
From Bakker’s equation (36), this means that the LV interfacial tension
$\gamma_{\mathrm{LV}}$ decreases with the temperature rise. Note also that the
bulk liquid density decreases with the the temperature rise, and this is
related to the increase of the average distance upon temperature rise
indicated in Fig. 4. In addition, the bulk stress decreases with the
temperature rise, and this means that the saturated vapor pressure
$p^{\mathrm{bulk}}_{\mathrm{V}}$ given by Eq. (32) increases with temperature.
Such basic features are indeed realized in the MD systems.

The right panel of Fig. 10 shows the LV interfacial tension
$\gamma_{\mathrm{LV}}$ and $p^{\mathrm{bulk}}$ at different temperature $T$.
As expected from the right panel, $\gamma_{\mathrm{LV}}$ decreases and
$p^{\mathrm{bulk}}$ increases with the temperature increase. The vapor
pressure is above the atmospheric pressure, and the extraordinary high
anisotropic local stress $\tau_{xx}$ gives the interfacial tension around 10
mN/m, and the resulting temperature dependence of about $-0.2$ mN/m$\cdot$ K.
Of course, the absolute value and the coefficient of the temperature
dependence depend on the molecular type and temperature range, but in general,
this non-negligible temperature gradient gives the local stress gradient
around the interface known as the temperature Marangoni effect, which drives
the liquid flow from a higher-temperature region to a lower-temperature
region.

##  5 Concluding remarks

Surface is physically as well as theoretically at the boundary of macroscopic
and microscopic physics, understanding of surface tension intrinsically needs
multiple knowledge of thermodynamics, statistical mechanics as well as
continuum mechanics. In this article, after introducing the interpretation of
the phase change, liquid-vapor coexistence and interfacial tension from these
points of view, molecular dynamics results were shown which indicate the
existence of microscopic stress anisotropy at the interface as the origin of
the macroscopic stress surface tension tangential to the interface. In
addition, a simple model to explain the stress anisotropy due to the density
anisotropy was provided. At present, molecular dynamics as a powerful tool is
available, and I hope that unsolved issues about this fascinating physics,
especially related to the Marangoni effect and wetting, will be elucidated.

###### Acknowledgements.

The simulation results in this article are provided by Minori Shintaku and
Yule Ding, former members of the author’s group. I sincerely appreciate their
cooperation.

## References

  * (1) Marchand, A., Weijs, J. H., Snoeijer, J. H., Andreotti, B.: Why is surface tension a force parallel to the interface?. Am. J. Phys. 79, 999–1008 (2011)
  * (2) de Gennes P.-G., Brochard-Wyart, F., Quere, D.: Capillarity and Wetting Phenomena: Drops, Bubbles, Pearls, Waves Springer, Berlin (2003)
  * (3) Gao, L., McCarthy, T. J.: Wetting 101∘. Langmuir 25, 14105-14115 (2009)
  * (4) Thompson, A. P., Aktulga, H. M. , Berger, R., Bolintineanu, D. S., Brown, W. M., Crozier, P. S., in’t Veld, P. J., Kohlmeyer, A., Moore, S. G., Nguyen, T. D., Shan, R., Stevens, M. J. , Tranchida, J., Trott, C., Plimpton, S. J.: LAMMPS - a flexible simulation tool for particle-based materials modeling at the atomic, meso, and continuum scales. Comp. Phys. Comm. 271, 108171 (2022).
  * (5) http://www.gromacs.org
  * (6) Rowlinson, J. S.: “Cohesion: A Scientific History of Intermolecular Forces. Cambridge University Press (2002)
  * (7) Ono, S.: Surface tension. Kyoritsu, Tokyo (1980) (in Japanese).
  * (8) Young, T.: III. An essay on the cohesion of fluids. Philosophical Transactions of the Royal Society of London 95, 65-87 (1805)
  * (9) Bernoulli, D.: Hydrodynamica, sive de viribus et motibus fluidorum commentarii : Opus Academicum.
  * (10) Fujiwara, K., Hyodo, T.: Introduction to Thermal Physics, From a Macroscopic to Microscopic Approach, University of Tokyo Press, Tokyo (1995) (in Japanese).
  * (11) van der Waals, J. D.: Thermodynamische theorie der capillariteit in de onderstelling van continue dichtheidsverandering. Verhandel. Konink. Akad. Weten. Amsterdam (Sect. 1), Vol. 1, No. 8 (1893)
  * (12) Rowlinson, J. S.: Translation of J. D. van der Waals’ “The Thermodynamic Theory of Capillarity Under the Hypothesis of a Continuous Variation of Density.” J. Stat. Phys. 20, 197-200 (1979)
  * (13) Ono, S., Kondo, S.: Molecular Theory of Surface Tension in Liquids. pp. 134-280. Springer, Berlin (1960)
  * (14) Rowlinson, J. S., Widom, B.: Molecular Theory of Capillarity. Dover, New York (1982)
  * (15) Wilhelmy, L.: Über die Abhängigkeit der Capillaritäts-Constanten des Alkohols von Substanz und Gestalt des benetzten festen Körpers, Ann. Phys. 195, 177-217 (1863)
  * (16) Thompson, S. M., Gubbins, K. E., Walton, J. P. R. B., Chantry, R. A. R., Rowlinson, J. S.: A molecular dynamics study of liquid drops, J. Chem. Phys. 81, 530–542 (1984)
  * (17) Yaguchi, H. Yano, T. Fujikawa, S.: Molecular Dynamics Study of Vapor-Liquid Equilibrium State of an Argon Nanodroplet and Its Vapor, J. Fluid. Sci. Tech. 5 180-191 (2010)
  * (18) Surblys, D., Yamaguchi, Y., , Kuroda, K., Kagawa, M., Nakajima, T., Fujimura, H.: Molecular dynamics analysis on wetting and interfacial properties of water-alcohol mixture droplets on a solid surface, J. Chem. Phys. 140, 034505 (2014)
  * (19) Yamaguchi, Y., Kusudo, H., Surblys, D., Omori, T., Kikugawa, G.: Interpretation of Young’s equation for a liquid droplet on a flat and smooth solid surface: Mechanical and thermodynamic routes with a simple Lennard-Jones liquid, J. Chem. Phys. 150, 044701 (2019)
  * (20) Kusudo, H., Omori, T., Yamaguchi, Y.: Local stress tensor calculation by the method-of-plane in microscopic systems with macroscopic flow: A formulation based on the velocity distribution function, J. Chem. Phys. 155, 184103 (2021)


Criticality of global monopole charges in diverse dimensions

Hong-Ming Cui and Zhong-Ying Fan

Department of Astrophysics, School of Physics and Electronic Engineering,

Guangzhou University, Guangzhou 510006, China

ABSTRACT

In this work, we construct charged AdS black holes with a global monopole
charge in diverse dimensions and study the thermodynamics. We find a critical
monopole charge below which the solution exhibits Van-der Waals like
behaviors. In the context of holography, this could be intepreted using the
boundary degrees of freedoms. As an example, we study the phase diagram in the
four dimensions analytically. We further analyze the microstructures of the
solution using Ruppeiner geometry. We find that repulsive interactions
dominates for black holes in a wide range of temperatures. However, around the
critical point, attractive interactions is dominant and the Ruppeiner scalar
curvature shows universal behaviors: it has a critical exponent $2$ and
coefficient $-1/8$ in diverse dimensions. Universality of the results is
interpreted from the scaling behavior of free energy near the critical point
for Van-der Waals like fluids.

###### Contents

  1. 1 Introduction
  2. 2 Charged black holes with a global monopole
  3. 3 Extended black hole thermodynamics
  4. 4 Monopole charge criticality in 4-dimensions
    1. 4.1 Coexistence curve and the phase diagram
    2. 4.2 Rupepeiner geometry and microstructures
  5. 5 Monopole charge criticality in higher dimensions
    1. 5.1 Critical exponents
  6. 6 Conclusion

##  1 Introduction

Einstein’s gravity with a negative cosmological constant admits black hole
solutions in asymptotically anti-de Sitter (AdS) spacetimes. Generally the
solutions have richer phase structures relative to the asymptotically flat
counterparts. In the past decade, people was interested in studying extended
thermodynamics of AdS black holes, in which coupling constants in the
Lagrangian were taken as thermodynamic variables. In [1], by treating the
cosmological constant as a thermodynamic variable, referred to as
thermodynamic pressure $P$, it was established that charged AdS black holes
exhibit $P-V$ criticality, similar to the liquid-gas transition of Van-der
Waals fluids. More interestingly, Ruppeiner geometry associated to
thermodynamic fluctuations of the volume (or pressure) partly illustrates the
microstructures of the black holes [2, 3]. Recently the thermodynamics of
charged AdS black holes was further extended by introducing central charges in
the boundary CFTs [4, 5]. This corresponds to varying both the cosmological
constant and the Newton’s gravitational constant. An advantage of this is the
Van-der Waals like behaviors of AdS black holes can be interpreted in the
context of holography. Generalization to higher curvature theories can be
found in [6, 7].

In this work, we would like to construct global monopoles swallowed by charged
AdS black holes in diverse dimensions and study the thermodynamics.
Previously, the topic was partly studied in [8, 9] in the four dimensions but
using a ghost-like Goldstone field. This is however less interesting in
physics since the solution is generally unstable. We resolve this issue by
using a canonical multiplet and obtain the correct solution. Similar to the
asymptotically flat case [10], the global monopole produces a strong
gravitational field but has little contributions to the black hole mass.
However, the metric becomes globally singular because of a deficit solid
angle. We define the conical deficit (up to a constant factor) to be the
global monopole charge. We find that the solution exhibits new Van-der Waals
like behaviors associated to this new charge and its conjugate. The critical
point depends on the ratio of the electric charge $Q$ to (certain power of)
the central charge $C$ in the boundary. As a consequence, all dual fluids with
a same ratio will have the same critical point, leading to a critical line on
the $Q-C$ plane.

The paper is organized as follows. In section 2, we consider Einstein-Maxwell
theory with a cosmological constant, which is minimally coupled to a multiplet
scalar field and construct charged AdS black holes with a global monopole
charge in diverse dimensions. In section 3, we derive the first law for the
solutions and extend the thermodynamics appropriately. In section 4, we study
criticality of global monopole charges in the four dimensions in details. We
derive the coexistence curve and analyze the phase diagram analytically. We
further study the microstructures of the solution using the Ruppeiner
geometry. In section 5, we generalize the results to higher dimensions and
compute various critical exponents. We briefly conclude in section 6.

##  2 Charged black holes with a global monopole

Global monopoles form as a result of breaking of a global gauge symmetry. The
simplest model producing a global monopole is described by a multiplet scalar
field (2). To find charged AdS black hole solutions outside a global monopole,
we consider Einstein-Maxwell theory with a cosmological constant, minimally
coupled to a multiplet $\phi^{a}$ in diverse dimensions

$\displaystyle I=\int d^{D}x\,\sqrt{-g}\,\left({\frac{1}{16\pi
G}}(R-2\Lambda)-{\frac{1}{16\pi}}F^{2}+\mathcal{L}[\phi^{a}]\right)\,,$ (1)

where $a=1\,,2\,,\cdots\,,D-1$, $\Lambda=-{\frac{(D-1)(D-2)}{2L^{2}}}$ is the
cosmological constant and $L$ is the AdS radius. The Lagrangian of the scalar
field is given by

$\displaystyle\mathcal{L}=-{\frac{1}{2}}(\partial\phi^{a})^{2}-V(\phi^{a})\,,\quad
V(\phi^{a})={\frac{\lambda}{4}}(\phi^{a}\phi^{a}-\eta^{2})^{2}\,,$ (2)

where $\eta$ is the scale of gauge symmetry breaking. The model has a global
O($D-1$) symmetry, which will be spontaneously broken to O($D-2$). For a
typical unification scale $\eta\sim 10^{16}\mathrm{Gev}\ll\ell_{P}^{-1}$ ,
where $\ell_{P}$ is the Plank length. This is the four dimensional case.
Generalization to diverse dimensions, the scale $\eta$ should be bounded by
$\ell_{P}^{-{\frac{D-2}{2}}}$ on dimensional grounds.

The Einstein field equations is given by

$\displaystyle R_{\mu\nu}-{\frac{1}{2}}Rg_{\mu\nu}+\Lambda g_{\mu\nu}=8\pi
G\left(T_{\mu\nu}^{\mathrm{EM}}+T_{\mu\nu}^{GM}\right)\,,$ (3)

where

$\displaystyle
T_{\mu\nu}^{\mathrm{EM}}={\frac{1}{8\pi}}\left(F_{\mu\alpha}F^{\alpha}_{\,\,\,\nu}-{\frac{1}{4}}g_{\mu\nu}F^{2}\right)\,,$
$\displaystyle
T_{\mu\nu}^{GM}=\partial_{\mu}\phi^{a}\partial_{\nu}\phi^{a}-g_{\mu\nu}\left[{\frac{1}{2}}(\partial\phi^{a})^{2}+V(\phi^{a})\right]\,.$
(4)

The equation of motion of the scalar field reads

$\square\phi^{a}={\frac{\partial V}{\partial\phi^{a}}}\,.$ (5)

The generally static spherically symmetric charged solution takes the form of

$\displaystyle
ds^{2}=-\tilde{f}(r)dt^{2}+f(r)^{-1}dr^{2}+r^{2}d\Omega_{D-2}^{2}\,,$
$\displaystyle A=a(r)\,dt\,,\qquad\phi^{a}=\eta h(r){\frac{x^{a}}{r}}\,,$ (6)

where $d\Omega_{D-2}^{2}$ is the metric of a unit $\mathrm{S}^{D-2}$ and
$x^{a}x^{a}=r^{2}$. Generally speaking, full solutions to $h(r)$ and the
metric functions $\tilde{f}(r)\,,f(r)$ can only be found numerically. However,
since most of the monopole energy is concentrated in a small region near the
core, whose size is essentially determined by the scale of gauge symmetry
breaking $\delta\sim\eta^{-{\frac{2}{D-2}}}$, we may take the monopole core
$\delta\ll r_{+}$ in diverse dimensions, where $r_{+}$ denotes the black hole
horizon radius. Under this approximation, $h(r)=1$ and the potential
$V(\phi^{a})$ vanishes outside the horizon. The stress tensor of the multiplet
greatly simplifies to

$\displaystyle
T^{0\mathrm{(GM)}}_{\,\,\,0}=T^{1\mathrm{(GM)}}_{\,\,\,1}=-{\frac{(D-2)\eta^{2}}{2r^{2}}}\,,\qquad
T^{i\mathrm{(GM)}}_{\,\,\,i}=-{\frac{(D-4)\eta^{2}}{2r^{2}}}\,,$ (7)

where $i=1\,,2\,,\cdots\,,D$. Consequently, the metric function
$\tilde{f}(r)=f(r)$ so that the electrostatic potential can be solved
immediately from Maxwell’s equation

$a(r)=-{\frac{4\pi Q}{(D-3)\Omega_{D-2}\,r^{D-3}}}\,,$ (8)

where $\Omega_{D-2}=2\pi^{(D-1)/2}/\Gamma\left({\frac{D-1}{2}}\right)$ is the
volume of a unit $\mathrm{S}^{D-2}$ and $Q$ is the electric charge carried by
the black hole. Finally, the metric function can be solved as

$f(r)=1-2\zeta-{\frac{16\pi
GM}{(D-2)\Omega_{D-2}\,r^{D-3}}}+{\frac{32\pi^{2}GQ^{2}}{(D-1)(D-3)\Omega_{D-2}^{2}\,r^{2D-6}}}+{\frac{r^{2}}{L^{2}}}\,,$
(9)

where

$\displaystyle\zeta={\frac{4\pi G\eta^{2}}{D-3}}\,,$ (10)

is referred to as the global monopole charge. Clearly the presence of a global
monopole charge does not contribute to the black hole mass111Precisely
speaking, the monopole has energy
$M_{\mathrm{core}}\sim\eta^{{\frac{2}{D-2}}}$ which is much less than the
black hole mass.. All its physical effect is inducing a conical deficit at
asymptotical infinity. To see this, let $t\rightarrow t/\sqrt{1-2\zeta}$,
$r\rightarrow r\sqrt{1-2\zeta}$, one finds asymptotically

$ds^{2}\Big{|}_{r\rightarrow\infty}=-\left({\frac{r^{2}}{L^{2}}}+1\right)dt^{2}+{\frac{dr^{2}}{{\frac{r^{2}}{L^{2}}}+1}}+(1-2\zeta)\,r^{2}d\Omega_{D-2}^{2}\,.$
(11)

Clearly there exists a conical deficit for the solid angle
$\Delta=\Omega_{D-2}\,\zeta$. Besides, the monopole charge should be upper
bounded as $\zeta<1/2$, which translates into
$\eta<\sqrt{{\frac{D-3}{8\pi}}}\,\ell_{P}^{-{\frac{D-2}{2}}}$.

##  3 Extended black hole thermodynamics

We proceed to study thermodynamics of the charged AdS black holes with a
global monopole. Since the monopole charge is also an integration constant of
the solution, the first law of thermodynamics acquires an extra contribution
from a new pair of conjugates $(U\,,\zeta)$

$\displaystyle dM=TdS+\Phi dQ+Ud\zeta\,,$ (12)

where $U$ is the chemical potential dual to $\zeta$. The various quantities
are given by

$\displaystyle M={\frac{(D-2)\Omega_{D-2}\,r_{+}^{D-3}}{16\pi
G}}\left(Z+{\frac{r_{+}^{2}}{L^{2}}}\right)+{\frac{2\pi
Q^{2}}{(D-3)\Omega_{D-2}\,r_{+}^{D-3}}}\,,$ $\displaystyle T={\frac{1}{4\pi
r_{+}}}\left[(D-3)Z+(D-1){\frac{r_{+}^{2}}{L^{2}}}\right]-{\frac{8\pi
GQ^{2}}{(D-2)\Omega_{D-2}^{2}r_{+}^{2D-5}}}\,,$ $\displaystyle
S={\frac{\Omega_{D-2}\,r_{+}^{D-2}}{4G}}\,,\quad\Phi={\frac{4\pi
Q}{(D-3)\Omega_{D-2}\,r_{+}^{D-3}}}\,,\quad U=-{\frac{(D-2)\Omega_{D-2}}{8\pi
Gr_{+}^{D-3}}}\,,$ (13)

where $Z\equiv 1-2\zeta$ is a useful combination here and after. It turns out
that there exists a critical monopole charge below which the solution shows
new Van-der Waals like behaviors in diverse dimensions.

However before studying this in details, we would like to extend the black
hole thermodynamics by introducing the thermodynamic pressure [1] and the
central charge [4] dual to the boundary

$\displaystyle P=-{\frac{\Lambda}{8\pi G}}\,,$ $\displaystyle
C={\frac{kL^{D-2}}{16\pi G}}\,,$ (14)

where $k$ is a constant depending on details of holographic systems. The first
law is extended to

$\displaystyle dM=TdS+\Phi dQ+Ud\zeta+V_{C}dP+\mu dC\,,$ (15)

where the thermodynamic volume $V_{C}$ and the chemical potential $\mu$ dual
to the central charge are given by

$\displaystyle
V_{C}={\frac{\Omega_{D-2}r_{+}^{D-1}}{D(D-1)}}+{\frac{\Omega_{D-2}L^{2}Zr_{+}^{D-3}}{D(D-1)}}+{\frac{32\pi^{2}GL^{2}Q^{2}}{D(D-1)(D-2)\Omega_{D-2}r_{+}^{D-3}}}\,,$
$\displaystyle\mu=-{\frac{2(D-1)\Omega_{D-2}r_{+}^{D-1}}{DkL^{D}}}+{\frac{2\Omega_{D-2}Zr_{+}^{D-3}}{DkL^{D-2}}}+{\frac{64\pi^{2}GQ^{2}}{D(D-2)\Omega_{D-2}kL^{D-2}r_{+}^{D-3}}}\,.$
(16)

Here the thermodynamic volume $V_{C}$ should not be confused with the original
one defined at constant $G$ [11]. The two volumes are related by

$V=V_{C}+\mu{\frac{\partial C}{\partial P}}\Big{|}_{G}=V_{C}-{\frac{(D-2)\mu
C}{2P}}\,.$ (17)

Though it is widely believed that $V$ obeys the reverse isoperimetric
inequality, this new thermodynamic volume $V_{C}$ generally does not [4]. The
Smarr relation is given by

$M={\frac{D-2}{D-3}}\big{(}TS+\mu C\big{)}+\Phi Q-{\frac{2}{D-3}}PV_{C}\,,$
(18)

where the monopole charge does not have any contributions (since $\zeta$ is
dimensionless).

Dependence of the various thermodynamic quantities on the pressure $P$ and the
central charge $C$ can be seen by using the inverse formula

$\displaystyle
G={\frac{1}{16\pi}}\left({\frac{k}{C}}\right)^{2/D}\left[{\frac{(D-1)(D-2)}{P}}\right]^{1-2/D}\,,\quad
L=\left[{\frac{(D-1)(D-2)C}{kP}}\right]^{1/D}\,.$ (19)

Since we are interested in $U-\zeta$ criticality, we will work in canonical
ensemble which has fixed charges $(Q\,,C\,,P\,,Z)$ or equivalently
$(Q\,,G\,,L\,,\zeta)$. The equation of state of the dual fluid can be formally
written as $T=T(U\,,\zeta)$. Explicit expressions for the equation of state
will be presented later when we focus on certain dimensions. Here we shall
point out that the existence of a critical point and the Van-der Waals-like
behaviors associated to the monopole charge does not rely on the extended
thermodynamics introduced above. All we have done in this part is trying to
interpret the results in the context of holography, as will be shown later.

Figure 1: The chemical potential $U$ is displayed as a function of $T$ for
fixed monopole charge. The various curves correspond to: $\zeta<\zeta_{c}$
(solid blue), $\zeta=\zeta_{c}$ (dashed red) and $\zeta>\zeta_{c}$ (dotted
brown). Here we have set $G=1$, $L=1$ and $Q=0.1$.

##  4 Monopole charge criticality in 4-dimensions

To study the $U-\zeta$ criticality in details, we consider the four dimensions
at first. The black hole solution reads

$\displaystyle ds^{2}=-f(r)dt^{2}+{\frac{1}{f(r)}}dr^{2}+r^{2}d\Omega^{2}\,,$
$\displaystyle
f(r)=1-2\zeta-{\frac{2GM}{r}}+{\frac{GQ^{2}}{r^{2}}}+{\frac{r^{2}}{L^{2}}}\,.$
(20)

The relevant thermodynamic quantities are given by

$\displaystyle
M={\frac{r_{+}^{3}}{2GL^{2}}}+{\frac{Z}{2G}}r_{+}+{\frac{Q}{2r_{+}}}\,,$
$\displaystyle T={\frac{3r_{+}}{4\pi L^{2}}}+{\frac{Z}{4\pi
r_{+}}}-{\frac{GQ^{2}}{4\pi r_{+}^{3}}}\,,$ $\displaystyle S={\frac{\pi
r_{+}^{2}}{G}}\,,\quad U=-{\frac{r_{+}}{G}}\,.$ (21)

Furthermore, the equation of state $T(U\,,\zeta)$ can be written explicitly as

$T=-{\frac{2Z}{U}}\sqrt{{\frac{2PC}{3k}}}-{\frac{3kU}{64\pi^{2}C}}+{\frac{32\pi
PCQ^{2}}{3kU^{3}}}\,,$ (22)

where we have adopted the relation $U=-8\pi r_{+}\sqrt{{\frac{2PC}{3k}}}$.
Using inflection point condition $0=\partial T/\partial
r_{+}=\partial^{2}T/\partial r_{+}^{2}$, we obtain the critical point

$\displaystyle r_{c}={\frac{3^{1/4}\sqrt{Q}}{2^{3/4}\pi^{1/4}P^{1/4}}}\,,$
$\displaystyle Z_{c}={\frac{3Q}{2}}\sqrt{{\frac{k}{\pi C}}}\,,$ $\displaystyle
T_{c}={\frac{P^{1/4}}{6^{1/4}\pi^{5/4}}}\sqrt{{\frac{kQ}{C}}}\,.$ (23)

It is immediately seen that the critical monopole charge depends on the ratio
of electric charge to (square root of) the central charge in the boundary. In
fact, this gives a critical line on the $Q-C$ plane: all the dual fluids with
a same ratio $Q/\sqrt{C}$ will have the same critical point. Interestingly,
the product

${\frac{T_{c}r_{c}}{Z_{c}}}={\frac{1}{3\pi}}\,,$ (24)

is a constant, independent of the charges on the boundary (but unfortunately
we do not find a proper interpretation for this). Define dimensionless
variables $t=T/T_{c}$, $u=U/U_{c}$, $z=Z/Z_{c}$. The equation of state
translates into the “law of corresponding state”

$\displaystyle t=\alpha u+{\frac{\beta}{u}}+{\frac{\gamma}{u^{3}}}\,,$ (25)

where the various parameters are given by

$\alpha={\frac{3}{8}}\,,\quad\beta={\frac{3z}{4}}\,,\quad\gamma=-{\frac{1}{8}}\,.$
(26)

This is a very useful relation in solving the coexistence curve as well as
computing various critical exponents.

Figure 2: $F-T$ diagram is displayed as a function of $\zeta$ for fixed $G=1$,
$L=1$ and $Q=0.1$. The various curves correspond to: $\zeta<\zeta_{c}$ (solid
blue), $\zeta=\zeta_{c}$ (dashed red) and $\zeta>\zeta_{c}$ (dotted brown). We
have set $k=16\pi$.

The $U-T$ diagram for a given monopole charge is depicted in Fig. 1. It
follows that for $\zeta<\zeta_{c}$, the chemical potential $U$ has multi-
values for a given temperature. This indicates a first order transition
between a pair of black holes, across which $U$ changes discontinuously. To
ensure this, we plot the $F-T$ diagram in Fig. 2, where $F=M-TS$ is the free
energy, given by

$F=-{\frac{2}{3}}\pi Pr_{+}^{3}+2\pi
Zr_{+}\sqrt{{\frac{2CP}{3k}}}+{\frac{3Q^{2}}{4r_{+}}}\,.$ (27)

It is easily seen that for a smaller monopole charge $\zeta<\zeta_{c}$, the
swallowtail behavior emerges. It implies a first order phase transition for a
small monopole charge but no such transition for a large monopole charge. The
situation is somehow similar to the electric charge criticality [12], in which
the transition between multi-phases is only possible when the number of
charged carriers is sufficiently small, though microstructures of the two
cases are quite different.

###  4.1 Coexistence curve and the phase diagram

A remarkable feature of the four dimensional solution is the coexistence curve
can be derived analytically. We will generalize the method developed in [13]
and adopt the Maxwell’s area law on the $T-S$ plane.

Figure 3: Maxwell’s area law on the $T-S$ graph. The oscillatory region is
unphysical and should be replaced by an isotherm $t^{*}$, according to the
area law. The rectangular region below the isotherm is the coexistence region.
We have set $z={\frac{9}{8}}$, that is $\zeta={\frac{18\zeta_{c}-1}{16}}$ and
$t^{*}\approx 1.09$.

For this purpose, it will be convenient for us to work with the law of
corresponding state (25). However, to keep the discussions below as general as
possible, we will not specify the parameters $(\alpha\,,\beta\,,\gamma)$ until
the end of this subsection. The equal area law gives

$\displaystyle
t^{*}\left(u_{g}^{2}-u_{l}^{2}\right)=\int^{u^{2}_{g}}_{u^{2}_{l}}t(u)du^{2}\,,$
(28)

where $u_{l}$ and $u_{g}$ ($u_{l}<u_{g}$) denotes the normalized chemical
potential of the small and the large black hole on the the isotherm $t=t^{*}$,
respectively (see Fig. 3 for an illustration). Integrating the law of
corresponding state yields

$t^{*}={\frac{2}{u_{g}+u_{l}}}\left[{\frac{\alpha}{3}}(u^{2}_{g}+u^{2}_{l}+u_{g}\,u_{l})+\beta+{\frac{\gamma}{u_{g}\,u_{l}}}\right]\,.$
(29)

On the other hand, the two phases must have the same temperature

$\displaystyle t^{*}=t(u_{g})=t(u_{l})\,,$ $\displaystyle
t^{*}={\frac{t(u_{g})+t(u_{l})}{2}}\,,$ (30)

which leads to

$\displaystyle\left(u_{g}+u_{l}\right)^{2}-{\frac{\alpha}{\gamma}}(u_{g}u_{l})^{3}+{\frac{\beta}{\gamma}}(u_{g}u_{l})^{2}-u_{g}u_{l}=0\,,$
$\displaystyle
t^{*}={\frac{u_{g}+u_{l}}{2}}\left[\alpha+{\frac{\beta}{u_{g}\,u_{l}}}+\gamma{\frac{(u_{g}+u_{l})^{2}-u_{g}\,u_{l}}{(u_{g}\,u_{l})^{3}}}\right]\,.$
(31)

It follows that the equations (29) and (31) can be solved analytically by
taking $u_{g}u_{l}$ and $u_{g}+u_{l}$ as independent variables. One finds

$\displaystyle\left[\alpha(u_{g}u_{l})^{2}+3\gamma\right]\left[\beta
u_{g}\,u_{l}-\alpha(u_{g}u_{l})^{2}+3\gamma\right]=0\,,$ (32)

which gives two set solutions

$\displaystyle\alpha(u_{g}u_{l})^{2}+3\gamma=0\quad$
$\displaystyle\Rightarrow$ $\displaystyle\quad
u_{g}u_{l}=\sqrt{-{\frac{3\gamma}{\alpha}}}\,,$ (33) $\displaystyle\beta
u_{g}\,u_{l}-\alpha(u_{g}u_{l})^{2}+3\gamma=0\quad$ $\displaystyle\Rightarrow$
$\displaystyle\quad
u_{g}u_{l}={\frac{\beta\pm\sqrt{\beta^{2}+12\alpha\gamma}}{2\alpha}}\,\,.$
(34)

However, the second solution is unphysical since it will lead to
$u^{2}_{g}=u^{2}_{l}={\frac{\beta\pm\sqrt{\beta^{2}+12\alpha\gamma}}{2\alpha}}$.
From (33), we obtain the physical solutions

$\displaystyle
u^{2}_{g,\,l}={\frac{1}{2\alpha}}\left[3\beta-4\sqrt{-3\alpha\gamma}\pm\sqrt{3(3\beta-2\sqrt{-3\alpha\gamma})(\beta-2\sqrt{-3\alpha\gamma})}\,\right]\,,$
(35)

where the sign “$+$” (“$-$”) corresponds to the saturated large (small) black
hole, respectively. Substituing the parameters (26) into (35) and (29), we
finally obtain

$\displaystyle u_{g,l}=\sqrt{3z-2\pm\sqrt{(3z-1)(3z-3)}}\,,$ (36)

and the coexistence curve on the $t-z$ plane

$\displaystyle t^{*}(z)=\sqrt{{\frac{3z-1}{2}}}\,.$ (37)

Notice that $z>1/3$ is guaranteed by the upper bound of the monopole charge
$\zeta<1/2$.

Figure 4: The coexistence curve of for charged AdS black holes with a global
monople on the $t-z$ plane. Figure 5: Phase diagram of charged AdS black holes
with a global monopole on the $t-u$ plane. Blue and red solid curves are
separated by the critical point (black dot), corresponding to the saturated
small and large black holes respectively. The three phases are represented by
the shaded regions with different color.

The coexistence line is depicted in Fig. 4. Notice that the transition occurs
for $T>T_{c}$, which is equivalent to $\zeta<\zeta_{c}$. Furthermore, combing
(36) and (37), we deduce the coexistence curve on the $t-u$ plane

$t^{*}(u)={\frac{u^{2}+1}{2u}}\,.$ (38)

With these results in hand, we are ready to draw the phase diagram, as shown
in Fig. 5. It is easily seen that given a small monopole charge
$\zeta<\zeta_{c}$, there exists three phases for the solution:
$\bullet$ for $u<u_{l}$, there is a single small black hole in the “liquid-
phase”;
$\bullet$ for $u>u_{g}$, there is a single large black hole in the “gas-
phase”;
$\bullet$ for $u_{l}<u<u_{g}$, the black hole is in a coexisting “liquid-gas
phases”, transiting from one phase to another.
These results are standard for Van-der Waals-like fluids.

###  4.2 Rupepeiner geometry and microstructures

Figure 6: Density plot for the normalized scalar curvature $R_{N}$ in the
$D=4$ dimensions.

Having established the phase diagram, we move to study microstructures of the
dual fluids using Ruppeiner geometry. In thermodynamic fluctuation theory, the
distance between two neighboring fluctuation states $\Delta l^{2}$ is related
to the probability of finding the macroscopic system between the neighboring
states $P(x^{0}\,,\,\cdots,x^{N})\varpropto e^{-{\frac{1}{2}}\Delta l^{2}}$,
where

$\displaystyle\Delta l^{2}=-{\frac{\partial^{2}S}{\partial x^{\mu}\partial
x^{\nu}}}\Delta x^{\mu}\Delta x^{\nu}\,.$ (39)

We choose the temperature $T$ and the chemical potential $U$ to be the
independent fluctuation variables. The line element reads

$\displaystyle
dl^{2}={\frac{C_{U}}{T^{2}}}dT^{2}-{\frac{1}{T}}\left({\frac{\partial\zeta}{\partial
U}}\right)_{T}dU^{2}\,,$ (40)

where $C_{U}=T\left({\frac{\partial S}{\partial T}}\right)_{U}$ is the heat
capacity at constant $U$. However, since $U\propto\sqrt{S}$, the heat capacity
$C_{U}$ vanishes identically. To solve this issue, we adopt the trick in [3]:
taking $C_{U}$ to be a nonzero constant at first and then sending
$C_{U}\rightarrow 0$. In this limit, the Ruppeiner scalar curvature diverges
as $R\propto 1/C_{U}$. Then we define a normalized scalar curvature
$R_{N}\equiv RC_{U}$ and propose that it describes the microstructures of the
dual fluids: a positive (negative) $R_{N}$ indicating that repulsive
(attractive) interactions dominates.

Evaluation of the normalized scalar curvature yields

$\displaystyle
R_{N}=-{\frac{(3u^{4}+1)(8tu^{3}-3u^{4}-1)}{2(4tu^{3}-3u^{4}-1)^{2}}}\,.$ (41)

The behavior of $R_{N}$ is depicted in Fig. 6. Roughly speaking, $R_{N}$ is
very small in most of the parameters space. However, close to the curve

$\displaystyle t_{div}={\frac{3u^{4}+1}{4u^{3}}}\,,$ (42)

the curvature $R_{N}$ changes dramatically. This curve exactly corresponds to
a divergent heat capacity $C_{\zeta}$ defined at constant monopole charge,
signaling critical point of the transition. On this curve, $R_{N}$ goes to
negative infinity. Besides, there exists a sign-changing curve across which
$R_{N}$ changes the sign

$\displaystyle t_{0}={\frac{3u^{4}+1}{8u^{3}}}\,.$ (43)

It follows that for $t>t_{0}$, $R_{N}<0$, implying the attractive interactions
dominates in this region. Otherwise, $R_{N}$ will be positive and the
repulsive interactions dominates instead.

Figure 7: Characteristic curves for $D=4$ dimensions. The black solid line is
the coexistence curve, and the red dot-dashed line is sign-changing curve
corresponding to $t_{0}$, on which $R_{N}=0$. The blue dashed line corresponds
to the temperature $t_{div}$, on which $R_{N}\rightarrow-\infty$. In the
shaded region $R_{N}>0$ and otherwise $R_{N}<0$.

In Fig. 7, we illustrate the coexistence curve (black solid line), the sign-
changing curve (red dot-dashed line) as well as the diverging curve (blue
dashed line), respectively. The shaded (unshaded) regions correspond to a
positive (negative) curvature $R_{N}$, indicating that the repulsive
(attractive) interactions dominates. However, any features above the
coexistence curve are tentative since the equation of state we adopted is
inapplicable in the coexistence region. The results that we have confidence to
state are those below the coexistence curve. In region II, repulsive
interactions dominates in the microstructures, including most of the saturated
small black holes. However, around the critical point, attractive interactions
is dominant. To explain this, let us expand the scalar curvature in powers of
$\delta t=t-1$. We deduce

$\displaystyle R_{N}=-{\frac{1}{8\delta t^{2}}}\pm{\frac{1}{(2\delta
t)^{3/2}}}+\cdots\,,$ (44)

where “+/-” sign corresponds to the saturated small/large black holes,
respectively. To leading order $R_{N}$ has a critical exponent $2$ and the
coefficient

$\displaystyle\lim_{\delta t\rightarrow 0}R_{N}\delta t^{2}=-{\frac{1}{8}}\,.$
(45)

In fact, these results are universal to charged AdS black holes with a global
monopole in diverse dimensions and match with the Van-der Waals fluid
perfectly [3]. Universality of the result could be interpreted from the
scaling behavior of free energy near the critical point, see the next section.

Figure 8:  The F-T digram for $D=5$ (left) and $D=6$ (right). In both panels,
$z>1$ (solid blue), $z=1$ (dashed red) and $z<1$ (dotted brown).

##  5 Monopole charge criticality in higher dimensions

We continue studying the critical phenomenon of monopole charges in higher
dimensions. The equation of state $T=T(U\,,\zeta)$ is listed as follows

$\displaystyle T={\frac{1}{4\pi
r_{+}}}\left[(D-3)Z+(D-1){\frac{r_{+}^{2}}{L^{2}}}\right]-{\frac{8\pi
GQ^{2}}{(D-2)\Omega_{D-2}^{2}r_{+}^{2D-5}}}\,,$ (46)

where

$\displaystyle r_{+}^{D-3}=-{\frac{(D-2)\Omega_{D-2}}{8\pi GU}}\,.$ (47)

Evaluation of the critical point using the inflection point condition yields

$\displaystyle r_{c}$ $\displaystyle=$
$\displaystyle\left[{\frac{2(D-3)(2D-5)\pi
Q^{2}}{\Omega_{D-2}^{2}P}}\right]^{{\frac{1}{2D-4}}}\,,$ (48) $\displaystyle
Z_{c}$ $\displaystyle=$
$\displaystyle{\frac{2\big{[}(D-2)(D-1)\big{]}^{{\frac{D-2}{D}}}[2\pi(D-3)(2D-5)]^{{\frac{1}{D-2}}}}{(2D-5)\pi}}{\frac{k^{{\frac{2}{D}}}Q^{{\frac{2}{D-2}}}}{\Omega_{D-2}^{{\frac{1}{D-2}}}C^{{\frac{2}{D}}}P^{{\frac{D-4}{D(D-2)}}}}}\,,$
(49) $\displaystyle T_{c}$ $\displaystyle=$
$\displaystyle{\frac{\big{[}(D-2)(D-1)\big{]}^{{\frac{D-2}{D}}}[2\pi(D-3)(2D-5)]^{{\frac{1}{2D-4}}}}{(2D-5)\pi}}{\frac{k^{{\frac{2}{D}}}Q^{{\frac{1}{D-2}}}}{\Omega_{D-2}^{{\frac{1}{D-2}}}C^{{\frac{2}{D}}}P^{{\frac{3D-8}{2D(D-2)}}}}}\,.$
(50)

Again the critical monopole charge depends on the ratio of the electric charge
to (certain power of) the central charge: $Q/C^{{\frac{D-2}{D}}}$. In
addition, the product

$\displaystyle{\frac{T_{c}r_{c}}{Z_{c}}}={\frac{(D-3)^{2}}{(2D-5)\pi}}\,,$
(51)

turns out to be a constant in any given dimensions.

Working with dimensionless variables $t=T/T_{c}\,,u=U/U_{c}\,,z=Z/Z_{c}$, we
deduce the law of corresponding state

$\displaystyle
t={\frac{2D-5}{4(D-2)}}u^{{\frac{1}{D-3}}}+{\frac{2D-5}{4(D-3)}}zu^{-{\frac{1}{D-3}}}-{\frac{u^{-{\frac{2D-5}{D-3}}}}{4(D-2)(D-3)}}\,,$
(52)

and the normalized free energy $f=F/F_{c}$

$\displaystyle
f={\frac{D-1}{4}}z\,u+{\frac{D-1}{4(D-2)u}}-{\frac{(D-3)^{2}}{4(D-2)}}\,u^{{\frac{D-1}{D-3}}}\,.$
(53)

The $F-T$ digrams for higher dimensional solutions are depicted in Fig. 8. The
existence of multi-phases and the transitions between them (the coexistence
curve) for a small monopole charge $\zeta<\zeta_{c}$ as well as the full phase
digram are quite similar to the four dimensional case. We shall not discuss
them anymore.

To proceed, we study the microstructures for higher dimensional solutions
using the Ruppeiner geometry. We deduce

$\displaystyle
R_{N}=-{\frac{\left[(2D-5)\,u^{{\frac{2(D-2)}{D-3}}}+1\right]\left[4(D-2)t\,u^{{\frac{2D-5}{D-3}}}-(2D-5)\,u^{{\frac{2(D-2)}{D-3}}}-1\right]}{2\left[2(D-2)t\,u^{{\frac{2D-5}{D-3}}}-(2D-5)\,u^{{\frac{2(D-2)}{D-3}}}-1\right]^{2}}}\,.$
(54)

Again there exists a diverging curve

$\displaystyle
t_{div}={\frac{(2D-5)\,u^{{\frac{2(D-2)}{D-3}}}+1}{2(D-2)\,u^{{\frac{2D-5}{D-3}}}}}\,,$
(55)

on which $R_{N}$ goes to negative infinity and a sign-changing curve

$\displaystyle
t_{0}={\frac{(2D-5)\,u^{{\frac{2(D-2)}{D-3}}}+1}{4(D-2)\,u^{{\frac{2D-5}{D-3}}}}}\,,$
(56)

across which the attractive and the repulsive interactions changes from one
another. In Fig. 9, we depict these curves as well as the coexistence lines
(which are computed numerically) for $D=5$ and $D=6$ dimensions, respectively.
Similar to the four dimensional case, small black holes are dominated by
repulsive interactions whereas the large black holes (at higher temperatures
$t>t_{0}$) are dominated by attractive interactions. In particular,
attractions always dominates around the critical point. The leading order
behavior of the scalar curvature (45) near the critical point is universal to
diverse dimensions.

Figure 9:  Behavior of the normalized scalar curvature $R_{N}$ of $D=5$ (left)
and $D=6$ (right) dimensions. All characteristic curves have the same meaning
as those in Fig. 7. In the shaded region $R_{N}>0$ and otherwise $R_{N}<0$.

Universality of the result could be interpreted using the two-scale factor
universality hypothesis [14]. Near the critical point, the singular part of
the free energy connected to critical phenomena can be promoted to be a
generalized homogeneous function [15, 16]:

$\displaystyle F_{s}=c_{1}|\delta t|^{\beta(\delta+1)}Y\big{(}c_{2}\,v|\delta
t|^{-\beta\delta}\big{)}\,,$ (57)

where $\beta\,,\delta$ are the critical exponents defined as usual (see next
subsection), $v=u-1$ is the order parameter and $Y()$ is a smooth function,
obeying the conditions $Y(0)<0$ and $Y^{\prime}(0)=0$. Note that the constants
$c_{1}\,,c_{2}$ are two scale factors, which are material dependent and
differentiate between systems even in the same universality class. Actually,
in modern approach to critical phenomena, the behavior of free energy can be
elegantly explained using renormalization group analysis in quantum field
theory.

By straightforward calculations, we deduce

$\displaystyle
R={\frac{\beta(\beta-\alpha)T_{c}}{2(2-\alpha)(1-\alpha)c_{1}Y(0)}}\,\delta
t^{\alpha-2}\,,$ (58) $\displaystyle
C_{U}=-{\frac{(2-\alpha)(1-\alpha)c_{1}Y(0)}{T_{c}}}\,\delta t^{-\alpha}\,,$
(59)

where we have made use of the scaling equation $2-\alpha=\beta(\delta+1)$. It
is immediately seen that the critical behaviors of both $R$ and $C_{U}$ depend
on details of the fluids, through the scale factor $c_{1}$. However, the
leading order behavior of the normalized curvature $R_{N}$ is universal

$\displaystyle\lim_{\delta t\rightarrow 0}R_{N}\delta
t^{2}=-{\frac{\beta(\beta-\alpha)}{2}}\,,$ (60)

depending only on the critical exponents $\alpha\,,\beta$. For Van-der Waals
like fluids, $\alpha=0\,,\beta=1/2$, this reproduces (45) exactly. The result
is robust since neither a constant $C_{U}$ nor the limit $C_{U}\rightarrow 0$
plays an essential role.

###  5.1 Critical exponents

We move to compute the various critical exponents for $U-\zeta$ criticality in
diverse dimensions. The exponents $\alpha$, $\beta$, $\gamma$, $\delta$ are
defined as follows:
$\bullet$ Exponent $\alpha$ describes the behavior of the specific heat at
constant $U$:

$\displaystyle C_{U}=T\left({\frac{\partial S}{\partial
T}}\right)_{U}\propto\lvert\delta t\rvert^{-\alpha}\,.$ (61)

$\bullet$ Exponent $\beta$ describes the behavior of the order parameter
$u_{g}-u_{l}$ on the given isotherm

$\displaystyle u_{g}-u_{l}\propto\lvert\delta t\rvert^{\beta}\,.$ (62)

$\bullet$ Exponent $\gamma$ determines the behavior of the isothermal
compressibility $\kappa_{T}$

$\displaystyle\kappa_{T}=-{\frac{1}{U}}\left({\frac{\partial
U}{\partial\zeta}}\right)_{T}\propto\lvert\delta t\rvert^{-\gamma}\,.$ (63)

$\bullet$ Exponent $\delta$ describes the following behavior on the critical
isotherm $T=T_{c}$:

$\displaystyle\lvert\zeta-\zeta_{c}\rvert\propto\lvert
U-U_{c}\rvert^{\delta}\,.$ (64)

Firstly, since the entropy $S$ is a function of $U$, the specific heat $C_{U}$
vanishes identically. This implies $\alpha=0$.

To compute $\beta$, we first rewrite the law of corresponding state as

$\displaystyle z={\frac{4(D-3)}{2D-5}}(\delta
t+1)(\upsilon+1)^{{\frac{1}{D-3}}}-{\frac{D-3}{D-2}}(\upsilon+1)^{{\frac{2}{D-3}}}+{\frac{1}{(D-2)(2D-5)(\upsilon+1)^{2}}}\,.$
(65)

Expanding the equation to relevant orders in the vicinity of the critical
point yields

$\displaystyle z=1+{\frac{4D-12}{2D-5}}\delta t+{\frac{4}{2D-5}}\delta
t\upsilon-{\frac{2}{3(D-3)^{2}}}\upsilon^{3}+O(\delta
t\upsilon^{2}\,,\,\upsilon^{4})\,.$ (66)

During the phase transition, the temperature remains a constant, which gives

$\displaystyle z$ $\displaystyle=1+{\frac{4D-12}{2D-5}}\delta
t+{\frac{4}{2D-5}}\delta t\upsilon_{g}-{\frac{2}{3(D-3)^{2}}}\upsilon_{g}^{3}$
(67) $\displaystyle=1+{\frac{4D-12}{2D-5}}\delta t+{\frac{4}{2D-5}}\delta
t\upsilon_{l}-{\frac{2}{3(D-3)^{2}}}\upsilon_{l}^{3}\,,$

Besides, using Maxwell’s area law, one has

$\displaystyle 0=\int^{\upsilon_{g}}_{\upsilon_{l}}\upsilon
dz=\int^{\upsilon_{g}}_{\upsilon_{l}}\left[{\frac{4\delta
t}{2D-5}}-{\frac{2\upsilon^{2}}{3(D-3)^{2}}}\right]\upsilon d\upsilon=0\,.$
(68)

Combining the equations, one finds the solution is
$\upsilon_{g}=-\upsilon_{l}\propto\sqrt{-\delta t}$. This leads to
$\beta={\frac{1}{2}}$.

Differentiating (66), we get

$\displaystyle\left({\frac{\partial
U}{\partial\zeta}}\right)_{T}={\frac{2D-5}{4}}{\frac{\upsilon_{c}\,z}{z_{c}\,\delta
t}}\qquad\Rightarrow\qquad\gamma=1\,.$ (69)

Finally, according to (66), the critical isotherm $\delta t=0$ adheres to

$\displaystyle
z-1=-{\frac{2\upsilon^{3}}{3(D-3)^{2}}}\qquad\Rightarrow\qquad\delta=3\,.$
(70)

This completes our derivations.

##  6 Conclusion

In this work, we constructed charged AdS black holes with a global monopole in
diverse dimensions. The presence of a monopole charge $\zeta$ and its
thermodynamic conjugate $U$ contributes to the first law of thermodynamics. We
found that there exists a critical monopole charge $\zeta_{c}$ below which the
solution exhibits Van-der Waals like behaviors. In the context of holography,
this could be interpreted using the boundary degrees of freedoms. The critical
point depends on the ratio of the electric charge to (certain power of) the
central charge.

As an example, we derived the coexistence curve in the four dimensions and
studied the phase diagram analytically. The results of higher dimensions are
qualitatively similar. We also investigated the microstructures of the
solution using thermodynamic geometry. We defined a normalized scalar curve
$R_{N}$ and established that the leading order behavior (45) near the critical
point is universal to diverse dimensions. Universality of the results is
interpreted from the scaling behavior of free energy near the critical point
for Van-der Waals like fluids. All these results improve our understanding of
the critical phenomenon for charged AdS black holes with a global monopole
charge significantly.

## Acknowledgments

Z.Y. Fan was supported in part by the National Natural Science Foundations of
China with Grant No. 11805041 and No. 11873025 and also supported in part by
Guangzhou Science and Technology Project 2023A03J0016.

## References

  * [1] D. Kubiznak and R. B. Mann, P-V criticality of charged AdS black holes, JHEP 07, 033 (2012) [arXiv:1205.0559 [hep-th]].
  * [2] S. W. Wei and Y. X. Liu, Insight into the Microscopic Structure of an AdS Black Hole from a Thermodynamical Phase Transition, Phys. Rev. Lett. 115, no.11, 111302 (2015) [erratum: Phys. Rev. Lett. 116, no.16, 169903 (2016)].
  * [3] S. W. Wei, Y. X. Liu and R. B. Mann, Repulsive Interactions and Universal Properties of Charged Anti–de Sitter Black Hole Microstructures, Phys. Rev. Lett. 123, no.7, 071103 (2019).
  * [4] W. Cong, D. Kubiznak and R. B. Mann, Thermodynamics of AdS Black Holes: Critical Behavior of the Central Charge, Phys. Rev. Lett. 127, no.9, 091301 (2021) [arXiv:2105.02223 [hep-th]].
  * [5] W. Cong, D. Kubiznak, R. B. Mann and M. R. Visser, Holographic CFT phase transitions and criticality for charged AdS black holes, JHEP 08, 174 (2022) [arXiv:2112.14848 [hep-th]].
  * [6] H. M. Cui and Z. Y. Fan, Criticality of central charges for Gauss-Bonnet black holes, [arXiv:2404.05945 [hep-th]].
  * [7] N. Kumar, S. Sen and S. Gangopadhyay, Breaking of the universal nature of the central charge criticality in AdS black holes in Gauss-Bonnet gravity, Phys. Rev. D 107, no.4, 046005 (2023) [arXiv:2211.00925 [gr-qc]].
  * [8] A. K. Ahmed, U. Camci and M. Jamil, Accretion on Reissner–Nordström–(anti)-de Sitter black hole with global monopole, Class. Quant. Grav. 33, no.21, 215012 (2016) doi:10.1088/0264-9381/33/21/215012 [arXiv:1610.01129 [gr-qc]].
  * [9] H. Zheng, Y. Chen and J. Tang, Holographic thermodynamics of a charged AdS black holes with global monopole, [arXiv:2405.12227 [hep-th]].
  * [10] M. Barriola and A. Vilenkin, Gravitational Field of a Global Monopole, Phys. Rev. Lett. 63, 341 (1989) doi:10.1103/PhysRevLett.63.341
  * [11] M. Cvetic, G. W. Gibbons, D. Kubiznak and C. N. Pope, Black Hole Enthalpy and an Entropy Inequality for the Thermodynamic Volume, Phys. Rev. D 84, 024037 (2011) doi:10.1103/PhysRevD.84.024037 [arXiv:1012.2888 [hep-th]].
  * [12] A. Chamblin, R. Emparan, C. Johnson, and R. Myers, Charged AdS black holes and catastrophic holography, Phys.Rev. D 60, (1999) 064018, [hep-th/9902170].
  * [13] E. Spallucci and A. Smailagic, Maxwell’s equal area law for charged Anti-deSitter black holes, Phys. Lett. B 723,436, (2013) [arXiv:1305.3379 [hep-th]]
  * [14] D. Stauffer, M. Ferer, and M. Wortis, Universality of Second-Order Phase Transitions: The Scale Factor for the Correlation Length, Phys. Rev. Lett. 29, 345.
  * [15] B. Widom, Equation of state in the neighborhood of the critical point, J. Chem, Phys. 43, 3898.
  * [16] H. E. Stanley, Introduction to Phase Transitions and Critical Phenomena, Oxford University, New York.