Localized density matrix minimization and linear scaling algorithms

We propose a convex variational approach to compute localized density matrices for both zero temperature and finite temperature cases, by adding an entry-wise $\ell_1$ regularization to the free energy of the quantum system. Based on the fact that the density matrix decays exponential away from the diagonal for insulating system or system at finite temperature, the proposed $\ell_1$ regularized variational method provides a nice way to approximate the original quantum system. We provide theoretical analysis of the approximation behavior and also design convergence guaranteed numerical algorithms based on Bregman iteration. More importantly, the $\ell_1$ regularized system naturally leads to localized density matrices with banded structure, which enables us to develop approximating algorithms to find the localized density matrices with computation cost linearly dependent on the problem size.

💡 Research Summary

The paper introduces a convex variational framework for computing localized density matrices (LDMs) in electronic structure calculations, applicable to both zero‑temperature (ground‑state) and finite‑temperature regimes. The key idea is to augment the standard free‑energy functional with an entry‑wise ℓ₁ penalty term. For the zero‑temperature case the functional reduces to the trace of the Hamiltonian times the density matrix, while for finite temperature the authors use the Mermin free energy containing the Fermi‑Dirac entropy ϕ(P)=P ln P+(1−P) ln (1−P). The constraints tr P = N, P = Pᵀ, and 0 ≼ P ≼ I remain convex, so the addition of the ℓ₁ term preserves overall convexity.

Theoretical analysis (Theorems 1 and 2) shows that as the regularization parameter η grows, the ℓ₁‑regularized solution P_{β,η} converges to the exact density matrix P_{β,∞} in Frobenius norm at a rate O(1/η). Moreover, the energy gap between the regularized and exact solutions is bounded by (1/η)‖P_{β,∞}‖₁. These results rely on the physical fact that for insulators or systems at finite temperature the density matrix decays exponentially away from the diagonal, which justifies the sparsity induced by the ℓ₁ term.

To solve the resulting optimization problems, the authors adopt Bregman iteration (equivalently ADMM). By introducing auxiliary variables Q and R, the problem is split into three sub‑problems per iteration: (i) a P‑update that involves the Hamiltonian and the entropy term, (ii) a Q‑update that handles the ℓ₁ penalty via soft‑thresholding (shrinkage) while enforcing the trace constraint, and (iii) an R‑update that projects onto the feasible set 0 ≼ R ≼ I via eigenvalue clipping. For the finite‑temperature case the P‑update leads to a nonlinear fixed‑point equation; the authors propose an inner iteration Z_{l+1}=1/(1+exp(βY_l)) with Y_l containing the current estimate and the Bregman multipliers. Under the condition β(λ+r)<4, this inner loop converges exponentially.

A major contribution is the development of linear‑scaling algorithms. By exploiting the exponential decay of the true density matrix, the authors restrict the search space to banded matrices with bandwidth b. All Bregman sub‑steps can be performed using only the banded entries, yielding O(nb) computational cost and memory, where n is the matrix dimension. Unlike traditional density‑matrix minimization (DMM) or purification methods, the banded formulation remains convex, guaranteeing global optimality and avoiding spurious local minima.

Numerical experiments on one‑dimensional chains and two‑dimensional lattices validate the theory. Varying η and the bandwidth demonstrates that high sparsity (up to 99.9 % zero entries) can be achieved while keeping the energy error below 10⁻⁶. Timing results show near‑linear growth of CPU time with system size up to 10⁵ degrees of freedom, outperforming standard purification schemes which either lose convergence or incur higher errors.

In conclusion, the paper delivers a mathematically rigorous, convex, ℓ₁‑regularized variational model for localized density matrices, an efficient Bregman‑based solver with provable convergence, and a practical linear‑scaling implementation. This combination addresses the long‑standing challenge of achieving both sparsity and guaranteed convergence in large‑scale electronic‑structure calculations, and opens avenues for extensions to more complex systems such as disordered materials, spin‑polarized calculations, and time‑dependent density functional theory.