# Spin-Orbit coupling¶

- Author
Ramón Cuadrado (University of Southampton) & Nils Wittemeier (ICN2)

The main goal of this tutorial is to cover the basic notions to perform calculations using the spin and spin-orbit coupling (SOC) in Siesta.

Note

When performing a SOC calculation it is crucial to:

carefully converge the k–points sampling

carefully converge the mesh cutoff

converge the density matrix with a tolerance below 10^-5

use relativistic pseudopotentials (PPs) and use the

*non–linear core corrections*when the PPs are built.

All the information regarding the basis set and all the parameters used in a common SIESTA calculation are valid for the SO.

In the next two examples the above mentioned
parameters (mesh cutoff, etc) are **not** converged because we want to speed
up the calculations to show an initial (and dirty) results, but
remember, for a real calculation those values have to be
converged.

An additional advice: remember to read the SOC section in the main Siesta manual!

## Performing calculations with different spin configuration.¶

Hint

Please move to directory Example-1.

Note

The original fdf for this example can be found in example_1.fdf.original. We have made some changes (can you see them?) to speed the calculation while maintaining reasonable quality results.

In this example you will calculate the band structure, density-of-states (DOS), and Mulliken charges of bulk FePt using different spin options.

### Getting started¶

To get started inspect the input file *example_1.fdf* file, and make sure that you
are familiar with all options in the input file. Refer to the user manual to clarify
any doubts. Alongside the input file you will find two pseudopotential files one for
Fe and one for Pt (*Fe_fept_SOC.psf*, *Pt_fept_SOC.psf*)

In addition, to the usual flags we have included the following options:

to calculate and store and band structure information

BandLinesScale ReciprocalLatticeVectors %block BandLines 1 0.0000 0.0000 0.0000 Gamma 20 0.5000 0.0000 0.0000 X 20 0.5000 0.5000 0.0000 M 20 0.0000 0.0000 0.0000 Gamma 20 0.0000 0.0000 0.5000 Z 20 0.0000 0.5000 0.5000 R 20 0.5000 0.5000 0.5000 A 20 0.0000 0.0000 0.5000 Z %endblock BandLines

to calculate (projected) density-of-states

%block PDOS.kgrid_Monkhorst_Pack 30 0 0 0.0 0 55 0 0.0 0 0 55 0.0 %endblock PDOS.kgrid_Monkhorst_Pack %block ProjectedDensityOfStates -10.00 5.00 0.02 1000 eV %endblock ProjectedDensityOfStates

to calculate and print mulliken charges

WriteMullikenPop 1

### Simulations¶

Run one SIESTA calculations for each possible options for the `Spin`

flag: `none`

,
`polarized`

, `non-colinear`

, `spin-orbit`

.
Please, change the `SystemLabel`

for each calculation to avoid overwriting
output files or perform each calculation in a separate subdirectory.

### Analysis¶

Visualize the bands structure, and compare the results. Do any of the calculations give the same band structure? Why?

Templates for plotting the bands/DOS with

`gnuplot`

can be found in the .gnp files.Note

The default version of gnuplot on MareNostrum is quite old (v4.6). If you plan on running

`gnuplot`

on MareNostrum, we recommend loading a recent version first:module load gnuplot/5.2.3

Note

If you prefer performing the analysis with

`sisl`

in a`jupyter`

notebook and you have both installed, please feel free to give it go. To get started please, have a look at example_1.ipynb. Unfortunately, we can not support running the`jupyter`

notebooks on MareNostrum.Visualize the density-of-states (

*.DOS*file).Note

For

`Spin none`

the*.DOS*file contains two columns: energy and DOS.For

`Spin polarized`

the*.DOS*file contains three columns: energy, DOS*↑*, DOS*↓*(DOS for spin-up/spin-down electrons). The integral of DOS*↑*(DOS*↓*) up to E_{F}will give the number of electrons with spin-up (spin-down). The total density of states is the sum of DOS*↑*and DOS*↓*.For

`Spin non-colinear`

or`Spin spin-orbit`

the*.DOS*file contains five columns: energy, DOS*↑*, DOS*↓*, Re{DOS*↑↓*} and Im{DOS*↑↓*}.

The integral of Re{DOS

*↑↓*} up to E_{F}will give the magnetization along*X*axis, M*x*, and the integral of Im DOS*↑↓*up to E_{F}, M*y*. The total density of states is the sum of DOS*↑*and DOS*↓*.Inspect the Mulliken charges printed in the output file (search for ‘Mulliken’). Are all atoms/orbitals equally spin-polarized? If not, electrons in which orbital do primarily contribute the magnetic properties?

*More advanced*: Visualize the**projected**density-of-states (pDOS)Hint

In the interest of time, you might want to go through the examples below first before doing this.

Note

The

*.PDOS*file can be processed as discussed in this how-to.You can plot the PDOS

*↑*, PDOS*↓*, PDOS*↑↓*, etc. for each atomic species or the total DOS (the sum of the relevant columns). It is worth noting that in general, for the transition metals the most relevant orbitals will be the*d*ones. You can play with the DOS/pDOS analysis programs to choose other PDOS curves.

## Calculation of the magnetic anisotropy of Pt_{2}.¶

Hint

Please move to directory Example-2

Note

The original fdf for this example can be found in example_2z.fdf.original. We have made some changes (can you see them?) to speed the calculation while maintaining reasonable quality results.

In this example, you will calculate the magnetic anisotropy of a Pt_{2}
dimer. This means you will have to run multiple SIESTA calculations each with
magnetic moments aligned along different directions and extract the total energy.

To get started inspect the input file *example_2z.fdf* file. Use this input file
and the pseudopotential of Pt atom, *Pt_pt2_SOC.psf* (note that the pseudo for
Pt is different from the one used in the first example), and calculate the total
self-consistent energy E*tot* for three highly symmetric orientations of the
initial magnetization, namely, X, Y and Z, keeping fixed the physical dimer
(the (x,y,z) coordinates of each atom).

Each calculation has to be initialized modifying the `DM.InitSpin`

, updating
the orientation of each atomic magnetization. As a result one has to have three
different energy values for each (*θ*,*ϕ*) pair.

Remember to rename each one of the *.fdf* files, change the SystemLabel
flag inside for each magnetic orientation and do not re-use the density
matrix (*.DM* file) obtained in other orientations. Reading the density matrix
will cause SIESTA to ignore the `DM.InitSpin`

. To calculate the total energy
for each specific magnetic configuration you need to start from scratch.

Note

You can also try to change the orientation of the dimer (say put it along the Z or Y axes) and compute again the energy values associated to different spin orientations and compare the results. Some of the total energies (which?) have to be the same.

Note

What would happen if you used `Spin non-colinear`

instead of
`Spin spin-orbit`

? Would you obtain the same curve?

## Calculation of the magnetic anisotropy of FePt bulk¶

Hint

Please move to directory Example-3

Note

The original fdf for this example can be found in example_3z.fdf.original. We have made some changes (can you see them?) to speed the calculation while maintaining reasonable quality results.

In this example, you will calculate the magnetic anisotropy bulk FePt.

You will find the *example_3z.fdf* file and the PPs of Fe and Pt: *Fe_fept_SOC.psf*
and *Pt_fept_SOC.psf*. Each calculation has to be initialized modifying
the `DM.InitSpin`

block, as above, updating the orientation of each atomic
magnetization. As a result one has to have different energy values
for each (*θ*,*ϕ*) pair.

A) Obtain the E

*tot*for three highly symmetric orientations of the initial magnetization, namely, X, Y and Z, keeping fixed the atoms in their unit cell. Again do not use the same density matrix for same calculations. Allow SIESTA to calculate it from scratch.B) Create a magnetization curve changing the angles at intervals of 20 degrees from the Z-axis initial orientation, to finalize along the X axis. Plot the resulting energy as a function of the varied angle. Check the orientation of the final magnetic moment in each output file. Does it coincide with the initial one?

If this set of calculations takes too much time to be performed completely, you can try to run only a few steps and see how the total energy evolves.

Note

If you do plot the magnetization curve, you might see that it has a minimum at around 60 degrees, which does not make too much sense, and it is quite shallow. In this case we went too far in the simplification of the parameters to gain speed. In particular, removing the polarization orbitals has a very significant effect. At the cost of a noticeably longer simulation time, you can restore those orbitals by using the block:

```
%Block PAO.Basis
Fe_fept_SOC 2
n=4 0 2 P
0.0 0.0
n=3 2 2
0.0 0.0
Pt_fept_SOC 2
n=6 0 2 P
0.00000 0.00000
n=5 2 2
0.00000 0.00000
%EndBlock PAO.Basis
```

You should now get a curve similar to that in Fig. 7: