Wannier Functions
Objectives:
- Learning how to set up FLEUR to calculate Wannier functions
- Studying of the input/output and the connection to Wannier90
- Doing a band interpolation
Introduction
In this tutorial, we will learn how to use FLEUR together with the WANNIER90 code.
As an example, we will obtain Wannier functions (WFs) for PbTiO3, one of the famous insulating ferroelectric materials. In general, the procedure consists of several steps.
!!! info 1. Converging the charge density and analyzing the electronic structure (from FLEUR) 2. Preparing the files necessary for WANNIER90 (from FLEUR) 3. Obtaining WFs (from WANNIER90)
Once you obtain WFs, we can proceed to post-processings. In this example, we will do the band interpolation.
1. Example: PbTiO3
Converging the charge density and analyzing the electronic structure
The first thing you need to do is to converge the charge density as usual FLEUR calculations. You can find the corresponding input in the PbTiO3 directory. Let's have a look:
cd PbTiO3; ls
You will find the usual input for FLEUR. To ease your job, we have provided the inp.xml and the corresponding kpts.xml and sym.xml already. You might want to check this input (you find the input used for inpgen at the end of the inp.xml) and make sure that you understand the setup. To save time, we also provide:
- a cdn.hdf file containing a self-consistent charge density. You may want to check this by running
fleur_MPIand thereby confirming that this charge is actually converged. - in the directory
precalcwe also provided ados_plot.pngand abandstructure.pngfile used to analyse the electronic structure as discussed next. Remember that we have parts of the tutorial dedicated to generating such plots. So if you want to reproduce these, please revisit these parts of the tutorial.
Electronic structure analysis
Before proceeding further, it is important to analyze the band structure and their orbital characters because from this information we decide which initial projections to use for the WFs. This is especially important when you wannierize a new class of materials for the first time. Once you know the recipe, you can apply same/similar projections to materials in the same class.
So we now will open the dos_plot.png and the bandstructure.png files:
display <precalc/dos_plot.png
display <precalc/bandstructure.png
The following main features of the electronic structure are clearly visible:
- A band located at ~ -7 eV comes from Pb s-orbital
- A group of bands from -5 eV to 0 originates in O p-orbitals
- A group of bands from 1 eV to 8 eV mainly comes from Ti d-orbitals and O p-orbitals.
What we are going to do is construct WFs that describe the second group by using p-like WFs located near the O atoms. There are in total 9 bands, and we will use px, py, pz orbitals for three O atoms (3x3=9).
It is important to understand that such an understanding of the electronic structure is important to be able to identify the correct states to use in the construction of Wannier functions. A sitation like encountered here, in which the states of interest form well separated bands is of course ideal for the transformation of Bloch states into Wannier states.
Preparing files necessary for Wannierization
To geneate files necessary for the wannierization from FLEUR, we need to start from a converged calculation, i.e. we need the files: cdn.hdf, inp.xml, sym.xml, and kpts.xml.
kpoints
To prepare for Wannier functions, we use an equally-distant kpoint mesh that spans the full Brillouin zone. This does not necessarily need to be same as the kpts used in the self-consistency calculation. In most of the cases, 8x8x8 k-mesh is sufficient as long as the final WFs do not have significant overlap over 8 unit cell distance. An additional kPointList can be generated by using inpgen:
inpgen -inp.xml -kpt wannier#gamma@grid=8,8,8 -noKsym
When you open kpts.xml, you can find kPointList whose name is "wannier"
<kPointList name="wannier" count="512" nx="8" ny="8" nz="8" type="mesh">
<kPoint weight=" 1.0000000000000"> 0.00/8.00 0.00/8.00 0.00/8.00</kPoint>
<kPoint weight=" 1.0000000000000"> 0.00/8.00 0.00/8.00 1.00/8.00</kPoint>
<kPoint weight=" 1.0000000000000"> 0.00/8.00 0.00/8.00 2.00/8.00</kPoint>
<kPoint weight=" 1.0000000000000"> 0.00/8.00 0.00/8.00 3.00/8.00</kPoint>
<kPoint weight=" 1.0000000000000"> 0.00/8.00 0.00/8.00 4.00/8.00</kPoint>
<kPoint weight=" 1.0000000000000"> 0.00/8.00 0.00/8.00 5.00/8.00</kPoint>
<kPoint weight=" 1.0000000000000"> 0.00/8.00 0.00/8.00 6.00/8.00</kPoint>
<kPoint weight=" 1.0000000000000"> 0.00/8.00 0.00/8.00 7.00/8.00</kPoint>
<kPoint weight=" 1.0000000000000"> 0.00/8.00 1.00/8.00 0.00/8.00</kPoint>
<kPoint weight=" 1.0000000000000"> 0.00/8.00 1.00/8.00 1.00/8.00</kPoint>
...
Open ```inp.xml``` and set kPointList ```<kPointListSelection listName="wannier"/>``` so that FLEUR uses the k-points we just created.
Important switches for Wannier functions in FLEUR
We will now discuss the most important ingredients to prepare the input for Wannier90 by FLEUR. These are connected to changes in the inp.xml file in a new tag in the output section. In this section, we have to add a new attribute wannier="T" and add the new section called wannier.
<output dos="F" band="F" slice="F" wannier="T">
<wannier>
<jobList> .... </jobList>
</wannier>
Most steps related to the generation of input for Wannier90 and further processing using Wannier funnctions are specified by keywords in the joblist tag.
projgen
Defining initial projections
The iterative approach to define maximally localized Wannier functions as implemented in Wannier90 needs a starting point which is given by projecting the calculated Bloch states onto some set of initial localized functions. These functions are typically defined in terms of atomic orbitals. While a more detailed setup using an proj file exists, in this tutorial, we will use a simpler steup in which this file is generated automatically.
This is done by providing FLEUR with a projgen_inp file. To define px, py, pz orbitals for O atoms, make a new file and name it as projgen_inp in the directory wann:
For our example this file will look like this:
O 1 1 0
O 1 2 0
O 1 3 0
The first column specifies the atom, the second and third column are for l and m values (FLEUR follows the convention from WANNIER90. For the definition of l and m, please check the Table of the official documentation for WANNIER90). The last column specifies the radial dependence, where 0 chooses FLEUR radial function corresponding to the chosen angular momentum.
We now create such a file:
cat >projgen_inp <<EOF
O 1 1 0
O 1 2 0
O 1 3 0
EOF
in addition to providing this file, we will have to specify the keyword projgen in the jobList tag of the wannier section in inp.xml. We will
summarize all changes later on.
Selection of states to use
The tag bandSelection within the wannier section in inp.xml is used to specify the Bloch states that we are going to use to construct WFs. minSpinUp and maxSpinUp define indices of the Bloch states, which are labeled in the increasing order in energy. In a non-spin-polarized calculation, the corresponding attributes minSpinDown and maxSpinDown can also be used, by default FLEUR assumes min/maxSpinDown=min/maxSpinUp.
<output dos="F" band="F" slice="F" wannier="T">
<wannier>
<bandSelection minSpinUp="14" maxSpinUp="22"/>
<jobList> projgen ... </jobList>
</wannier>
</output>
````
```<bandSelection minSpinUp="14" maxSpinUp="22"/>``` means that we use from 14th Bloch state to 22nd Bloch state (in total 9 Bloch states). This is because
1. We defined local orbitals, 5d states for Pb and 3s3p states for Ti, 5+4=9 (Please check ``` <lo type="SCLO" ... > ``` in ```inp.xml```.). So, we have 9 semicore states far below the Fermi energy.
2. There are 3 O s-states (one s-state from each O atoms).
3. We do not include Pb s-state, which is located at ~ -9 eV.
### Input for WANNIER90
FLERU can also be used to generate further required input for Wannier90:
* `bkpts`: contains vectors connecting neighboring kpoints, which are necessary for evaluating finite-difference expressions in $M_{mn}^{(\mathbf{k},\mathbf{b})}$.
* `WF1.win`: the input file for the WANNIER90 program.
These files can be geneated by specifying ```prepwan90``` in `jobList`.
### $A_{mn}(\mathbf{k})$ and $M_{mn}^{(\mathbf{k},\mathbf{b})}$
To run WANNIER90 we need to provide $A_{mn}(\mathbf{k})$ (initial projection) and $M_{mn}^{(\mathbf{k},\mathbf{b})}$ (overlap of the Bloch states between neighboring kpoints).
These are stored in the files:
* `WF1.amn` (in case of collinear magnetic setups the WF2
* `WF1.mmn`
(in case of collinear magnetic setups the WF2.xxx files contain the data for the second spin.)
and can be generated by specifying ```matrixamn``` and ```matrixmmn``` in ```jobList```.
FELUR will also generate ```WF1.eig``` that contains energy eigenvalues of the Bloch states.
### Summary of changes in the inp.xml
So to perform all these steps at once, we specify in `inp.xml`
Please modify the `inp.xml` accordingly.
Then run FLEUR to generate the file for Wannier90.
fleur_MPI
Running WANNIER90
Now, we are finally ready to run wannier90.x. Before running the program, let's quickly check WF1.win, which specifies options and parameters for obtaining WFs. Then we can run
wannier90.x WF1
Post-processing: Band interpolation
Post-processings in general consist of three separate steps. Let's assume that we want to calculate an observable . What we have to do is:
- Evaluation of in the Bloch basis.
- Obtaining in the Wannier basis by Fourier transformation together with .
- Interpolation in k-space by Fourier-transforming .
FLEUR provides several routines for the steps 1 and 2. WANNIER90 also has an user-friendly interface for this. Especially, the band interpolation can be done very easily within WANNIER90. For this, the WF1.eig file was already created by FLEUR. To obtain the band structrue from the WFs, you need to open WF1.win and change a few keys.
- First, uncomment
restart=plotandbands_plot=true. - Second, we also need to provide high-symmetry lines at the end of
WF1.win.
For example, WF1.win should look like this
length_unit=Bohr
num_wann= 9
num_iter= 100
num_bands= 9
search_shells= 200
!optional parameters for wannierization
!num_cg_steps=
!trial_step=
!fixed_step=
!restart=wannierise
! optional parameters for plotting
spin=up
restart=plot
!wannier_plot=true
!wannier_plot_supercell=3
bands_plot=true
!fermi_surface_plot=true
...
begin kpoint_path
G 0.0 0.0 0.0 X 0.5 0.0 0.0
X 0.5 0.0 0.0 M 0.5 0.5 0.0
M 0.5 0.5 0.0 G 0.0 0.0 0.0
G 0.0 0.0 0.0 Z 0.0 0.0 0.5
Z 0.0 0.0 0.5 R 0.5 0.0 0.5
R 0.5 0.0 0.5 A 0.5 0.5 0.5
A 0.5 0.5 0.5 Z 0.0 0.0 0.5
end kpoint_path
Please modify the `WF1.win` file accordingly.
wannier90.x WF1
Now, the band structure is written in WF1_band.dat.
head -10 WF1_band.dat
The first column is for the kpoint and the second column is for the energy eigenvalues. In order to visualize this interpolated bandstructure, you can use a plotting program of your choice. Here, for your conviencence, we provide a small python script that generates a plot with data from the original DFT calculation and the interpolated band from wannier90.
python ../compare.py
display <compare.png
Of course the interpolated wannier bands will only cover the O-p bands we constructed the wannier functions from. And hopefully you find perfect agreement in your plot here.
Postprocessing: calculation of the Ferroelectric polarization
Wannier functions can be used to calculate the ferroelectric polarization of a materials. In detail there are several contributions to this:
- Ionic part: The ions are considered to carry some nominal ionic charges. In general, the polarization is calculated from the difference between the polarization term of the present system and the polarization term of a reference system which does not carry an electronic polarization.
- Electronic part: The files
proj.1/proj.2/projdefine the wannier function host sites. The centers of the wannier functions are shifted away from these host sites, whereby dipoles are created. These dipole moments are added up to give the electronic contribution to the polarization.
The charges of the ions can be given by providing a file called IONS. In our example it should look like this:
cat >IONS <<EOF
3
Ti 4.0
Pb 2.0
O -2.0
EOF
Now we can run FLEUR again with the following keywords in the jobList section:
* dipole3 to specify that we want to calculate the dipole moment from the polarization.
* stopopt to prevent FLEUR to continue after this initial step, i.e. to stop before calculating DFT eigenfunctions.
<wannier>
<bandSelection minSpinUp="14" maxSpinUp="22"/>
<jobList> dipole3 stopopt</jobList>
</wannier>
Please modify `inp.xml` accordingly.
fleur_MPI
Now FLEUR will output the polarization splitted in different contributions. We should find a total polarization of about .
Please understand, that this is not the correct and full value since we included only the contribution of the Oxygen-p orbitals. In order to get the full polarization you would have to include more orbitals in your construction of Wannier functions.
2. Example: Fe
We now switch to the second example in this tutorial.
cd ../Fe
For this second example, we again provided you with some input and some precalculated quantities:
Check the inp.xml file to verify:
- We calculate bulk Fe, with a single atom per unit cell. We treat the 3s and 3p states with local orbitals, so we have 8 electrons in semi-core bands.
- We already generated a k-point set suitable for Wannier calculations using the
inpgencommand as in the first example. The corresponding name of the k-point set was already set ininp.xml.
We also already provide a cdn.hdf file with a converged density.
Setting up the input for Wannier calculations in FLEUR
Again, we first need to define initial projections. For bcc Fe, it works fine to use "usual" spd projections. Hence, we create a file
cat >projgen_inp <<EOF
Fe 0 0 0
Fe 1 0 0
Fe 2 0 0
EOF
This should define total 9 WFs (s, px, py, pz, dxy, dyz, dzx, dz2, dx2y2). Note that in the third column, which is supposed to indicate "m" for the angular part of the WFs, is "0". Specifying "0" chooses all the orbital having a value of "l" For example,
Fe 1 0 0
is equivalent to
Fe 1 1 0
Fe 1 2 0
Fe 1 3 0
We defined bandSelection from 5 to 13 (9 Bloch states) because
there are 4 semicore states far below the valence states (we defined 3s3p states as local orbitals).
Also, set numbands="20" so that FLEUR calculates enough states. By default, FLEUR will not calculate too many unoccupied states, and thus you might need to set this parameter to obtain all the Bloch states needed.
Summarizing, we need to modify the inp.xml at two positions:
<calculationSetup>
<cutoffs Kmax="3.50000000" Gmax="10.50000000" GmaxXC="10.50000000" numbands="20"/>
....
</calculationSetup>
<output dos="F" band="F" slice="F" wannier="T">
<wannier>
<bandSelection minSpinUp="5" maxSpinUp="13"/>
<jobList> projgen prepwan90 matrixamn matrixmmn </jobList>
</wannier>
...
</output>
Please modify `inp.xml` by including the `numbands="20"` attribute, the `wannier="T"` attribute and the `wannier` section
Now FLEUR can be executed
fleur_MPI
Running WANNIER90
Now, we are finally ready to run wannier90.x. We will only use the first spin in this example, the second spin data is also available in the files starting with WF2.
wannier90.x WF1
In this example, we will not calculate the bandstructure, but a Fermi surface. To achieve this, we first have to identify the Fermi energy.
This can be done by greping for it in the out-file (if you have a out file of a SCF calculation available). You will see something like this:
--> new fermi energy : 0.361941 htr
Take this number and convert it to eV ()
Now we will calculate a Fermi surface. Within this tutorial we have only limited visualization capabilities, so we will only make a 2D slice through the Fermi surface which we can visualize using gnuplot.
Open `WF1.win` and insert the following three lines:
kslice = true
kslice_task = fermi_lines
fermi_energy = 9.85
Then run the postw90.x postprocessing tool of wannier90.
postw90.x WF1
This will generate the data needed for producing the plotting
gnuplot -e "set term pdf; set output 'fermi.pdf'" -c WF*.gnu
Now you will see a file fermi.pdf in your list of files, which you can inspect. Wannier90 can also generate 3D data which can be visualized e.g. with the Fermisurfer to generate an image like this
Learn more
This tutorial only shows very basic functionality of Wannier90, but focusses on the FLEUR interface. Hence, many subjects like disentangling, convergence of the wannier90 localization process etc are not discussed at all. We refer to the Wannier90 documentation for these.
A good starting point for further exploration is: * The Wannier90 docu * The FLEUR docu on the wannier interface
(Note: this tutorial is based on work of D.Go and F.Freimuth).