🔬 materials science

CP2K on the road to exascale

This paper presents the CP2K atomistic simulation package, highlighting its second-generation Car-Parrinello molecular dynamics capabilities and detailing development strategies for massively parallel low-scaling solvers and approximate computing on low-precision hardware to prepare for the exascale era.

Original authors: Thomas D. Kühne, Christian Plessl, Robert Schade, Ole Schütt

Published 2026-09-29

📖 5 min read🧠 Deep dive

CC BY 4.0

Original authors: Thomas D. Kühne, Christian Plessl, Robert Schade, Ole Schütt

Original paper licensed under CC BY 4.0 (http://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Imagine trying to understand how a drop of water behaves, not just as a static liquid, but as a living, breathing system where every single atom is constantly moving, colliding, and reacting. To do this, scientists use a powerful type of computer simulation that treats matter according to the strict laws of quantum mechanics, the physics that governs the very small. These simulations are incredibly demanding because they must track the behavior of electrons and nuclei simultaneously, a task that requires immense computing power. For decades, researchers have relied on supercomputers to run these calculations, but as the systems they want to study grow larger—moving from tiny clusters to vast materials with millions of atoms—the time required to get an answer has become a bottleneck. The goal has always been to reach a point where computers can solve these problems so fast that scientists can watch complex chemical processes unfold in real time, simulating entire worlds of atoms rather than just a few.

A team of researchers behind the CP2K software package, a widely used tool for these kinds of atomic simulations, has been working to bridge the gap between current capabilities and the next generation of supercomputers, known as exascale machines. These new computers are designed to perform a quintillion calculations per second, offering a massive leap in speed. However, simply having faster hardware is not enough; the software must be rewritten to handle the unique way these machines operate. The team's work focuses on refining the core algorithms that drive the simulations, ensuring they can scale up to use millions of processor cores simultaneously without getting bogged down. They are particularly interested in a specific method called second-generation Car-Parrinello molecular dynamics, which allows the computer to simulate the movement of atoms over very long periods, such as nanoseconds, while keeping the electrons in their most stable state. This approach is crucial because it allows for the study of realistic, finite-temperature systems, like liquids or solids, rather than just frozen snapshots of matter.

The researchers have identified three major trends in modern computing that their software must adapt to: the need for massive parallel processing, the widespread use of graphics processing units (GPUs) for calculation, and the availability of hardware that performs calculations with varying levels of precision. To meet these challenges, the team has restructured how their code handles data. They have developed new ways to divide large mathematical problems into smaller, independent chunks that different parts of a computer can solve at the same time without needing to constantly stop and check in with each other. This is achieved by organizing data into fixed blocks assigned to specific threads of processing, which eliminates the delays that usually occur when many processors try to access the same information simultaneously. They have also integrated support for GPUs, the powerful chips originally designed for video games but now essential for scientific computing, allowing the software to offload the heaviest calculations to these specialized accelerators.

A significant portion of their effort involves tackling the mathematical complexity of electron interactions, which are the most expensive part of these simulations. Traditionally, calculating how electrons repel each other in a large system required a number of steps that grew explosively as the system got bigger, making large-scale simulations nearly impossible. The team has introduced new techniques to screen out negligible interactions, effectively ignoring the tiny forces that do not matter, which reduces the computational cost dramatically. They have also developed methods to approximate the density of electrons using smaller, easier-to-calculate matrices, correcting the results later with a more precise but cheaper calculation. This allows them to maintain high accuracy while drastically cutting down the time needed. Furthermore, they are exploring the use of "mixed precision," where the computer performs the bulk of its work using a lower level of numerical accuracy that is faster to compute, and then uses a specialized mathematical correction to ensure the final result remains exact.

The results of these developments are already impressive. In simulations of bulk water, the researchers demonstrated that their optimized code could handle systems with up to one million atoms on a supercomputer with tens of thousands of cores. More recently, by combining their new submatrix method with the purification scheme for electron density, they achieved a record-breaking simulation of a system containing more than 100 million atoms. This massive calculation ran on a supercomputer with a sustained performance of 324 petaflops, utilizing more than 67 percent of the machine's total theoretical power. This efficiency proves that the software can effectively harness the raw speed of modern hardware to solve problems that were previously out of reach.

Looking ahead, the authors see a future where these classical simulations are augmented by other emerging technologies. They suggest that artificial intelligence, specifically deep neural networks, could be used to make decisions during a simulation, such as determining whether a complex calculation needs to be done from scratch or if a trained network can provide a sufficiently accurate answer. They also envision hybrid approaches where a quantum computer, which operates on the principles of quantum mechanics, could work alongside the classical supercomputer to solve specific parts of the problem, such as calculating the behavior of electrons in a way that is currently too difficult for classical machines. While these future applications are still in the planning stages, the current work on CP2K lays the essential groundwork, ensuring that the software is ready to evolve alongside the hardware, turning the promise of exascale computing into a practical reality for understanding the material world.

Technical Summary: CP2K on the Road to Exascale

Problem Statement
The CP2K program package, a versatile tool for atomistic simulations ranging from condensed matter physics to quantum chemistry, faces the challenge of scaling to the exascale era. While CP2K has historically excelled in parallel scalability, the transition to exascale computing requires addressing three specific trends: the necessity for massive parallelism (increasing node and thread counts), the widespread adoption of GPUs, and the emergence of mixed-precision hardware accelerators. The core computational bottlenecks involve the efficient treatment of large-scale orbital-free and Kohn–Sham density functional theory (DFT) calculations, as well as wavefunction-based correlation methods (such as Hartree-Fock, MP2, RPA, and GW). Specifically, the calculation of four-center two-electron repulsion integrals (ERIs) and the solution of eigenvalue problems for systems containing millions of atoms demand algorithms that maintain linear scaling with system size while minimizing synchronization overhead and leveraging low-precision arithmetic without sacrificing accuracy.

Methodology
The paper outlines a multi-faceted development strategy to optimize CP2K for future architectures:

Algorithmic Optimization for Post-Hartree-Fock Methods:
- Integral Screening: To reduce the scaling of ERIs from $O(N^4)$ to $O(N^2)$ and eventually linear scaling, the authors employ Schwarz inequality screening and density matrix screening. The latter utilizes the largest density matrix element ( $P_{max}$ ) to filter integrals, ideally using the propagated density matrix from the previous timestep in second-generation Car-Parrinello molecular dynamics (AIMD).
- Approximation Schemes: The implementation relies on the Auxiliary Density Matrix Method (ADMM) and the Resolution of Identity (RI) approach. ADMM approximates the exact Hartree-Fock exchange energy using a smaller, auxiliary density matrix, correcting the difference at the semi-local DFT level. RI factorizes four-center integrals into 2- and 3-center integrals using an auxiliary Gaussian basis set, computed via Cholesky decomposition of the Coulomb metric.
- Linear-Scaling Solvers: The code utilizes the matrix sign function, $sign(A) = A(A^2)^{-1/2}$ , to compute inverses and square roots of large sparse matrices. This function is also used to purify the Kohn–Sham matrix into the density matrix $P$ , enabling linear-scaling electronic structure computations.
Architectural Adaptation for Exascale:
- Parallelism and Threading: Building on its MPI foundation, CP2K integrates OpenMP for multi-threading. To prevent bottlenecks from reductions and atomic access collisions, the authors developed the Distributed Block-sparse Matrices (DBM) library, which uses fixed assignments of matrix block rows to threads to eliminate synchronization needs.
- GPU Acceleration: Support for GPUs is implemented at the library level (e.g., DBM, COSMA, SpFFT, DBCSR). The authors note a shift in bottlenecks from register file constraints to host-to-device communication, addressing this via GPU-aware MPI.
- Mixed-Precision and the Submatrix Method: To exploit hardware accelerators capable of high FLOPs in reduced precision (e.g., FP16 with FP32 accumulation), the authors introduce the submatrix method. This approach converts the evaluation of a matrix function on a large distributed sparse matrix into a parallel evaluation on many dense, smaller submatrices. When combined with the purification scheme, this allows for the use of approximate computing paradigms where numerical inaccuracies are compensated by a modified Langevin-type equation.

Key Contributions and Results

Implementation of Advanced Methods: The paper details the implementation of analytic nuclear gradients for arbitrary boundary conditions within post-HF methods (HF, MP2, RPA, hybrid/double-hybrid DFT, GW), which are prerequisites for accurate finite-temperature simulations.
Scalability Demonstrations: Using the linear-scaling approach with a double-ζ valence polarization basis set, the authors demonstrated the ability to simulate systems with up to one million atoms on a Cray XT5.
Record-Breaking Performance: By combining the submatrix method with the purification scheme, the team achieved record-sized linear-scaling electronic structure computations on systems with more than 100 million atoms. This effort yielded a sustained performance of 324 PFLOPs with an efficiency exceeding 67%.
Software Ecosystem: The modular structure of CP2K, utilizing modern Fortran 2008, has facilitated the extraction of stand-alone libraries (e.g., DBCSR) and the integration of GPU-accelerated libraries, enhancing the code's adaptability to diverse hardware.

Significance
The paper positions CP2K as a critical tool for the "immanent exascale era." Its significance lies not only in its ability to handle increasingly complex and accurate simulations but also in its strategic adaptation to hardware trends. By successfully integrating the submatrix method with approximate computing, CP2K addresses the challenge of utilizing low-precision hardware architectures effectively. Furthermore, the authors outline a roadmap for future integration of deep neural networks (DNN) for "learning on-the-fly" schemes and hybrid quantum-classical algorithms for exact configuration interaction and reduced density-matrix functional theory, suggesting that CP2K will remain at the forefront of atomistic simulation methodologies as computing paradigms evolve.

Technical Summary: CP2K on the Road to Exascale

More like this