How to use analytical continuation in GW? [Edit:Bug Found]]

[KalaShayminS]

Hello,
I am using abinit 7.2.1. From literature it said analytical continuation is implented in abinit, but I cannot find how to use it.
In test v67mbpt, only seperated screen and Sigma calculations are found. When I try to combine them, it prompts error
"Frequencies in the SCR file are not compatible with the analytic continuation. Verify the frequencies in the SCR file."
It seems that adjusting nfreqim does not help.

Is there any examples how to run this kind of calculation?
Thanks very much.

[gabriel.antonius]

What do you mean, "when you try to combine them" ?
The proper way of doing it is to compute the screening in one dataset, and Sigma in the next one. Like in "tests/v67mbpt/Input/t02.in". This is as "combined" as it gets!
Post your input if that doesn't help.

[KalaShayminS]

What do you mean, "when you try to combine them" ?
The proper way of doing it is to compute the screening in one dataset, and Sigma in the next one. Like in "tests/v67mbpt/Input/t02.in". This is as "combined" as it gets!
Post your input if that doesn't help.

Thank you very much for the reply.
From documentation it says "nomegasi" and "omegasimax " only works if optdriver=4 and gwcalctyp=x1. However in t02.in, "gwcalctyp" is not set, which should means plasmon-pole model by default. Is analytical continuation only used for spectral function in this example?
If I set gwcalctyp = 1 in t02.in then the error message in the previous thread appears.

[gabriel.antonius]

Just to make sure, we are talking about analytical continuation of Sigma to the complex plane, as in the example v67mbpt/t02.in, and not analytical continuation of the screening, for a contour deformation G0W0 calculation, as in example "v5/t71.in".

Indeed, in the example, a plasmon-pole model is used. With gwcalctyp=1, I manage to suppress the error message by setting nfreqim2 = 6
.
A higher value triggers the error message. I think this is a bug. Could someone comment on this?

[KalaShayminS]

Just to make sure, we are talking about analytical continuation of Sigma to the complex plane, as in the example v67mbpt/t02.in, and not analytical continuation of the screening, for a contour deformation G0W0 calculation, as in example "v5/t71.in".

Indeed, in the example, a plasmon-pole model is used. With gwcalctyp=1, I manage to suppress the error message by setting nfreqim2 = 6
.
A higher value triggers the error message. I think this is a bug. Could someone comment on this?

I think I have managed to figure out where the bug is. In 70_gw/m_screening.F90:1154, there is a frequency filter which ignore all frequency with Im(omega)<0.001*Ha_eV. However, omega itself is stored in unit Ha. When nfreqim2 <=6, the minimum Gauss-Lengendre knot is larger than 0.027 Hartree so there is no problem. But if nfreqim2>6, the minimum knot is filtered and is not counted in Sigma calculation, thus make the program stop.

[gabriel.antonius]

Yes, you're right!

As a matter of fact, this is also the cause of an other bug I had notice, but did not report (shame on me!) In contour deformation g0w0, using nfreqim > 29 causes an imaginary frequency to be dropped, and the result is erroneous.

Anyway, here is a quick fix. In 70_gw/m_screening.F90:1143

replace the old

Code: Select all

 if (Er%nomega==2) then
   Er%nomega_r=1
   Er%nomega_i=1
 else ! Real frequencies are packed in the first locations.
   Er%nomega_r=1
   do iomega=1,Er%nomega
     if ((REAL(Er%omega(iomega))>0.001*Ha_eV).AND.&
&     (AIMAG(Er%omega(iomega))<0.001*Ha_eV)) Er%nomega_r=iomega
   end do
   Er%nomega_i=0
   do iomega=Er%nomega_r+1,Er%nomega
     if ((REAL(Er%omega(iomega))<0.001*Ha_eV).AND.&
&     (AIMAG(Er%omega(iomega))>0.001*Ha_eV)) Er%nomega_i=Er%nomega_i+1
   end do
   Er%nomega_c=Er%nomega-Er%nomega_r-Er%nomega_i
 end if

with the new

Code: Select all

 if (Er%nomega==2) then
   Er%nomega_r=1
   Er%nomega_i=1
 else ! Real frequencies are packed in the first locations.
   Er%nomega_r=1
   do iomega=1,Er%nomega
     if ((REAL(Er%omega(iomega))>0.00001*Ha_eV).AND.&
&       (AIMAG(Er%omega(iomega))<0.00001*Ha_eV)) Er%nomega_r=iomega
   end do
   Er%nomega_i=0
   Er%nomega_c=0
   do iomega=Er%nomega_r+1,Er%nomega
     if ((REAL(Er%omega(iomega))<0.00001*Ha_eV).AND.&
&       (AIMAG(Er%omega(iomega))>0.00001*Ha_eV)) then
       Er%nomega_i=Er%nomega_i+1
     else if ((REAL(Er%omega(iomega))<0.00001*Ha_eV).AND.&
&            (AIMAG(Er%omega(iomega))<0.00001*Ha_eV)) then
       Er%nomega_c=Er%nomega_c+1
     end if
   end do
   if (Er%nomega_c > 0) then
     write(msg,'(3a,i6)')&
&      '  Some complex frequencies are too small to qualify as real or imaginary.',ch10,&
&      '  Number of unidentified frequencies = ', Er%nomega_c
     MSG_WARNING(msg)
   end if
   Er%nomega_c=Er%nomega-Er%nomega_r-Er%nomega_i
 end if

The tolerance criterion is a bit more stringent, and at least it gives a warning when frequencies are dropped.

I think it solves the problem, does it?
Thanks