3.8 Kinetic energy, scalar relativity, spin-orbit coupling, and full relativity

While the non-relativistic level of theory is exactly defined and will be the same in any first-principles implementation (at a complete basis set, all-electron level anyway), there are many different versions of scalar-relativistic approximations which can yield considerably different total energies for different systems. Their unifying feature is that any two scalar-relativistic methods should still yield the same energy differences for properties that concern valence electrons: Binding energies, valence eigenvalues, etc.

The recommended level of scalar relativity in FHI-aims is the so-called “atomic ZORA” approximation, as defined specifically in Equations (55) and (56) of Ref. [36]. It is important to refer to this specific definition since there are other variants of ZORA (“zero-order regular approximation”) in the literature and in other codes, including variants also called “atomic ZORA” but following a different mathematical definition.

The keyword needed to use this level of theory is

relativistic atomic_zora scalar

That’s it.

The “atomic ZORA” level of theory as implemented in FHI-aims has held up extremely well in large, high-accuracy benchmarks of scalar-relativistic total-energy based properties [195] as well as energy band structures [148]. It works for the right mathematical reasons. It can, in principle, be used across the entire periodic table (there should be no need to resort to non-relativistic calculations except for benchmarking purposes).

In addition, a non-selfconsistent treatment of spin-orbit coupling for band structures, densities of states, absorption properties and for the independent-particle dielectric response is also available and can be used in addition to (on top of) scalar relativistic calculations using the atomic ZORA. This is described in detail in Ref. [148], including a simple discussion of relativistic treatments in general and of how the “atomic ZORA” and the spin-orbit coupling formalism on top of it are related.

The keyword to add post-scf spin-orbit coupling is

include_spin_orbit

That’s it. Note that this keyword can be used as a followup to both non-spinpolarized and spin-polarized scalar-relativistic calculations.

An important fact to keep in mind is that a scalar-relativistic calculation yields two distinct spin sets of spin states, one for each spin channel. However, after the perturbative spin-orbit coupling treatment, only a single set of states emerges as output, since spin-orbit coupling mixes the scalar-relativistic spin states and the spin channels are no longer distinct. Thus, the output files for any quantities derived from spin-orbit coupled calculations (densities of states, band structures, etc.) are not and cannot be printed as separate spin channels - only one set of files is written that includes states derived fron both former spin channels.

More details follow below, but here are three additional important point:

• •

  Never mix results from different scalar relativistic treatments in total-energy differences. Absolute total-energy differences between different relativistic treatments can be very large because the deep core state energies change.

• •

  The absolute core level energies in the “atomic ZORA” approximation are far away from measured core level energies that would appear in experiment or in the actual Dirac equation. However, the relative core level shifts (differences) between different chemical systems are still reliable.

• •

  FHI-aims also includes another relativistic treatment called “scaled ZORA” but this seems to be slightly less accurate and does not have support for forces or stresses or any other use cases. We do not recommend to use “scaled ZORA” any more (“atomic ZORA” simply seems to do the better job).

Spin-orbit coupling (SOC) is a simple consequence of transforming Dirac’s equation to a (two-component) Schrödinger-like form. This leads to an approximate “spin-orbit-coupled” Hamiltonian of the form

$$\widehat{H}={\widehat{t}}_{SR}+\widehat{v}+{\widehat{v}}_{SOC},$$

(3.36)

where ${\widehat{t}}_{SR}$ is the usual scalar relativistic kinetic energy operator (e.g., atomic ZORA), $\widehat{v}$ is the local or non-local potential, and ${\widehat{v}}_{SOC}$ is the spin-orbit coupling operator,

$$v_{SOC}=\frac{i}{4c^2}𝝈\cdot 𝒑v\times 𝒑.$$

(3.37)

FHI-aims currently implements a treatment of spin-orbit coupling which adds spin-orbit coupling corrections to the Kohn-Sham eigenvalues, band structures, and densities of states in a single evaluation after the scalar-relativistic s.c.f. cycle has converged. This means that it is a post-processed implementation of spin-orbit coupling. It works in the Hilbert space of calculated scalar-relativistic eigenstates, as opposed to the “full” space spanned by the computational basis set, to dramatically reduce the problem size. This is known as the “second-variational” method. It only calculates and diagonalizes the spin-orbit-coupled Hamiltonian once; therefore, the resulting spin-orbit-coupled eigenstates are non-self-consistent.

Full details on the implementation of spin-orbit coupling in FHI-aims, as well as a derivation of the spin-orbit-coupled Hamiltonian from the Dirac equation and a detailed benchmark of the effect of spin-orbit coupling on band structures, may be found in Ref. [148] When publishing results using spin-orbit coupling in FHI-aims, please remember to cite this reference.

Applying the SOC operator as a correction to scalar-relativistic eigenvectors is quantitatively accurate (to a few 0.01 eV for valence band structures) for elements below Xe ($Z$=54) when combined with atomic ZORA. For heavy elements (approximately Au and beyond) this level of theory is only qualitatively accurate. It captures the majority of the SOC effect, but quantitative deviations above 0.1 eV for band structures must be expected. Similarly, the corrections for any core levels would require one to go beyond non-self-consistent SOC.

Since this implementation of spin-orbit coupling operates in the Hilbert space spanned by the calculated scalar-relativistic eigenvectors, for accurate high-lying bands one must include sufficiently many unoccupied states. This may be done by setting empty_states to a higher value or, if you are feeling particularly paranoid, setting the
calculate_all_eigenstates keyword to include all possible eigenstates. It is the opinion of the authors that this is only relevant for materials containing Au and heavier elements.

As of FHI-aims version 20231212, empty_states are automatically increased by:

n_empty_states = n_empty_states + 10*n_atoms

by default for spin-orbit coupled calculations (i.e. when include_spin_orbit is enabled). As SOC is evaluated in the basis of the available scalar-relativistic KS states, the number of states included in the scalar-relativistic SCF sycle is smaller than the full basis size leading to inadequate numerical convergence. This effect is more noticeable in hybrid systems. This does make SOC calculations more expensive, however, it will lead to better numerically converged results. If a user wishes to override this default setting, they could manually set empty_states as desired.

Note: While FHI-aims also prints out a corrected total-energy expression based on the SOC-corrected eigenvalues, do not use this value. It is experimental.

Note: This will set the empty_states keyword to a higher value automatically.

This flag is presently kept for test purposes only (the electronic kinetic energy is separately computed and printed for each scf iteration anyway) but may be useful for some future modifications.

For example, this will allow you to run a physically incorrect calculation of heavy elements (think Au) with Schrödinger’s expression for the kinetic energy, instead of a scalar relativistic treatment. The results will be wrong, so this flag should only be set for test purposes. When set, the code assumes that the user must know what they are doing.

Detailed expressions for the scalar relativistic treatments available here are given in Ref. [36]. We here only repeat the salient options and expressions. Possible options for r-type are:

• •

  none : Non-relativistic kinetic energy. In this case, s-type and threshold need not be provided.

• •

  atomic_zora : Atomic ZORA approximation as described in Ref. [36]. threshold need not be provided. This is the currently recommended option for energy differences and valence and unoccupied eigenvalues.

• •

  zora The ZORA approximation is used throughout the s.c.f. cycle, followed by a “scaled ZORA” [303] post-processing step (rescaling of all eigenvalues). WARNING: Do not rely on intermediate, simple ZORA total energies, but only on the final, rescaled total energies instead! ZORA (unscaled) is not the same as “atomic ZORA” and cannot be trusted. We also no longer recommend scaled ZORA values since there is no clear advantage. Just use atomic_zora unless there is need to do otherwise.

Forces are only provided for none and atomic_zora. The default is atomic_zora scalar.

Remember to never take energy differences between calculations performed with different “relativistic” settings.

We recommend to simply use atomic_zora for all calculations, unless there is a particular need to stay with relativistic none.

If you really do want to use scaled ZORA (the case of zora keyword), the threshold option is required. It specifies the threshold value above which the difference between the sum-of-free-atoms ZORA expression and that for the actual potential during the s.c.f. cycle will be calculated. In areas of shallow potentials, where both expressions are substantially similar, this saves the extra integration effort associated with ZORA. For (very!) safe settings, threshold may be set to 10^-12; in our experience, also 10^-9 does not lead to any noticeable accuracy loss.

For 1120, the setting none will be accepted, but a warning will be issued. For $Z$>20, choosing none will cause the code to stop with an error message in order to avoid accidental calculations with incorrect relativity. If a non-relativistic calculation is still desired, for example for test purposes, this “stop” can be disabled by setting the flag override_relativity—but use this only if you know what you are doing.

This flag is presently based on the perturbative spin-orbit coupling method. Therefore, you should in the same time set:

relativistic atomic_zora scalar

include_spin_orbit

Since the use_spin_texture keyword is computed on bands, you should also set:

output band

or

output band_mulliken

If you wish to compute the spin texture of the S.C.F. k-grid, please see the keyword use_spin_texture_scf.

Spin texture is defined as the expectation value of the vector of Pauli matrices. The complete formulae we use, with examples, are given in the methods section and supporting information of Ref. [154]. By setting this flag, the code will calculate the expectation values for three components (${\sigma }_x$, ${\sigma }_y$, and ${\sigma }_z$), which can help understand the spin polarization behaviour of individual bands. The values are printed out in separate files named spin_texture.out and spin_texture.dat which contain the same kind of data but simply in different formats.

Please note that for degenerate eigenstates on symmetric k-points, the value of the spin texture will be meaningless, since for these states the eigenvectors are arbitrary. Because the eigenvectors are arbitrary, different eigensolvers will output different eigenvectors, and the eigenvectors will also be different on different hardware or for calculations with different numbers of tasks. As of version 250304 the spin_texture.dat additionally contains the k-points in relative coordinates of the reciprocal basis vectors (i.e. how they are specified in an output band command) to allow more straightforward determination of the symmetric k-points. The symmetric k-points are those which correspond to relative fractional k-vectors of the form [0 0 $\alpha$], [0 $\alpha$ $\alpha$], or [$\alpha$ $\alpha$ $\alpha$] for some $\alpha$.

As of version 250224 use_spin_texture can be used in combination with the ScaLAPACK eigensolver (KS_method parallel). Previously, only use of the LAPACK eigensolver (KS_method serial) in combination with use_spin_texture was possible.

Like the use_spin_texture keyword, use_spin_texture_scf is presently based on the perturbative spin-orbit coupling method. Therefore, you should in the same time set:

relativistic atomic_zora scalar

include_spin_orbit

Also analogous to the use_spin_texture keyword, the three components (${\sigma }_x$, ${\sigma }_y$, and ${\sigma }_z$), are printed out in separate files named spin_texture_scf.out and spin_texture_scf.dat.

As of version 250304, use_spin_texture_scf can be used in combination with the ScaLAPACK eigensolver (KS_method parallel).