Ab Initio Many Body Quantum Embedding and Local Correlation in Crystalline Materials using Interpolative Separable Density Fitting

This paper presents an efficient, linear-scaling implementation of ab initio many-body quantum embedding and local correlation methods for infinite periodic systems using interpolative separable density fitting, enabling accurate thermodynamic limit estimates of coupled cluster ground-state energies for both weakly and strongly correlated solids.

The Gist

Imagine you are trying to understand how a massive, infinite city behaves by studying just a single house. In the world of quantum chemistry, this "city" is a crystal (like diamond or a metal oxide), and the "house" is a tiny repeating unit called a unit cell. Scientists want to know the exact energy of this infinite city to predict its properties, but calculating the interactions between every single electron in an infinite grid is like trying to count every grain of sand on every beach on Earth simultaneously. It's computationally impossible with traditional methods because the work grows too fast as the grid gets bigger.

This paper introduces a clever new "shortcut" called Interpolative Separable Density Fitting (ISDF), combined with a technique called FFTISDF, to solve this problem. Here is how it works, using simple analogies:

The Problem: The "Too Many Neighbors" Issue

In a crystal, electrons don't just interact with their immediate neighbors; they feel the pull of electrons far away. To get an accurate answer, you need to sample the crystal with a grid of points (called k-points).

• The Old Way: Imagine trying to calculate the noise level in a stadium by asking every single person what they are hearing from every other person. If you double the number of people (k-points), the number of conversations you have to track explodes. This is why previous methods hit a wall when trying to simulate large, infinite crystals.
• The Goal: The authors wanted to calculate the energy of these infinite crystals using up to 1,000 k-points (a very dense grid) to get a result that represents the "thermodynamic limit" (the true, infinite size of the material).

The Solution: The "Smart Summarizer"

The authors developed a method that acts like a smart summarizer or a translator.

1. The Interpolation Points (The "Key Witnesses"):
   Instead of asking every single electron in the crystal about its interactions, the method picks a small, strategic set of "witnesses" (called interpolation points). Think of these as key reporters in a newsroom. Instead of interviewing every single citizen in a city to understand the mood, you interview a carefully selected group of 100 people who can accurately represent the feelings of the whole city.

   • The paper shows that by using these "witnesses," they can reconstruct the behavior of the entire electron cloud with high accuracy, but with a fraction of the work.

2. The Linear Scaling (The "Magic Elevator"):
   In the old methods, if you doubled the size of your simulation (more k-points), the time it took to run the calculation would quadruple or even increase much faster (like climbing a steep, endless hill).

   • With this new method, the time it takes to calculate grows linearly. If you double the number of k-points, it only takes double the time. It's like having a magic elevator that lets you go up the mountain without getting tired, no matter how high you go. This allowed them to run simulations with up to 1,000 k-points, which was previously impossible.

The Tools: "Embedding" and "Local Correlation"

To get the most accurate energy numbers, the paper uses two specific strategies:

• Density Matrix Embedding (The "Focus Group"): This is like taking a small, representative group of people (a fragment of the crystal) and studying them in deep detail, while treating the rest of the city as a simplified background. This allows for a very precise look at the "local" interactions.
• Local Natural Orbital Correlation (The "Efficient Sorting"): This method sorts the electrons so that only the ones that really matter for a specific interaction are calculated in high detail, ignoring the ones that are too far away to matter.

What They Tested

The team tested this new "smart summarizer" on four different types of materials:

1. Diamond: A hard, wide-gap semiconductor.
2. Carbon Dioxide (CO2): A molecular crystal (like dry ice).
3. Nickel Oxide (NiO): A material where electrons are "strongly correlated" (they act like a chaotic crowd rather than independent individuals).
4. CaCuO2: A cuprate superconductor with a layered structure.

The Results

• Accuracy: They showed that their method could predict the energy of these materials with extreme precision, matching the results of much slower, older methods but doing it in a fraction of the time.
• The "Thermodynamic Limit": By using up to 1,000 k-points and then mathematically "extrapolating" (predicting the trend to infinity), they were able to give the most accurate estimates yet for the ground-state energy of these infinite crystals.
• Magnetic Properties: For Nickel Oxide and CaCuO2, they calculated how the atoms interact magnetically (specifically the "exchange constants"). Their results were much closer to real-world experimental values than previous calculations, proving that including these "strong correlations" is vital for understanding these materials.

The Bottom Line

This paper presents a new computational engine that makes it possible to simulate infinite crystals with the same level of detail previously reserved for tiny molecules. By using a "smart summarizer" (ISDF) to reduce the complexity of electron interactions, they turned a task that was computationally impossible into one that is now efficient and scalable. This allows scientists to finally get reliable answers about the true, infinite nature of solid materials without needing a supercomputer the size of a planet.

Technical Summary

Technical Summary: Ab Initio Many Body Quantum Embedding and Local Correlation in Crystalline Materials using Interpolative Separable Density Fitting

Problem Statement
Correlated wavefunction methods offer a pathway to controllable accuracy in computational materials science through systematically improvable many-body expansions. However, their application to infinite periodic systems is hindered by steep computational scaling with system size, creating a significant obstacle to obtaining thermodynamic limit (TDL) results. While existing strategies utilize translational symmetry (Bloch orbitals) to reduce costs by factors of $N_k^{-1}$ to $N_k^{-2}$, and employ locality or low-rank approximations (such as Cholesky decomposition or resolution-of-the-identity) to manage correlation, a unified approach that achieves linear scaling with the number of k-points ($N_k$) for high-level many-body methods remains elusive. Specifically, standard density fitting techniques often suffer from quadratic scaling or prohibitive storage requirements when applied to large k-point meshes necessary for TDL convergence.

Methodology
The authors present an efficient implementation of ab initio many-body quantum embedding and local correlation methods for infinite periodic systems. The core of the methodology is the integration of Translational Symmetry Adapted Interpolative Separable Density Fitting (FFTISDF) within the Gaussian-Plane-Wave framework (implemented in PySCF).

Key technical components include:

1. FFTISDF Representation: The Coulomb interaction is represented via a tensor hypercontraction (THC) form. Products of basis functions on a real-space grid are approximated using values from a select set of interpolating points (IPs) determined by a pivoted Cholesky decomposition. This allows the two-electron integrals to be expressed in a factorized form:
   $$(\mu k_\mu \nu k_\nu | \lambda k_\lambda \sigma k_\sigma) = \sum_{IJ} X^{k_\mu}_{I\mu} X^{-k_\nu}_{I\nu} X^{k_\lambda}_{J\lambda} X^{-k_\sigma}_{J\sigma} W^q_{IJ}$$
   where $X$ are interpolating vectors and $W$ is the Coulomb kernel.
2. Linear Scaling Algorithms: By leveraging the convolution theorem and Fast Fourier Transforms (FFT), the evaluation of metric tensors, Coulomb kernels, and exchange matrices is accelerated. This reduces the scaling of integral evaluation, Coulomb matrix construction, and exchange matrix construction to strictly linear with respect to $N_k$.
3. Quantum Embedding and Local Correlation: The FFTISDF formalism is applied to:
   • Density Matrix Embedding Theory (DMET): Accelerating the construction of the fragment Hamiltonian by transforming the k-space crystal orbital Hamiltonian to the fragment space using the factorized integrals.
   • Local Natural Orbital (LNO) Correlation: Facilitating global low-level correlation calculations (e.g., k-SOS-MP2, k-dRPA) and local high-level corrections (LNO-CCSD, LNO-CCSD(T)).
4. Extrapolation Strategies: The authors employ two extrapolation techniques to reach the TDL:
   • k-point Extrapolation: Linear fitting of energies against $N_k^{-1}$ for k-point sampling.
   • Correlation Domain Extrapolation: Using k-SOS-MP2 as a baseline to extrapolate LNO-CCSD and LNO-CCSD(T) energies to the full correlation domain by varying truncation thresholds.

Key Contributions

• Algorithmic Implementation: A robust implementation of FFTISDF for periodic systems that enables linear scaling with $N_k$ for many-body calculations, overcoming the quadratic scaling bottlenecks of traditional Gaussian Density Fitting (GDF) and FFT-based Density Fitting (FFTDF) in exchange and Hamiltonian construction.
• Scalability Demonstration: Successful execution of correlated ground-state calculations (DMET and LNO-based CCSD/CCSD(T)) on crystalline systems using up to 1000 k-points.
• Storage Efficiency: A dramatic reduction in storage requirements, scaling as $N_k N_{IP}^2 + N_k N_{IP} N_{AO}$, compared to the $N_k^2 N_{AO}^2 N_{aux}$ scaling of GDF. For a 1000 k-point diamond calculation, storage was reduced from ~881 GB to ~2.1 GB.
• Thermodynamic Limit Estimates: The ability to extrapolate reliable CCSD(T) ground-state energies to the thermodynamic limit for both weakly and strongly correlated solids.

Results
The methodology was benchmarked on four representative systems:

1. Diamond (Wide-bandgap semiconductor):
   • FFTISDF with $c_{IP}=14$ achieved per-atom energy errors below $10^{-4}$ a.u. relative to high-cutoff FFTDF references.
   • Wall-clock time analysis confirmed linear scaling for FFTISDF across all steps (integrals, Coulomb, exchange, Hamiltonian), whereas GDF and FFTDF exhibited quadratic scaling for exchange and Hamiltonian construction.
   • TDL extrapolation yielded converged CCSD(T) energies, demonstrating the method's ability to sample the energy-increasing regime of k-point convergence that is inaccessible to standard methods due to cost.
2. Carbon Dioxide (Molecular crystal):
   • Validated an extrapolation procedure where LNO energies (truncated domains) were linearly fitted against the difference from global k-SOS-MP2 energies ($\Delta_{LNO}$).
   • Extrapolated TDL energies for MP2, CCSD, and CCSD(T) showed excellent agreement with global reference values (errors $\sim 10^{-5}$ to $10^{-4}$ a.u. per atom), confirming the efficacy of the domain extrapolation strategy.
3. Nickel Monoxide (NiO) and CaCuO2 (Strongly correlated cuprates):
   • Applied FFTISDF-DMET with a CCSD solver to determine magnetic exchange constants ($J$).
   • For NiO, DMET yielded a next-nearest-neighbor exchange constant $J_2 = -14.05$ meV, an improvement over the underestimation by k-HF ($-7.59$ meV), though still slightly below experimental ranges.
   • For CaCuO2, DMET yielded $J = -165.44$ meV, significantly larger in magnitude than the k-HF value ($-38.97$ meV) and in good agreement with experimental values derived from spin-wave dispersion ($-142$ to $-158$ meV).

Significance and Claims
The paper claims that this work establishes a robust computational framework capable of obtaining converged ground-state energies from high-level correlated methods (DMET, LNO-MP2, LNO-CCSD, LNO-CCSD(T)) for crystals in the thermodynamic limit. The primary significance lies in the linear scaling with respect to the number of k-points, which enables the practical extraction of TDL results for systems that were previously computationally prohibitive.

The authors note that while the method is successful, it currently requires a relatively large number of interpolation points (typically 10 times the number of atomic orbitals) to achieve desired accuracy, making the determination of these points computationally demanding for large unit cells. Additionally, the real-space grid size remains large for systems with transition metals. The authors conclude that while the framework is a significant step forward, future work should focus on developing more efficient protocols for generating compact interpolation points and exploring alternative correlating spaces for LNO methods.