comp-chem/cp2k
1
0
Fork 0
mirror of https://github.com/cp2k/cp2k.git synced 2026-07-28 22:25:32 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions

Expose LBFGS print control as input keyword, increase default verbosity, and expand document for optimization (#5274)

Browse source
This commit is contained in:
HE Zilong 2026-05-26 17:05:23 +08:00 • committed by GitHub
parent a7cc81a68e
commit 3274ff38b0
No known key found for this signature in database
GPG key ID: B5690EEEBB952194
5 changed files with 303 additions and 83 deletions
Download patch file Download diff file Expand all files Collapse all files

305
docs/methods/optimization/geometry_and_cell_opt.md
View file

@ -32,10 +32,20 @@ only internal atomic coordinates should relax.
Use `CELL_OPT` when the equilibrium cell volume, shape, pressure, or residual stress is part of the
question. It is most common for bulk crystals, strained solids, two-dimensional materials with
relaxed in-plane lattice constants, and structures prepared for fixed-cell molecular dynamics. A
cell optimization is not a replacement for finite-temperature pressure sampling; if the desired
quantity is a thermal average at finite temperature, use an appropriate NPT molecular-dynamics
workflow instead.
relaxed in-plane lattice constants, and structures prepared for fixed-cell molecular dynamics.
```{note}
An optimization is not a replacement for finite-temperature pressure sampling; if the desired
quantity is a thermal average at finite temperature, use a molecular-dynamics workflow with
appropriate ensemble instead.
For the usual electronic-structure methods based on the Born-Oppenheimer approximation
(which, by neglecting nuclear motion, provides the concept of "potential-energy surface"
in the first place), there is no involvement of the temperature (kinetic energy of the
nuclei) in an optimization, and the energy that is being optimized does not include
components like zero-point energy and thermal corrections (harmonic or anharmonic).
The notion "optimization at 0 K" is really a misnomer in this regard.
```
## Basic setup
@ -58,7 +68,7 @@ block:
&GEO_OPT
TYPE MINIMIZATION
OPTIMIZER BFGS
MAX_ITER 200
MAX_ITER 200
MAX_DR 3.0E-03
RMS_DR 1.5E-03
MAX_FORCE 4.5E-04
@ -87,7 +97,7 @@ sense as well:
&MOTION
&CELL_OPT
OPTIMIZER BFGS
MAX_ITER 200
MAX_ITER 200
EXTERNAL_PRESSURE 1.01325E+00
PRESSURE_TOLERANCE 100.0
MAX_DR 3.0E-03
@ -98,6 +108,14 @@ sense as well:
&END MOTION
```
```{note}
The availability of analytical stress tensor depends on the method for evaluating energy and force.
Some electronic-structure methods may only support numerical stress tensor from finite differences,
which would take quite more time at each iteration of optimization. This can be confirmed by setting
[FORCE_EVAL/PRINT/STRESS_TENSOR](#CP2K_INPUT.FORCE_EVAL.PRINT.STRESS_TENSOR) and seeing if numerical
stress is mentioned in the output.
```
For `GEO_OPT`, [TYPE](#CP2K_INPUT.MOTION.GEO_OPT.TYPE) selects the optimization target.
`MINIMIZATION` searches for a local minimum and is the normal choice. `TRANSITION_STATE` activates
transition-state search machinery; it requires additional method-specific settings in
@ -108,30 +126,72 @@ treated as an ordinary minimization.
optimization iteration can require more than one force evaluation, depending on the optimizer, line
search, and the selected `CELL_OPT` mode.
An optimization may terminate if one of the following events takes place:
- The convergence criteria (see below) are satisfied before reaching `MAX_ITER`;
- The `MAX_ITER` is reached, regardless of whether the convergence criteria are satisfied or not;
- Some external control for the program is activated, such as hitting the global walltime as
specified by [WALLTIME](#CP2K_INPUT.GLOBAL.WALLTIME), or discovering a file in the working
directory with filename `EXIT` or `EXIT_GEO` that is created by user in case of emergency.
## Starting structure and cell
Start from a chemically and physically sensible structure. Severe close contacts, unrealistic
coordination environments, or a badly prepared cell can lead to unstable forces, SCF failures,
unphysical atom displacements, or optimizer steps that are difficult to recover from.
unphysical atom displacements, or optimizer steps that are difficult to recover from. For instance,
"amorphous cell" structures generated by randomly packing molecules and/or atoms in a box can be
extremely challenging for optimization, especially if the chemical bonds between atoms are
significantly scrambled. Same for snapshots from a high-temperature molecular-dynamics simulation in
vaporized or molten conditions with lots of broken chemical bonds.
For `CELL_OPT`, the initial cell matters as much as the initial coordinates. The starting volume and
shape should be close enough to the expected structure that the pressure and stress are not
dominated by preparation artefacts. For slabs, two-dimensional materials, or systems with vacuum, do
not relax the vacuum direction unless this is physically intended; constrain the appropriate cell
shape should be close enough to the expected structure and density, such that the pressure and
stress are not dominated by preparation artefacts.
A visualization of the structure and cell in modelling programs, with the box and the neighboring
periodic images displayed, will be very helpful; neither large vacuous gaps nor crowded cluster of
atoms should occur near the boundary of the box on the directions consistent with the periodicity,
and conversely, sufficient vacuum space on the non-periodic directions is crucial for eliminating
unwanted interactions across the boundary.
For surface slabs, two- or one-dimensional materials, or systems with vacuum, do not relax the
vacuum direction unless physically intended; it is better to constrain the appropriate cell
components or use a fixed-cell optimization after choosing the desired cell.
For variable-cell optimizations, a conventional cell representation with A along the X-axis and B on
the XY plane is expected, while non-standard cell orientations can usually still work.
```{warning}
An experimental `cif` should not be used as the initial structure without validation; refinement and
careful revision is necessary to make it "computation-ready" during pre-processing:
**Do not use experimental structure blindly.**
- Watch out for crystallographic disorder and fractional occupation of atoms: using the
superposition of all atoms as if every occupancy is 1.00 is a frequent pitfall as it introduces
contacting or even overlapping atoms.
- Beware of composition: the atomic structure may not match the intended charge-neutral chemical
formula, especially hydrogen atoms, which may be missing or duplicated due to disorder.
If available, **computational** materials databases are the most recommended avenue
to obtain structures that are "computation-ready", or even better, already optimized
with some electronic-structure methods. On the other hand, structures that are from
**experimental** characterization are frequently not "computation-ready", and thus
should not be subject to optimization without careful validation in pre-processing.
This can be prominent for `cif` and `pdb` structures determined by powder or single-
crystal XRD which can be affected by sample quality and thermal motion.
- Watch out for crystallographic disorder and atoms with low resolution or fractional
occupation: using the superposition of all atoms as if every occupancy is 1.00 is
highly likely to introduce contacting or even overlapping atoms.
- Beware of composition: the atomic structure may not match the intended macroscopic,
charge-neutral chemical formula, owing to missing or duplicated hydrogen atoms,
small counter ions, solvent or ligand molecules.
Possible resolutions vary from simple manual editing in the modelling stage, to
utilization of supercells and enumeration of special quasirandom structures (common
for materials with dopants), and to more rigorous XRD refinement and application of
quantum crystallography methods.
```
```{note}
For variable-cell optimizations, the [CELL_OPT](#CP2K_INPUT.MOTION.CELL_OPT) section
clarifies that a cell representation with the first vector A along the X-axis and
the second vector B on the XY plane is expected. This is to ensure that the cell
vectors are updated correctly in response to the external pressure in the current
implementation; even when an unconventional vector definition could sometimes feel
easy (for instance the primitive rhombohedral cell of the face-centered cubic lattice),
an appropriate linear transformation to the cell and the atomic coordinates should
still be applied to create the model for input. For more details, see online
[github discussion](https://github.com/cp2k/cp2k/issues/3384).
```
## Choosing an optimizer
@ -184,19 +244,17 @@ preliminary optimization with looser settings before tightening the convergence
The most important atomic convergence criteria are
[MAX_FORCE](#CP2K_INPUT.MOTION.GEO_OPT.MAX_FORCE),
[RMS_FORCE](#CP2K_INPUT.MOTION.GEO_OPT.RMS_FORCE), [MAX_DR](#CP2K_INPUT.MOTION.GEO_OPT.MAX_DR), and
[RMS_DR](#CP2K_INPUT.MOTION.GEO_OPT.RMS_DR). Force criteria are given in Hartree/Bohr, while
displacement criteria are given in Bohr. For `CELL_OPT`, the pressure target and pressure
convergence are controlled by [EXTERNAL_PRESSURE](#CP2K_INPUT.MOTION.CELL_OPT.EXTERNAL_PRESSURE) and
[PRESSURE_TOLERANCE](#CP2K_INPUT.MOTION.CELL_OPT.PRESSURE_TOLERANCE). The optimization is considered
converged only when all active criteria are satisfied.
[RMS_DR](#CP2K_INPUT.MOTION.GEO_OPT.RMS_DR). Force criteria are given in unit of Hartree/Bohr, while
displacement criteria are given in Bohr. In addition, for `CELL_OPT`, the pressure target and
pressure convergence in bar are specified by
[EXTERNAL_PRESSURE](#CP2K_INPUT.MOTION.CELL_OPT.EXTERNAL_PRESSURE) and
[PRESSURE_TOLERANCE](#CP2K_INPUT.MOTION.CELL_OPT.PRESSURE_TOLERANCE) respectively.
```{note}
The convergence criteria of the LBFGS optimizer are not the same with that of the BFGS and CG
optimizers, so you may see a converged report without all criteria listed above meeting the
expectation.
```
A typical convergence report for a variable-cell optimization looks like:
A typical convergence report for a variable-cell optimization, as controlled by
[PROGRAM_RUN_INFO](#CP2K_INPUT.MOTION.CELL_OPT.PRINT.PROGRAM_RUN_INFO), looks like the following
where forces are reported as gradients. For a fixed-cell optimization, the
[PROGRAM_RUN_INFO](#CP2K_INPUT.MOTION.GEO_OPT.PRINT.PROGRAM_RUN_INFO) does not mention the pressure
but is otherwise the same.
```none
OPT| **************************************************************************
@ -232,7 +290,44 @@ A typical convergence report for a variable-cell optimization looks like:
OPT| **************************************************************************
```
In the output, forces are reported as gradients.
Satisfying all of these criteria within `MAX_ITER` steps is a sufficient condition for the
convergence of optimization. With the `BFGS` and `CG` optimizers, this is also the necessary
condition. With the `LBFGS` optimizer, however, this is not a necessary condition because `LBFGS`
has its own extra pair of criteria named
[WANTED_PROJ_GRADIENT](#CP2K_INPUT.MOTION.GEO_OPT.LBFGS.WANTED_PROJ_GRADIENT) and
[WANTED_REL_F_ERROR](#CP2K_INPUT.MOTION.GEO_OPT.LBFGS.WANTED_REL_F_ERROR) which may override the
aforementioned general ones. It is possible to see an optimization task with the `LBFGS` optimizer
finish with the message below, even when one or more of the general criteria have not been met yet;
try restarting the optimization task with tightened criteria or alternative optimizers until
satisfactory convergence.
```none
***********************************************
* Specific L-BFGS convergence criteria
* WANTED_PROJ_GRADIENT and WANTED_REL_F_ERROR
* satisfied .... run CONVERGED!
***********************************************
```
The `LBFGS` optimizer has its own [PRINT_LEVEL](#CP2K_INPUT.MOTION.GEO_OPT.LBFGS.PRINT_LEVEL) which
controls the verbosity of details printed to the main output. In addition, a separate `iterate.dat`
file is also created by `LBFGS` optimizer where the columns `projg` and `f` can be checked against
the `WANTED_PROJ_GRADIENT` and `WANTED_REL_F_ERROR` criteria. Omitting some headers, the main body
of an `iterate.dat` excerpt is reproduced below.
```none
it nf nseg nact sub itls stepl tstep projg f
0 1 - - - - - - 4.632D-02 -1.082D+00
1 2 1 0 --- 0 1.0D+00 6.6D-02 5.016D-02 -1.087D+00
2 4 1 0 --- 1 1.9D+00 1.3D-01 5.707D-02 -1.097D+00
3 6 1 0 --- 1 1.7D+00 1.3D-01 6.212D-02 -1.108D+00
4 8 1 0 --- 1 1.5D+00 1.3D-01 6.413D-02 -1.120D+00
5 10 1 0 --- 1 1.5D+00 1.3D-01 6.112D-02 -1.132D+00
6 11 3 2 con 0 1.0D+00 1.3D-01 4.979D-02 -1.143D+00
7 12 3 2 con 0 1.0D+00 1.3D-01 2.440D-02 -1.150D+00
8 14 1 0 con 1 5.1D-01 6.6D-02 3.685D-03 -1.152D+00
9 15 1 0 con 0 1.0D+00 1.2D-02 7.569D-04 -1.152D+00
```
Geometry and cell optimization normally try to lower the energy, but convergence is determined by
the active force, displacement, and pressure criteria rather than by the total energy alone.
@ -241,6 +336,15 @@ is usually acceptable if all active convergence criteria are satisfied. In some
may indicate that the structure is not a real minimum, which can be checked by vibrational analysis
where imaginary frequencies reveal unstable modes.
Generally speaking, the usual convergence criteria may be satisfied within dozens or a few hundreds
of iterations. The default `MAX_ITER` value of 200 is mostly sufficient and may be raised to about
300 ~ 400 to deal with difficult cases like massive flexible structures. Running a single-point
calculation with `RUN_TYPE ENERGY_FORCE` beforehand, and then multiplying the elapsed time by
expected iterations, provides an estimate for the total time cost of an optimization task. To
monitor a long optimization, do not focus on these numeric criteria alone; it is recommended to
download the geometry (trajectory) file and observe the evolution of the structure and cell in a
visualizer program occasionally, and if anything goes wrong, halt the program in time.
## Constraints, cell degrees of freedom, and symmetry
Atomic constraints are defined in [MOTION/CONSTRAINT](#CP2K_INPUT.MOTION.CONSTRAINT). A common
@ -264,7 +368,9 @@ contains atom indices in the order used in [COORD](#CP2K_INPUT.FORCE_EVAL.SUBSYS
are useful, for example, for fixing the bottom layers of a slab, holding part of a support
structure, or imposing a chemically intended restriction. Check them carefully before production
calculations: accidental constraints are a common reason for apparently converged but physically
wrong structures.
wrong structures. For example, in a model of surface or micropore adsorption (physisorption or
chemisorption), the substrate should be allowed to relax locally rather than fully fixed to account
for the realistic interaction between the substrate (adsorbent) and the adsorbate.
For `CELL_OPT`, also decide which cell degrees of freedom should be optimized. Useful keywords
include [CONSTRAINT](#CP2K_INPUT.MOTION.CELL_OPT.CONSTRAINT) for fixing selected cell components,
@ -280,25 +386,58 @@ symmetry; otherwise they can prevent the structure from relaxing to a lower-symm
periodic fixed-cell `GEO_OPT`, see [KEEP_SPACE_GROUP](#CP2K_INPUT.MOTION.GEO_OPT.KEEP_SPACE_GROUP)
and related symmetry keywords.
```{note}
The `KEEP_SYMMETRY` keyword should always be used in conjunction with the cell symmetry
specified by [FORCE_EVAL/SUBSYS/CELL/SYMMETRY](#CP2K_INPUT.FORCE_EVAL.SUBSYS.CELL.SYMMETRY).
These keywords act only on the cell vectors and their gradients in a variable-cell
optimization, and they both assume a conventional cell where A is along the X-axis and
B is on the XY plane as noted in "Starting structure and cell" above.
To allow for the use of `KEEP_SPACE_GROUP`, the `spglib` library should be installed
and linked to the CP2K build in order to detect and preserve the space group. Usually
`KEEP_SPACE_GROUP` is used together with `KEEP_SYMMETRY`.
```
## Electronic-structure settings
Density functional theory (DFT) is an electronic-structure (wavefunction) method routinely used for
optimization. This section elaborate on the relevant aspects, assuming basic knowledge about the
method which can be found at [](../dft/index.md).
### SCF quality
For ab-initio geometry or cell optimization, every optimization step requires an
electronic-structure calculation. The force quality depends on the SCF convergence, grid settings,
basis set, pseudopotentials, and method-specific thresholds. For `CELL_OPT`, the stress tensor and
pressure must also be accurate enough for reliable cell updates. Important Quickstep settings
include [EPS_DEFAULT](#CP2K_INPUT.FORCE_EVAL.DFT.QS.EPS_DEFAULT),
[CUTOFF](#CP2K_INPUT.FORCE_EVAL.DFT.MGRID.CUTOFF),
[REL_CUTOFF](#CP2K_INPUT.FORCE_EVAL.DFT.MGRID.REL_CUTOFF), and
[EPS_SCF](#CP2K_INPUT.FORCE_EVAL.DFT.SCF.EPS_SCF).
Every optimization step requires an electronic-structure calculation, and the force quality depends
on the SCF convergence, grid settings, basis set, pseudopotentials, and method-specific thresholds.
For `CELL_OPT`, the stress tensor and pressure must also be accurate enough for reliable cell
updates. (The SCF convergence for each iteration is not to be confused with the convergence of
geometry in terms of stepsize, gradient and pressure across iterations in the optimization process.)
Optimization is driven by forces, and `CELL_OPT` is additionally driven by stress, rather than by
total energies alone. Therefore, simply choosing `EPS_SCF`, `CUTOFF`, `REL_CUTOFF`, or `EPS_DEFAULT`
from a final static-energy test may not be reliable. These settings usually need to be at least as
strict as, and often stricter than, those used for routine single-point energy calculations, because
noisy or poorly converged forces and stresses can slow down the optimization, cause oscillations,
produce line-search failures, or lead to an unreliable structure or cell.
total energies alone. Therefore, simply adjusting settings and tuning parameters from a final
static-energy test may not be sufficient. They usually need to be at least as strict as, and often
stricter than, those used for routine single-point energy calculations, because noisy or poorly
converged forces and stresses can slow down the optimization, cause oscillations, produce
line-search failures, or lead to an unreliable structure or cell.
Important Quickstep settings include [CUTOFF](#CP2K_INPUT.FORCE_EVAL.DFT.MGRID.CUTOFF),
[REL_CUTOFF](#CP2K_INPUT.FORCE_EVAL.DFT.MGRID.REL_CUTOFF),
[EPS_DEFAULT](#CP2K_INPUT.FORCE_EVAL.DFT.QS.EPS_DEFAULT), and
[EPS_SCF](#CP2K_INPUT.FORCE_EVAL.DFT.SCF.EPS_SCF).
```{note}
The mesh spacing of the plane-wave and real-space grid for electronic-structure
calculations is determined both by the `CUTOFF` setting and by the dimensions of
the cell. A variable-cell optimization may witness significant variations of the
cell, which will also affect the number of grid points and introduce artificial
discontinuities to the energy, force, stress and other properties even if the
change in structure is small.
It is recommended for better stability to set up a reference cell in the section
[FORCE_EVAL/SUBSYS/CELL/CELL_REF](#CP2K_INPUT.FORCE_EVAL.SUBSYS.CELL.CELL_REF),
with the dimension of the reference cell larger than the expected upper bound within
reach during the whole variable-cell optimization process, which will be kept constant
so that the number of grid points is held fixed and the change in energy more smooth.
```
### Extrapolation
@ -328,31 +467,71 @@ this warning appears together with abnormal changes in the geometry, cell, total
stress, or SCF behaviour, stop using that density-matrix-based method and switch to a safer
alternative.
```{note}
The applicability of extrapolation schemes depends on versions. To exemplify, extrapolation
was formerly available for gamma-only calculations, but the latest versions of CP2K has
started supporting extrapolation with full k-point sampling. However, when k-points are
requested to be reduced by symmetry, the dynamic refreshing of symmetry of k-points alongside
geometry updates could make extrapolation unavailable and automatically fall back to `USE_GUESS`.
```
```{danger}
**Be responsible, and do not ignore SCF convergence failure blindly.**
By default, a failure in SCF convergence aborts the program with a message like:
`SCF run NOT converged. To continue the calculation regardless, please set the keyword
IGNORE_CONVERGENCE_FAILURE.` Setting
[IGNORE_CONVERGENCE_FAILURE](#CP2K_INPUT.FORCE_EVAL.DFT.SCF.IGNORE_CONVERGENCE_FAILURE)
to `.TRUE.` turns it into a warning `SCF run NOT converged` that allows for optimization
to continue. However, bad SCF convergence leads to unreliable energy, force, and stress
on the current step, introducing error to the updated structure and the extrapolated
wavefunction on the next step. An optimization process with most or all steps failing
to converge SCF cycles does not yield trustworthy results in the end.
Therefore, `IGNORE_CONVERGENCE_FAILURE` should only be considered as a **last resort**
out of desperation, rather than a universal remedy used on a regular basis or even as
the default. There are dozens of measures available towards SCF convergence on top of
a realistic, chemically sensible model structure and appropriate net charge, spin
multiplicity, atomic magnetization, and k-point sampling; please try achieving SCF
convergence on the starting structure in a single-point calculation in the first place.
```
## Output files and restarts
The main output file contains the optimization progress and convergence checks. CP2K also writes
geometry and restart files whose names are based on [PROJECT_NAME](#CP2K_INPUT.GLOBAL.PROJECT_NAME)
(or `PROJECT`):
The main output file contains the optimization progress and convergence checks. Geometry and restart
files are also written with names based on the global
[PROJECT_NAME](#CP2K_INPUT.GLOBAL.PROJECT_NAME) (or `PROJECT`) as well as relevant printkeys under
the [MOTION/PRINT](#CP2K_INPUT.MOTION.PRINT) section.
- `project-pos-1.xyz` contains the sequence of geometries visited during the optimization. The last
frame is the latest optimized geometry.
- Cell information and stress output, when requested, are controlled by print keys such as
[MOTION/PRINT/CELL](#CP2K_INPUT.MOTION.PRINT.CELL) and
[MOTION/PRINT/STRESS](#CP2K_INPUT.MOTION.PRINT.STRESS).
- `project-FINAL-1_{cycles}.cif` is a CIF file containing the cell information of the final
structure.
- `project-FINAL-1_{cycles}.xyz` is an extented XYZ file containing the cell information of final
structure
- `project-pos-1.xyz` contains the sequence of geometries visited during the optimization. This is
controlled by the [TRAJECTORY](#CP2K_INPUT.MOTION.PRINT.TRAJECTORY) section with the default
format of XYZ (or XMOL). Note that the original XYZ specification does not convey periodicity, and
not all of the visualizers can parse the additional cell information in the comment line of XYZ
(e.g. in extended XYZ specification); switching to other formats using the `FORMAT` keyword may be
needed for visualization.
- `project-FINAL-1_{iter}.cif` and `project-FINAL-1_{iter}.xyz` is a CIF and an extended XYZ file
for the final structure, where `{iter}` is the number of iterations (steps) during optimization
until reaching convergence or hitting `MAX_ITER`. These are controlled by
[FINAL_STRUCTURE](#CP2K_INPUT.MOTION.PRINT.FINAL_STRUCTURE).
- `project-1.cell` contains the cell vectors during optimization, as controlled by
[CELL](#CP2K_INPUT.MOTION.PRINT.CELL).
- `project-1.stress` contains the stress tensors during optimization, as controlled by
[STRESS](#CP2K_INPUT.MOTION.PRINT.STRESS).
- `project-1.restart` is a CP2K restart input containing the latest geometry, cell, and relevant
settings.
- `project-1.restart.bak-*` are backup restart files from previous optimization steps.
settings, as controlled by [RESTART](#CP2K_INPUT.MOTION.PRINT.RESTART). When `BACKUP_COPIES` under
this section is non-zero, backup restart files from previous optimization steps will also be
preserved as `project-1.restart.bak-*`.
If a job stops before convergence, continue from the latest restart file, for example:
The restart file can be used directly as a new input file like
```none
cp2k -i project-1.restart -o project-restart.out
```
Alternatively, the [EXT_RESTART](#CP2K_INPUT.EXT_RESTART) section can be used for finer restart
controls. These will be convenient in case an optimization does not converge and a new optimization
starting from the latest structure is needed.
For expensive electronic-structure calculations, wavefunction restart files can reduce the cost of
continuation or follow-up calculations; see
[FORCE_EVAL/DFT/SCF/PRINT/RESTART](#CP2K_INPUT.FORCE_EVAL.DFT.SCF.PRINT.RESTART). For BFGS
@ -361,7 +540,7 @@ optimization driver.
## Practical checklist
Before trusting an optimized structure or cell, it's suggested to check that:
Before trusting an optimized structure or cell, it is suggested to check that:
- All active force and displacement convergence criteria are satisfied;
- For `CELL_OPT`, the pressure criterion is satisfied and the final cell is physically sensible;