Nuclear–electronic orbital method at scale: A two-stage hybrid-layout distributed Cholesky decomposition approach

Cholesky decomposition (CD) is widely used to reduce the cost of storing two-particle integrals and has been extended to multicomponent nuclear–electronic orbital (NEO) methods. Although NEO-CD lowers memory demands, single-node implementations can still face memory bottlenecks due to the size of the mixed asymmetric nuclear–electronic tensors. In this study, we present a two-stage hybrid-layout distributed NEO-CD framework suitable for large NEO calculations at scale. Cholesky vectors are first generated and orthonormalized using an atomic orbital (AO)-index distribution. During Fock matrix construction, mixed nuclear–electronic vectors retain the AO-index layout, while pure electronic and pure nuclear vectors are redistributed over the auxiliary index. This design avoids costly inter-rank communication and provides efficient construction of both Coulomb and exchange matrices. Benchmark calculations on water clusters with all protons treated quantum mechanically show excellent scaling and load balance. Using only 104 MPI ranks, a full NEO-Fock build for (H2O)250 with 500 quantum protons, the largest NEO calculation to date, completes in 96.2 s per Fock matrix construction. The distributed NEO-CD framework enables NEO-DFT calculations for liquid-like water clusters up to (H2O)100. Extrapolation of the binding energies yields a bulk-limit value of ΔEn→∞ = −10.35 ± 0.19 kcal/mol, in excellent agreement with the experimental vaporization enthalpy of water. These results demonstrate that fully quantum mechanical simulations of large molecular systems are now practical on modern computing architectures.