Spin-spiral calculations

Spin-spiral calculations model continuously rotating magnetic structures without the need for large supercells. This page provides step-by-step instructions for setting up, running, and analyzing spin-spiral calculations. For the underlying formalism (the generalized Bloch condition, the modified Hamiltonian, and the basis-set requirements) see Spin spirals.

Prerequisites

Spin spirals are a noncollinear-magnetism feature. You need the vasp_ncl executable that supports noncollinear magnetism (LNONCOLLINEAR = True), and it is helpful to be familiar with collinear magnetic calculations, the MAGMOM tag, and Setting up an electronic minimization before starting.

Step-by-step instructions

Step 1: Start from a converged magnetic state

Begin with a converged collinear or noncollinear ground-state calculation of the primitive cell, following Setting up an electronic minimization, and keep the resulting CHGCAR as a starting point for the spin-spiral run. The total energy and on-site magnetic moments (LORBIT) should be converged with respect to k-mesh density and energy cutoff (ENCUT).

The introduction of a spin spiral will lower the symmetry of the system. At present, VASP cannot correctly account for the presence of a spin spiral in its symmetry analysis. Therefore, the use of symmetry must be switched off in the subsequent steps by setting ISYM = -1. Generally, the number of irreducible k points changes for different symmetry settings (ISYM and thus restarting from CHGCAR should be preferred over restarting from WAVECAR. You may set LWAVE = False to avoid writing the WAVECAR. It is not advantageous to switch off symmetry from the get-go, as using symmetry reduces the computational effort and hence saves substential computational cost during the course of convergence studies.

Magnetism usually arises in d-electron and f-electron bands. To write the augmentation charges of those bands, set LMAXMIX = 4 (d electrons) or 6.

Choose the remaining electronic-minimization tags (ALGO, EDIFF, NELM, ISMEAR) as described in Setting up an electronic minimization; frustrated or incommensurate magnets often need a larger NELM and ALGO = Conjugate is often most reliable.

Initialize the magnetic order with the MAGMOM tag. In the noncollinear case, this takes three components per atom. If no initial moments are given, the initial order is ferromagnetic, which may trap the system in a local minimum even when the ground state is not ferromagnetic.

Step 2: Set up the spin-spiral INCAR

Add the following tags to the INCAR:

LNONCOLLINEAR = .TRUE.
LSPIRAL       = .TRUE.
QSPIRAL       = q1 q2 q3
ISYM          = -1
ENINI         = 300
ENCUT         = 400

Here, q1 q2 q3 are the components of the spin-spiral propagation vector in direct (fractional) coordinates of the reciprocal lattice. In a spin-spiral calculation, the magnetic configuration is set both by the initial moments within the cell using the MAGMOM tag, in the noncollinear (three-component-per-atom) format, and by how the magnetization rotates between cells (the QSPIRAL vector). Worked examples are given in the Examples section below. Choose the remaining electronic-minimization tags as in step 1.

Mind: Symmetry must be switched off completely (ISYM = -1). VASP cannot account for the symmetry reduction introduced by a spin spiral, so leaving symmetry on produces incorrect results.

Tip: Set ENCUT (or equivalently ENMAX) at least 100 eV above ENINI. VASP prints a runtime warning if the value is too small for the chosen q vector, e.g.:

 ----------------------------------------------------------------------------- 
|                                                                             |
|           W    W    AA    RRRRR   N    N  II  N    N   GGGG   !!!           |
|           W    W   A  A   R    R  NN   N  II  NN   N  G    G  !!!           |
|           W    W  A    A  R    R  N N  N  II  N N  N  G       !!!           |
|           W WW W  AAAAAA  RRRRR   N  N N  II  N  N N  G  GGG   !            |
|           WW  WW  A    A  R   R   N   NN  II  N   NN  G    G                |
|           W    W  A    A  R    R  N    N  II  N    N   GGGG   !!!           |
|                                                                             |
|      To represent the spin spiral you requested, with a kinetic             |
|      energy cutoff of ENINI=  300.00 eV, choose ENMAX >  331.21 eV          |
|      Currently ENMAX=  400.00 eV                                            |
|                                                                             |
 -----------------------------------------------------------------------------

Step 3: Constrain to a planar spiral (optional)

To prevent the magnetization density from developing a z-component, set:

LZEROZ = .TRUE.

This forces [math]\displaystyle{ m_z({\bf r})=0 }[/math] at each step of the electronic minimization.

Step 4: Extract local magnetic moments (optional)

Analyzing site-resolved local moments is less straightforward than usual, because the spin-spiral period is generally incommensurate with the unit cell, so the magnetization density is not cell-periodic:

The standard analysis via LORBIT (output in the PROCAR file and at the end of the OUTCAR file) does not account for this. As a workaround, use the constrained-moment infrastructure with a zero penalty potential:

I_CONSTRAINED_M = 1
LAMBDA          = 0.0
RWIGS           = <one radius per species>

This switches on the constrained-magnetic-moment approach but sets the penalty potential to zero (LAMBDA = 0.0). The magnetization density is then integrated inside site-centered spheres of radius RWIGS, and the resulting local moments are written under M_int in the OSZICAR file, e.g.:

 E_p =  0.00000E+00  lambda =  0.000E+00
<lVp>=  0.00000E+00
 DBL =  0.00000E+00
 ion        MW_int                 M_int
  1  1.178  0.000  0.000    1.573  0.000  0.000
RMM:   8    -0.819213822792E+01    0.53417E-07   -0.43965E-08  2542   0.310E-03

Here, the local moment on ion 1 (after iteration 8) is [math]\displaystyle{ M=1.573\,\hat{x}\;\mu_{\rm B} }[/math].

Mind: Do not forget to set RWIGS for every species; without it the integration spheres are undefined.

Step 5: Run the calculation

Run VASP as usual and monitor convergence in the OSZICAR file. A spin-spiral calculation has approximately the same cost as a standard noncollinear calculation of the primitive cell.

Examples

Example 1: Initializing the magnetic configuration

The initial moments within the cell, together with the QSPIRAL vector, determine which magnetic configuration the calculation starts from.

Double-layer antiferromagnet. Two magnetic atoms with initial moments M along y, and [math]\displaystyle{ {\bf q}=(0,0,\tfrac{1}{2}) }[/math]:

MAGMOM  = 0 M 0  0 M 0
QSPIRAL = 0.0 0.0 0.5

Flat spin spiral. Two magnetic atoms with initial moments M along y and x, respectively, and [math]\displaystyle{ {\bf q}=(0,0,\tfrac{1}{2}) }[/math]:

MAGMOM  = 0 M 0  M 0 0
QSPIRAL = 0.0 0.0 0.5

Tip: Both configurations obey the same generalized Bloch condition, [math]\displaystyle{ {\bf q}=(0,0,0.5) }[/math], and during the electronic minimization, one may transform into the other if that lowers the total energy. The Bloch condition fixes the change in magnetization density from one cell to the next, but does not constrain the magnetic order within a cell.

Example 2: Spin-spiral energy of a NiI₂ monolayer

This example uses a spin-spiral calculation to obtain the total energy of a flat spin spiral in a NiI₂ monolayer. The spin spiral is scanned along the QSPIRAL path Γ–M–K–Γ to find the magnetic ground state. The triangular Ni sublattice is magnetically frustrated, which stabilizes an incommensurate spiral.

Structure: A 1T (CdI₂-type) NiI₂ monolayer, hexagonal a = 3.97 Å with ≈8 Å of vacuum (POSCAR):

NiI2 1T monolayer
1.0
   3.9697   0.0000    0.0000
  -1.9849   3.4379    0.0000
   0.0000   0.0000   11.0271
Ni I
1 2
Direct
  0.0000 0.0000 0.5000
  0.3333 0.6667 0.6373
  0.6667 0.3333 0.3627

INCAR: One self-consistent calculation per q, changing only QSPIRAL. The Ni moment is initialized in the xy-plane and iodine is nonmagnetic:

LNONCOLLINEAR = .TRUE.
LSPIRAL       = .TRUE.
QSPIRAL       = q1 q2 0.0
ISYM          = -1
ENINI         = 600
ENCUT         = 700
ISMEAR        = 0
SIGMA         = 0.05
MAGMOM        = 2 0 0  0 0 0  0 0 0

A Γ-centered 12×12×1 k-mesh and Ni_pv and I PAW potentials were used. The total energy E(q) is read from the OUTCAR.

Result — spin-spiral energy relative to the ferromagnetic state (q = Γ):

q (fractional)	point / direction	E(q) − E(Γ) (meV)
(0, 0, 0)	Γ (ferromagnetic)	0.0
(0.143, 0, 0)	Γ→M	−22.7
(0.214, 0, 0)	Γ→M (minimum)	−29.7
(0.5, 0, 0)	M	+36.3
(1/3, 1/3, 0)	K	−12.5
(0.133, 0.133, 0)	Γ→K (minimum)	−32.9

The energy is lowest for an incommensurate spin spiral, not for the ferromagnet (Γ) or the M-point antiferromagnet. Two shallow minima appear along inequivalent in-plane directions — near q ≈ (0.21, 0, 0) along Γ–M and q ≈ (0.13, 0.13, 0) along Γ–K — separated by a barrier of only ≈3 meV, and the spiral is stabilized by ≈30 meV per formula unit relative to the ferromagnet. This reproduces the qualitative result of the GPAW tutorial: a frustration-driven incommensurate ground state with two near-degenerate minima.

Tip: The exact position of the incommensurate minimum depends on the functional, PAW potentials, k-mesh, and cutoff, so the numbers above reproduce the physics of the GPAW tutorial rather than identical values. Determining the orientation of the spiral plane additionally requires a non-self-consistent spin–orbit-coupling step, which is not covered here.

Magnon dispersion and exchange interactions

Scanning the spin-spiral energy E(q) as in Example 2 is the basis of the frozen-magnon method for extracting magnetic interactions. Mapping the computed energies onto a classical Heisenberg model,

[math]\displaystyle{ E({\bf q}) = E_0 - \sum_{\bf R} J({\bf R})\, \cos({\bf q}\cdot{\bf R}), }[/math]

gives the interatomic exchange constants [math]\displaystyle{ J({\bf R}) }[/math] as a Fourier transform of E(q). If the ferromagnetic state is the ground state, the adiabatic magnon dispersion follows from the same energies, [math]\displaystyle{ \hbar\omega({\bf q}) \propto [E({\bf q}) - E(0)]/M }[/math], with M the local moment. If instead the minimum of E(q) lies at a finite q (as for NiI₂) the system is predicted to order as an incommensurate spin spiral with that propagation vector.

Practical notes:

• Use a fixed k-mesh and the same ENINI/ENMAX for every q, since the relevant energy differences are only a few meV.
• Restart each point from the same converged charge density (ICHARG = 1) for consistency across the scan.
• The mapping assumes rigid local moments; verify that the local moment is roughly q-independent (in Example 2 it varies by only a few percent).

Mind: The spin-spiral approach is incompatible with spin–orbit coupling (LSORBIT = TRUE), because spin–orbit coupling breaks the spin-rotational symmetry that the generalized Bloch condition relies on.

Related tags and articles

Spin spirals, Setting up an electronic minimization

LNONCOLLINEAR, MAGMOM, ISYM, I_CONSTRAINED_M, LAMBDA, RWIGS, LORBIT

LSPIRAL, QSPIRAL, LZEROZ, , ENINI