Examples

Hard rod fluid

Hard rods in 1D can be solved exactly (up to numerical accuracy) using the PercusFunctional.

using Plots
using ClassicalDFT

L = 10
dx = 0.01
geom = CartesianGeometry{1}(L, dx; boundary=:walls)
functional = PercusFunctional(geom)
T, μ = 1.0, 2.0
ρ = solve(functional; T, μ).ρ
plot(geom, ρ, label="ρ")
[ Info: starting iteration
[ Info: converged after 486 iterations (‖ΔρEL‖ = 9.870572026393631e-9 < 1.0e-8)

Hard sphere fluid

For hard spheres in 3D, FMT functionals such as the RosenfeldFunctional can be employed.

Planar hard-wall confinement

A flat-wall setup in PlanarGeometry is considered in the following, such that quantities are resolved only along the $z$ coordinate and remain translationally invariant in the $x$ and $y$ directions.

using Plots
using ClassicalDFT

L = 10
dz = 0.01
geom = PlanarGeometry(L, dz; boundary=:walls)
functional = RosenfeldFunctional(geom)
T, μ = 1.0, 3.0
ρ = solve(functional; T, μ).ρ
plot(geom, ρ, label="ρ")
[ Info: starting iteration
[ Info: converged after 388 iterations (‖ΔρEL‖ = 9.988906479740933e-9 < 1.0e-8)

Sedimentation

In the following, a linearly decreasing external potential is imposed in order to model gravitational sedimentation.

using Plots
using ClassicalDFT

L = 100
dz = 0.01
geom = PlanarGeometry(L, dz; boundary=:walls)
functional = RosenfeldFunctional(geom)
T, μ = 1.0, 3.0
Vext(z) = 0.05 * z
ρ = solve(functional; T, μ, Vext).ρ
plot(geom, ρ, label="ρ")
[ Info: starting iteration
[ Info: converged after 388 iterations (‖ΔρEL‖ = 9.754028473807352e-9 < 1.0e-8)

Lennard-Jones fluid

The (truncated) Lennard-Jones fluid is addressed with a hard-sphere functional plus the LennardJonesMeanFieldFunctional, which treats the attractive part of the interparticle pair potential.

Planar hard-wall confinement

In the following, we utilize the AndersonSolver, which usually converges much more rapidly compared to the standard PicardSolver. The progress of the iterations is monitored with PrintCallback.

using Plots
using ClassicalDFT

L = 10
dz = 0.01
geom = PlanarGeometry(L, dz; boundary=:walls)
functional = CompositeFunctional(geom, RosenfeldFunctional, LennardJonesMeanFieldFunctional)
T, μ = 1.5, 1.0
solver = AndersonSolver(functional)
callbacks = [PrintCallback(["iteration", "‖ΔρEL‖"]; every=10)]
ρ = solve(functional; T, μ, solver, callbacks).ρ
plot(geom, ρ, label="ρ")
[ Info: starting iteration
iteration      ‖ΔρEL‖
0              1.94773
10             0.0750246
20             0.00287115
30             1.02017e-05
40             1.06138e-07
[ Info: converged after 47 iterations (‖ΔρEL‖ = 8.865916667488705e-9 < 1.0e-8)

Liquid-gas coexistence

Instead of the chemical potential, it is also possible to fix the mean density ρ̄. This is used in the following to predict liquid-gas coexistence, whereby the iteration is started with a step-like initial density profile ρ_init.

using Plots
using ClassicalDFT

L = 50
dx = 0.01
geom = PlanarGeometry(L, dx; boundary=:periodic)
T = 0.9
ρ̄ = 0.3
ρ_init(x) = x < L/2 ? 1.0 : 0.0
functional = CompositeFunctional(geom, RosenfeldFunctional, LennardJonesMeanFieldFunctional)
ρ = solve(functional; T, ρ̄, ρ_init).ρ
plot(geom, ρ, label="ρ")
[ Info: starting iteration
[ Info: converged after 972 iterations (‖ΔρEL‖ = 9.936846456781723e-9 < 1.0e-8)

Spherical geometries

SphericalGeometry and SphericalShellGeometry describe spherically symmetric setups, with profiles being resolved along the radial coordinate $r$. The difference between both geometries is how the radial $r$ axis is constrained: In SphericalGeometry, the $r$ axis starts at the origin and reaches up to a radius $R$, while in SphericalShellGeometry, only a shell between an inner radius $R_\mathrm{inner}$ and an outer radius $R_\mathrm{outer}$ is considered. Of course, for equivalent setups, both geometries yield equivalent results, as demonstrated in the following:

using Plots
using ClassicalDFT

R_inner = 5
R_outer = 15
dr = 0.01
geom = SphericalShellGeometry(R_inner, R_outer, dr; boundary=:walls)
functional = RosenfeldFunctional(geom)
T, μ = 1.0, 3.0
ρ = solve(functional; T, μ).ρ
plot(geom, ρ, label="ρ")
[ Info: starting iteration
[ Info: converged after 394 iterations (‖ΔρEL‖ = 9.922726085243028e-9 < 1.0e-8)

using Plots
using ClassicalDFT

R = 15
dr = 0.01
geom = SphericalGeometry(R, dr; boundary=:wall)
functional = RosenfeldFunctional(geom)
T, μ = 1.0, 3.0
Vext(r) = r < 5 ? Inf : 0.0
ρ = solve(functional; T, μ, Vext).ρ
plot(geom, ρ, label="ρ")
[ Info: starting iteration
[ Info: converged after 394 iterations (‖ΔρEL‖ = 9.900106456850466e-9 < 1.0e-8)

Test particle

The spherical geometries are suitable for test particle calculations, where one sets the external potential to the interparticle pair potential, $V_\mathrm{ext}(r) = \phi(r)$. The idea goes back to Percus and amounts to fixing a particle at the origin. The resulting density profile $\rho(r)$ is, up to normalization by the bulk density, the radial distribution function $g(r)$.

using Plots
using ClassicalDFT

R = 20
dr = 0.01
geom = SphericalGeometry(R, dr; boundary=:bulk)
functional = RosenfeldFunctional(geom)
T, μ = 1.0, 3.0
Vext(r) = r < 1.0 ? Inf : 0.0
ρ = solve(functional; T, μ, Vext).ρ
plot(geom, ρ, label="ρ")
[ Info: starting iteration
[ Info: converged after 363 iterations (‖ΔρEL‖ = 9.753315932670148e-9 < 1.0e-8)