1. Noncollinear magnetism and force theorem

For many materials the magnetic unit cell is neither ferromagnetic nor antiferromagnetic but features a more complex magnetic structure in which the orientation of the magnetic moments related to the different atoms changes within the unit cell. Such magnetic structures are called noncollinear and may be caused by different circumstances. One possible cause is magnetic frustration, i.e., a preferred alignment of the magnetic moments at the different atoms is in contradiction with the geometry of the system. Another cause may be spin-orbit coupling. In this tutorial we will calculate several noncollinear magnetic structures for a magnetically frustrated system.

In Fleur noncollinearity is treated such that within the MT sphere of a certain atom a single magnetization direction is assumed in the setup of the Hamiltonian but this direction can differ between the atoms. The magnetization directions in the MT spheres are given in terms of the two angles $\alpha$ and $\beta$ analogous to the $\phi$ and $\theta$ angles used for the specification of the SQA in SOC calculations. In the interstitial region between the MT spheres the orientation of the magnetization continuously changes.

1.1. Noncollinear magnetism in a Cr Monolayer with hexagonal lattice geometry

NOTE: While performing the calculations for this example you should download the visualization tool XCrysDen and get it to run if it doesn't already run on your computer. You will need it. Alternatively you can also work with VESTA.

NOTE 2: This tutorial example is motivated by the article Ph Kurz, G. Bihlmayer, and S. Blügel, Noncollinear magnetism of Cr and Mn monolayers on Cu(111), Journal of Applied Physics 87, 6101 (2000).

For Cr atoms the preferred magnetic ordering is antiferromagnetic. In a triangular lattice this is not possible and such a setup therefore is an example for a prototypical magnetically frustrated system. To obtain such a lattice for Cr one can grow a single layer of Cr atoms on a Cu(111) surface. In such a sample the Cr atoms will occupy the positions of Cu atoms in the respective layer. As an approximation to such a setup we calculate an unsupported Cr monolayer with atoms at such positions. Assume an fcc Cu crystal with a lattice constant of $6.65~a_0$. This translates into a nearest-neighbor distance between the atoms of $4.70~a_0$. In this tutorial example we consider a unit cell of two atoms, and in the exercises a unit cell containing three atoms is investigated. The setup for this example is sketched in the following figure.

Calculate the lattice parameters for the sketched unit cell and set up an inpgen input for such a film system with two inequivalent Cr atoms.

Calculations with noncollinear magnetism require parameters that are not provided within the default Fleur input template generated by the input generator. We have to use inpgen with the -noco command line option to obtain a template containing such parameters. Generate a Fleur input file in this way.

To actually perform the calculations we also have to make some adjustments to the Fleur input file:

1. Set the number of iterations to 90: This system does not show trivial convergance behavior. Note that the 90 iterations might not be enough for some of the calculations. For such a case run Fleur multiple times and remove the file "mixing_history" between the runs.

2. Set /calculationSetup/scfLoop/@maxIterBroyd to 10. This forces the deletion of the charge density mixing history after every 10 iterations of the SCF loop. This is a measure to obtain a more benign convergence behavior.

3. Set /calculationSetup/coreElectrons/@ctail to "F". Noncollinear features are incompatibel with this feature. The switch turns on a certain treatment of the part of the core electron density that reaches out of the MT spheres.

4. Set /calculationSetup/magnetism/@l_noco to "T". This actually switches on the noncollinear calculation mode.

5. Set /calculationSetup/nocoParams/@l_mperp to "T". This activates the calculation of the offdiagonal (Spin-Up - Spin-Down) part of the density in the MT spheres. We only need this for visualization purposes.

6. Reduce all MT radii to 90% and to keep the product of the MT radii and $K_\text{max}$ constant increase $K_\text{max}$ accordingly.

7. Change the XC functional to the LSDA formulation of Moruzzi, Janak, and Williams (choose "mjw" for xcFunctional/@name). This is just done in accordance with the article cited above.

8. Increase the initial magnetic moments to $3.0~\mu_\text{B}$. This is done because lower initial magnetic moments may lead to the calculation of a nonmagnetic self-consistent solution. All magnetic moments are supposed to be positive.

We will perform calculations with different in-plane magnetic setups. To restrict the SQA to the plane of the film set /atomGroups/atomGroup/nocoParams/@beta to "Pi/2.0" for both atoms.

For one of the atoms keep the angle $\alpha$ (/atomGroups/atomGroup/nocoParams/@alpha) constant at $0.0$ but vary it for the other one in steps of "Pi/6.0" from "2.0*Pi/6.0" (this is the hardest calculation) to "Pi" in different subdirectories. Calculate the self-consistent solutions for each of these setups.

Plot the total energy and the magnetic moment in dependence of $\alpha$. Your results should be similar to the following plot.

We now want to visualize for the $\alpha=4\pi/6$ case the spatially-dependent orientation of the magnetization within the unit cell. For this copy everything from this calculation to a new directory and in this directory change in the inp.xml file /output/plotting/@iplot to "1". If you invoke fleur for this input it will complain that you don't provide a plot_inp file. But on the other hand it will generate a template for such a file. In the plot_inp file one specifies which plots are to be generated. By default 2 plots are generated, one 2D plot and one 3D plot. We slightly adapt this template:

1. We extend the first line with ",polar=t". Note that you should not introduce any unwanted spaces (no spaces in front of the comma) at this place. This line is sensitve to the column in which a variable is set.

2. for the 3D plot modify the grid such that it becomes a 35 x 60 x 35 grid (the 60 should be in the elongated dimension of the unit cell) and set the offset (zero) to (-0.25,-0.25,0.0).

Invoke Fleur again. Now the following XCrysDen structure files are generated:

• cdn_pl.xsf: a file containing the charge density (this is also produced for calculations that are not noncollinear)

• mdnx_pl.xsf: magnetization density in x direction

• mdny_pl.xsf: magnetization density in y direction

• mdnz_pl.xsf: magnetization density in z direction

if polar=t you also obtain the following files:

• mabs_pl.xsf: absolute value of the magnetization density

• mphi_pl.xsf: phi angle for the magnetization direction at every point of the magnetization density

• mtha_pl.xsf: theta angle for the magnetization direction at every point of the magnetization density

Run XCrysden and open the file mphi_pl.xsf that contains the magnetization directions in the plane. Go to Tools->Data grid and select "plot2". You can then define different parametrizations for the visualization. Play with it, use a rainbow color coding together with a linear scale. If you open one of the other files not related to angles you will realize that most of the density and magnetization density is found very near to the atom positions.

For the magnetization directions you should be able to obtain a color-encoded visualization similar to the following figure.

If you have time you can also go to earlier tutorial examples from which you have self-consistent solutions and generate charge density plots for those. For example Si should have a large density between the atoms due to the covalent bonding. Is this really the case?

2. Exercises

2.1 Three-atom unit cell of the Cr film

Extracting parameters for a Heisenberg model from DFT calculations with simple magnetic configurations of a Cr film and performing simulations on the basis of this model predicts that the ground state for this system features an in-plane rotation of 120° between the magnetic moments of neighboring Cr atoms. It turns out that DFT calculations for such a magnetic configuration and similar configurations agree with this prediction. Note that such an agreement between Heisenberg model predictions and more accurate DFT calculations is not for all systems the case.

Here we want to perform the DFT calculation for the magnetic ground state structure.

Print out or draw a reasonable region of a triangular lattice (you can use the sketch above) and select a triangle of neighboring atoms. For this triangle draw arrows at the atoms pointing from the center of the triangle outwards. Assume that these arrows indicate the magnetic moments. Draw the magnetic moments for the other atoms according to the 120° rotations from atom to atom.

On the basis of this drawing select a magnetic unit cell containing three inequivalent atoms. Verify that the repetitions of this unit cell feature the same magnetization. Set up an inpgen input with this Cr film unit cell and the nearest-neighbor distance provided in task 1.1. Generate a Fleur input with parametrization for noncollinear setups and adapt the parameters as good as possible to the Fleur input files from task 1.1. Set the SQAs at the different atoms according to the construction of the magnetic unit cell.

Perform a Fleur calculation to obtain the ground-state density for this setup. Compare the total energy per atom to the energetically most likely configuration from task 1.1. Plot the magnetization direction at each point in the unit cell in the film plane.

Results to be delivered: 1. The comparison of the total energy per atom to task 1.1. 2. A plot of the magnetization directions at every point in the unit cell in the film plane.

2.2. Magnetocrystalline anisotropy in Co revisited with the magnetic force theorem

Last week the magnetocrystalline anisotropy in Co was determined by performing self-consistent calculations for different spin-quantization axes (SQAs). There is also a faster way to do this. According to the magnetic force theorem already single-shot calculations for different SQAs (only determination of the eigenvalue sum for a fixed input density) on top of a representative self-consistent density for one SQA should provide reasonable approximate results. To do this we have to add a special XML element to one of the inp.xml files for Co from last week. Choose the one with $\theta=0.0$ and $\phi=0.0$ and add

   <forceTheorem>
      <MAE theta="0.0 0.5*Pi 0.25*Pi 0.5*Pi" phi="0.0 0.0 0.0 0.5*Pi" />
   </forceTheorem>

directly after the output section of the inp.xml file. In the forceTheorem XML element we provide the parametrization for such a calculation. It can be used to determine the magnetocrystalline anisotropy energy (MAE), for other quantities like the $\mathbf{D}_{ij}$ parameters of the Dzyaloshinskii-Moriya interaction (DMI), or for the determination of spin-spiral dispersion relations.

For the MAE calculation a MAE xml element has to be provided that features the two xml attributes theta and phi, both representing a list of respective angles. The i-th angles from both lists define the i-th SQA. The example above provides exactly those angles that were also used last week.

Run fleur with this input to first calculate a self-consistent density. The force theorem calculation starts if the iteration number of the current run equals the maximum iteration number specified in the inp.xml file or the convergence criterion is reached. Note that the usage of the force theorem replaces the 4 self-consistency calculations from last week by a single self-consistency calculation. We thus nearly obtain a factor 4 speedup.

In the out.xml file you should now have a "Forcetheorem_MAE" section in which a list with the different angles together with respective eigenvalue sums is provided. Determine the "K_i" parameters from last weeks uniaxial MAE model and compare them to the parameters determined with self-consistent calculations.

Results to be delivered: The eigenvalue sums for each combination of angles and the parameters $K_1$ and $K_2$ for the MAE model obtained from these eigenvalue sums.