Dear Sitangshu,
Thank you for pointing this out. We have tested many times the units of the electron-phonon matrix elements across different systems and we are pretty sure they are correct. The full |g| matrix element has energy units, and since |g| = |dvscf|/sqrt(E_ph), that means that |dvscf| has units of energy^{3/2}. In Yambo, the standard energy units are hartree so the |dvscf| is given in hartree^{3/2}, and in order to obtain the correct |g| in hartree we need to temporarily convert back the phonon energies, since they by default they are given in eV like any other energy eigenvalue in yambopy.
Now, looking at Figure 6(a) of the paper you cite, it looks like they plot the value of |g| in eV. Indeed, I made a test and tried to plot that quantity by converting the gkkp array in eV (parameters: 12x12x1 k/q-grid plus spin-orbit so I have i_k=143, then for the TA mode I took i_nu=1 and for bottom conduction band i_c1=26):
Code: Select all
g_of_q = np.abs(yelph.gkkp[:,i_k,i_nu,i_c1,i_c1])*ha2ev
By doing this I get values pretty similar to the 2013 paper (see figure below).
As for the kpoint units: you are right, they are different from the quantum espresso "cart. coord. in units 2pi/alat". I also used to get confused all the time, so there is a function in YamboLatticeDB called get_units_info() as a reminder about reciprocal-space unit conventions between QE and Yambo. It always prints this output:
Yambo cartesian units [cc in yambo]:
:: self.car_kpoints*2.*pi
QE cartesian unists [cart. coord. in units 2pi/alat in QE:
:: self.car_kpoints*self.alat[0]
Internal yambo units [iku]:
:: self.iku_kpoints
Reduced coordinates [rlu in yambo, cryst. coord. in QE]:
:: self.red_kpoints
From which you can see that in Yambo, the Cartesian coordinates are different from QE because they are not in units of 2pi/alat. In yambopy, the Cartesian coordinates are in units of 2pi but not in units of 1/alat. So, in order to get the same as they use in the 2013 paper, we need to rescale like this: ylat.car_kpoints*ylat.alat[0].
Below an example on how one could do the rescaling while still using the function plot_elph from YamboElectronPhononDB (but of course it may be more convenient to do your own easily editable plotting script, this is taken from tutorial/databases_yambopy/elph_plot.py):
Code: Select all
if Qspace_Plot:
from yambopy.plot.plotting import BZ_hexagon # We need to rescale the BZ border to QE cc in 2pi/a
QE_cc_kpoints = ylat.car_kpoints*ylat.alat[0] # Rescale kpoints to QE cc in 2pi/a
yelph.rlat = yelph.rlat*ylat.alat[0] # Rescale rec. lattice vectors to QE cc in 2pi/a
yelph.plot_elph(g_of_q,kcoords=QE_cc_kpoints,s=100,plt_cbar=True,marker='H',cmap='hot') # Give explicitly rescaled kpts
yelph.ax.patches[-1].remove() # Remove default BZ borders in previous units
scaled_BZ = BZ_hexagon(yelph.rlat,center=(0.,0.),orientation=np.radians(30),color='white',linewidth=2) # New BZ borders
yelph.ax.add_patch(scaled_BZ) # Add new BZ borders
yelph.ax.set_title('|g(q)| (eV)')
yelph.ax.set_xlabel('q_x')
yelph.ax.set_ylabel('q_y')
plt.show()
plt.savefig('test_g_q_sitangshu_b%d.png'%i_c1,dpi=150)
Below is what I get in my test with updated units for |g| and qpoints. This is just a test with a coarse grid and I'm using SOC which probably wasn't used in the 2013 paper (plus different pseudopotential, lattice parameter optimization, etc). All in all, it looks pretty similar to their plot to me. But maybe I missed something, in that case please don't hesitate to let me know!
test_g_q_sitangshu_b2.png
Best,
Fulvio
You do not have the required permissions to view the files attached to this post.
Dr. Fulvio Paleari
S3-CNR Institute of Nanoscience and MaX Center
Modena, Italy