cp2k-workflow

CP2K Workflow

Description

This skill covers end-to-end DFT, AIMD, and atomistic simulation workflows using CP2K, from input file construction through convergence testing, geometry optimization, and ab initio molecular dynamics. CP2K uses the Gaussian and Plane Waves (GPW) method, which combines a Gaussian-type orbital (GTO) basis with an auxiliary plane-wave density grid, making it efficient for large systems, liquid-phase simulations, and condensed-phase AIMD. Invoke this skill when setting up CP2K calculations, selecting basis sets and GTH pseudopotentials, converging CUTOFF and REL_CUTOFF parameters, diagnosing SCF convergence failures, generating AIMD trajectories, or producing DFT reference data for ML potential training with an open-source code.

Domain Context

CP2K solves the Kohn-Sham DFT equations using the GPW (Gaussian and Plane Waves) or GAPW (Gaussian Augmented Plane Waves) method. The key idea is a dual representation of the electron density: the wavefunction and density matrix are represented in a localized Gaussian basis, while the Hartree potential is evaluated on an auxiliary plane-wave grid using fast Fourier transforms. This avoids the O(N³) cost of dielectric screening in pure plane-wave codes for very large systems, and enables linear-scaling DFT (LS-DFT) via the linear-scaling module.

Physical approximations and structural concepts that must be understood before using this workflow:

• GPW method: The total electron density is expanded as a sum of atom-centered Gaussian functions, then mapped onto a regular real-space grid at a grid spacing controlled by CUTOFF (in Ry). The Hartree energy is evaluated in reciprocal space via FFT. The accuracy of this mapping depends on both CUTOFF and REL_CUTOFF — both must be converged explicitly. Unlike VASP or QE, a high CUTOFF alone is insufficient if REL_CUTOFF is too low.
• GAPW method: Augments GPW with atom-centered hard and soft auxiliary density expansions near nuclei (analogous to PAW), allowing all-electron calculations or harder core regions without increasing the plane-wave grid cost. Used for core-level spectroscopy (XPS, XANES) and when soft pseudopotentials are insufficient. Higher cost than GPW for the same basis.
• GTH pseudopotentials: Goedecker-Teter-Hutter norm-conserving pseudopotentials are the standard for CP2K. They are available for LDA, PBE, BLYP, PBE0, and several other functionals. The pseudopotential and basis set must be generated for and matched to the same XC functional — using a GTH-PBE pseudopotential with a BLYP functional is inconsistent and produces wrong results.
• Gaussian basis sets: CP2K ships MOLOPT-optimized basis sets (short-range SR variants for solids, standard for molecules). Quality tiers: SZV (minimal, screening only), DZVP (double-ζ + polarization, standard production), TZV2P (triple-ζ, higher accuracy), TZV2PX (extended). Short-range SR basis sets are preferred for periodic solids because they minimize linear dependence and basis set superposition error in condensed phases.
• CUTOFF (Ry): Controls the plane-wave grid density for the auxiliary density. Analogous to ecutrho in QE but expressed per-calculation, not per-element. Higher CUTOFF increases grid density and accuracy of the Hartree term. Typical production values: 400–800 Ry for most GTH basis sets.
• REL_CUTOFF (Ry): Controls the accuracy of the Gaussian-to-grid mapping for each atom type. Each Gaussian is mapped onto a subgrid whose fineness is determined by REL_CUTOFF. If REL_CUTOFF is too small, Gaussian functions are mapped poorly regardless of CUTOFF. Default: 60 Ry; recommended production: 60–100 Ry. Converge by scanning REL_CUTOFF at fixed CUTOFF.
• OT (Orbital Transformation) method: An iterative direct minimization approach for finding the ground-state wavefunction without diagonalizing the Kohn-Sham matrix. Efficient and stable for insulators and semiconductors with a band gap. Requires a good preconditioner. Cannot handle fractional occupations and must not be used for metallic systems. OUTER_SCF is required alongside OT to converge the total energy, not just the orbital gradient.
• Diagonalization + smearing for metals: For metallic systems, OT must be replaced with standard diagonalization (&DIAGONALIZATION ON) plus ADDED_MOS (extra molecular orbitals above the HOMO) and Fermi-Dirac or Gaussian smearing. Broyden or Pulay mixing is then required for charge density convergence.
• Periodic boundary conditions and Poisson solver: CP2K supports fully periodic (3D), slab (2D), wire (1D), and isolated (0D) boundary conditions via the POISSON section. The correct Poisson solver must match the periodicity: PERIODIC for bulk, Martyna-Tuckerman (MT) for isolated molecules or charged slabs, and wavelet-based WAVELET for general non-orthorhombic cells with mixed periodicity.
• AIMD integration: CP2K uses the Verlet algorithm for NVE, and thermostat-augmented integrators for NVT (CSVR, NOSE) and NPT (BAROSTAT + thermostat). The timestep must be ≤ 1/10 of the period of the fastest vibrational mode — typically 0.5 fs for hydrogen-containing systems, 1–2 fs for heavy-element systems.