Cp2k

CP2K Input Generation Skill

Reference for generating, validating, and debugging CP2K input files with DFT-aware checks.

Quick Start

Validate an Input File

# Input file validation (section structure, keyword checks, parameter sanity)
python scripts/validate_input.py path/to/input.inp

# Parse CP2K output for energies, SCF convergence, timings
python scripts/parse_output.py path/to/output.out

All script paths are relative to the skill directory (.claude/skills/cp2k/). Run commands from that directory, or use full paths from the repo root (e.g., python .claude/skills/cp2k/scripts/validate_input.py input.inp).

Important Constraints

DO NOT generate files in the examples/ directory. Always output to the user's working directory or a user-specified path.

Environment-agnostic: this skill assumes cp2k (or cp2k.psmp / cp2k.ssmp) is available on $PATH but does not assume a specific installation path.

CP2K_DATA_DIR must point to the CP2K data directory containing BASIS_MOLOPT, GTH_POTENTIALS, etc. Alternatively, place copies of the required basis set and pseudopotential files in the working directory and reference them by relative path in the input.

---

Input File Structure

CP2K uses a hierarchical, section-based input format. Sections open with &SECTION_NAME and close with &END SECTION_NAME. Keywords are set inside sections. The key nesting:

&GLOBAL                          # Project name, run type, print level
&FORCE_EVAL                      # Method selection (Quickstep, FIST, etc.)
  &DFT                           # All electronic structure settings
    &MGRID                       #   Planewave grid: CUTOFF, REL_CUTOFF, NGRIDS
    &QS                          #   Quickstep method: EPS_DEFAULT, METHOD (GPW/GAPW)
    &SCF                         #   Self-consistent field convergence
      &OT                        #     Orbital transformation (insulators)
      &DIAGONALIZATION            #     Standard diag (metals)
      &OUTER_SCF                 #     Outer loop for tight convergence
      &MIXING                    #     Density mixing for diagonalization
    &XC                          #   Exchange-correlation functional
      &XC_FUNCTIONAL             #     PBE, BLYP, B3LYP, PBE0, etc.
      &VDW_POTENTIAL             #     Dispersion corrections (DFT-D3, rVV10)
  &SUBSYS                        # Atomic structure definition
    &CELL                        #   ABC vectors, PERIODIC XYZ/XY/NONE
    &COORD                       #   Atomic positions (Cartesian or scaled)
    &KIND               #   Per-element basis set and pseudopotential
  &MM                            # Classical force field (FIST method only)
&MOTION                          # Geometry opt, MD, cell opt, NEB
  &GEO_OPT                      #   Geometry optimizer settings
  &CELL_OPT                     #   Cell vector optimization
  &MD                            #   Molecular dynamics settings
  &BAND                          #   NEB settings
  &PRINT                         #   Trajectory, restart, energy output

Minimal Example: Water Single-Point Energy

&GLOBAL
  PROJECT water
  RUN_TYPE ENERGY
  PRINT_LEVEL MEDIUM
&END GLOBAL

&FORCE_EVAL
  METHOD Quickstep
  &DFT
    BASIS_SET_FILE_NAME BASIS_MOLOPT
    POTENTIAL_FILE_NAME GTH_POTENTIALS
    &MGRID
      CUTOFF 400
      REL_CUTOFF 60
    &END MGRID
    &QS
      EPS_DEFAULT 1.0E-12
    &END QS
    &SCF
      EPS_SCF 1.0E-6
      MAX_SCF 50
      &OT
        MINIMIZER DIIS
        PRECONDITIONER FULL_ALL
      &END OT
    &END SCF
    &XC
      &XC_FUNCTIONAL PBE
      &END XC_FUNCTIONAL
    &END XC
  &END DFT
  &SUBSYS
    &CELL
      ABC 10.0 10.0 10.0
      PERIODIC NONE
    &END CELL
    &COORD
      O   0.000  0.000  0.000
      H   0.757  0.586  0.000
      H  -0.757  0.586  0.000
    &END COORD
    &KIND O
      BASIS_SET DZVP-MOLOPT-GTH
      POTENTIAL GTH-PBE-q6
    &END KIND
    &KIND H
      BASIS_SET DZVP-MOLOPT-GTH
      POTENTIAL GTH-PBE-q1
    &END KIND
  &END SUBSYS
&END FORCE_EVAL

---

Method Selection Decision Tree

digraph method_selection {
    rankdir=TB;
    node [shape=diamond];

    calc [label="What calculation?" shape=ellipse];
    sys [label="System type?"];
    metal [label="Metallic?"];

    node [shape=box];
    sp [label="RUN_TYPE ENERGY\n(single point)"];
    geo [label="RUN_TYPE GEO_OPT\n&MOTION > &GEO_OPT"];
    md [label="RUN_TYPE MD\n&MOTION > &MD"];
    cell [label="RUN_TYPE CELL_OPT\n&MOTION > &CELL_OPT\nSTRESS_TENSOR ANALYTICAL"];
    neb [label="RUN_TYPE BAND\n&MOTION > &BAND"];
    vib [label="RUN_TYPE VIBRATIONAL_ANALYSIS\n&VIBRATIONAL_ANALYSIS"];

    diag [label="Diagonalization\n+ SMEARING\n+ MIXING (Broyden)"];
    ot [label="OT minimizer\n+ FULL_ALL preconditioner\n+ OUTER_SCF"];

    calc -> sp [label="energy"];
    calc -> geo [label="geometry opt"];
    calc -> md [label="AIMD"];
    calc -> cell [label="cell opt"];
    calc -> neb [label="NEB"];
    calc -> vib [label="frequencies"];

    sp -> sys; geo -> sys; md -> sys; cell -> sys; neb -> sys; vib -> sys;
    sys -> metal;
    metal -> diag [label="yes (metal/\nsmall gap)"];
    metal -> ot [label="no (insulator/\nmolecular)"];
}

Rule of thumb: OT is faster and more robust for systems with a band gap. Diagonalization is required when fractional occupations exist (metals, radicals with near-degenerate states).

---

Basis Set & Pseudopotential Quick Reference

| Quality | Basis Set | Use Case | |---------|-----------|----------| | Minimum | SZV-MOLOPT-GTH | Quick tests only | | Standard | DZVP-MOLOPT-GTH | Production (most cases) | | Standard (condensed) | DZVP-MOLOPT-SR-GTH | Condensed phase, saves cost | | High | TZVP-MOLOPT-GTH | Benchmark, high accuracy | | Very High | TZV2P-MOLOPT-GTH | Reference calculations |

Functional / Pseudopotential Matching Rules

The pseudopotential family must match the XC functional used to generate it:

| XC Functional | Pseudopotential | Notes | |---------------|-----------------|-------| | PBE | GTH-PBE | Most common choice | | BLYP | GTH-BLYP | Also used for B3LYP | | B3LYP | GTH-BLYP | No GTH-B3LYP exists | | PBE0 | GTH-PBE | PBE-derived hybrid | | SCAN | GTH-PBE | Use PBE pseudo (standard practice) |

Valence electron count must be consistent between basis set and pseudopotential. The -qN suffix on the pseudopotential denotes N valence electrons. Verify: oxygen uses -q6, carbon -q4, hydrogen -q1, nitrogen -q5, silicon -q4, etc.

Data Files

| File | Contents | |------|----------| | BASIS_MOLOPT | MOLOPT basis sets for standard DFT | | BASIS_MOLOPT_UZH | Extended basis sets, needed for hybrid functionals (HFX) | | GTH_POTENTIALS | GTH pseudopotentials for all functional families | | BASIS_ADMM | Auxiliary basis sets for ADMM (hybrid functional acceleration) | | POTENTIAL | All-electron potentials (rarely needed) |

---

Cutoff Selection

CUTOFF controls the planewave expansion of the electron density. REL_CUTOFF controls the finest grid level.

| Parameter | Typical Range | Default Recommendation | |-----------|--------------|----------------------| | CUTOFF | 300-600 Ry | 400 Ry (safe starting point) | | REL_CUTOFF | 50-60 Ry | 60 Ry | | NGRIDS | 4-5 | 5 (default, rarely changed) |

Convergence Testing

Run single-point calculations at increasing CUTOFF (200, 300, 400, 500, 600 Ry). Plot total energy vs cutoff. Converged when energy change 500 atoms) | | FULLSINGLEINVERSE | Medium | Compromise between FULLALL and FULLKINETIC |

Diagonalization -- Metals and Small-Gap Systems

Required when fractional occupations are needed.

&SCF
  EPS_SCF 1.0E-6
  MAX_SCF 100
  ADDED_MOS 20            # Extra empty states (increase for metals)
  &DIAGONALIZATION
  &END DIAGONALIZATION
  &SMEARING
    METHOD FERMI_DIRAC
    ELECTRONIC_TEMPERATURE [K] 300
  &END SMEARING
  &MIXING
    METHOD BROYDEN_MIXING
    ALPHA 0.2              # Mixing parameter (lower = more stable, slower)
  &END MIXING
&END SCF

EPS_SCF Guidelines

| Calculation Type | EPS_SCF | Rationale | |-----------------|---------|-----------| | MD (AIMD) | 1E-6 | Sufficient for forces | | Geometry optimization | 1E-8 | Tight convergence needed for gradients | | Cell optimization | 1E-8 | Stress tensor requires tight SCF | | Vibrational analysis | 1E-8 | Numerical derivatives amplify noise |

Rule: sqrt(EPS_DEFAULT) should be smaller than EPS_SCF. Default EPS_DEFAULT of 1E-12 pairs with EPS_SCF of 1E-6.

---

Simulation Workflow

0. Structure Preparation

Obtain coordinates from .xyz, .cif, databases (Materials Project, COD), or build manually. Set up &CELL (lattice vectors or ABC + angles) and &COORD (Cartesian or SCALED).

For isolated molecules: use PERIODIC NONE in &CELL and POISSON_SOLVER MT (or WAVELET) in &POISSON.

For periodic systems: ensure cell vectors are consistent with crystal structure. Use PERIODIC XYZ (3D), PERIODIC XY (slab), or PERIODIC X (wire).

1. Geometry Optimization

&GLOBAL
  RUN_TYPE GEO_OPT
&END GLOBAL
&MOTION
  &GEO_OPT
    OPTIMIZER BFGS          # BFGS (default, fast) or LBFGS (large systems) or CG
    MAX_ITER 500
    MAX_DR 3.0E-3           # Max displacement [bohr]
    MAX_FORCE 4.5E-4        # Max force [bohr/hartree]
    RMS_DR 1.5E-3
    RMS_FORCE 3.0E-4
  &END GEO_OPT
&END MOTION

Use tight EPS_SCF (1E-8). For cells under pressure or with soft modes, also optimize cell vectors (use RUN_TYPE CELL_OPT instead).

2. Equilibration (AIMD)

Short run to thermalize the system. Use a robust thermostat (CSVR recommended for equilibration).

&GLOBAL
  RUN_TYPE MD
&END GLOBAL
&MOTION
  &MD
    ENSEMBLE NVT
    STEPS 5000
    TIMESTEP 0.5            # fs
    TEMPERATURE 300
    &THERMOSTAT
      TYPE CSVR
      &CSVR
        TIMECON 50           # fs, coupling constant
      &END CSVR
    &END THERMOSTAT
  &END MD
&END MOTION

Add wavefunction extrapolation to reduce SCF iterations per MD step:

&DFT
  WFN_RESTART_FILE_NAME RESTART.wfn   # restart from previous wavefunction
  &SCF
    SCF_GUESS RESTART
  &END SCF
  &XC
    ...
  &END XC
&END DFT

In &FORCE_EVAL > &DFT:

  EXTRAPOLATION PS
  EXTRAPOLATION_ORDER 3

3. Production (AIMD)

Continue from equilibrated state. Longer run with trajectory output.

&MOTION
  &MD
    ENSEMBLE NVT
    STEPS 50000
    TIMESTEP 0.5
    TEMPERATURE 300
    &THERMOSTAT
      TYPE NOSE
      REGION MASSIVE
      &NOSE
        LENGTH 3
        TIMECON 100          # fs
      &END NOSE
    &END THERMOSTAT
  &END MD
  &PRINT
    &TRAJECTORY
      &EACH
        MD 10
      &END EACH
    &END TRAJECTORY
    &VELOCITIES
      &EACH
        MD 10
      &END EACH
    &END VELOCITIES
    &FORCES
      &EACH
        MD 100
      &END EACH
    &END FORCES
    &RESTART
      &EACH
        MD 500
      &END EACH
    &END RESTART
    &RESTART_HISTORY
      &EACH
        MD 5000
      &END EACH
    &END RESTART_HISTORY
  &END PRINT
&END MOTION

Thermostat Comparison

| Thermostat | CP2K Keyword | Use For | Notes | |------------|-------------|---------|-------| | Nose-Hoover | NOSE | General NVT production | Correct canonical; use REGION MASSIVE | | CSVR | CSVR | Equilibration, robust | Fast thermalization, correct canonical ensemble | | GLE | GLE | Path integrals, NQE | Generalized Langevin, advanced use |

Ensemble Selection

| Ensemble | CP2K Keyword | Use For | |----------|-------------|---------| | NVE | NVE | Microcanonical, transport properties, energy conservation tests | | NVT | NVT | Most simulations, thermalized sampling | | NPTI | NPTI | Isotropic pressure (liquids, amorphous) | | NPTF | NPTF | Anisotropic pressure (crystals, anisotropic cells) |

AIMD Timestep

• 0.5 fs -- safe default for systems with hydrogen
• 1.0 fs -- acceptable for heavy-element-only systems (no H, no Li)
• Reduce to 0.25 fs for high-temperature simulations (> 1000 K) or light-element-rich systems

Wavefunction Extrapolation

Always enable for MD. Extrapolation predicts the initial wavefunction guess from previous steps, reducing SCF iterations by 30-50%.

EXTRAPOLATION PS
EXTRAPOLATION_ORDER 3

Use PS (polynomial scheme). Order 3 is the recommended default. Higher orders can introduce instabilities.

---

Output Files

| File | Content | |------|---------| | PROJECT.out | Main log: SCF convergence, energies, forces, timings | | PROJECT-pos-1.xyz | Trajectory (atomic positions) | | PROJECT-vel-1.xyz | Velocities | | PROJECT-frc-1.xyz | Forces on atoms | | PROJECT-1.ener | Step, time, kinetic E, temperature, potential E, conserved quantity | | PROJECT-1.cell | Cell parameters per step (NPT/CELLOPT) | | PROJECT-1.restart | Human-readable restart file | | PROJECT-RESTART.wfn | Binary wavefunction restart (for SCFGUESS RESTART) | | PROJECT-RESTART.wfn.bak-1 | Previous wavefunction backup |

Key Output Patterns to Check

• SCF convergence: grep for SCF run converged or SCF run NOT converged
• Total energy: grep for ENERGY| Total FORCE_EVAL
• Forces: grep for SUM OF ATOMIC FORCES (should be near zero for equilibrium)
• Geometry opt: grep for OPTIMIZATION COMPLETED or MAXIMUM NUMBER OF OPTIMIZATION STEPS REACHED

---

Validation

validate_input.py Checks

| Check | What It Catches | |-------|-----------------| | Section nesting | Unclosed sections, mismatched &END tags | | Required sections | Missing &GLOBAL, &FORCE_EVAL, &SUBSYS | | Basis/pseudo consistency | Functional-pseudopotential family mismatch | | Valence electrons | -qN suffix inconsistency with element | | CUTOFF range | Values below 200 Ry or above 1200 Ry flagged | | RELCUTOFF range | Values below 40 Ry flagged | | EPSSCF vs run type | Loose EPSSCF with geometry/cell optimization | | OT on metals | OT without band gap warning | | CELLOPT stress | Missing STRESS_TENSOR ANALYTICAL | | Periodicity/Poisson | PERIODIC NONE without appropriate Poisson solver | | Wavefunction extrapolation | MD without EXTRAPOLATION flagged | | Thermostat region | NVT without REGION MASSIVE flagged |

---

Common Mistakes

1. Basis/pseudopotential functional mismatch. GTH-PBE pseudopotentials with BLYP functional. Match the family: PBE functional uses GTH-PBE, BLYP uses GTH-BLYP.

1. Valence electron count mismatch. Using -q4 pseudopotential where -q6 is needed (e.g., oxygen). Check element-specific valence counts in GTH_POTENTIALS.

1. CUTOFF too low. Values below 200 Ry produce inaccurate forces and energies. Start at 400 Ry and converge systematically.

1. OT on metallic systems. OT assumes integer occupations. Metals require diagonalization with SMEARING. Symptom: SCF oscillates and never converges.

1. Missing STRESSTENSOR ANALYTICAL for CELLOPT. Cell optimization needs stress tensors. Without this keyword, CP2K uses numerical stress (slow and less accurate) or fails.

1. EPS_SCF too loose for geometry/cell optimization. MD tolerates 1E-6; geometry and cell optimization need 1E-8. Loose SCF produces noisy gradients that prevent convergence.

1. No wavefunction extrapolation for MD. Without EXTRAPOLATION PS, each MD step starts SCF from atomic guess, wasting 2-5x SCF iterations.

1. Missing OUTERSCF with OT. A single OT SCF loop may not converge tightly enough. OUTERSCF provides an additional convergence layer that re-evaluates the Kohn-Sham matrix.

1. Wrong periodicity for isolated molecules. Use PERIODIC NONE in &CELL and set &POISSON solver to MT (Martyna-Tuckerman) or WAVELET. Default periodic boundary conditions introduce spurious interactions with images.

1. Forgetting REGION MASSIVE for NVT thermostat. Without REGION MASSIVE, Nose-Hoover couples to global kinetic energy only, leading to poor temperature equipartition. Always set REGION MASSIVE

…

Source & license

This open-source skill is cataloged on AgentStack and links to its original source — we do not rehost the code.

• Author: WuGroup-XJTLU
• Source: WuGroup-XJTLU/cc-skills-ZhenghaoWu-Group
• License: Apache-2.0

Install and usage instructions live in the source repository linked above.