Jellium Problem in Electron Gas Theory

• The jellium problem is defined as the study of charged particles immersed in a uniform neutralizing background, modeling systems like electron gases and plasmas.
• It provides a fundamental framework for analyzing energy asymptotics, mean-field limits, and exchange–correlation effects in condensed matter physics.
• The model enables rigorous mathematical and numerical investigations into screening, phase transitions, and renormalized energies across different dimensions.

The jellium problem is the study of interacting charged particles embedded in a uniform neutralizing background, with the background introduced to enforce global neutrality and to regularize long-range Coulomb energies. In condensed-matter usage, it is the paradigmatic model of the uniform electron gas; in statistical mechanics it is the one-component plasma; and in two dimensions it also appears as the logarithmic Coulomb gas that underlies several exact and asymptotic results. Across these settings, the problem serves as a reference point for simple metals, warm dense matter, density-functional theory, Coulomb-gas asymptotics, and the plasma analogy for fractional quantum Hall states (Fantoni, 2021, Rougerie, 2022, Dufty, 2017).

1. Canonical formulations

In three dimensions, jellium is most commonly written as a gas of electrons of charge $e$ moving in a compensating positive background. A compact Hamiltonian is

$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$

where $V_{\mathrm{bg}}$ collects the electron–background and background–background terms needed for neutrality and a finite thermodynamic limit. In the units used for finite-temperature simulations, lengths are scaled by $a=(4\pi n/3)^{-1/3}$, energies by the Rydberg, and $r_s=a/a_B$ parametrizes density; small $r_s$ is the high-density regime and large $r_s$ the strongly correlated regime (Fantoni, 2021).

In two dimensions, the Coulomb kernel is logarithmic,

$G(z)=-\log|z|,\qquad -\Delta G = 2\pi \delta_0,$

and the classical one-component plasma on a bounded domain $D\subset \mathbb{R}^2$ with background density $\rho_b$ has energy

$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$0

Neutrality is imposed by $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$1. Hard-wall confinement or a smooth confining potential growing at infinity are standard choices (Rougerie, 2022).

A broader formulation replaces the Coulomb law by a Riesz kernel $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$2, $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$3. In dimension $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$4, Coulomb corresponds to $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$5, while the logarithmic case is $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$6 in $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$7. In this language, the neutralized jellium energy in a domain $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$8, $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$9, is

$V_{\mathrm{bg}}$0

This Riesz form is the setting for several sharp asymptotic results and for the comparison between jellium and the uniform electron gas in density-functional theory (Cotar et al., 2017).

2. Thermodynamic limit, mean-field theory, and density functionals

At finite temperature, the classical $V_{\mathrm{bg}}$1D jellium Gibbs measure is

$V_{\mathrm{bg}}$2

Under the mean-field scaling $V_{\mathrm{bg}}$3, with $V_{\mathrm{bg}}$4 smooth and confining, the leading equilibrium is governed by

$V_{\mathrm{bg}}$5

minimized over probability measures $V_{\mathrm{bg}}$6. If $V_{\mathrm{bg}}$7 is the unique minimizer, then the $V_{\mathrm{bg}}$8-point equilibrium densities converge weakly after the natural rescaling $V_{\mathrm{bg}}$9, and the free energy has the local-density expansion

$a=(4\pi n/3)^{-1/3}$0

where $a=(4\pi n/3)^{-1/3}$1 is defined by the homogeneous thermodynamic limit of jellium (Rougerie, 2022).

For homogeneous jellium, Lieb–Narnhofer-type thermodynamic limits are central. In $a=(4\pi n/3)^{-1/3}$2D, if $a=(4\pi n/3)^{-1/3}$3 is a box of side $a=(4\pi n/3)^{-1/3}$4, $a=(4\pi n/3)^{-1/3}$5, and $a=(4\pi n/3)^{-1/3}$6 is the background-generated confinement with $a=(4\pi n/3)^{-1/3}$7 outside $a=(4\pi n/3)^{-1/3}$8, then

$a=(4\pi n/3)^{-1/3}$9

exists and is independent of the shape of $r_s=a/a_B$0 in $r_s=a/a_B$1D and, more generally, in $r_s=a/a_B$2. This homogeneous free-energy density is the local input of the local-density approximation (Rougerie, 2022).

In density-functional theory, the exact relationship between the inhomogeneous electron gas and inhomogeneous jellium is unusually sharp. If $r_s=a/a_B$3 and $r_s=a/a_B$4 are the Mermin intrinsic free-energy functionals of the electron gas and jellium, then

$r_s=a/a_B$5

with $r_s=a/a_B$6 the electron–background potential and $r_s=a/a_B$7 the background self-energy. Since the Hartree sector exactly absorbs the background dependence, the exchange–correlation functionals coincide: $r_s=a/a_B$8 This identity justifies the use of jellium-based exchange–correlation information in local and gradient approximations for real inhomogeneous electronic systems (Dufty, 2017).

3. Next-order asymptotics, renormalized energy, and equivalent definitions

A major modern formulation of the jellium problem concerns the next-order term beyond the mean-field $r_s=a/a_B$9 energy. For Riesz interactions $r_s$0, the minimal confined energy has the asymptotic form

$r_s$1

where $r_s$2 is defined via a renormalized infinite-volume energy $r_s$3. In the Caffarelli–Silvestre extended formulation,

$r_s$4

and

$r_s$5

For the many-marginals optimal-transport formulation of the uniform electron gas, one similarly gets

$r_s$6

The main theorem is

$r_s$7

for $r_s$8 and $r_s$9, and for $r_s$0 and $r_s$1. In the Coulomb case $r_s$2, this verifies the long-standing conjecture that the second-order constants of jellium and the uniform electron gas coincide (Cotar et al., 2017).

At the microscopic scale, the jellium energy equidistributes. After blow-up at the interparticle spacing $r_s$3, the local renormalized energy density in a microscopic cube converges to the minimum jellium energy corresponding to the local macroscopic density. In the Coulomb case this leads to surface-order discrepancy bounds for point counts,

$r_s$4

in the regime described in the theorem. This is a rigorous form of microscopic rigidity that does not amount to a crystallization proof but does identify the local energy density with the jellium minimum (Petrache et al., 2016).

A complementary question is whether apparently different definitions of the jellium ground-state energy give the same thermodynamic constant. In $r_s$5, one has

$r_s$6

where $r_s$7 is the classical background-neutralized energy density, $r_s$8 is the periodic Coulomb or Ewald-type definition, and $r_s$9 is the indirect-energy definition with prescribed constant one-body density. In $G(z)=-\log|z|,\qquad -\Delta G = 2\pi \delta_0,$0D, a modified floating crystal trial state with exact characteristic-function density extends this equivalence and connects $G(z)=-\log|z|,\qquad -\Delta G = 2\pi \delta_0,$1 to Serfaty’s renormalized energy and to the order-$G(z)=-\log|z|,\qquad -\Delta G = 2\pi \delta_0,$2 term of logarithmic energy on the sphere (Lauritsen, 2021). In $G(z)=-\log|z|,\qquad -\Delta G = 2\pi \delta_0,$3D Coulomb jellium, a related floating Wigner crystal construction suppresses the boundary charge fluctuations that otherwise lead to a macroscopic energy increase, and proves the coincidence of the background, constant-density, and periodic definitions (Lewin et al., 2019).

4. Plasma analogy, screening, and incompressibility in two dimensions

In fractional quantum Hall theory, Laughlin’s wavefunction at filling $G(z)=-\log|z|,\qquad -\Delta G = 2\pi \delta_0,$4,

$G(z)=-\log|z|,\qquad -\Delta G = 2\pi \delta_0,$5

has squared modulus

$G(z)=-\log|z|,\qquad -\Delta G = 2\pi \delta_0,$6

with $G(z)=-\log|z|,\qquad -\Delta G = 2\pi \delta_0,$7. Thus Laughlin’s state is the Gibbs measure of a $G(z)=-\log|z|,\qquad -\Delta G = 2\pi \delta_0,$8D one-component plasma with quadratic confinement. Analytic quasi-hole factors $G(z)=-\log|z|,\qquad -\Delta G = 2\pi \delta_0,$9 insert additional positive charges, because the terms $D\subset \mathbb{R}^2$0 are superharmonic in each variable (Rougerie, 2022).

The central structural result is an incompressibility estimate. In physical units,

$D\subset \mathbb{R}^2$1

In the paper’s rescaled units, for minimizing configurations of

$D\subset \mathbb{R}^2$2

with $D\subset \mathbb{R}^2$3 superharmonic in each variable, every disk $D\subset \mathbb{R}^2$4 satisfies

$D\subset \mathbb{R}^2$5

A mesoscopic version for Laughlin-like states states that for disks of radius $D\subset \mathbb{R}^2$6, $D\subset \mathbb{R}^2$7,

$D\subset \mathbb{R}^2$8

and the probability of violating this bound has an exponential large-deviation estimate.

The mechanism is a screening construction. For any finite set of point charges $D\subset \mathbb{R}^2$9, one builds a screening region $\rho_b$0 of constant negative charge density such that

$\rho_b$1

satisfies

$\rho_b$2

This is equivalent to a Thomas–Fermi-like obstacle problem with constraint $\rho_b$3, $\rho_b$4. The resulting exclusion rule yields the density upper bound, and it feeds directly into a stability theorem for the Laughlin phase: for smooth external potential $\rho_b$5 and weak smooth interaction $\rho_b$6, there exists $\rho_b$7 such that

$\rho_b$8

meaning that among all holomorphic perturbations of Laughlin’s state, independent quasi-hole states of product form are asymptotically optimal (Rougerie, 2022).

5. Broken symmetries, Hartree–Fock pathologies, and screened approximations

Within Hartree–Fock, the uniform Fermi gas is not the full story. In three-dimensional jellium, broken-symmetry states are energetically favored at any density against the homogeneous Fermi gas with isotropic Fermi surface. At high density, Hartree–Fock yields metallic spin-unpolarized incommensurate crystals in which the charge and spin densities form an incommensurate crystal with more crystal sites than electrons; as $\rho_b$9, these states approach pure spin-density waves with modulation wavevector $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$00. At lower density, the commensurate Wigner crystal is favored for $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$01, the system undergoes several structural phase transitions, and the polarization transition occurs around $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$02 (Baguet et al., 2013).

A distinct Hartree–Fock issue is the unscreened single-particle dispersion. In $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$03D jellium, the Hartree–Fock spectrum

$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$04

with

$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$05

has a logarithmically divergent derivative at $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$06. Consequently, the density of states

$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$07

vanishes at the Fermi energy. Slater’s hyper-Hartree–Fock equations do not remove this anomaly; instead they shift it from $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$08 to a tunable boundary $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$09 between optimized and unoptimized orbitals. This supports the interpretation that the pathology is an artifact of the nonlocal exchange operator rather than exclusively a failure of screening (Blair et al., 2015).

Self-consistent $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$10 addresses different aspects of the same problem. In $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$11D jellium, the fully self-consistent approximation is built from

$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$12

with dielectric function

$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$13

Its most serious deficiency is an incorrect dielectric response, traceable to violation of particle-number conservation in $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$14. Enforcing

$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$15

restores the physically required $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$16 behavior. At $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$17, $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$18, and $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$19, the corrected dielectric function exhibits a plasmon zero close to $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$20, with

$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$21

and the same calculations yield benchmark $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$22 values for $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$23, the quasiparticle residue $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$24, and $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$25 for $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$26 (Houcke et al., 2016).

A recent extension shows how quantum geometry modifies the classical phase logic. In $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$27D $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$28-jellium, the lower band has Berry curvature

$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$29

while the dispersion remains quadratic. Self-consistent Hartree–Fock then yields two Fermi liquids, two Wigner crystals, and an anomalous Hall crystal; the anomalous Hall crystal occupies a large region of the $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$30 phase diagram, and a continuous transition separates it from a halo Wigner crystal (Soejima et al., 17 Mar 2025).

6. Numerical methods, finite temperature, and open directions

At finite temperature, fermionic path-integral Monte Carlo for jellium is obstructed by the sign problem. A restricted-path, fixed-node formulation combined with a canonical worm algorithm gives a practical route through this obstruction. In the implementation for fully polarized $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$31 electrons with periodic boundary conditions, long-range Coulomb interactions are handled through the Fraser correction,

$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$32

which is equivalent in spirit to Ewald summation. The restricted path uses free-fermion trial nodes

$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$33

and two canonical variants were studied: Algorithm A, which under-samples exchange at low temperature, and Algorithm B, which restores exchange sampling by a tailored $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$34-sector. For $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$35 and $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$36, the method reproduces canonical restricted PIMC energies and pair correlations with good agreement at high density, and substantially improves low-$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$37 results when permutations are restored (Fantoni, 2021).

The same model continues to motivate questions outside the standard metallic-fluid regime. In the dielectric Kirzhnits–Maksimov–Khomskii formalism for the electron–proton jellium model, the screened interaction

$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$38

leads to a linearized integral equation for $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$39. Solving that equation directly gives a dome-shaped $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$40 with maximum

$$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$41

within Thomas–Fermi–RPA screening, while static Hubbard local-field corrections raise $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$42 only several-fold and still keep it below $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$43 K across the explored densities. This supports the conclusion that superconductivity in the jellium model remains weak-coupling over the accessible density range (Sadovskii, 28 Jun 2026).

Several mathematical and physical questions remain open. In the $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$44D Laughlin/jellium regime, the universal incompressibility bounds are presently established on disks of radius $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$45 with $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$46, and extension to all mesoscopic scales $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$47 is identified as a natural open problem; stronger or singular interactions, realistic short-scale disorder, and generalizations beyond independent quasi-holes are likewise unresolved (Rougerie, 2022). In the Riesz/Coulomb asymptotic theory, the logarithmic $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$48 case remains exceptional, and the crystallization hypotheses behind the conjectured optimal constants—triangular in $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$49D logarithmic systems, bcc or fcc in $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$50D Riesz systems depending on $$H = \sum_{i=1}^N \frac{p_i^2}{2m} + \frac{1}{2} \sum_{i\neq j} \frac{e^2}{|r_i-r_j|} + V_{\mathrm{bg}},$$51—remain open (Cotar et al., 2017). These unresolved points help explain why the jellium problem remains a live interface between Coulomb-gas analysis, many-body theory, and computational condensed matter.