Requests for technical support from the VASP group should be posted in the VASP-forum.

# LCHIMAG

LCHIMAG = .TRUE. | .FALSE.
Default: LCHIMAG = .FALSE.

Description: calculate the chemical shifts by means of linear response.

For LCHIMAG=.TRUE., VASP calculates the chemical shift tensors.

The chemical shielding tensor is defined as:

${\displaystyle \sigma (\mathbf {R} )_{ij}=-{\frac {\partial B^{\mathrm {ind} }(\mathbf {R} )_{i}}{\partial B_{j}^{\mathrm {ext} }}}}$

Here ${\displaystyle {\mathbf {R}}}$ denotes the atomic nuclear site, ${\displaystyle i}$ and ${\displaystyle j}$ denote cartesian indices, ${\displaystyle B^{\mathrm {ext} }}$ an applied DC external magnetic field and ${\displaystyle B^{\mathrm {ind} }(\mathbf {R} )}$ the induced magnetic field at the nucleus. NMR experiments yield information on the symmetric part of the tensor. Typical NMR experiments yield information on the shielding relative to that of a reference compound:

${\displaystyle \delta (\mathbf {R} )_{ij}=\sigma _{ij}^{\mathrm {ref} }-\sigma (\mathbf {R} )_{ij}}$

In this (approximate) formula ${\displaystyle \sigma _{ij}^{\mathrm {ref} }}$ is the isotropic shielding of the nucleous in the reference compound. ${\displaystyle \delta (\mathbf {R} )_{ij}}$ is the chemical shift tensor.

In VASP the chemical "shift" tensor is calculated as:

${\displaystyle \delta (\mathbf {R} )_{ij}\mathrm {[VASP]} ={\frac {\partial B^{\mathrm {ind} }(\mathbf {R} )_{i}}{\partial B_{j}^{\mathrm {ext} }}}}$

This is minus the shielding tensor. It is not the true chemical shift tensor. To convert it to the real shift tensor one should add the reference shielding:

${\displaystyle \delta (\mathbf {R} )_{ij}=\sigma _{ij}^{\mathrm {ref} }+\delta (\mathbf {R} )_{ij}\mathrm {[VASP]} }$

VASP can calculate chemical "shifts" for non-metallic crystalline systems using the linear response method of Yates, Pickard and Mauri.[1][2]

INPUT

A typical INCAR could look like this:

PREC = A               # nice
ENCUT = 600.0          # typically higher cutoffs than usual are needed
EDIFF = 1E-8           # you need much smaller EDIFFs than normal.
ISMEAR = 0; SIGMA = 0.1 # no fancy smearings, SIGMA sufficiently small
LREAL = A              # helps for speed for large systems, not necessary per se

LCHIMAG = .TRUE.       # to switch on linear response for chemical shifts
DQ = 0.001             # often the default is sufficient
ICHIBARE = 1           # often the default is sufficient
LNMR_SYM_RED = .TRUE.  # be on the safe side
NSLPLINE = .TRUE.      # only needed if LREAL is NOT set.


The first block of tags in the INCAR above expresses the fact that the calculations of chemical shifts by means of linear response often require a high accuracy (PREC=A, EDIFF≤1E-8, high ENCUT).

The chemical shifts are calculated from the orbital magnetic response under the assumption that the system is an insulator. It makes no sense to use smearing schemes intended for metals, indeed, doing so can generate nonsense. It is safe to use ISMEAR=0 and make SIGMA so small that states have no fractional occupancies.

The seçond block of tags switches on the calculation of the chemical shifts (LCHIMAG=.TRUE.), and controls several aspects of the finite difference k-space derivatives (Eqs. 38, 40, and 47 in the work of Yates et al.[2]):

• DQ is the step size for the finite difference k-space derivative. Typical values are in the range [0.001 - 0.003]. The default (DQ=0.001) is often sufficient.
• ICHIBARE is the order of the finite difference stencil used to calculate the magnetic susceptibility (second order derivative in Eq. 47 of Yates et al.[2]). ICHIBARE may be set to 1, 2, or 3. Often the default (ICHIBARE=1) is sufficient. A higher ICHIBARE results in a substantial increase of the computational load.
• For NLSPLINE=.TRUE., the PAW projectors in reciprocal space (LREAL=.FALSE.) are set up using a spline interpolation so that they are k-differentiable. This improves the susceptibility contribution to the shifts (via the aforementioned Eq. 47). It only slightly affects the other contributions to the shifts (Eqs. 38 and 40). It is advised to set NLSPLINE=.TRUE. if PAW projectors are applied in reciprocal space, but only in case of calculation of chemical shifts. As this option also gives slightly different total energies, it is advised to use the default NLSPLINE=.FALSE. in all other calculations for reasons of compatibility. Real space projectors are k-differentiable by construction, hence do not require to set NLSPLINE=.TRUE.
• The star on which the k-space derivative is calculated is oriented along the cartesian directions in k-space. If the symmetry operations in k-space do not map this star onto itself, erroneous results can be obtained. To have VASP check for such operations, set LNMR_SYM_RED=.TRUE., and such operations will be discarded, resulting in a larger IBZ. In case of any doubt set LNMR_SYM_RED=.TRUE. Beware: It matters how the real space lattice vectors are set up relative to the cartesian coordinates in POSCAR. It determines the orientation of the k-space star and hence can affect the efficiency via the number of k-points in the IBZ.

The chemical shift is calculated via the induced current.[1][2] It has contributions from the plane wave grid and one-center contributions (the induced field at the center of a PAW sphere due to the augmentation current inside that sphere). Two-center contributions (induced fields due to augmentation currents in other PAW spheres) are standard neglected. These contributions can be switched on using LLRAUG.

For very high accuracy calculations use LASPH.

No special POTCAR files are necessary. The GIPAW is applied using the projectors functions and partial waves that are stored in the regular POTCAR files. A few remarks, however, on accuracy in relation to the different POTCAR flavours:

• Results sensitively depend on the quality, i.e., completeness of the partial wave/projector function set in the energy range needed for good chemical transferability. Result obtained with different POTCAR flavours can differ a few ppm for first and second row sp-bonded elements (except for H).
• Use POTCAR files generated with a consistent exchange-correlation functional. The PAW reconstruction with AE partial waves is crucial as the field on the nucleus needs to be calculated. So avoid, if possible, overriding LEXCH from POTCAR with an explicit GGA-tag in the INCAR.

OUTPUT

At the end of the OUTCAR file, VASP prints the chemical shift tensors both before and after space group symmetrization. These are the absolute tensors for the infinite lattice, excluding core contributions. Next lines Q=0 CONTRIBUTION TO CHEMICAL SHIFT are printed. This is a shift tensor arising solely from the ${\displaystyle \mathbf {G=0} }$ component of the induced field. This component is related to the shape of the sample and depends only on the induced macroscopic surface currents (via the orbital magnetic susceptibility). It is printed for a spherical sample (for which is it nucleus independent), and calculated according to Eqs. 46-48 of Yates et al.[2], i.e., using the so-called pGv-approximation to the magnetic susceptibility. To obtain the full absolute tensor the contribution for ${\displaystyle \mathbf {G=0} }$ has to be added to the nuclear shifts. The approximate susceptibility itself is also printed.

Finally the tensor is processed and its (CSA) characteristics are printed on OUTCAR. The tensor is symmetrized (${\displaystyle \sigma _{ij}=\sigma _{ji}}$ is enforced) and diagonalized. From the three diagonal values the isotropic chemical "shift" ${\displaystyle \delta _{\mathrm {iso} }\mathrm {[VASP]} }$, span ${\displaystyle \Omega }$ and skew ${\displaystyle \kappa }$ are calculated and printed.[3] Note that ${\displaystyle \kappa }$ is ill-defined if ${\displaystyle \Omega =0}$. Note that the isotropic chemical shift ${\displaystyle \delta _{\mathrm {iso} }\mathrm {[VASP]} }$ (ISO_SHIFT) as printed is actually minus the isotropic shielding. To make it a real shift one should add the reference shielding. Also note that ${\displaystyle \Omega }$ (SPAN) and ${\displaystyle \kappa }$ (SKEW) are unambiguously defined.[3] Units are ppm, except for the skew. This typically looks like:

  ---------------------------------------------------------------------------------
CSA tensor (J. Mason, Solid State Nucl. Magn. Reson. 2, 285 (1993))
---------------------------------------------------------------------------------
EXCLUDING G=0 CONTRIBUTION             INCLUDING G=0 CONTRIBUTION
-----------------------------------   -----------------------------------
ATOM    ISO_SHIFT        SPAN        SKEW     ISO_SHIFT        SPAN        SKEW
---------------------------------------------------------------------------------
(absolute, valence only)
1    4598.8125      0.0000      0.0000     4589.9696      0.0000      0.0000
2     291.5486      0.0000      0.0000      282.7058      0.0000      0.0000
3     736.5979    344.8803      1.0000      727.7550    344.8803      1.0000
4     736.5979    344.8803      1.0000      727.7550    344.8803      1.0000
5     736.5979    344.8803      1.0000      727.7550    344.8803      1.0000
---------------------------------------------------------------------------------
(absolute, valence and core)
1   -6536.1417      0.0000      0.0000    -6547.9848      0.0000      0.0000
2   -5706.3864      0.0000      0.0000    -5718.2296      0.0000      0.0000
3   -2369.4015    344.8803      1.0000    -2381.2446    344.8803      1.0000
4   -2369.4015    344.8803      1.0000    -2381.2446    344.8803      1.0000
5   -2369.4015    344.8803      1.0000    -2381.2446    344.8803      1.0000
---------------------------------------------------------------------------------
IF SPAN.EQ.0, THEN SKEW IS ILL-DEFINED
---------------------------------------------------------------------------------



The columns excluding the ${\displaystyle \mathbf {G=0} }$ contribution are useful for supercell calculations on molecules. The columns including the ${\displaystyle \mathbf {G=0} }$ contribution are for crystals. The upper block gives the shielding due to only the electrons included in the SCF calculation. The lower block has the contributions due to the frozen PAW cores added. These core contributions are rigid.[4] They depend on POTCAR and are isotropic, i.e. affect neither SPAN nor SKEW.

Beware: the treatment of the orbital magnetism is non-relativistic. This is fine for light nuclei. The standard POTCARs are scalar-relativistic and account partially for relativistic effects. The accuracy can be improved using LBONE, that restores the small B-component of the wave function inside the PAW spheres. Spin-orbit coupling is not implemented for chemical shift calculations.

What to do in case of insufficient memory? VASP trades off memory savings against speed, opting for the latter. The response calculation is inherently parallel over k-points. This can be used to economize on memory: First do a regular self-consistent calculation at high accuracy for the full k-point mesh. Save the CHGCAR file. Next do a chemical shift calculation for each k-point in the IBZ separately, starting from CHGCAR, i.e., using ICHARG=11. Finally calculate the shifts as a k-point weighted average of the symmetrized shifts of the individual k-points.