TDepES Legacy code question(s)

mjocic · August 30, 2021, 1:32pm

Hello,

I am currently trying to perform TDepES calculation on a polar material with SOC.
I’m following the example provided in tutorial EPH-Legacy/tdepes ("/tutorial/tdepes/") that uses the legacy code.

Following TDepES (/topics/TDepES/) the restoration of charge neutrality,
as seen in Ponce2015, can be invoked in code by performing an ddk calculation, followed
by electric field perturbation for q==[0,0,0]. These DDBs are then used for el-ph calculations
for q!=[0,0,0] using previously obtained wfk and wfk+q, as shown in examples:
example59: /tests/v7/Input/t59.abi
and
example58: /tests/v7/Input/t58.abi
using getddb or irdddb.

So, is there any significant difference between getddb and irdddb and which one would you recommend using? I parallelize over q-points: one job submition per q point then I merge all DDBs with mrgddb tool, so would I benefit to use getddb_filepath instead?
Why there isn’t any getddb/irdddb in tutorial example EPH-Legacy/tdepes (/tutorial/tdepes/)?
Should the input for el-ph dataset be different for q==[0,0,0] from q!=[0,0,0]?
Do I just set nspinor=2 and double the number of el. bands for SOC?
In TDepES (/topics/TDepES/) there is mention of temperature_para.py script, is this the same as temperature_final.py found in the tutorial EPH-Legacy/tdepes (/tutorial/tdepes/), only updated?
What would be the best way to get temperature depenence of electronic band(s) for k!=[0,0,0]?
Is SOC supported in temperature_final.py?

Kind regards,
Milan Jocic

PS Because I can’t post more than 2 links, in parentheses there should be a https://docs.abinit.org preceding the part of every link given…

sponce · September 2, 2021, 11:58am

Hello Milan Jocic,

I’m afraid I cannot easily answer to all your questions since it was a long time ago and the code might have changed. Also note that it is now considered as a legacy way of doing ELPH in Abinit and the recommended way relies on interpolation of the perturbed potential.

However what I can say is that the data used for the Ponce2015 paper can be accessed here https://archive.materialscloud.org/record/2021.137

The SOC has not been explicitly tested and is therefore not officially supported in
temperature_final.py but it might well work. You can try it out but use with care.

For 5., yes temperature_final.py is an update of temperature_para.py
For 6, if you follow the EPH-Legacy tutorial, it should show how to do this Tdepes - abinit

Hope this helps,
Samuel

mjocic · September 2, 2021, 2:09pm

Hello Ponce,

Thank you very much for the response. I understand that this topic isn’t recent but I would still appreciate if you could consider the following:

Considering temperature_final.py, I don’t think the SOC would not work by default:
This segment in rf_final.py code:

    occtmp = EIGR2D.nsppol*EIGR2D.occ[:,:,:]/2 # jband # should be 1 !
    delta_E_ddw = N.einsum('lij,k->lijk',eig0[:,:,:].real,N.ones(EIGR2D.nband)) - \
              N.einsum('lij,k->likj',eig0[:,:,:].real,N.ones(EIGR2D.nband)) - \
              N.einsum('i,ljk->ljik',N.ones((EIGR2D.nband)),(2*occtmp-1))*smearing*1j # spin,ikpt,iband,jband

    ddw_tmp = N.einsum('ijkln,lm->mijkn',ddw_addQ,2*bose+1.0) # itemp,ikpt,iband,jband,ispin
    ddw_add = N.einsum('ijklm,mjkl->imjk',ddw_tmp,1.0/delta_E_ddw) # temp,spin,ikpt,iband
    delta_E = N.einsum('lij,k->lijk',eig0[:,:,:].real,N.ones(EIGR2D.nband)) - \
              N.einsum('lij,k->likj',eigq.EIG[:,:,:].real,N.ones(EIGR2D.nband)) # spin,ikpt,iband,jband
    delta_E_sm = N.einsum('i,ljk->ljik',N.ones((EIGR2D.nband)),(2*occtmp-1))*smearing*1j # spin,ikpt,iband,jband
    num1 = N.einsum('ij,mkl->mkijl',bose,N.ones((EIGR2D.nsppol,EIGR2D.nkpt,EIGR2D.nband))) +1.0 \
          - N.einsum('ij,mkl->mkijl',N.ones((3*EIGR2D.natom,ntemp)),occtmp) # spin,k,mod,temp,band # bef was (imode,tmp,band)
    deno1 = N.einsum('mijk,l->mijkl',delta_E,N.ones(3*EIGR2D.natom),dtype=N.complex)

Has a problem with SOC, since occupations are 1 per band (instead of 2) so this array occtmp gets filled with 1/2 which results division with 0 in ddw_add i.e. (2*occtmp-1) in delta_E_ddw becomes zero.

I suspect that line is related to eq (15) and (16) in Ponce2015. There the 1j*smearing is the ad-hoc parameter of numerical nature i*delta inferred as finite lifetime. So, the occtmp variable is an array over states that takes 1 for occupied and 0 for unoccupied in nonSOC case (0.5 and 0 for SOC) so the delta_E_ddw would suggest that i*delta changes sign in equations (15) and (16) from + to - for unoccupied states. I am confused why is that the case?

Anyhow, changing (2*occtmp-1) to N.ones(occtmp.shape) in delta_E_ddw and other places where it occurs can remove this division with 0, when dealing with SOC case. Since I am only interested in temperature dependence, that takes only the real part of temperature dependent renormalization delta Epsilon for adiabatic/dynamic case of R.I.A. , I guess the sign of i*delta shouldn’t be an issue?

There is a part of rf_final.py the code that calculates average renormalization for degenerate states def make_average that for SOC should consider 4-fold degenerate states, that part I edited too but the results of my renormalization do not look promising.

Another question that I have is about III B section of Ponce2015, Restoration of charge neutrality:
Since polar materials already exibit non-zero Born effective charges, will the restoration of charge neutrality help the convergence over k-points for them like it does for non-polar materials? Is it even correct to impose charge neutrality in that case? If so, why? Does the code add the long-range contribution to GKK when getddb is activated for el-ph interaction ( phonon=1 and ieig2rf=5 ).
I use the Legacy code since I was under the impression that SOC is not supported in the new EPH code? Is that still the case or has the support been added, because I couldn’t find it in release notes? In the new EPH tutorial, eph-intro, long range dipolar fields are adressed as in Verdi2015, so does the code know when it is dealing with polar or non-polar materials?

Kind regards,
Milan Jocic

sponce · September 10, 2021, 6:44am

Hello Milan Jocic,

Great that you adapted the code to make it work with SOC.

The Born effective charges indeed have to fulfil a sum rule (both for IR active and non-IR materials).
However in the case of polar materials, it is not important to restaure charge neutrality and will make no difference. Note that if you are using the non-adiabatic formulation for non-polar materials, it will also probably not be very important. The issue was simply that there was residual charges and therefore the adiabatic formulation was diverging even though it should not for diamond.
You are correct that SOC is not supported at present in the new EPH code but I think it is under development and might be soon ?
Yes the new code is dealing properly with polar and non-polar materials when doing the Fourier interpolation of the perturbed potential and can include long range dipoles as well as quadrupoles.

Best,
Samuel

gmatteo · September 11, 2021, 9:59pm

You are correct that SOC is not supported at present in the new EPH code but I think it is under development and might be soon

The next version adds supports for SOC in the EPH code but only in the non-magnetic case (nspinor 2 nspden 1). The non-collinear magnetic case (nspinor 2 nspden 4) is not yet supported.
Also note that the computation of the dynamical quadrupoles does not support nspinor 2.

mjocic · September 13, 2021, 9:44am

The next version adds supports for SOC in the EPH code but only in the non-magnetic case (nspinor 2 nspden 1).

Thank you. When will the next version be released?

jmbgmx · September 13, 2021, 10:15am

Expected before the end of the month…

mjocic · November 17, 2021, 11:58am

Hello,

I have a follow up question: temperature_final.py and rf_final.py compute lifetimes, they are contained in total_corr[3,:,:,:,:] array. I’m interested is this the value of 1/tau or 1/(2*tau) from equation (18) in Ponce2015?

Thanks,
Milan Jocic