DIY tutorial
As we learned from the basic tutorial, in WannierBerri calculations are done by `Calculator
<https://wannier-berri.org/docs/calculators.html#wannierberri.calculators.Calculator>`__s. In this tutorial we will se how to create our own calculators.
As an example let’s consider anomalous Nernst conductivity (Xiao et al. 2006 ). In the zero-temperature limit \(\alpha_{\alpha\beta}^{\rm ANE}\) may be obtained from AHC \(\sigma_{\alpha\beta}(\epsilon)^{\rm AHE}\) evaluated over a dense grid of Fermi levels \(\epsilon\)
where \(f(\varepsilon)=1/\left(1+e^\frac{\varepsilon-\mu}{k_{\rm B}T}\right)\) and the AHC is defined as a Fermi-sea integral of Berry curvature
\begin{equation} \sigma^{\rm AHC}_{xy} = -\frac{e^2}{\hbar} \sum_n \int\frac{d\mathbf{k}}{(2\pi)^3} \Omega^n_\gamma f(\epsilon_n) \label{eq:ahc}\tag{2} \end{equation}
In the zero-temperature limit it reduced to
Where we omit the dimensional factor for simplicity.
Thus, there are two ways of calculating ANE :
via eq (\ref{eq:ane-ahc})
via eq (\ref{eq:ane})
let’s try eq (\ref{eq:ane-ahc}) first, because AHC is already implemented
Warning : this is an advanced tutorial. You something may not work from first try, you need to work with documentation to resolve the problem, but do not hesitate to ask the TA if you cannot.
[ ]:
# Preliminary
# Set environment variables - not mandatory but recommended if you use Parallel()
#import os
#os.environ['OPENBLAS_NUM_THREADS'] = '1'
#os.environ['MKL_NUM_THREADS'] = '1'
import wannierberri as wberri
import wannierberri.models
print (f"Using WannierBerri version {wberri.__version__}")
import numpy as np
import scipy
import matplotlib.pyplot as plt
from termcolor import cprint
# This block is needed if you run this cell for a second time
# because one cannot initiate two parallel environments at a time
try:
parallel.shutdown()
except NameError:
pass
# Chiose one of the options:
parallel = wberri.Parallel(num_cpus=2)
#parallel = wberri.Parallel() # automatic detection
#parallel = wberri.Serial()
Using WannierBerri version 1.2.1
initializing ray with {'num_cpus': 2}
2025-05-14 11:55:56,422 INFO worker.py:1879 -- Started a local Ray instance. View the dashboard at http://127.0.0.1:8265
(raylet) [2025-05-14 15:04:56,665 E 18313 18313] (raylet) node_manager.cc:3287: 2 Workers (tasks / actors) killed due to memory pressure (OOM), 0 Workers crashed due to other reasons at node (ID: 41f5f8eb248e8ca703d3f988466034c2421af859bb67197d42528a9a, IP: 128.179.188.185) over the last time period. To see more information about the Workers killed on this node, use `ray logs raylet.out -ip 128.179.188.185`
(raylet)
(raylet) Refer to the documentation on how to address the out of memory issue: https://docs.ray.io/en/latest/ray-core/scheduling/ray-oom-prevention.html. Consider provisioning more memory on this node or reducing task parallelism by requesting more CPUs per task. To adjust the kill threshold, set the environment variable `RAY_memory_usage_threshold` when starting Ray. To disable worker killing, set the environment variable `RAY_memory_monitor_refresh_ms` to zero.
[2]:
# Instead of Wannier functions, we work with a Hlldane model here https://wannier-berri.org/docs/models.html
model=wberri.models.Haldane_ptb()
system = wberri.System_PythTB(model)
number of wannier functions: 2
R=0 found at position(s) [[3]]
shape of Ham_R = (7, 2, 2)
Real-space lattice:
[[1. 0. 0. ]
[0.5 0.8660254 0. ]
[0. 0. 1. ]]
Number of wannier functions: 2
Number of R points: 7
Recommended size of FFT grid [3 3 1]
Reading the system from PythTB finished successfully
anomalous Hall conductivity
Now, let’s evaluate the AHC on a grid of Ef-points, and then take a derivtive
[3]:
calculators = {}
Efermi = np.linspace(-3,3,21)
# Set a grid
grid = wberri.Grid(system, length=30 ) # length [ Ang] ~= 2*pi/dk
calculators ["ahc"] = wberri.calculators.static.AHC(Efermi=Efermi,
tetra=True,
constant_factor=1,
kwargs_formula={"external_terms":False})
result_run_ahc = wberri.run(system,
grid=grid,
calculators = calculators,
parallel=parallel,
adpt_num_iter=0,
fout_name='Fe',
restart=False,
file_Klist="Klist_ahc.pickle", # needed to restart a calculation in future
print_Kpoints=False
)
Minimal symmetric FFT grid : [3 3 1]
Starting run()
Using the follwing calculators :
############################################################
'ahc' : <wannierberri.calculators.static.AHC object at 0x74552017eb70> : 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
ahc <wannierberri.calculators.static.AHC object at 0x74552017eb70>
All calculators are compatible
Grid is regular
The set of k points is a Grid() with NKdiv=[7 7 1], NKFFT=[5 5 1], NKtot=[35 35 1]
generating K_list
Done in 0.00014543533325195312 s
excluding symmetry-equivalent K-points from initial grid
Done in 0.0006546974182128906 s
K_list contains 49 Irreducible points(100.0%) out of initial 7x7x1=49 grid
Done, sum of weights:1.0000000000000007
processing 49 K points : using 2 processes.
# K-points calculated Wall time (sec) Est. remaining (sec)
time for processing 49 K-points on 2 processes: 2.4343 ; per K-point 0.0497 ; proc-sec per K-point 0.0994
time1 = 0.0007510185241699219
Totally processed 49 K-points
run() finished
[4]:
!ls
Fe-ahc_iter-0000.dat Fe-ane_iter-0000.npz tutorial-wb-DIY-solution.ipynb
Fe-ahc_iter-0000.npz Klist_ahc.pickle
Fe-ane_iter-0000.dat tutorial-wb-DIY.ipynb
[5]:
# taking derivative and plotting
ef = result_run_ahc.results["ahc"].Energies[0]
ahc_z = result_run_ahc.results["ahc"].data[:,2]
# take the derivative
d_ef=ef[1]-ef[0]
d_ahc_z = (ahc_z[1:]-ahc_z[:-1])/d_ef
efnew = (ef[1:]+ef[:-1])/2
plt.plot(efnew, d_ahc_z)
[5]:
[<matplotlib.lines.Line2D at 0x7454fe78dee0>]

Look into documentation
Look into documentation to see how the AHC calculator is defined : look for “Calculator” on https://docs.wannier-berri.org and press the “source” link.
One can see that it is based on the StaticCalculator, and only redefines the __init__
method. Namely, it redefines formula that stands under the integral (Berry curvature) and fder=0
meaning that the formula is weighted by the 0th derivative of the Fermi distribution.
copy the definition of AHC calculator and redefine the and below and define the init class
(hint : leave the factor the same, for simplicity) (hint : formula shouldbe imported from wannierberri.formula.covariant.DerOmega )
[6]:
from wannierberri.calculators.static import StaticCalculator
from wannierberri.formula.covariant import Omega
from wannierberri import factors
class ANC(StaticCalculator):
r"""Anomalous Nernst conductivity
| Output: :math:`O = \int [dk] \Omega f' `
"""
def __init__(self,
constant_factor=1, # Figure out the pre-factor yourself
**kwargs):
"""describe input parameters here"""
self.Formula = Omega
self.fder = 1
super().__init__(constant_factor=constant_factor, **kwargs)
Evaluate ANE using the new calculator and plot the result
[7]:
# insert the needed code below
calculators = {}
# Set a grid
grid = wberri.Grid(system, length=50 ) # length [ Ang] ~= 2*pi/dk
calculators ["ane"] = ANC(Efermi=Efermi,tetra=True)
result_run_ane = wberri.run(system,
grid=grid,
calculators = calculators,
parallel=parallel,
adpt_num_iter=0,
fout_name='Fe',
restart=False,
# file_Klist="Klist_ahc.pickle" # needed to restart a calculation in future
)
Minimal symmetric FFT grid : [3 3 1]
Starting run()
Using the follwing calculators :
############################################################
'ane' : <__main__.ANC object at 0x7455200d9fa0> : Anomalous Nernst conductivity
| Output: :math:`O = \int [dk] \Omega f' `
############################################################
Calculation on grid - checking calculators for compatibility
ane <__main__.ANC object at 0x7455200d9fa0>
All calculators are compatible
Grid is regular
The set of k points is a Grid() with NKdiv=[19 19 1], NKFFT=[3 3 1], NKtot=[57 57 1]
generating K_list
Done in 0.0009851455688476562 s
excluding symmetry-equivalent K-points from initial grid
Done in 0.003620147705078125 s
K_list contains 361 Irreducible points(100.0%) out of initial 19x19x1=361 grid
Done, sum of weights:0.999999999999993
processing 361 K points : using 2 processes.
# K-points calculated Wall time (sec) Est. remaining (sec)
/home/stepan/github/wannier-berri/wannierberri/grid/grid.py:223: UserWarning: the requested k-grid [58 58 50] was adjusted to [57 57 50].
warnings.warn(f" the requested k-grid {NK} was adjusted to {NKFFT * NKdiv}. ")
time for processing 361 K-points on 2 processes: 0.6829 ; per K-point 0.0019 ; proc-sec per K-point 0.0038
time1 = 0.0052797794342041016
Totally processed 361 K-points
run() finished
[8]:
ef = result_run_ane.results["ane"].Energies[0]
ane_z = result_run_ane.results["ane"].dataSmooth[:,2]
plt.scatter(efnew,d_ahc_z)
plt.plot(ef, ane_z)
plt.show()

Questions:
Make the k-points or Efermi denser. Will the agreement improve?
we calculated the AHC and ANE with zero refinement iterations, and the results matched well (or at least they should). If we calculate them separately, but with some refinement iterations, will the results match?
What if we use refinement itrations and run the two calculations together, in one
run()
call ?
FormulaProduct
Look at the definition of erryDipole_FermiSurf
https://wannier-berri.org/_modules/wannierberri/calculators/static.html
It evaluates
\begin{equation} D_{ab} =\sum_n \int\frac{d\mathbf{k}}{(2\pi)^3} v_a\Omega^n_b f(\epsilon_n) \label{eq:ahc}\tag{2} \end{equation}
also look at definition of the Feormula https://wannier-berri.org/docs/formula.html
you can see that it is base on a FormulaProduct. Essentially, by analogy you may define any tensor, e.g.
i.e. product of several quantities weighted by rge second derivative of Fermi distribution. You may write any analogous formula, if you want. Try to define corresponding FormulaProduct and StaticCalculator classes, and try to evaluate them.
[ ]:
[ ]:
Problem :
try to prove (analytically and numerically) that
and tell me what is the factor missing?
** Have fun ! **
[ ]: