howto:cdft
Differences
This shows you the differences between two versions of the page.
Both sides previous revisionPrevious revisionNext revision | Previous revisionNext revisionBoth sides next revision | ||
howto:cdft [2018/07/24 05:38] – Update to CP2K 6.1 nholmber | howto:cdft [2018/11/02 13:11] – [Structure of input file] nholmber | ||
---|---|---|---|
Line 21: | Line 21: | ||
A more exhaustive list of potential applications has been presented in this [[doi> | A more exhaustive list of potential applications has been presented in this [[doi> | ||
- | The charge and spin localized states are created by enforcing electron and spin density localization within atom centered regions of space. The relevant theory has been derived by Wu and Van Voorhis in a series of key papers: [[doi> | + | The charge and spin localized states are created by enforcing electron and spin density localization within atom centered regions of space. The relevant theory has been derived by Wu and Van Voorhis in a series of key papers: [[doi> |
- | In this tutorial, only the main results | + | In this tutorial, only the main theoretical aspects |
\begin{equation} | \begin{equation} | ||
Line 35: | Line 35: | ||
\end{equation} | \end{equation} | ||
- | where $c_j$ are atomic coefficients which determine how each atom is included in the constraint (more on this later), $P_j$ is the so-called cell function which determines the volume occupied by atom $j$ according to some population analysis method, and $\mathcal{N}$ is the set of all atoms in a system. Currently, only the Becke partitioning scheme is fully supported in CP2K, which will be elaborated in the following section. Different types of constraints can be constructed by modifying the weight function according to the following conventions | + | where $c_j$ are atomic coefficients which determine how each atom is included in the constraint (more on this later), $P_j$ is the so-called cell function which determines the volume occupied by atom $j$ according to some population analysis method, and $\mathcal{N}$ is the set of all atoms in a system. Different types of constraints can be constructed by modifying the weight function according to the following conventions |
* charge density constraint ($\rho^\uparrow + \rho^\downarrow$): | * charge density constraint ($\rho^\uparrow + \rho^\downarrow$): | ||
* magnetization density constraint ($\rho^\uparrow - \rho^\downarrow$): | * magnetization density constraint ($\rho^\uparrow - \rho^\downarrow$): | ||
* spin specific constraint ($\rho^{\uparrow/ | * spin specific constraint ($\rho^{\uparrow/ | ||
- | When CDFT is used in a molecular dynamics or a geometry optimization simulation, additional force terms arising from the constraints are calculated | + | The Becke and Hirshfeld space partitioning schemes can be used as constraint weight functions in CP2K. The main differences between these two constraints will be explained in a subsequent section. Please note that Becke constraints have been tested much more extensively. |
+ | |||
+ | When CDFT is used in a molecular dynamics or a geometry optimization simulation, additional force terms arising from the constraints are calculated | ||
\begin{equation} | \begin{equation} | ||
Line 75: | Line 77: | ||
===== Using the CDFT module ===== | ===== Using the CDFT module ===== | ||
- | The input sections | + | The input section |
==== Defining CDFT SCF parameters | ==== Defining CDFT SCF parameters | ||
Line 81: | Line 83: | ||
Settings for the CDFT SCF loop are controlled by the input section [[inp> | Settings for the CDFT SCF loop are controlled by the input section [[inp> | ||
- | <note important> | + | <note important> |
<code cp2k> | <code cp2k> | ||
&QS | &QS | ||
... | ... | ||
- | & | ||
- | ... | ||
- | &END BECKE_CONSTRAINT | ||
! CDFT loop settings | ! CDFT loop settings | ||
+ | ! Please note that prior to CP2K version 7.0, | ||
+ | ! Becke constraints were separate from the CDFT section | ||
&CDFT | &CDFT | ||
TYPE_OF_CONSTRAINT BECKE | TYPE_OF_CONSTRAINT BECKE | ||
+ | ! Compute CDFT charges? | ||
+ | ATOMIC_CHARGES | ||
+ | ! Constraint strength and target values | ||
+ | ! Give one value per constraint | ||
+ | STRENGTH | ||
+ | TARGET | ||
+ | ! Constraint definitions, | ||
+ | & | ||
+ | ATOMS 1 | ||
+ | COEFF 1 | ||
+ | CONSTRAINT_TYPE CHARGE | ||
+ | &END ATOM_GROUP | ||
+ | ! No constraint applied but calculate charges | ||
+ | & | ||
+ | ATOMS 2 | ||
+ | &END DUMMY_ATOMS | ||
+ | ! CDFT convergence and optimizer settings | ||
& | & | ||
- | TYPE BECKE_CONSTRAINT | + | TYPE CDFT_CONSTRAINT |
EXTRAPOLATION_ORDER 2 | EXTRAPOLATION_ORDER 2 | ||
MAX_SCF 10 | MAX_SCF 10 | ||
! Convergence threshold | ! Convergence threshold | ||
EPS_SCF 1.0E-3 | EPS_SCF 1.0E-3 | ||
- | ! Optimizer selection: | + | ! Optimizer selection: |
+ | ! Now Newton' | ||
OPTIMIZER NEWTON_LS | OPTIMIZER NEWTON_LS | ||
- | ! Optimizer step size | + | ! Optimizer |
STEP_SIZE -1.0 | STEP_SIZE -1.0 | ||
! Note that the section CDFT_OPT exists in CP2K version >= 6.1 | ! Note that the section CDFT_OPT exists in CP2K version >= 6.1 | ||
Line 109: | Line 128: | ||
CONTINUE_LS | CONTINUE_LS | ||
FACTOR_LS 0.5 | FACTOR_LS 0.5 | ||
- | ! Finite difference settings for calculation of Jacobian matrix | + | ! Finite difference settings for Jacobian matrix |
JACOBIAN_STEP 1.0E-2 | JACOBIAN_STEP 1.0E-2 | ||
JACOBIAN_FREQ 1 1 | JACOBIAN_FREQ 1 1 | ||
Line 115: | Line 134: | ||
JACOBIAN_RESTART FALSE | JACOBIAN_RESTART FALSE | ||
&END CDFT_OPT | &END CDFT_OPT | ||
- | & | + | &END |
+ | ! Settigs specific to Becke constraints | ||
+ | & | ||
+ | ... | ||
+ | &END BECKE_CONSTRAINT | ||
+ | ! Print information about CDFT calculation | ||
+ | & | ||
+ | &EACH | ||
+ | QS_SCF 1 | ||
+ | &END EACH | ||
+ | COMMON_ITERATION_LEVELS 2 | ||
+ | ADD_LAST NUMERIC | ||
+ | FILENAME ./${NAME} | ||
+ | &END PROGRAM_RUN_INFO | ||
&END CDFT | &END CDFT | ||
&END QS | &END QS | ||
</ | </ | ||
- | The structure of this input section is quite straightforward. The keyword [[inp> | + | The structure of this input section is quite straightforward |
+ | - Constraint definitions (type, which atoms to include, constraint target, etc) | ||
+ | - CDFT SCF loop settings (solver, convergence criterion, etc) | ||
+ | - Constraint weight function specific settings (Becke/ | ||
+ | |||
+ | In the above example, a Becke constraint is selected using the keyword [[inp> | ||
+ | |||
+ | {{ howto: | ||
+ | **Figure 2.** Using a fragment based CDFT constraint. The system is first divided into two fragments with atomic positions fixed in the same configuration as in the full system. The electron and spin densities of the fragment systems are then saved to cube files and subsequently used as input files for the CDFT calculation, | ||
+ | |||
+ | The OUTER_SCF section within the CDFT section defines settings for the CDFT SCF loop. The keyword [[inp> | ||
Above, for instance, the Jacobian is explicitly calculated every CDFT SCF iteration and MD step by perturbing each constraint Lagragian using a first order forward difference stencil with a step size of $10^{-2}$. The Newton step size is optimized with backtracking line search using the update formula $\alpha_n = 0.5*\alpha_{n-1}$ for a maximum of 5 steps as long as the CDFT constraint error decreases. | Above, for instance, the Jacobian is explicitly calculated every CDFT SCF iteration and MD step by perturbing each constraint Lagragian using a first order forward difference stencil with a step size of $10^{-2}$. The Newton step size is optimized with backtracking line search using the update formula $\alpha_n = 0.5*\alpha_{n-1}$ for a maximum of 5 steps as long as the CDFT constraint error decreases. | ||
- | ==== Defining | + | ==== Available |
- | As explained in the previous section, the CDFT module in CP2K currently | + | The CDFT module in CP2K currently supports using Becke or Hirshfeld (no forces) |
- | {{ howto: | ||
- | **Figure 2.** Comparison of the Voronoi (lines) and Becke partitioning (contours) schemes. At left, the Becke partitioning is performed without atomic size information. At right, the size of the red atom is 30 % larger than the black atoms, and the contours of the red atom extend farther than without atomic size adjustments. | ||
- | The Voronoi and, by extension, the Becke partitioning methods treat each element equally. This leads to unphysical partial charges in most systems. For example, the Becke scheme predicts a positive charge on oxygen and a negative charge on hydrogen in water (see examples for input files). This problem can be remedied by accounting for atomic radii during the partitioning. This behavior is activated by the keyword [[inp> | + | == Becke constraints == |
- | The algorithmic implementation of the Becke density partitioning method | + | The Becke density partitioning method |
- | The actual constraints are defined in the section [[inp> | + | {{ howto: |
+ | **Figure 3.** Comparison | ||
- | {{ howto: | + | The Voronoi and, by extension, the Becke partitioning methods treat each element equally. This leads to unphysical partial charges in most systems. For example, the Becke scheme predicts a positive charge on oxygen and a negative charge on hydrogen in water (see examples for input files). This problem can be remedied by accounting for atomic radii during the partitioning. This behavior |
- | **Figure 3.** Using a fragment based Becke constraint. The system | + | |
+ | The algorithmic implementation | ||
An example of a Becke constraint input section is given below. This choice of parameters should be reasonable for most systems, ignoring the atomic radii and constraint definitions which are system dependent. Decreasing the partitioning cutoff might be useful for solvated system MD simulations, | An example of a Becke constraint input section is given below. This choice of parameters should be reasonable for most systems, ignoring the atomic radii and constraint definitions which are system dependent. Decreasing the partitioning cutoff might be useful for solvated system MD simulations, | ||
<code cp2k> | <code cp2k> | ||
- | &QS | + | &CDFT |
... | ... | ||
& | & | ||
Line 149: | Line 191: | ||
ADJUST_SIZE | ADJUST_SIZE | ||
ATOMIC_RADII | ATOMIC_RADII | ||
- | ! Compute Becke charges? | ||
- | ATOMIC_CHARGES | ||
- | ! Constraint strength and target values | ||
- | ! Give one value per constraint | ||
- | STRENGTH | ||
- | TARGET | ||
! Cutoff scheme | ! Cutoff scheme | ||
CUTOFF_TYPE | CUTOFF_TYPE | ||
Line 166: | Line 202: | ||
IN_MEMORY | IN_MEMORY | ||
SHOULD_SKIP | SHOULD_SKIP | ||
- | ! Constraint definitions, | + | &END CDFT |
- | & | + | |
- | ATOMS 1 | + | |
- | COEFF 1 | + | |
- | CONSTRAINT_TYPE CHARGE | + | |
- | | + | |
- | ! No constraint applied but calculate charges | + | |
- | & | + | |
- | ATOMS 2 | + | |
- | &END DUMMY_ATOMS | + | |
- | ! Print information about CDFT calculation | + | |
- | & | + | |
- | &EACH | + | |
- | QS_SCF 1 | + | |
- | &END EACH | + | |
- | COMMON_ITERATION_LEVELS 2 | + | |
- | ADD_LAST NUMERIC | + | |
- | FILENAME ./${NAME} | + | |
- | &END PROGRAM_RUN_INFO | + | |
- | &END BECKE_CONSTRAINT | + | |
- | &END QS | + | |
</ | </ | ||
+ | |||
+ | == Hirshfeld constraints == | ||
+ | |||
+ | Hirshfeld constraints are cheaper to construct than Becke constraints in large systems because Hirshfeld constraints are essentially just weighted sums of spherical Gaussian functions. The keywords [[inp> | ||
+ | * The first choice implies that the CDFT weight function for each atom is a single Gaussian function whose radius is controlled by the GAUSSIAN_SHAPE keyword. By default, tabulated covalent radii are used as the radii of the Gaussian, but it is also possible to select van der Waals radii or to define custom radii. | ||
+ | * The latter choice implies that the atomic weight function are constructed from isolated atomic densities which are expanded in terms of multiple spherical Gaussians. This choice is affected by the selected basis set and pseudopotential but not by the Gaussian shape keyword. | ||
==== Selected examples | ==== Selected examples | ||
Line 293: | Line 315: | ||
- How much charge is transferred between the two water molecules according to the different constraint methods? Note that net charges for fragment based constraints are not reported with respect to the core charge, but they can be recovered by post-processing the reported absolute populations. Instead, each atomic charge is referenced to the number of electrons per atom in the system where the isolated densities are superimposed. | - How much charge is transferred between the two water molecules according to the different constraint methods? Note that net charges for fragment based constraints are not reported with respect to the core charge, but they can be recovered by post-processing the reported absolute populations. Instead, each atomic charge is referenced to the number of electrons per atom in the system where the isolated densities are superimposed. | ||
- Compare the calculated charge transfer energies to the [[doi> | - Compare the calculated charge transfer energies to the [[doi> | ||
+ | |||
+ | |||
+ | === Zn dimer cation with Hirshfeld constraints === | ||
+ | |||
+ | <note important> | ||
+ | |||
+ | This tutorial is exactly the same as the Zn dimer example above but using Hirshfeld partitioning based constraints instead of Becke constraints. You can find the input files {{: | ||
+ | |||
+ | It might be instructive to visualize how the Becke and Hirshfeld weight function schemes differ, in particular, how the methods assign a volume to each atom in the system. You can activate the section [[inp> | ||
Line 403: | Line 434: | ||
In the above example input file, a common file '' | In the above example input file, a common file '' | ||
- | The keyword [[inp> | + | The keyword [[inp> |
+ | |||
+ | The keyword [[inp> | ||
+ | |||
+ | A special run type is available for '' | ||
==== Example: Electronic coupling of Zn cation dimer ==== | ==== Example: Electronic coupling of Zn cation dimer ==== | ||
Line 450: | Line 485: | ||
Other files created during the execution are related to the individual CDFT states. The '' | Other files created during the execution are related to the individual CDFT states. The '' | ||
+ | |||
+ | |||
+ | ==== Example: Configuration interaction calculations with CDFT (CDFT-CI): The case of $\mathrm{H}_2^+$ ==== | ||
+ | |||
+ | DFT exchange-correlation functionals suffer from varying degrees of self-interaction error (SIE). A pathological example of a system where SIE leads to unphysical results with DFT is the simple dissociation reaction | ||
+ | |||
+ | $\mathrm{H}_2^+ \rightarrow \mathrm{H}^+ + \mathrm{H}$. | ||
+ | |||
+ | Even though this system contains only 1 electron, the dissociation profile obtained with PBE notably deviates from the exact Hartree-Fock profile as shown in Figure 4 below. | ||
+ | |||
+ | {{ howto: | ||
+ | ** Figure 4. ** Illustration of DFT self-interaction error for the reaction $\mathrm{H}_2^+ \rightarrow \mathrm{H}^+ + \mathrm{H}$. PBE notably deviates from the exact Hartree-Fock dissociation profile. The correct profile can be recovered with constrained DFT configuration interaction (CDFT-CI) using fragment constraint states $|\mathrm{H^+H}> | ||
+ | |||
+ | We can use CDFT states as the basis of a configuration interaction (CI) simulation to correct for SIE in this system. As the figure above shows, CDFT-CI using the PBE functional is able to reproduce the exact dissociation profile. You can read up on the theory behind CDFT-CI simulations from the references given at the start of this tutorial. Very briefly, CDFT-CI simulations involve representing the system' | ||
+ | |||
+ | In this tutorial, you will reproduce the DFT and CDFT results from the figure above. You can find the input files {{: | ||
+ | |||
+ | - Start by examining the simulation script '' | ||
+ | - While the CDFT simulations are running, look at the results from the DFT simulations with PBE. Can you figure out the reason why PBE predicts an unphysical dissociation profile? (Hint. Compute the partial charges). | ||
+ | - Inspect the output files produced by the CDFT-CI simulations once they are done. Find the output from the CDFT-CI module in the main output files. Look at the CI expansion coefficients in terms of the CDFT states. How would you characterize the CI wavefunction? | ||
+ | - Plot the CDFT-CI and DFT dissociation profiles with your favorite plotting tool. Use the Hartree-Fock data from the provided data file as a reference. You can produce a similar data file from your simulations by calling '' |
howto/cdft.txt · Last modified: 2024/01/03 13:20 by oschuett