JuliaCon 2025

Calculating electronic properties of materials using the density matrix formalism enables
linear scaling with system size, compared to the cubic scaling of direct diagonalization. To
obtain density matrices, one typically applies the Fermi–Dirac operator to the Hamiltonian
of a given system. The second-order spectral projection (SP2) scheme has proven effective
across a broad range of temperatures. We generalize this SP2 method, drawing inspiration
from machine learning models, to recursively expand the Fermi–Dirac operator. This approach
accurately captures properties such as fractional occupation numbers and electronic entropy
while significantly reducing computational cost by replacing highly nonlinear matrix
functions with iterative matrix multiplications. Performance further improves by
running these calculations on GPUs.

In this presentation, we show how Julia enables rapid prototyping of this algorithm and the
subsequent design of a user-friendly API. We share strategies for optimizing performance
(e.g., removing type instabilities, reducing allocations, and avoiding code duplication)
with the help of robust benchmarking tools from the Julia community. We have been pleasantly
surprised to see that writing GPU-compatible Julia code adds minimal overhead compared to
calling CUDA C++ directly, while substantially improving readability—particularly when
leveraging tensor cores. We also take advantage of automatic differentiation packages to
generalize our approach to differentiating physically meaningful functions. Test cases for
this model use a rich ecosystem of random distribution generation code, showcasing the
flexibility and power of Julia for computational chemistry applications.