Evaluating the oscillators using simmetries

[marco.govoni]

Dear Yambo developers,
we are studying the Yambo code. In particular we are interested in how Yambo calculates the oscillators using symmetries.
According to eq. (39) of ref. ( A. Marini et a., yambo: An ab initio tool for excited state calculations, Computer Physics Communications (2009) ) we would like to understand the following code-line (written in the file src/wf_and_fft/scatter_Bamp.F).

Code: Select all

 isc%rhotw(1:isc%ngrho)=  rho_tw_rs(fft_g_table(g_rot(qs,1:isc%ngrho),iGo),1)

As far as we know, g_rot(qs, 1:isc%ngrho) simply applies the simmetry qs to the G vectors, rotating them. We would like to know what fft_g_table (contained in src/wf_and_fft/fft_setup.F) does.
Thank you,
Sincerely
Marco Govoni and Ivan Marri

[myrta gruning]

Dear Marco Govoni and Ivan Marri,

Well I am not really into these technical details (Andrea is but he cannot answer at this moment), but I will try to answer, I hope it helps.

The subroutine calculates
rhotw(G)=<ib ik_is i_sp| e^[I(Rs(G-G0)).r] | ob ok_os o_sp>

to do this, first it multiplies the two wfs in direct space, then fft transforms this product (resulting in rho_tw_rs) and finally (in the line you pointed out) applies the rotation and subtracts G0 using two mappings:
1. Given the indexes of the symmetry operation and a G vector g_rot gives the index corresponding to the rotated Rs G vector.
2. Then, my understanding is that fft_g_table maps the index corresponding to the (rotated) G vector and the index corresponding to the G0 vector into the index corresponding to the G-G0 vector.

Cheers,

m

[andrea marini]

Let me only add that when performing and FFT from real to reciprocal space, the G components are stored in the array used by the FFT using an order different from the one used in the code. The table

Code: Select all

fft_g_table(iG,iGo)

then, gives the index in the FFT ordering of the component relative to the RL vector (G_{iG}-G_{iGo}). In this way the oscillator that is returned by scatter_Bamp is ordered in the correct way.

[marco.govoni]

Dear developers,

we are studying how the oscillators strengths are implemented in the yambo code.
As pointed out in your Com.Phys.Commun. (A2 appendix), the matrix elements are defined as:

rho_{nm} (k,q,G) = Int_dr psi*_{nk} ( r) e^{i (q+G).r } psi_{mk'} ( r)

where k' = k - q - G_0, k' and k are in the 1BZ, q is in 1QBZ. Int_dr denote an integration in the real space over the unit cell. psi are the Bloch wavefunctions:

psi_{nk} ( r) = e^{i k.r} u_{nk} ( r)

so that

rho_{nm} (k,q,G) = Int_dr u*_{nk} ( r) e^{i (G-G_0).r } u_{mk'} ( r)

by applying symmetries on k and k' : k= R k_{IBZ} and k' = R' k'_{IBZ}
the periodic function u transforms as : u_{nk} ( r) = u_{n R k_{IBZ}} ( r) = u_{n k_{IBZ}} ( R^{-1} r)

as a consequence we have

(*) rho_{nm} (k,q,G) = Int_dr u*_{n k_{IBZ} } ( R^{-1} r) e^{i (G-G_0).r } u_{m k'_{IBZ} } ( R'^{-1} r)

----------------------

The procedure is implemented in Scatter_Bamp.F as follows:

Step #1 :
Load the wavefunctions in real space: u_{n k_{IBZ} } ( r) and u_{m k'_{IBZ} } ( r)

Step #2 :
Use symmetries in real space to evaluate the product :
u*_{n k_{IBZ} } ( R^{-1} r) u_{m k'_{IBZ} } ( R'^{-1} r)
where complex conjugate operations can be removed or added depending on the presence of Time-Reversal (litr and lotr variables)

Step #3 :
Perform a Fourier transformation of this product

Step #4 :
Rotate the G-space using symmetry qs and shift the whole rotated space of a G_0 vector.

Code: Select all

       isc%rhotw(1:isc%ngrho)=  rho_tw_rs(fft_g_table(g_rot(qs,1:isc%ngrho),iGo),1)

Comparing this procedure with eq. (*), we don't understand the meaning of rotation qs, which as far as we know is defined as:

q_{QBZ} = R_{qs} q_{QIBZ}

Are we omitting something?
Thank you in advance.

Sincerely,

Marco Govoni and Ivan Marri

[andrea marini]

we are studying how the oscillators strengths are implemented in the yambo code.
As pointed out in your Com.Phys.Commun. (A2 appendix), the matrix elements are defined as:

rho_{nm} (k,q,G) = Int_dr psi*_{nk} ( r) e^{i (q+G).r } psi_{mk'} ( r)

In the appendix it is explained how to move the symmetry acting on k on the right, in order to reduce the number of operations to do in the evaluation of the BS kernel.

This procedure, however, is implemented in the routines the call scatter_Bamp.F and not in the routine itself, where it is calculated the most general oscillator of the type

Code: Select all

 !
 ! rhotw(G)=<ib ik_is i_sp| e^[I(RsG-G0).r] | ob ok_os o_sp>   
 !
 !  is(3)=(ib,ik,is,i_sp) --<--:--<-- os(3)=(ob,ok,os,o_sp)
 !                             :
 !                            /:\ iqs=(ig0,iq,qs)
 !                             :
 !
 !                     \:/
 ! REMEMBER that     ->-:->-  =  conjg [-<-:-<-]
 !                                        /:\
 !  iq_is = ik_is-ok_os + Go

as a consequence we have
(*) rho_{nm} (k,q,G) = Int_dr u*_{n k_{IBZ} } ( R^{-1} r) e^{i (G-G_0).r } u_{m k'_{IBZ} } ( R'^{-1} r)

Comparing this procedure with eq. (*), we don't understand the meaning of rotation qs, which as far as we know is defined as:

q_{QBZ} = R_{qs} q_{QIBZ}

Are we omitting something?

The operation of symmetry acting on G is due to the fact that all summations on G are rotated by the inverse of the symmetry R_{qs}. For example

Code: Select all

\sum_R \sum_q \sum_G <exp(i(Rq+G))> v(Rq+G) = \sum_R \sum_q \sum_G' <exp(iR(q+G'))> v(R(q+G'))