Wannierisation (including SAWF)
This tutorial shows how to construct (Symmetry adapted) Wannier functions with WannierBerri, with magnetic symmetries. We will use the example of bcc Fe.
0. Compute QuantumEspresso files
in tthe tutorial repository only the input files for Quantum ESPRESSO are provided to obtain the necessary files one needs to run
pw.x < Fe_pw_scf_in > Fe_pw_scf_out
pw.x < Fe_pw_nscf_in > Fe_pw_nscf_out
wannier90.x -pp Fe
pw2wannier90.x < Fe_pw2wan_in > Fe_pw2wan_out
1. Setup
First import modules and set up the parallel environment
[1]:
import ray
# Initialize Ray with 8 CPU cores
# Do this only once at the beginning of your script. Initializing multiple times will lead to errors.
ray.init(num_cpus=8)
# If needed, you can do ray.shutdown() at the end of your script to clean up resources.
/home/stepan/github/WannierBerri-tutorial/.conda/lib/python3.12/site-packages/tqdm/auto.py:21: TqdmWarning: IProgress not found. Please update jupyter and ipywidgets. See https://ipywidgets.readthedocs.io/en/stable/user_install.html
from .autonotebook import tqdm as notebook_tqdm
2026-02-22 03:16:16,460 INFO util.py:154 -- Missing packages: ['ipywidgets']. Run `pip install -U ipywidgets`, then restart the notebook server for rich notebook output.
2026-02-22 03:16:18,264 INFO worker.py:1918 -- Started a local Ray instance. View the dashboard at http://127.0.0.1:8266
[1]:
| Python version: | 3.12.12 |
| Ray version: | 2.48.0 |
| Dashboard: | http://127.0.0.1:8266 |
[2]:
import os
from matplotlib import pyplot as plt
import scipy
import wannierberri as wb
import numpy as np
path_data = "./pwscf/" # adjust path if needed to point to the data in the tests fo wannier-berri repository
assert os.path.exists(path_data), f"Path {path_data} does not exist"
import irrep, spglib
print (f"using wannier-berri version: {wb.__version__}")
print (f"using irrep version: {irrep.__version__}")
print (f"using spglib version: {spglib.__version__}")
using wannier-berri version: 1.8.0
using irrep version: 2.6.1
using spglib version: 2.6.0
2. Read the bandstructure from Quantum ESPRESSO
we use the Bnadstructure object form irrep to read the bandstructure from Quantum ESPRESSO. It also can be used to read from VASP, ABINIT, GPAW etc. , see documentation of irrep for more details.
[3]:
from irrep.bandstructure import BandStructure
assert os.path.exists(path_data), f"Path {path_data} does not exist"
bandstructure = BandStructure(code='espresso', # to work with VASP or abinit please refer to the documentation of irrep
prefix=os.path.join(path_data, "Fe"),
magmom=[[0,0,1]], # set the magnetic moments for a magnetic system (units do not matter)
include_TR=True) # set include_TR=False if you do not want to include the symmetries involving time reversal (magnetic symmetries)
spacegroup = bandstructure.spacegroup
# spacegroup.show() # uncomment to see the detected symmetries of the system
3. Create the “w90files”
Here we use the Wannier90data object, to create the eig, symmetrizer, and mmn file (see below)
[4]:
w90data = wb.w90files.Wannier90data.from_bandstructure(
bandstructure=bandstructure,
files = ["eig","symmetrizer", "mmn"] # remove "mmn" if you do the next step (see below)
)
got irreducible=None, mp_grid=None, seedname=wannier90, files=['eig', 'symmetrizer', 'mmn'], projections=None, unk_grid=None, normalize=True
mpgrid = [4 4 4], 64
detected grid=(np.int64(4), np.int64(4), np.int64(4)), selected_kpoints=[ 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47
48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63]
self.irreducible=False
mpgrid = [4 4 4], 64
Shells found with weights [0.41728408] and tolerance 7.693398823301265e-16
NK= 64, selected_kpoints = [ 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47
48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63], kptirr = [ 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47
48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63]
3.1 Read the mmn file from pw2wannier90 (if needed)
irrep only reads the plane wave coefficients of the wavefunctions, and not the projections (EXCEPT for GPAW). This is totally fine for the symmetrizer, and ok for the amn file, but is very inaccurate for the mmn file.
This step is necessary if :
QE/VAST/Abinit with ultrasoft or PAW pseudopotentials.
Not needed if:
GPAW (projections are read by irrep)
QE/VAST/Abinit with norm-conserving pseudopotentials (plane wave coefficients are sufficient to construct the mmn file in this case)
[5]:
mmn = wb.w90files.mmn.MMN.from_w90_file(os.path.join(path_data,"Fe"), bkvec=w90data.bkvec)
w90data.set_file("mmn", mmn, overwrite=True)
3.2 Save the w90files with to npz for later use (optional)
You can save the created w90files to npz format for later use, to save time. Later thay can be loaded
[6]:
w90data.to_npz(os.path.join(path_data,"Fe_w90data"))
w90data_loaded = wb.w90files.Wannier90data.from_npz(os.path.join(path_data,"Fe_w90data"))
saving to ./pwscf/Fe_w90data.sawf.npz :
saving to ./pwscf/Fe_w90data.bkvec.npz :
saving to ./pwscf/Fe_w90data.chk.npz :
saving to ./pwscf/Fe_w90data.eig.npz :
saving to ./pwscf/Fe_w90data.mmn.npz :
files = ['mmn', 'eig', 'amn', 'bkvec']
Trying to read file mmn from npz ./pwscf/Fe_w90data.mmn.npz
setting file mmn from npz ./pwscf/Fe_w90data.mmn.npz as <wannierberri.w90files.mmn.MMN object at 0x7b09bd0dd5b0>
Trying to read file eig from npz ./pwscf/Fe_w90data.eig.npz
setting file eig from npz ./pwscf/Fe_w90data.eig.npz as <wannierberri.w90files.eig.EIG object at 0x7b0c7c14f080>
Trying to read file amn from npz ./pwscf/Fe_w90data.amn.npz
Trying to read file bkvec from npz ./pwscf/Fe_w90data.bkvec.npz
setting file bkvec from npz ./pwscf/Fe_w90data.bkvec.npz as <wannierberri.w90files.bkvectors.BKVectors object at 0x7b09b37abbf0>
/home/stepan/github/wannier-berri/wannierberri/w90files/w90data.py:371: UserWarning: file ./pwscf/Fe_w90data.amn.npz not found, cannot read amn file ([Errno 2] No such file or directory: './pwscf/Fe_w90data.amn.npz').
Set it manually, if needed
warnings.warn(f"file {seedname}.{f}.npz not found, cannot read {f} file ({e}).\n Set it manually, if needed")
4 Choose projections
This is similar to wannier90. However, in this case the projections already include the symmetry information.
‘position_num’ may contain only one position of the ornit, in this case the others will be generated by symmetry. However, this should be used with caution, as the generated positions may differ by a lattice vector, which may cause problems in some cases.
4.1 set the amn and symmetry of future Wannier functions
The last line serves two purposes:
creates and sets the amn file - initial guess for the Wannier functions
sets how the Wannier functions should transform under the symmetry operations.
Both are based on the projections_set object.
[7]:
from wannierberri.symmetry.projections import Projection, ProjectionsSet
# now set the transformations of WFs. Make sure, the projections are consistent with the amn file
proj_s = Projection(position_num = [[0,0,0]], orbital='s', spacegroup=spacegroup)
proj_p = Projection(position_num = [[0,0,0]], orbital='p', spacegroup=spacegroup)
proj_d = Projection(position_num = [[0,0,0]], orbital='d', spacegroup=spacegroup)
projections_set = ProjectionsSet(projections=[proj_s, proj_p, proj_d])
w90data.set_projections(projections_set, bandstructure=bandstructure)
finding num points from 3 projections
Creating amn. Using projections_set
ProjectionsSet with 9 Wannier functions and 0 free variables
Projection 0, 0, 0:['s'] with 1 Wannier functions on 1 points (1 per site)
Projection 0, 0, 0:['p'] with 3 Wannier functions on 1 points (3 per site)
Projection 0, 0, 0:['d'] with 5 Wannier functions on 1 points (5 per site)
5. Wannierise
Note, that in wannierberri there is not separation into “disentangle” and “wannierise” steps, as in wannier90. Both procedures are done together at every iteration. If you want to do only disentanglement, you can set localise=False. then num_iter will mean the number of disentanglement iterations.
[8]:
froz_max=25
_ = wb.wannierise.wannierise(
w90data=w90data,
froz_min=4,
froz_max=froz_max,
outer_min=0,
outer_max=100,
print_progress_every=20,
num_iter=100,
conv_tol=1e-6,
localise=True,
sitesym=True,
)
####################################################################################################
starting WFs
----------------------------------------------------------------------------------------------------
wannier centers and spreads
----------------------------------------------------------------------------------------------------
0.000000000000 0.000000000000 0.000000000000 | 20.981979926676
0.000000000000 0.000000000000 0.000000000000 | 17.481272996893
0.000000000000 0.000000000000 0.000000000000 | 21.205837708017
0.000000000000 0.000000000000 0.000000000000 | 21.342570639751
0.000000000000 0.000000000000 0.000000000000 | 23.929723510145
0.000000000000 0.000000000000 0.000000000000 | 20.588261495264
0.000000000000 0.000000000000 0.000000000000 | 23.929723510145
0.000000000000 0.000000000000 0.000000000000 | 20.588261495264
0.000000000000 0.000000000000 0.000000000000 | 22.270813393077
0.000000000000 0.000000000000 0.000000000000 | 15.155350709404
0.000000000000 0.000000000000 0.000000000000 | 18.964394055881
0.000000000000 0.000000000000 0.000000000000 | 22.497006646651
0.000000000000 0.000000000000 0.000000000000 | 18.964394055881
0.000000000000 0.000000000000 0.000000000000 | 22.497006646651
0.000000000000 0.000000000000 0.000000000000 | 20.688871040090
0.000000000000 0.000000000000 0.000000000000 | 21.018422522393
0.000000000000 0.000000000000 0.000000000000 | 21.975981977632
0.000000000000 0.000000000000 0.000000000000 | 23.737355990674
----------------------------------------------------------------------------------------------------
0.000000000000 0.000000000000 0.000000000000 | 377.817228320486 <- sum
maximal spread = 23.929723510145
####################################################################################################
####################################################################################################
Iteration 0 (from wannierizer)
----------------------------------------------------------------------------------------------------
wannier centers and spreads
----------------------------------------------------------------------------------------------------
0.000000000000 0.000000000000 0.000000000000 | 9.695443700873
0.000000000000 0.000000000000 0.000000000000 | 11.023837184599
0.000000000000 0.000000000000 0.000000000000 | 11.178168524439
0.000000000000 0.000000000000 0.000000000000 | 14.557934158409
0.000000000000 0.000000000000 0.000000000000 | 12.745155904959
0.000000000000 0.000000000000 0.000000000000 | 11.642554627883
0.000000000000 0.000000000000 0.000000000000 | 12.745155904959
0.000000000000 0.000000000000 0.000000000000 | 11.642554627883
0.000000000000 0.000000000000 0.000000000000 | 12.644254725857
0.000000000000 0.000000000000 0.000000000000 | 7.105420992061
0.000000000000 0.000000000000 0.000000000000 | 9.752314791331
0.000000000000 0.000000000000 0.000000000000 | 11.344262807596
0.000000000000 0.000000000000 0.000000000000 | 9.752314791331
0.000000000000 0.000000000000 0.000000000000 | 11.344262807596
0.000000000000 0.000000000000 0.000000000000 | 11.086601936003
0.000000000000 0.000000000000 0.000000000000 | 11.713154450500
0.000000000000 0.000000000000 0.000000000000 | 11.354228517606
0.000000000000 0.000000000000 0.000000000000 | 11.755896123511
----------------------------------------------------------------------------------------------------
0.000000000000 0.000000000000 0.000000000000 | 203.083516577397 <- sum
maximal spread = 14.557934158409
standard deviation = 0.0
####################################################################################################
####################################################################################################
Iteration 20 (from wannierizer)
----------------------------------------------------------------------------------------------------
wannier centers and spreads
----------------------------------------------------------------------------------------------------
0.000000000000 0.000000000000 0.000000000000 | 1.348582299678
0.000000000000 0.000000000000 0.000000000000 | 1.389795515256
0.000000000000 0.000000000000 0.000000000000 | 1.422163757020
0.000000000000 0.000000000000 0.000000000000 | 1.497286103259
0.000000000000 0.000000000000 0.000000000000 | 1.393750879851
0.000000000000 0.000000000000 0.000000000000 | 1.491353763982
0.000000000000 0.000000000000 0.000000000000 | 1.393750879851
0.000000000000 0.000000000000 0.000000000000 | 1.491353763982
0.000000000000 0.000000000000 0.000000000000 | 0.479241892284
0.000000000000 0.000000000000 0.000000000000 | 0.435110992388
0.000000000000 0.000000000000 0.000000000000 | 0.432289142928
0.000000000000 0.000000000000 0.000000000000 | 0.406130651457
0.000000000000 0.000000000000 0.000000000000 | 0.432289142928
0.000000000000 0.000000000000 0.000000000000 | 0.406130651457
0.000000000000 0.000000000000 0.000000000000 | 0.472048782391
0.000000000000 0.000000000000 0.000000000000 | 0.434139865644
0.000000000000 0.000000000000 0.000000000000 | 0.432848818697
0.000000000000 0.000000000000 0.000000000000 | 0.399576584524
----------------------------------------------------------------------------------------------------
0.000000000000 0.000000000000 0.000000000000 | 15.757843487578 <- sum
maximal spread = 1.497286103259
standard deviation = 0.00027294468598274105
####################################################################################################
####################################################################################################
Iteration 40 (from wannierizer)
----------------------------------------------------------------------------------------------------
wannier centers and spreads
----------------------------------------------------------------------------------------------------
0.000000000000 0.000000000000 0.000000000000 | 1.347224358570
0.000000000000 0.000000000000 0.000000000000 | 1.389634727117
0.000000000000 0.000000000000 0.000000000000 | 1.422282922495
0.000000000000 0.000000000000 0.000000000000 | 1.497610058586
0.000000000000 0.000000000000 0.000000000000 | 1.393357997897
0.000000000000 0.000000000000 0.000000000000 | 1.491125696332
0.000000000000 0.000000000000 0.000000000000 | 1.393357997897
0.000000000000 0.000000000000 0.000000000000 | 1.491125696332
0.000000000000 0.000000000000 0.000000000000 | 0.479205014200
0.000000000000 0.000000000000 0.000000000000 | 0.435118689102
0.000000000000 0.000000000000 0.000000000000 | 0.432544241226
0.000000000000 0.000000000000 0.000000000000 | 0.406142082296
0.000000000000 0.000000000000 0.000000000000 | 0.432544241226
0.000000000000 0.000000000000 0.000000000000 | 0.406142082296
0.000000000000 0.000000000000 0.000000000000 | 0.472040986643
0.000000000000 0.000000000000 0.000000000000 | 0.434088025033
0.000000000000 0.000000000000 0.000000000000 | 0.432671725092
0.000000000000 0.000000000000 0.000000000000 | 0.399594227660
----------------------------------------------------------------------------------------------------
0.000000000000 0.000000000000 0.000000000000 | 15.755810769997 <- sum
maximal spread = 1.497610058586
standard deviation = 1.0875991549849182e-05
####################################################################################################
####################################################################################################
Iteration 60 (from wannierizer)
----------------------------------------------------------------------------------------------------
wannier centers and spreads
----------------------------------------------------------------------------------------------------
0.000000000000 0.000000000000 0.000000000000 | 1.347211202246
0.000000000000 0.000000000000 0.000000000000 | 1.389639352551
0.000000000000 0.000000000000 0.000000000000 | 1.422446446208
0.000000000000 0.000000000000 0.000000000000 | 1.497638425144
0.000000000000 0.000000000000 0.000000000000 | 1.393259887296
0.000000000000 0.000000000000 0.000000000000 | 1.491124710911
0.000000000000 0.000000000000 0.000000000000 | 1.393259887296
0.000000000000 0.000000000000 0.000000000000 | 1.491124710911
0.000000000000 0.000000000000 0.000000000000 | 0.479173634118
0.000000000000 0.000000000000 0.000000000000 | 0.435118560379
0.000000000000 0.000000000000 0.000000000000 | 0.432568899662
0.000000000000 0.000000000000 0.000000000000 | 0.406149193813
0.000000000000 0.000000000000 0.000000000000 | 0.432568899662
0.000000000000 0.000000000000 0.000000000000 | 0.406149193813
0.000000000000 0.000000000000 0.000000000000 | 0.472028930898
0.000000000000 0.000000000000 0.000000000000 | 0.434078133912
0.000000000000 0.000000000000 0.000000000000 | 0.432637407642
0.000000000000 0.000000000000 0.000000000000 | 0.399598398625
----------------------------------------------------------------------------------------------------
0.000000000000 0.000000000000 0.000000000000 | 15.755775875088 <- sum
maximal spread = 1.497638425144
standard deviation = 3.749179382775483e-06
####################################################################################################
####################################################################################################
Iteration 80 (from wannierizer)
----------------------------------------------------------------------------------------------------
wannier centers and spreads
----------------------------------------------------------------------------------------------------
0.000000000000 0.000000000000 0.000000000000 | 1.347214947968
0.000000000000 0.000000000000 0.000000000000 | 1.389640569526
0.000000000000 0.000000000000 0.000000000000 | 1.422494060479
0.000000000000 0.000000000000 0.000000000000 | 1.497642247863
0.000000000000 0.000000000000 0.000000000000 | 1.393233075489
0.000000000000 0.000000000000 0.000000000000 | 1.491127311604
0.000000000000 0.000000000000 0.000000000000 | 1.393233075489
0.000000000000 0.000000000000 0.000000000000 | 1.491127311604
0.000000000000 0.000000000000 0.000000000000 | 0.479165289696
0.000000000000 0.000000000000 0.000000000000 | 0.435118524666
0.000000000000 0.000000000000 0.000000000000 | 0.432573306432
0.000000000000 0.000000000000 0.000000000000 | 0.406150849447
0.000000000000 0.000000000000 0.000000000000 | 0.432573306432
0.000000000000 0.000000000000 0.000000000000 | 0.406150849447
0.000000000000 0.000000000000 0.000000000000 | 0.472026578114
0.000000000000 0.000000000000 0.000000000000 | 0.434076349307
0.000000000000 0.000000000000 0.000000000000 | 0.432630352201
0.000000000000 0.000000000000 0.000000000000 | 0.399599188666
----------------------------------------------------------------------------------------------------
0.000000000000 0.000000000000 0.000000000000 | 15.755777194429 <- sum
maximal spread = 1.497642247863
standard deviation = 9.778359094176371e-07
####################################################################################################
Converged after 80 iterations
####################################################################################################
Final state (from wannierizer)
----------------------------------------------------------------------------------------------------
wannier centers and spreads
----------------------------------------------------------------------------------------------------
0.000000000000 0.000000000000 0.000000000000 | 1.347214947968
0.000000000000 0.000000000000 0.000000000000 | 1.389640569526
0.000000000000 0.000000000000 0.000000000000 | 1.422494060479
0.000000000000 0.000000000000 0.000000000000 | 1.497642247863
0.000000000000 0.000000000000 0.000000000000 | 1.393233075489
0.000000000000 0.000000000000 0.000000000000 | 1.491127311604
0.000000000000 0.000000000000 0.000000000000 | 1.393233075489
0.000000000000 0.000000000000 0.000000000000 | 1.491127311604
0.000000000000 0.000000000000 0.000000000000 | 0.479165289696
0.000000000000 0.000000000000 0.000000000000 | 0.435118524666
0.000000000000 0.000000000000 0.000000000000 | 0.432573306432
0.000000000000 0.000000000000 0.000000000000 | 0.406150849447
0.000000000000 0.000000000000 0.000000000000 | 0.432573306432
0.000000000000 0.000000000000 0.000000000000 | 0.406150849447
0.000000000000 0.000000000000 0.000000000000 | 0.472026578114
0.000000000000 0.000000000000 0.000000000000 | 0.434076349307
0.000000000000 0.000000000000 0.000000000000 | 0.432630352201
0.000000000000 0.000000000000 0.000000000000 | 0.399599188666
----------------------------------------------------------------------------------------------------
0.000000000000 0.000000000000 0.000000000000 | 15.755777194429 <- sum
maximal spread = 1.497642247863
standard deviation = 9.778359094176371e-07
####################################################################################################
####################################################################################################
Final state (from chk)
----------------------------------------------------------------------------------------------------
wannier centers and spreads
----------------------------------------------------------------------------------------------------
-0.000000000046 0.000000000392 -0.000000000180 | 1.347214998306
0.000000000047 0.000000000226 -0.000000000011 | 1.389640582026
-0.000000000266 0.000000003457 -0.000000001014 | 1.422494615901
-0.000000000106 0.000000002668 0.000000000016 | 1.497642281188
-0.000000000074 0.000000001104 -0.000000001468 | 1.393232763254
0.000000000173 0.000000002392 0.000000000035 | 1.491127351215
0.000000000235 0.000000000220 0.000000000509 | 1.393232762304
0.000000000480 0.000000000559 0.000000000805 | 1.491127350685
-0.000000000010 0.000000000061 -0.000000000058 | 0.479165193349
-0.000000000006 0.000000000041 0.000000000006 | 0.435118524347
-0.000000000025 0.000000000202 0.000000000008 | 0.432573352941
-0.000000000011 0.000000000237 -0.000000000032 | 0.406150866697
-0.000000000022 -0.000000000005 -0.000000000001 | 0.432573352918
-0.000000000029 0.000000000000 -0.000000000012 | 0.406150866674
0.000000000012 0.000000000115 -0.000000000010 | 0.472026554856
0.000000000036 0.000000000240 0.000000000004 | 0.434076332439
0.000000000043 0.000000000049 -0.000000000009 | 0.432630276020
-0.000000000003 -0.000000000004 -0.000000000019 | 0.399599196509
----------------------------------------------------------------------------------------------------
0.000000000428 0.000000011954 -0.000000001429 | 15.755777221627 <- sum
maximal spread = 1.497642281188
####################################################################################################
time for creating wannierizer 0.2476518154144287
time for iterations 2.3662729263305664
time for updating 1.9084126949310303
total time for wannierization 5.285660028457642
saving to wannier90.chk :
6. Create System_w90 object
[9]:
system = wb.system.System_w90(w90data= w90data, berry=True, transl_inv_JM=True,
symmetrize=True)
# optionally - save it for later use
system.save_npz("Fe_system")
irreducible : False, symmetrize set to True
setting Rvec
setting AA..
setting AA - OK
recentering JM - OK
Real-space lattice:
[[ 1.4349963 1.4349963 1.4349963]
[-1.4349963 1.4349963 1.4349963]
[-1.4349963 -1.4349963 1.4349963]]
Number of wannier functions: 18
Number of R points: 89
Recommended size of FFT grid [4 4 4]
num_blocks_left = 3, num_blocks_right = 3
number o R-vectors before symmetrization: 89
number o R-vectors after symmetrization: 89
saving system of class System_R to Fe_system
properties: ['num_wann', 'real_lattice', 'iRvec', 'periodic', 'is_phonon', 'wannier_centers_cart', 'pointgroup']
saving num_wann
saving num_wann to Fe_system/num_wann.npz
- Ok!
saving real_lattice
saving real_lattice to Fe_system/real_lattice.npz
- Ok!
saving iRvec
saving iRvec to Fe_system/iRvec.npz
- Ok!
saving periodic
saving periodic to Fe_system/periodic.npz
- Ok!
saving is_phonon
saving is_phonon to Fe_system/is_phonon.npz
- Ok!
saving wannier_centers_cart
saving wannier_centers_cart to Fe_system/wannier_centers_cart.npz
- Ok!
saving pointgroup
saving pointgroup to Fe_system/pointgroup.npz
- Ok!
saving Ham - Ok!
saving AA - Ok!
7. Bands along path
7.1 calculate bands
[10]:
# all kpoints given in reduced coordinates
path=wb.Path.from_nodes(system,
nodes=[
[0.0000, 0.0000, 0.0000 ], # G
[0.500 ,-0.5000, -0.5000], # H
[0.7500, 0.2500, -0.2500], # P
[0.5000, 0.0000, -0.5000], # N
[0.0000, 0.0000, 0.000 ]
] , # G
labels=["G","H","P","N","G"],
length=200 ) # length [ Ang] ~= 2*pi/dk
bands_path=wb.evaluate_k_path(system, path=path)
Starting run()
Using the follwing calculators :
############################################################
'tabulate' : <wannierberri.calculators.tabulate.TabulatorAll object at 0x7b0bd8ae73e0> :
TabulatorAll - a pack of all k-resolved calculators (Tabulators)
Includes the following tabulators :
--------------------------------------------------
"Energy" : <wannierberri.calculators.tabulate.Energy object at 0x7b0bd8e80fb0> : calculator not described
--------------------------------------------------
############################################################
Calculation along a path - checking calculators for compatibility
tabulate <wannierberri.calculators.tabulate.TabulatorAll object at 0x7b0bd8ae73e0>
All calculators are compatible
Symmetrization switched off for Path
Grid is regular
The set of k points is a Path() with 215 points and labels {0: 'G', 70: 'H', 130: 'P', 165: 'N', 214: 'G'}
generating K_list
Done
Done, sum of weights:215.0
############################################################
Iteration 0 out of 0
processing 215 K points : using 8 processes.
# K-points calculated Wall time (sec) Est. remaining (sec) Est. total (sec)
/home/stepan/github/wannier-berri/wannierberri/grid/path.py:272: UserWarning: symmetry is not used for a tabulation along path
warnings.warn("symmetry is not used for a tabulation along path")
time for processing 215 K-points on 8 processes: 1.9090 ; per K-point 0.0089 ; proc-sec per K-point 0.0710
time1 = 4.76837158203125e-07
Totally processed 215 K-points
run() finished
7.2 plot bands
[11]:
# plot the bands and compare with pw
# EF = 12
A = np.loadtxt("./pwscf/Fe_bands_pw.dat")
bohr_ang = scipy.constants.physical_constants['Bohr radius'][0] / 1e-10
alatt = 5.4235* bohr_ang
A[:,0]*= 2*np.pi/alatt
A[:,1] = A[:,1]
plt.scatter(A[:,0], A[:,1], c="black", s=5)
bands_path.plot_path_fat(path,
quantity=None,
# save_file="Fe_bands.pdf",
Eshift=0,
Emin=-10, Emax=50,
iband=None,
mode="fatband",
fatfactor=20,
cut_k=False,
linecolor="red",
close_fig=False,
show_fig=False,
label=f"WB"
)
plt.ylim(4, 40)
plt.hlines(froz_max, 0, A[-1,0], linestyles="dashed", label="frozen window", color="black")
plt.hlines(12.6, 0, A[-1,0], linestyles="dashed", label="Fermi level", color="green")
plt.legend()
plt.savefig("Fe_bands.pdf")
8. AHC and Ohmic conductivity
8.1 calculate
[13]:
results_grid = {}
efermi = np.linspace(12.4,12.8,1001)
param = dict(Efermi=efermi)
calculators_grid = {
"CDOS": wb.calculators.static.CumDOS(**param),
"ohmic": wb.calculators.static.Ohmic_FermiSea(**param),
"ahc_internal": wb.calculators.static.AHC(kwargs_formula={"external_terms":False}, **param),
"ahc_external": wb.calculators.static.AHC(kwargs_formula={"internal_terms":False}, **param ),
}
grid = wb.Grid(system, NKFFT=6, NK=48)
result_grid = wb.run(system,
grid,
calculators_grid,
fout_name="Fe_grid",
adpt_num_iter=0,
symmetrize=False, # we do not symmetrize here so that we can chaeck how symmetric are the WFs
use_irred_kpt=False,
)
# plot the bands to compare with pw
Starting run()
Using the follwing calculators :
############################################################
'CDOS' : <wannierberri.calculators.static.CumDOS object at 0x7b0bd863e060> : Cumulative density of states
'ohmic' : <wannierberri.calculators.static.Ohmic_FermiSea object at 0x7b0bd863f410> : Ohmic conductivity (:math:`S/m`)
| With Fermi sea integral. Eq(31) in `Ref <https://www.nature.com/articles/s41524-021-00498-5>`__
| Output: :math:`\sigma_{\alpha\beta} = e^2/\hbar \tau \int [dk] \partial_\beta v_\alpha f`for \tau=1fs| Instruction: :math:`j_\alpha = \sigma_{\alpha\beta} E_\beta`
'ahc_internal' : <wannierberri.calculators.static.AHC object at 0x7b0bd8721ca0> : Anomalous Hall conductivity (:math:`s^3 \cdot A^2 / (kg \cdot m^3) = S/m`)
| With Fermi sea integral Eq(11) in `Ref <https://www.nature.com/articles/s41524-021-00498-5>`__
| Output: :math:`O = -e^2/\hbar \int [dk] \Omega f`
| Instruction: :math:`j_\alpha = \sigma_{\alpha\beta} E_\beta = \epsilon_{\alpha\beta\delta} O_\delta E_\beta`
'ahc_external' : <wannierberri.calculators.static.AHC object at 0x7b0bd8722600> : Anomalous Hall conductivity (:math:`s^3 \cdot A^2 / (kg \cdot m^3) = S/m`)
| With Fermi sea integral Eq(11) in `Ref <https://www.nature.com/articles/s41524-021-00498-5>`__
| Output: :math:`O = -e^2/\hbar \int [dk] \Omega f`
| Instruction: :math:`j_\alpha = \sigma_{\alpha\beta} E_\beta = \epsilon_{\alpha\beta\delta} O_\delta E_\beta`
############################################################
Calculation on grid - checking calculators for compatibility
CDOS <wannierberri.calculators.static.CumDOS object at 0x7b0bd863e060>
ohmic <wannierberri.calculators.static.Ohmic_FermiSea object at 0x7b0bd863f410>
ahc_internal <wannierberri.calculators.static.AHC object at 0x7b0bd8721ca0>
ahc_external <wannierberri.calculators.static.AHC object at 0x7b0bd8722600>
All calculators are compatible
Grid is regular
The set of k points is a Grid() with NKdiv=[8 8 8], NKFFT=[6 6 6], NKtot=[48 48 48]
generating K_list
Done in 0.0022728443145751953 s
K_list contains 512 Irreducible points(100.0%) out of initial 8x8x8=512 grid
Done, sum of weights:1.0
############################################################
Iteration 0 out of 0
processing 512 K points : using 8 processes.
# K-points calculated Wall time (sec) Est. remaining (sec) Est. total (sec)
80 5.2 28.2 33.5
160 10.4 22.9 33.4
240 15.5 17.6 33.1
320 20.7 12.4 33.1
400 26.0 7.3 33.3
480 31.4 2.1 33.5
time for processing 512 K-points on 8 processes: 33.4957 ; per K-point 0.0654 ; proc-sec per K-point 0.5234
time1 = 4.76837158203125e-07
Totally processed 512 K-points
run() finished
8.2 Plot
[14]:
def plotxyz(axes, x, data, pre = "", label="",ls="-"):
for i in range(3):
ax =axes[i]
ax.plot(x, data[:,i], ls, label=label)
ax.set_title(f"{pre}{'xyz'[i]}")
quantities = ["ahc_internal", "ahc_external","ohmic"]
nfig = len(quantities)
for quantity in quantities:
fig = None
res = result_grid.results[quantity]
data = res.data
E = res.Energies[0]
if fig is None:
if data.ndim == 2:
nfigx = 1
nfigy = data.shape[1]
elif data.ndim == 3:
nfigx = data.shape[1]
nfigy = data.shape[2]
fig, axes = plt.subplots(nfigx, nfigy, figsize=(5*nfigy,5*nfigx))
if nfigx ==1:
plotxyz(axes, E, data)
else:
for i in range(nfigx):
plotxyz(axes[i], E, data[:,i], pre="xyz"[i], )
for ax in axes.flat:
ax.legend()
fig.suptitle(quantity)
plt.show()
plt.close()
/tmp/ipykernel_65542/1474735278.py:31: UserWarning: No artists with labels found to put in legend. Note that artists whose label start with an underscore are ignored when legend() is called with no argument.
ax.legend()
[ ]: