Geometry Optimization in Quantum ESPRESSO: rela…

Before you compute any property — band structure, phonons, formation energies — you need a structure that sits at a genuine minimum of the potential energy surface. Geometry optimization relaxes atomic positions and, optionally, the unit cell until forces and stresses vanish. This tutorial explains the difference between relax and vc-relax, the namelists that control them, and the pitfalls that trap most newcomers.

Why Relax a Structure at All

Experimental coordinates carry thermal expansion, measurement error, and finite-temperature effects, while your DFT calculation is athermal and functional-specific. A PBE calculation has its own equilibrium geometry that differs from experiment by a percent or two. If you compute forces on an unrelaxed structure, they will be nonzero, and any derived quantity — vibrational frequencies especially — becomes meaningless. Relaxation finds the geometry consistent with your chosen functional and pseudopotentials.

There are two flavors:

relax: optimize atomic positions only; the cell is fixed.
vc-relax: variable-cell relaxation; optimize both atomic positions and the lattice vectors.

flowchart TD
  A[Starting geometry] --> B{Cell known / fixed?}
  B -->|Yes| C["calculation = 'relax'"]
  B -->|No| D["calculation = 'vc-relax'"]
  C --> E["&IONS: ion_dynamics = bfgs"]
  D --> F["&IONS + &CELL<br/>cell_dynamics = bfgs"]
  E --> G[Minimize forces F_I]
  F --> H[Minimize forces F_I<br/>and stress σ]
  G --> I{forc_conv_thr met?}
  H --> J{forc + press thresholds met?}
  I -->|No| G
  J -->|No| H
  I -->|Yes| K[Fixed-cell minimum]
  J -->|Yes| L[Restart vc-relax once<br/>to reduce Pulay stress]
  L --> M[Production geometry]
  K --> M

Ionic forces follow from the Hellmann–Feynman theorem:

$$ \mathbf{F}_{I} = -\frac{\partial E}{\partial\mathbf{R}_{I}} $$

and cell relaxation is driven by the stress tensor $\boldsymbol{\sigma}$.

Fixed-Cell Relaxation with calculation = ‘relax’

Use relax when the cell is known or fixed by symmetry — for example, relaxing atoms in a slab, a molecule in a box, or a defect in a host lattice. Quantum ESPRESSO moves atoms along the calculated Hellmann–Feynman forces using, by default, the BFGS quasi-Newton algorithm.

&CONTROL
    calculation   = 'relax'
    prefix        = 'batio3'
    outdir        = './out'
    pseudo_dir    = './pseudo'
    forc_conv_thr = 1.0d-4
    etot_conv_thr = 1.0d-5
    nstep         = 100
/
&SYSTEM
    ibrav     = 1
    celldm(1) = 7.57
    nat       = 5
    ntyp      = 3
    ecutwfc   = 50.0
    ecutrho   = 400.0
/
&ELECTRONS
    conv_thr    = 1.0d-9
    mixing_beta = 0.7
/
&IONS
    ion_dynamics = 'bfgs'
/
ATOMIC_SPECIES
 Ba  137.327  Ba.pbe-spn-kjpaw_psl.1.0.0.UPF
 Ti   47.867  Ti.pbe-spn-kjpaw_psl.1.0.0.UPF
 O    15.999  O.pbe-n-kjpaw_psl.1.0.0.UPF
ATOMIC_POSITIONS (crystal)
 Ba  0.00  0.00  0.00
 Ti  0.50  0.50  0.50
 O   0.50  0.50  0.00
 O   0.50  0.00  0.50
 O   0.00  0.50  0.50
K_POINTS (automatic)
 6 6 6 0 0 0

pw.x -in batio3.relax.in > batio3.relax.out

The &IONS namelist

The &IONS namelist is required for any relax, vc-relax, or molecular dynamics run. Key options:

ion_dynamics: 'bfgs' (default, quasi-Newton, best for local minima) or 'damp' (damped dynamics).
bfgs_ndim: number of previous steps used to build the Hessian; increasing to 3-4 can speed up difficult relaxations.
pot_extrapolation and wfc_extrapolation: reuse the previous step’s potential and wavefunctions to accelerate the SCF at each new geometry.

Convergence Criteria: forc_conv_thr and etot_conv_thr

A relaxation stops when both conditions are satisfied simultaneously:

Total energy change between steps is below etot_conv_thr (Ry): $\lvert E_{i+1}-E_{i}\rvert < \varepsilon_{E}$.
The largest force component is below forc_conv_thr (Ry/bohr): $\max_{I,\alpha}\lvert F_{I\alpha}\rvert < \varepsilon_{F}$.

Parameter	Namelist	Default	Typical tight value
`etot_conv_thr`	`&CONTROL`	1.0d-4 Ry	1.0d-5 Ry
`forc_conv_thr`	`&CONTROL`	1.0d-3 Ry/bohr	1.0d-4 Ry/bohr
`press_conv_thr`	`&CELL`	0.5 kbar	0.1 kbar
`conv_thr`	`&ELECTRONS`	none	1.0d-9

A critical rule: the electronic conv_thr must be much tighter than the ionic thresholds. Forces are derivatives of the energy, so a loosely converged SCF produces noisy forces that prevent the ionic loop from ever converging. A good heuristic is conv_thr around 1.0d-9 when forc_conv_thr is 1.0d-4.

Variable-Cell Relaxation with calculation = ‘vc-relax’

When you do not know the equilibrium lattice parameters — a new compound, a high-pressure phase, or any bulk optimization — use vc-relax. It adds cell degrees of freedom driven by the calculated stress tensor and requires the &CELL namelist.

&CONTROL
    calculation    = 'vc-relax'
    prefix         = 'batio3'
    outdir         = './out'
    pseudo_dir     = './pseudo'
    forc_conv_thr  = 1.0d-4
    etot_conv_thr  = 1.0d-5
    nstep          = 200
/
&SYSTEM
    ibrav     = 1
    celldm(1) = 7.57
    nat       = 5
    ntyp      = 3
    ecutwfc   = 50.0
    ecutrho   = 400.0
/
&ELECTRONS
    conv_thr    = 1.0d-9
    mixing_beta = 0.7
/
&IONS
    ion_dynamics = 'bfgs'
/
&CELL
    cell_dynamics  = 'bfgs'
    press          = 0.0
    press_conv_thr = 0.1
    cell_dofree    = 'all'
/
ATOMIC_SPECIES
 Ba  137.327  Ba.pbe-spn-kjpaw_psl.1.0.0.UPF
 Ti   47.867  Ti.pbe-spn-kjpaw_psl.1.0.0.UPF
 O    15.999  O.pbe-n-kjpaw_psl.1.0.0.UPF
ATOMIC_POSITIONS (crystal)
 Ba  0.00  0.00  0.00
 Ti  0.50  0.50  0.50
 O   0.50  0.50  0.00
 O   0.50  0.00  0.50
 O   0.00  0.50  0.50
K_POINTS (automatic)
 6 6 6 0 0 0

pw.x -in batio3.vcrelax.in > batio3.vcrelax.out

The &CELL namelist

cell_dynamics: 'bfgs' for local minimization; 'sd' (steepest descent) rarely needed.
press: target external pressure in kbar. Set to 0.0 for ambient-pressure equilibrium; set positive for high-pressure phases.
press_conv_thr: stress convergence threshold in kbar.
cell_dofree: constrains which cell parameters move. 'all' relaxes everything; 'ibrav' preserves the Bravais symmetry; 'volume', '2Dxy', 'z', and others fix specific degrees of freedom. Use this to enforce known symmetry — for example, keep a cubic cell cubic.

Reading Forces and Stress from the Output

Quantum ESPRESSO prints forces and stress at every ionic step. To monitor progress:

grep "Total force" batio3.vcrelax.out
grep -A5 "total   stress" batio3.vcrelax.out

The total force should decrease toward forc_conv_thr. The stress tensor is printed in both Ry/bohr^3 and kbar; the reported scalar pressure P= should approach your target press. When you see bfgs converged in N scf cycles, the run succeeded and the final block prints the optimized CELL_PARAMETERS and ATOMIC_POSITIONS — copy these into your production input.

Cell Optimization Pitfalls

Variable-cell runs have failure modes that fixed-cell runs do not, and nearly all of them trace back to the plane-wave basis.

Pulay stress and the ecutwfc trap

The plane-wave basis set is defined by ecutwfc at a fixed cell volume. When vc-relax changes the volume, the number of plane waves would ideally change too, but the basis is frozen at the initial cell. This mismatch produces a spurious Pulay stress that biases the equilibrium volume — typically making the relaxed cell too small at low cutoffs. Two defenses:

Use a generous ecutwfc, 30-40% above your converged value, so the residual Pulay stress is negligible.
Two-step protocol: after the first vc-relax converges, restart a second vc-relax from the new geometry so the basis is re-generated at the updated cell. Repeat until the volume stops changing between restarts.

Other common traps

Too few k-points: stress is even more sensitive to k-sampling than energy. Under-converged grids give wrong lattice constants.
Loose SCF convergence: noisy stress prevents cell convergence; keep conv_thr tight.
Broken symmetry: without cell_dofree, a cubic crystal can drift into a slightly distorted cell due to numerical noise. Constrain it if the phase is known.
Getting stuck in metastable minima: BFGS finds the nearest local minimum, not the global one. For structure prediction, sample multiple starting configurations.

A Practical Recipe

Converge ecutwfc, ecutrho, and the k-grid on total energy for a fixed cell.
Run vc-relax with a generous cutoff and tight conv_thr.
Restart a second vc-relax from the relaxed cell to remove Pulay stress.
Confirm forces are below forc_conv_thr and pressure below press_conv_thr.
Use the final geometry for production property calculations.

Monitoring and Restarting a Relaxation

Long relaxations sometimes hit the nstep limit before converging, or a walltime limit kills the job. Quantum ESPRESSO writes restart information continuously, so you can resume rather than start over. To restart, set restart_mode = 'restart' in &CONTROL and resubmit with the same input — the code picks up from the last completed ionic step using the saved wavefunctions and charge density in outdir.

&CONTROL
    calculation   = 'vc-relax'
    restart_mode  = 'restart'
    prefix        = 'batio3'
    outdir        = './out'
    max_seconds   = 82800
/

Setting max_seconds slightly below your scheduler walltime lets pw.x stop cleanly and flush a restartable checkpoint instead of being terminated mid-step. To watch convergence live, track the printed energy and force at each step:

grep "!    total energy" batio3.vcrelax.out | tail
grep "Total force" batio3.vcrelax.out | tail

A healthy relaxation shows the total energy decreasing monotonically and the total force shrinking toward the threshold. Oscillating energies usually mean mixing_beta is too aggressive or the SCF is under-converged.

Constraining and Fixing Atoms

Not every atom should move. In slab calculations you often fix the bottom layers to mimic the bulk while relaxing only the surface, and in adsorption studies you may freeze the substrate. Quantum ESPRESSO supports this directly in ATOMIC_POSITIONS by appending three integer flags (0 or 1) after each coordinate line: 1 allows motion, 0 freezes that Cartesian component.

ATOMIC_POSITIONS (angstrom)
 O   0.000  0.000  8.000  1 1 1
 Ti  0.000  0.000  4.000  1 1 1
 Ti  0.000  0.000  0.000  0 0 0

Here the bottom titanium is completely fixed while the upper atoms relax freely. This selective relaxation dramatically reduces the degrees of freedom and keeps the optimization anchored to a physically sensible reference, which is essential for surfaces, interfaces, and defect calculations.

For related workflows, browse the blog and the platform’s computational engines; the full input reference lives in the official Quantum ESPRESSO documentation.

Run it on Simatra

Relaxations are iterative — each ionic step is a full SCF, and vc-relax restarts multiply the cost, especially for large or low-symmetry cells. Simatra runs Quantum ESPRESSO on GPU-accelerated clusters (instance GPU-Opt-V2) with up to 5x faster convergence and support for supercells up to ~2,000 atoms, turning multi-day relaxations into hours. Start a free trial with $100 in credits at app.simatra.io, and explore the computational engines — Quantum ESPRESSO plus our native CUDA/HIP/Metal KRONOS engine — that power the platform.