Approaching the thermodynamic limit of a bounded one-component plasma

TL;DR

This work establishes the thermodynamic limit for a bounded one-component plasma (BOCP) by MD, avoiding periodic boundaries and Ewald summation. It introduces convergent excess energies and the ionic compressibility factor, derives a wide-range ionic EOS, and demonstrates a consistent melting point near $Γ_m \approx 174$ with smooth extrapolation to $N \to \infty$ via $u(Γ,N) = u_\infty(Γ) + c_1(Γ)/(\ln N + c_2(Γ))$. The authors also connect BOCP results with conventional OCP simulations by calibrating the LAMMPS cutoff $r_c(Γ)$, showing its impact on interfacial properties and metastability, and providing practical guidance for accurate long-range Coulomb simulations. Overall, the paper delivers benchmark energies, a validated EOS, and a methodology to reconcile bounded-cluster and periodic-boundary approaches, with implications for transport, phase behavior, and simulations of strongly coupled plasmas.

Abstract

The classical one-component plasma (OCP) bounded by a spherical surface reflecting ions (BOCP) is studied using molecular dynamics (MD). Simulations performed for a series of sufficiently large BOCP's make it possible to establish the size dependencies for the investigated quantities and extrapolate them to the thermodynamic limit. In particular, the total electrostatic energy per ion is estimated in the limit of infinite BOCP in a wide range of the Coulomb coupling parameter $Γ$ from 0.03 to 1000 with the relative error of the order 0.1%. Calculated energies are by about 0.5% lower as compared to the modern Monte Carlo (MC) simulation data obtained by different authors at $Γ<30$ and almost coincide with the MC results at $Γ>175$. We introduce two more converging characteristic energies, the excess interatomic electrostatic energy and the excess ion-background electrostatic energy, which enable us to calculate the ionic compressibility factor inaccessible in conventional MC and MD simulation of the OCP with periodic boundary conditions. The derived wide-range ionic equation of state can be recommended for testing OCP simulations with various effective interaction potentials. Based on this equation, we propose an improved cutoff radius for the interionic forces implemented in LAMMPS and perform MD simulation of the OCP to demonstrate that location of the metastable region of the fluid-solid phase transition depends sensitively on this radius.

Figures (12)

• Figure 1: Ion density distribution $f_c(r)$ (\ref{['seq:refEquation7']}) and the radial distribution function $f_p(r)$ (\ref{['seq:refEquation6']}) of BOCP at different $\Gamma$: (a) and (b) $\Gamma=0.1$, 30, and 200 (red, green, and blue lines, respectively), $N=5000$; (c) $f_c(r)$ (red line) and $f_p(r)$ (blue line) at $\Gamma=1000$ and $N=14000$.
• Figure 2: Size dependence of the BOCP total electrostatic energy from MD simulation (dots) and its fit by Eq. (\ref{['seq:refEquation18']}) (lines): (a) diamonds and solid lines correspond to $\Gamma=0.1$ and circles and dashed line, to $\Gamma=30$; (b) circles indicate MD data for $\Gamma=3$, solid line is their fit, and dashed line marks the thermodynamic limit $u_\infty(3)$.
• Figure 3: Thermodynamic limit of the BOCP electrostatic energy compared with the results of MC simulations for the OCP with PBC. MD results of this work (circles) and their approximation by the functions (\ref{['seq:refEquation37']}) and (\ref{['seq:refEquation38']}) (solid line), MC results by Brush et al.Brush_1966 (squares), Hansen Hansen_1973 (diamonds), Caillol and Gilles Caillol_2010 (stars), Demyanov and Levashov Demyanov_2022 (triangles). Approximation of the results available in studies by Slattery et al.Slattery_1982 and Caillol (MC within hyperspherical boundary conditions) Caillol_1999_3 are shown by dotted and dashed-dotted lines, respectively. Dashed-dotted-dotted line indicates the Debye--Hückel approximation and dashed line, the Madelung constant for bcc lattice. Inset shows a magnified fragment of this figure.
• Figure 4: Thermal fraction of the BOCP electrostatic energy as a function of $\Gamma$ for the BOCP size $N=5000$ (solid line), 10000 (dashed line), 30000 (dashed-dotted line), and 50000 (dotted line). Inset shows the size-dependent Coulomb coupling parameter $\Gamma_{\mathrm{m}}(N)$ (open circles); solid line indicates the curve fit and dashed line, the thermodynamic limit $\Gamma_{\mathrm{m}}(\infty)$.
• Figure 5: Thermal fraction of the electrostatic energy in the thermodynamic limit as a function of $\Gamma$. This work: MD results (dots) and their approximation by the function (\ref{['seq:refEquation37']}) (dashed line) and (\ref{['seq:refEquation38']}) (dashed-dotted line). Solid line indicates calculations based on the harmonic lattice approximation Chugunov_2005 and dotted line, MC results Slattery_1982.