BCD convergence with NK #499
Replies: 5 comments 27 replies
|
BCD is really quite a delicate thing to claculate, so maybe you actually need a denser grid. However, my advice is :
|
|
Wait, GeTe is a semiconductor, and you are placing the Fermi level in the gap? Than you should get zero in both cases, just because you are integrating a full derivative. |
|
Hi Stepan, Thanks for the reply. I was doping the system to get the BCD at the Fermi level. 
I was looking the range of 0.5 eV above and below the fermi level, sampling 300 points to do the BCD integration. From the below plots, the magnitude between the sea vs surface approach itself was fluctuating a lot. 
I will try implementing the SAWF and denser grids to see how it goes and shall let you know. Thanks |
|
Is there a systematic way to choose the energy windows for SAWF, or is it purely trial and error? I used the same windows that worked perfectly for my Maximum Localization (MLWF) run, but when I switched to Symmetry-Adapted Wannier Functions (SAWF), the interpolated bands failed to match the Quantum ESPRESSO results. 
This was the code I was using for the SAWF process: import os
import ray
import numpy as np
from matplotlib import pyplot as plt
import scipy
import wannierberri as wb
import irrep
import spglib
NUM_CPUS = 16
ray.init(num_cpus=NUM_CPUS)
print(f"Python version:\t\t{__import__('sys').version.split()[0]}")
print(f"WannierBerri version:\t{wb.__version__}")
print(f"irrep version:\t\t{irrep.__version__}")
print(f"spglib version:\t\t{spglib.__version__}")
QE_PREFIX = "gete" # prefix in scf/nscf input → looks for tmp/gete.save
SEEDNAME = "my" # wannier90 seedname (my.win, my.mmn, my.eig, …)
PATH_DATA = "./" # folder containing tmp/, my.mmn, my.eig, …
assert os.path.exists(PATH_DATA), f"PATH_DATA={PATH_DATA} not found"
assert os.path.exists(os.path.join(PATH_DATA, "tmp")), \
"QE outdir 'tmp/' not found; adjust PATH_DATA or QE_PREFIX"
from irrep.bandstructure import BandStructure
print("\n── 2. Reading bandstructure & detecting symmetries ──")
bandstructure = BandStructure(
code="espresso",
prefix=os.path.join(PATH_DATA, "tmp", QE_PREFIX),
magmom=None, # no magnetic moments → pure TR symmetry
include_TR=True, # include time-reversal (and its combinations with spatial ops)
)
spacegroup = bandstructure.spacegroup
print(f"Detected spacegroup number: {spacegroup.number}")
# spacegroup.show() # uncomment to list all detected symmetry operations
print("\n── 3. Creating w90data (eig + symmetrizer from irrep) ──")
w90data = wb.w90files.Wannier90data.from_bandstructure(
bandstructure=bandstructure,
files=["eig", "symmetrizer"], # MMN will be read separately
)
print(" Reading MMN from pw2wannier90 output …")
mmn_path = os.path.join(PATH_DATA, SEEDNAME) # → ./my (w/o extension)
mmn = wb.w90files.mmn.MMN.from_w90_file(mmn_path, bkvec=w90data.bkvec)
w90data.set_file("mmn", mmn, overwrite=True)
npz_path = os.path.join(PATH_DATA, f"{SEEDNAME}_w90data")
w90data.to_npz(npz_path)
print(f" Saved to {npz_path}.*.npz (load later with Wannier90data.from_npz)")
# Fractional coordinates taken directly from nscf input (ATOMIC_POSITIONS crystal)
from wannierberri.symmetry.projections import Projection, ProjectionsSet
pos_Te = [[0.5227847276, 0.5227847276, 0.5227847276]]
pos_Ge = [[0.9982152724, 0.9982152724, 0.9982152724]]
proj_Te_s = Projection(position_num=pos_Te, orbital="s", spacegroup=spacegroup)
proj_Te_p = Projection(position_num=pos_Te, orbital="p", spacegroup=spacegroup)
proj_Ge_s = Projection(position_num=pos_Ge, orbital="s", spacegroup=spacegroup)
proj_Ge_p = Projection(position_num=pos_Ge, orbital="p", spacegroup=spacegroup)
projections_set = ProjectionsSet(projections=[proj_Te_s, proj_Te_p, proj_Ge_s, proj_Ge_p])
w90data.set_projections(projections_set, bandstructure=bandstructure)
froz_max = 13.65
_ = wb.wannierise.wannierise(
w90data=w90data,
froz_min=-4.79,
froz_max=froz_max,
outer_min=-4.79,
outer_max=16.75,
print_progress_every=500,
num_iter=5000,
conv_tol=1e-8, # tighter than default; loosen to 1e-6 if too slow
localise=True, # True = disentangle + localise; False = only disentangle
sitesym=True, # enforces symmetry (SAWF)
)
system = wb.system.System_w90(w90data= w90data, berry=True, transl_inv_JM=True,
symmetrize=True)
# optionally - save it for later use
system.save_npz("my_system")
# K-path from your my.win (kpoint_path block), in reduced coordinates
path = wb.Path.from_nodes(
system,
nodes=[
[0.0000000000, 0.0000000000, 0.0000000000], # G
[0.0000000000, 0.5000000000, 0.0000000000], # L
[0.3714285714, 0.6285714286, 0.0000000000], # U
[0.5000000000, 0.5000000000, 0.0000000000], # X
[0.0000000000, 0.0000000000, 0.0000000000], # G
[0.5000000000, 0.5000000000, 0.5000000000], # Z
[0.6285714286, 0.2432432432, 0.6285714286], # U
],
labels=["G", "L", "U", "X", "G", "Z", "U"],
length=2000,
)
bands_path = wb.evaluate_k_path(system, path=path)
fig, ax = plt.subplots(figsize=(8, 6))
bands_path.plot_path_fat(
path,
quantity=None,
Eshift=0,
Emin=-5.0,
Emax=15.0,
iband=None,
mode="fatband",
fatfactor=20,
linecolor="red",
close_fig=False,
show_fig=False,
label="WannierBerri (SAWF)",
)Is there any error in the code for wannierisation? Will playing with the energy windows be enough to get the good overlap with the DFT bands? I got these warnings in the process as well:
wberri_3500095.txt - file with the wannierisation data |
2.1 what if you do not set any windows at all?
I assume you are using the latest pip version |
|
I am using the 1.8.0 version. Doing wannierberri with symmetrize=False, also gave bad wannierization. At the end of 5000 iterations, the spread was: The overlap plot with QE bands were still bad. I tried with few different values for the energy windows(shrinking frozen, extending outer window, not setting windows) and they did not make much difference in the band overlap plot(still wavy). Is the solution to play with the energy window values to get good SAWF plots? |
|
I am not sure what the problem is here. Please, upload the whole folder somewhere (including the gete.save folder and the w90 files .mmn, .eig, .amn .win) To reduce the data size, please reduce the nscf grid to 4x4x4, it should be enough to get reasonable results. I will investigate this, when I have time. |
|
I am attaching the files in my GDrive link. The calculations were done with a 4x4x4 grid. Thanks in advance for the help. |
|
My first suggestion - get more bands in the nscf run, and increase the upper edge of the outer window. Because, also with Wannier90, on the 0th iteration (right after the disentanglement), you get huge spreads for the Ge WFs This indicates that something is wrong, probably some of the Ge bands are actually higher than band 40. Also, you may try to first get a good wannierisation for scalar bands (no spin-orbit), and than proceed to the spinor case. Please, try the scalar case with a 4x4x4 grid and ~50 bands in the nscf run. If it still fails, upload the data. |
Can you mention how to use the updated branch? I was not getting the updates in the code with the above command. Thanks |
|
Such general questions may be figured out with AI chatbots, they suggest: But as far ass the PR is already merged, you can install from master If it has no effect, maybe first uninstall the previous installation, because they are under the same version number, they may look as "the same" I usually work with editable installs git clone https://github.com/USER/REPO.git
cd REPO
git checkout branchname
pip install -e . |
|
Thanks for the fix. I was to able to get good SAWF bands. Now with this, I was comparing BCD with the Fermi sea vs Fermi surface approach using the following code: Using a 300x300x300 k-grid, Fermi sea vs Surface approach for BCD converges closely (which might get even better with higher k-grid?) Also comparing the different NK for fermi sea, looks well converged (ignoring the spikes) The fermi surface approach seems to be better than the fermi sea approach in this case as the values are less fluctuating with the different grids. Thanks for the suggestions to fix this. |
|
great. You may also try "tetra=True" for both calculators, it should improve convergence, especially for Fermi surface calc. Next, try increasing adpt_num_iter (upto smth like 100), it should help also |
|
I applied tetrahedral smearing and increased the adapt_num_iter to 50, and the plots looks much cleaner. |
|
Hi, How do we ensure SOC is enabled in the SAWF process? I was using s and p orbitals of Ge and Te, so with SOC, num_wann = 16. In the tutorial, I couldn't find the part with SOC for wannierisation. The current code was giving |
|
Ok, but did the wannierisation give 16 WFs ? It is on, do ntot worry, if it was on the the dft part. Otherwise the script would crash. It was just a wrong behaviour of the num_wann property, which was recently fixed, see https://github.com/wannier-berri/wannier-berri/pull/511/changes 
|
|
The wannierisation was giving 16 as expected. I though we had to turn on with spinor=True command somewhere in the script |
|
I had a question with respect to the convergence of BCD with the Fermi sea vs surface methods. I did the SAWF with conv_tol = 1e-8 converged in 3000 iterations. And in the calculation of BCD, I used denser grids of 250 and 300, 
The Fermi surface seems to be giving a converged BCD value but the Fermi sea approach gives a different value. Could this be a issue with the choice of the energy window for the wannierisation? Is the solution in this case to increase the grid size/ increase the number of adaptive iterations from 50 to higher value like 100 to get more converged BCD value? |
|
What is worrying is that your spreads are quite large. From your data (without SOC) I was getting for pz +sp2 projections With SOC it should not be much different. |
|
First, make sure you use the latest version of WannierBerri. add wannierberri.welcome()to see versions of all essential packages |
|
For the files that I had posted earlier in the Drive folder, I got clean convergence and same value for the BCD components (with SOC) as I had posted earlier (#499 (reply in thread)). I was varying the doping amounts to see how the coefficient would vary specifically and for some of the doping concentrations, I was having some difficulty with the BCD convergence from Sea vs surface approaches. I think I might get good results if I keep playing with the froz_max values. Sorry for the confusion. |
|
Just to understand better what you are doing, how is doping introduced into the DFT calculation? |
|
Doping was done with |
Uh oh!
There was an error while loading. Please reload this page.
Hi everyone,
I am currently calculating the Berry Curvature Dipole (BCD) for GeTe. By symmetry, we expect$D_{xy} = -D_{yx}$ . However, I am running into severe convergence issues, both when scaling up my NK grid and when comparing the FermiSea vs. FermiSurface formulations. I saw a very similar discussion in #239, and I was wondering if there have been any recent recommended workarounds for this feature?
To ensure a strict mathematical comparison, I evaluated both methods at strictly$T=0$ K
Calculators separated: FermiSea and FermiSurface evaluated in distinct run() calls to prevent confusing the adaptive mesh.
Method: tetra=True applied to both calculators.
Symmetry: system.set_symmetry(['TimeReversal', 'C3z', 'My']) and use_irred_kpt=True.
Adaptive Refinement: Heavy refinement focused on a narrow energy window (adpt_num_iter=75, adpt_fac=10, adpt_mesh=3).
Base Grids: NK = 100, 150, 200.
Despite these heavy optimizations, the BCD values are fluctuating wildly, particularly for the volume integration (FermiSea).
At E=E_F,
NK=100:
NK=150:
NK= 200:
As you can see, the FermiSurface and the FermiSea values are fluctuating.
Is the FermiSea volume integration mathematically unstable for BCD in materials like GeTe? Should I increase the adaptive iterations and the NK grid even further?
Are there any other suggestions or undocumented tweaks to improve the convergence of the tetrahedron method here?
Thanks
Dominic
All reactions