Finding initial projections based on symmetry indicators
This tutorial explains how to find suitable projections for a given DFT bandstructure and specified frozen window.
Theory and background
Automating Wannier Functions Inspired by Topological Quantum Chemistry
Wannier functions (WFs) have become a popular and powerful technique for studying the diverse electronic properties of materials. However, the construction of high-quality WFs remains a non-trivial task, often requiring manual intervention and expertise, which obstructs high-throughput calculations using WF-based methods. A crucial aspect of WF construction is the selection of initial projections, which serve as an initial guess for the WFs. Traditionally, these projections are chosen based on knowledge of the orbital character of the bands under study (within a frozen window), which requires visual inspection of the bands, and hence quite some expertise and time.
In this tutorial, I present an alternative approach for selecting initial projections based on the symmetry indicators (irreducible representations) of DFT Bloch bands. The method is inspired by Topological Quantum Chemistry (TQC) [Nature 547, 298 (2017)]. Briefly, from TCS, we know the following :
If we take an isolated set of bands, and compute the irreducible representations (irreps) of the bands at high-symmetry points in the Brillouin zone, we can determine the topological properties of the bands.
If the same set of irreps can be generated by placing Wannier functions at certain Wyckoff positions in the unit cell, then the bands are either topologically trivial (if all WFs are on atomic sites) or an obstructed atomic limit (if some WFs are on empty Wyckoff positions).
If the irreps cannot be generated by placing WFs at Wyckoff positions, then the bands are topologically non-trivial.
In practice we are interested to describe precisely the bands in a certain energy range (frozen window), and we may include some more states to the outer window to ensure that the WFs are well-localized. The method presented here is based on the following steps:
Compute the irreps of the bands in the frozen window at all irreducible k-points in the Brillouin zone.
Define a set of trial wyckoff positions (WP) and trial orbitals for the WFs. (as a maximum one can include all Wyckoff positions and all orbitals (s,p,d,f) for each WP, but to speed up the calculations, one can include only the most relevant ones). each pair (WP, orbital) defines a set of trial projections, or an elementary band representation (EBR).
Compute the irreps of the trial WFs at all irreducible k-points in the Brillouin zone. (the symmetry indicators of each EBR)
Find such integer (non-negative) coefficients for the EBRs that for each k-point the sum of the EBRs yields at least all the irreps of the bands in the frozen window, and at the same time the sum of the EBRs fits into irreps of the bands in the outer window.
If some of EBRs are generated by WPs with free parameters (like (x,x,1/2) or (x,y,z)), then the free parameters can be optimized to maximize the minimal distance between different WPs.
The resulting trial positions are printed in the order of increasing number of Wannier functions. The input for wannier90.win file is also printed for the user’s convenience.
Furhter user can try different sets of trial projections, generate the wannier90.amn file and construct the wannier functions, either using Wannier90, or WannierBerri. It is preferred to use Symmetry adapted wannier functions (SAWFs) [Sakuma, PRB 87, 235109 (2013)] in this case.
In complicated cases some of the projections sets may result WFs that are not well localized, and do not provide good wannierisation. But there will be a number of sets, so by try-and-error one can find the best set of projections. Further, this try-and-error can be automated (work in progress).
This method does not require manual inspection of the orbital character of the bands, which makes it suitable for high-throughput calculations.
Import some modules and classes
[11]:
import warnings
from irrep.bandstructure import BandStructure
from fractions import Fraction
import sympy
from wannierberri.w90files import DMN, EIG, WIN
import wannierberri as wb
tested_version = '1.0.1'
if wb.__version__ != tested_version:
warnings.warn(f'This tutorial was tested with version {tested_version} of wannierberri')
from wannierberri.wannierise.projections_searcher import EBRsearcher
from wannierberri.wannierise.projections import Projection, ProjectionsSet
generate the dmn file
We will need the dmn file, which can be generated by pw2wannier90.x
However, if you do not have one, e.g. if you are using a code different from Quantum ESPRESSO, you can generate the file with the following code, using the irrep interface.
Notes: * actually the method is not propelly tested with the original dmn files, so it is anyway better to generate the dmn file with the code below. * the dmn file created now does not contain any information on the Wannier functions, but at this moment we do not need it. Later, when we know the projections, we can generate the amn file with the correct information.
[2]:
print("calculating DMN")
path = "./wannier-berri/tests/data/diamond/" # path to the data
bandstructure = BandStructure(prefix=path + "di", code="espresso",
Ecut=100, include_TR=False)
spacegroup = bandstructure.spacegroup
# spacegroup.show()
try:
dmn = DMN("diamond-only-bands")
except FileNotFoundError:
dmn = DMN(empty=True)
dmn.from_irrep(bandstructure)
dmn.to_npz("diamond-only-bands.dmn")
calculating DMN
calling w90 file with diamond-only-bands, dmn, tags=['D_wann_block_indices', '_NB', 'kpt2kptirr', 'kptirr', 'kptirr2kpt', 'kpt2kptirr_sym', '_NK', 'num_wann', 'comment', 'NKirr', 'Nsym', 'time_reversals'], read_npz=True, write_npz=True, kwargs={}
calling w90 file with ./wannier-berri/tests/data/diamond/diamond, eig, tags=['data'], read_npz=True, write_npz=True, kwargs={}
We will also need the eig and win files, so let’s load them.
[3]:
prefix = path + "diamond"
eig = EIG(prefix)
win = WIN(prefix)
calling w90 file with ./wannier-berri/tests/data/diamond/diamond, eig, tags=['data'], read_npz=True, write_npz=True, kwargs={}
Define the trial projections
Of course we know which are the good projections for diamond (from examples of Wannier90), but let’s pretend we do not know, and we want to find them. Let’s put some trial orbitals on some trial Wyckoff positions. In the complete case we should use all Wyckoff positions and all orbitals, but for the sake of simplicity we will use only a few.
[14]:
trial_projections = ProjectionsSet()
x, y, z = sympy.symbols('x y z')
F12 = Fraction(1, 2)
F14 = Fraction(1, 4)
F18 = Fraction(1, 8)
WP = [[0, 0, 0], [x, 0, 0], [F12, F12, F12], [F14, F14, F14], [F18, F18, F18], [0, x, z]]
# in principle, those should be all wyckoff position for the spacegroup
# but we will only consider a few random positions
positions = [",".join(str(y) for y in x) for x in WP]
print(positions)
for p in positions:
for o in ['s','sp3']:
proj = Projection(position_sym=p, orbital=o, spacegroup=spacegroup)
trial_projections.add(proj)
print("trial_projections")
print(trial_projections.write_with_multiplicities(orbit=False))
['0,0,0', 'x,0,0', '1/2,1/2,1/2', '1/4,1/4,1/4', '1/8,1/8,1/8', '0,x,z']
trial_projections
--------------------------------------------------------------------------------
1 X | 0,0,0:['s']
1 X | 0,0,0:['sp3']
1 X | x,0,0:['s']
1 X | x,0,0:['sp3']
1 X | 1/2,1/2,1/2:['s']
1 X | 1/2,1/2,1/2:['sp3']
1 X | 1/4,1/4,1/4:['s']
1 X | 1/4,1/4,1/4:['sp3']
1 X | 1/8,1/8,1/8:['s']
1 X | 1/8,1/8,1/8:['sp3']
1 X | 0,x,z:['s']
1 X | 0,x,z:['sp3']
total number of Wannier functions = 450
--------------------------------------------------------------------------------
the “total number of Wannier functions” is the number of Wannier functions if we include all projections once (which we will not do), so this number is meaningless here.
Create the EBRsearcher object
this object combines the needed data (the dmn, eig, win files ) and the energies of the frozen and outer window)
In the example the frozen window covers the valence bands, and the outer window includes also some conduction bands. )
[15]:
ebrsearcher = EBRsearcher(
win=win, # the win file is used only to know the k-points.
dmn=dmn,
eig=eig,
spacegroup=spacegroup,
trial_projections=trial_projections,
froz_min=-10,
froz_max=20,
outer_min=-20,
outer_max=50,
debug=False # set to True to see more printed information
)
Find combinations of EBRs
the combinations are just arrays of integers, where each integer is the number of times the corresponding EBR is used in the combination.
[16]:
combinations = ebrsearcher.find_combinations(max_num_wann=20)
print (combinations)
[array([1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]), array([1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0]), array([0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0]), array([0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0])]
for more information let’s print them in a more readable way, including the input for wannier90.win file.
[17]:
for c in combinations:
print(("+" * 80 + "\n") * 2)
print(trial_projections.write_with_multiplicities(c))
newset = trial_projections.get_combination(c)
newset.join_same_wyckoff()
newset.maximize_distance()
print(newset.write_wannier90(mod1=True))
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
--------------------------------------------------------------------------------
1 X | 0,0,0:['s']
total number of Wannier functions = 4
--------------------------------------------------------------------------------
finding num points from 1 projections
finding num points from 1 projections
minimal distance 1.1412636225506292
positions
[[0. 0. 0. ]
[0.5 0. 0. ]
[0. 0. 0.5]
[0. 0.5 0. ]]
distances
[[2.28 1.14 1.14 1.14]
[1.14 2.28 1.14 1.14]
[1.14 1.14 2.28 1.14]
[1.14 1.14 1.14 2.28]]
num_wann = 4
begin projections
f=0.000000000000, 0.000000000000, 0.000000000000: s
f=0.500000000000, 0.000000000000, 0.000000000000: s
f=0.000000000000, 0.000000000000, 0.500000000000: s
f=0.000000000000, 0.500000000000, 0.000000000000: s
end projections
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
--------------------------------------------------------------------------------
1 X | 0,0,0:['s']
1 X | 1/8,1/8,1/8:['s']
total number of Wannier functions = 6
--------------------------------------------------------------------------------
finding num points from 2 projections
finding num points from 2 projections
minimal distance 0.6988783843123343
positions
[[ 0. 0. 0. ]
[ 0.5 0. 0. ]
[ 0. 0. 0.5 ]
[ 0. 0.5 0. ]
[ 0.125 0.125 0.125]
[ 0.875 -0.125 -0.125]]
distances
[[2.28 1.14 1.14 1.14 0.7 0.7 ]
[1.14 2.28 1.14 1.14 0.7 0.7 ]
[1.14 1.14 2.28 1.14 0.7 0.7 ]
[1.14 1.14 1.14 2.28 0.7 0.7 ]
[0.7 0.7 0.7 0.7 2.28 1.4 ]
[0.7 0.7 0.7 0.7 1.4 2.28]]
num_wann = 6
begin projections
f=0.000000000000, 0.000000000000, 0.000000000000: s
f=0.500000000000, 0.000000000000, 0.000000000000: s
f=0.000000000000, 0.000000000000, 0.500000000000: s
f=0.000000000000, 0.500000000000, 0.000000000000: s
f=0.125000000000, 0.125000000000, 0.125000000000: s
f=0.875000000000, 0.875000000000, 0.875000000000: s
end projections
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
--------------------------------------------------------------------------------
1 X | 1/8,1/8,1/8:['sp3']
total number of Wannier functions = 8
--------------------------------------------------------------------------------
finding num points from 1 projections
finding num points from 1 projections
minimal distance 1.3977567686246695
positions
[[ 0.125 0.125 0.125]
[ 0.875 -0.125 -0.125]]
distances
[[2.28 1.4 ]
[1.4 2.28]]
num_wann = 8
begin projections
f=0.125000000000, 0.125000000000, 0.125000000000: sp3
f=0.875000000000, 0.875000000000, 0.875000000000: sp3
end projections
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
--------------------------------------------------------------------------------
1 X | 1/4,1/4,1/4:['s']
total number of Wannier functions = 8
--------------------------------------------------------------------------------
finding num points from 1 projections
finding num points from 1 projections
minimal distance 1.141263622550628
positions
[[ 0.25 0.25 0.25]
[ 1.25 -0.25 -0.25]
[ 0.25 0.25 -0.25]
[-0.25 1.25 -0.25]
[ 0.25 -0.25 0.25]
[-0.25 -0.25 1.25]
[-0.25 0.25 0.25]
[-0.25 -0.25 -0.25]]
distances
[[2.28 1.14 1.14 1.14 1.14 1.14 1.14 1.61]
[1.14 2.28 1.14 1.14 1.14 1.14 1.61 1.14]
[1.14 1.14 2.28 1.14 1.14 1.61 1.14 1.14]
[1.14 1.14 1.14 2.28 1.61 1.14 1.14 1.14]
[1.14 1.14 1.14 1.61 2.28 1.14 1.14 1.14]
[1.14 1.14 1.61 1.14 1.14 2.28 1.14 1.14]
[1.14 1.61 1.14 1.14 1.14 1.14 2.28 1.14]
[1.61 1.14 1.14 1.14 1.14 1.14 1.14 2.28]]
num_wann = 8
begin projections
f=0.250000000000, 0.250000000000, 0.250000000000: s
f=0.250000000000, 0.750000000000, 0.750000000000: s
f=0.250000000000, 0.250000000000, 0.750000000000: s
f=0.750000000000, 0.250000000000, 0.750000000000: s
f=0.250000000000, 0.750000000000, 0.250000000000: s
f=0.750000000000, 0.750000000000, 0.250000000000: s
f=0.750000000000, 0.250000000000, 0.250000000000: s
f=0.750000000000, 0.750000000000, 0.750000000000: s
end projections
So, we get the known result that we can put s orbitals on the bonds to get the Wannier functions for valence bands only (4 WFs), and we can put sp3 orbitals on the atoms to the Wannier functions for the valence and conduction bands (8 WFs).
But we also get some other combinations, which are one also may try to use.