Computing magnetic interaction parameters with TB2J¶
- Author
Xu He (ICN2)
Introduction¶
To study the magnetic properties of materials from first principles, it is often required to map the full electronic interactions to effective Heisenberg spin Hamiltonian, which will help to understand the origin of the magnetic interaction and enables large-scale magnetic property simulation, e.g. with spin dynamics or Monte Carlo. One of the most successful method for doing this is to use the so-called [LKAG] method or the Liechtenstein formula, It takes the local rigid spin rotation as a perturbation and use the Green’s function method to compute the magnetic interaction parameters. The method is described in detail in reference [TB2J] and references therein. It requires only one self-consistent calculation (or three in the case of non-collinear spin) in a unitcell for the computation of the magnetic interactions up to a large distance given a dense enough k-point mesh.
TB2J is a python package which implements such method, which can make use of the SIESTA output. In this exercise you will calculate the magnetic interaction parameters with TB2J package. We’ll first learn how to install TB2J. Then we’ll do two exercises of using TB2J to calculate the exchange parameters in the collinear and non-collinear cases, with the examples of bcc-Fe and the 2D CrI3 mono-layer, respectively.
Note
This tutorial shows how to use the interface of SIESTA and TB2J in a minimalist manner. More information about the usage of the TB2J package can be found in the full TB2J documentation, and the background can be found in [TB2J] and references therein.
Installation¶
Note
These packages are already installed if you are using a school
virtual machine to run the exercise, type workon siesta_school
to activate the virtual environment where TB2J is installed.
TB2J is a python package, and requires Python (version>=3.7) to be installed. You can install TB2J with pip by the following command:
pip install TB2J
To use TB2J with SIESTA, sisl is also needed. It can be installed with:
pip install sisl
The Heisenberg model¶
Before you use the TB2J, please read very carefully this section. There are many conventions of the Heisenberg Hamiltonian. We strongly suggest that you clearly specify the convention you use in any published work. Here we describe the convention used in TB2J.
The Heisenberg Hamiltonian contains three different parts and reads as
where the first term represents the isotropic exchange, and the second term is the symmetric anisotropic exchange, where \(\mathbf{J}^{ani}\) is a \(3\times 3\) tensor with \(J^{ani}=J^{ani,T}\). The final term is the Dzyaloshinskii-Moriya interaction (DMI).
Note
More details of conventions in the Heisenberg Hamiltonian can be found in this page .
Exercise 1: isotropic exchange parameter of BCC-Fe.¶
In this exercise, we’ll compute the exchange parameters for bcc Fe. The files needed can be downloaded with the link below.
Assuming that we already have a relaxed structure, we do a self-consistent calculation with collinear spin first. We need to turn on the option to save the Hamiltonian and overlap matrices, which will be read by TB2J. Although the CDF.Save option is not necessary, we recommend the use of it for saving the matrices to the netcdf format.
SaveHS True
CDF.Save True
After the calculation, run the following command:
siesta2J.py --fdf_fname Fe.fdf --elements Fe --kmesh 7 7 7 --nz 80
where siesta2J.py is a script in TB2J, which compute the magnetic interaction parameters from SIESTA output. Use the –fdf_fname option to specify the name of the SIESTA input file, and the –elements Fe option to specify that the magnetic element in this structure is Fe; use the –kmesh 7 7 7 option to specify the sampling of the Brillouin zone within TB2J calculation; the –nz option specifies the number of poles used for the integration along the energy axis. The number of k-points and poles should be tested to make sure that the magnetic interactions are converged to the precision required.
There are a few more options, which can be viewed with the following command:
siesta2J.py --help
After running siesta2J.py, a folder called TB2J_results is generated, in which there is file named exchange.txt. We can open this file with a text editor. The description of the output file can be found in the TB2j documentation .
For each atom pair, the exchange parameter and the orbital decomposition is given. In this exercise, we can grep the exchange parameters linked with Fe by :
grep Fe1 exchange.txt
which yields the following contents, showing the 1st, 2nd, 3rd, and 4th nearest-neighbor exchange parameters as 15.84 meV, 11.02 meV, -0.89 meV, -0.48 meV, respectively. Further neighbors are also shown.
Fe1 Fe1 ( -1, 0, 0) 19.6542 (-2.467, 0.000, 0.000) 2.467
Fe1 Fe1 ( 0, -1, 0) 19.6563 ( 0.822, -2.325, 0.000) 2.467
Fe1 Fe1 ( 0, 0, -1) 19.6484 ( 0.822, 1.163, -2.014) 2.467
Fe1 Fe1 ( 0, 0, 1) 19.6484 (-0.822, -1.163, 2.014) 2.467
Fe1 Fe1 ( 0, 1, 0) 19.6563 (-0.822, 2.325, 0.000) 2.467
Fe1 Fe1 ( 1, 0, 0) 19.6542 ( 2.467, 0.000, 0.000) 2.467
Fe1 Fe1 ( -1, -1, -1) 19.6631 (-0.822, -1.163, -2.014) 2.467
Fe1 Fe1 ( 1, 1, 1) 19.6631 ( 0.822, 1.163, 2.014) 2.467
Fe1 Fe1 ( 0, -1, -1) 12.3918 ( 1.644, -1.163, -2.014) 2.848
Fe1 Fe1 ( 0, 1, 1) 12.3918 (-1.644, 1.163, 2.014) 2.848
Fe1 Fe1 ( -1, -1, 0) 12.3937 (-1.644, -2.325, 0.000) 2.848
Fe1 Fe1 ( -1, 0, -1) 12.3904 (-1.644, 1.163, -2.014) 2.848
Fe1 Fe1 ( 1, 0, 1) 12.3904 ( 1.644, -1.163, 2.014) 2.848
Fe1 Fe1 ( 1, 1, 0) 12.3937 ( 1.644, 2.325, 0.000) 2.848
Fe1 Fe1 ( -1, 0, 1) -2.6483 (-3.289, -1.163, 2.014) 4.028
Fe1 Fe1 ( -1, 1, 0) -2.6415 (-3.289, 2.325, 0.000) 4.028
Fe1 Fe1 ( 0, -1, 1) -2.6464 ( 0.000, -3.488, 2.014) 4.028
Fe1 Fe1 ( 0, 1, -1) -2.6464 ( 0.000, 3.488, -2.014) 4.028
Fe1 Fe1 ( 1, -1, 0) -2.6415 ( 3.289, -2.325, 0.000) 4.028
Fe1 Fe1 ( 1, 0, -1) -2.6483 ( 3.289, 1.163, -2.014) 4.028
Fe1 Fe1 ( -2, -1, -1) -2.6517 (-3.289, -1.163, -2.014) 4.028
Fe1 Fe1 ( -1, -2, -1) -2.6522 (-0.000, -3.488, -2.014) 4.028
Fe1 Fe1 ( -1, -1, -2) -2.6541 (-0.000, -0.000, -4.028) 4.028
Fe1 Fe1 ( 1, 1, 2) -2.6541 ( 0.000, 0.000, 4.028) 4.028
Fe1 Fe1 ( 1, 2, 1) -2.6522 ( 0.000, 3.488, 2.014) 4.028
Fe1 Fe1 ( 2, 1, 1) -2.6517 ( 3.289, 1.163, 2.014) 4.028
Fe1 Fe1 ( -1, -1, 1) 0.0218 (-2.467, -3.488, 2.014) 4.723
Fe1 Fe1 ( -1, 1, -1) 0.0246 (-2.467, 3.488, -2.014) 4.723
Fe1 Fe1 ( -1, 1, 1) 0.0236 (-4.111, 1.163, 2.014) 4.723
Fe1 Fe1 ( 1, -1, -1) 0.0236 ( 4.111, -1.163, -2.014) 4.723
Fe1 Fe1 ( 1, -1, 1) 0.0246 ( 2.467, -3.488, 2.014) 4.723
Fe1 Fe1 ( 1, 1, -1) 0.0218 ( 2.467, 3.488, -2.014) 4.723
Fe1 Fe1 ( -2, -1, 0) 0.0265 (-4.111, -2.325, 0.000) 4.723
Fe1 Fe1 ( -2, 0, -1) 0.0273 (-4.111, 1.163, -2.014) 4.723
.............
In the TB2J_results directory, several other formats of the exchange parameters are also provided, which can be used in spin dynamics code, e.g. MULTIBINIT.
Exercise 2: isotropic/anisotropic exchange and DMI of CrI3 monolayer¶
In this exercise, we’ll compute the isotropic exchange, anisotropic exchange, and the DMI in CrI3 monolayer. The anisotropic exchange and DMI are due to the spin-orbit coupling, which should be enabled in the SIESTA calculation (see this tutorial for more details).
The files needed for this exercise can be downloaded:
Three calculations of the cell with z-axis roated to x, y, and z axis need to be performed. The spins are set to be along the z-axis. Then the magnetic interactions from these three calculations can be then merged. This rotate-and-merge strategy is discussed in detail here .
Given a structure in a file, e.g. CrI3.xyz, one can use the following command to generate the rotated structures, which can then be used in the SIESTA SCF calculations:
TB2J_rotate.py CrI3.xyz --ftype xyz
In this exercise, the three rotated structure are already included in the siesta-x.fdf, siesta-y.fdf, and siesta-z.fdf files. We can perform the three scf calculations in three directories, say x, y, z. In each of them, we run the following command after the SIESTA SCF calculation.
siesta2J.py --fdf_fname siesta.fdf --elements Cr --kmesh 9 9 1 --nz 60
Then we can merge the results of the three calculations using:
TB2J_merge.py x y z --type structure
You’ll get the result in the TB2J_results directory. In the exchange.out file, the isotropic exchange (J_iso), anisotropic exchange tensor (J_ani), and the DMI could be found. For example, the magnetic interaction parameters for one of the first nearest neighbors are shown below:
Exchange:
i j R J_iso(meV) vector distance(A)
----------------------------------------------------------------------------------------
Cr1 Cr2 ( 0, 0, 0) 1.5450 (-3.893, -0.000, 0.000) 3.893
J_iso: 1.5450
[Testing!] Jprime: 10.394, B: -4.666
[Testing!] DMI: ( 0.0000 -0.0000 -0.0000)
[Testing!]J_ani:
[[-0.17 0. -0. ]
[ 0. -0.873 -0.203]
[-0. -0.203 -0.396]]
It can be found that the 1NN exchange is ferromagnetic, with a strong anisotropy. The DMI is 0.
The results for one of the second nearest neighbors show that the exchange is also ferromagnetic and there is a non-zero DMI:
Cr1 Cr1 ( -1, 1, 0) 0.9050 (-0.000, -6.745, -0.001) 6.745
J_iso: 0.9050
[Testing!] Jprime: 3.456, B: -1.316
[Testing!] DMI: ( 0.0000 0.0373 0.0197)
[Testing!]J_ani:
[[-0.037 0. 0. ]
[ 0. -0.129 -0.022]
[ 0. -0.022 -0.129]]