The ABINIT Users Mailing List ( CLOSED )

[Geoffrey Stenuit <geoffrey.stenuit@tyndall.ie>]

Hi Matthieu,

Thank you for your advices. I have tried to hack abinit (version - 
4.6.5) and in particular the following file ~/Src_8seqpar/cgwf.F90 (at 
line 435). See below to see the very few  modifications.
To my understand, this routine is called to compute the eigenvalue 
problem with the conjugate gradient technique. Since i just want to 
extract the the expectation value <psi|H|psi> with H being fixed by the 
provided density file (BSonlyi_DEN) and \psi> being taken from the 
BSonlyi_WFK file, before perfoming any diagonalisation, i decided to 
sent to the output the values of  lam0=<psi|H|psi> (for each band) 
before it will be modified.
I am not sure if my approach is correct ! Could you please comment it ?

Please find in attachment my abinit input file (which call for the 
reading of *_DEN and* _WFK files, consider non-scf loop () and number of 
step=1, since if nstep=0, then the density is fixed by the _WFK file).

Thank you very much for you comments and corrections,

Regards,

Joe

!Compute <g|H|c>
   sij_opt=0;if(prtvol<0) call status(0,filstat,iexit,level,'call 
getghc   ')
   call 
getghc(cwavef,dimffnl,ffnl,filstat,ghc,gsc_dummy,gs_hamk,gvnlc,kg_k,&
&    kinpw,lmnmax,matblk,mgfft,mpi_enreg,mpsang,mpssoang,natom,1,&
&    
npw,nspinor,ntypat,nvloc,n4,n5,n6,ph3d,prtvol,sij_opt,tim_getghc,vlocal)
   if(prtvol<0) call status(0,filstat,iexit,level,'after getghc  ')

!!!!!!!!!!!! MODIFICATIONS !!!!!!!!!!!!
if (prtvol == 10) then
   call dotprod_g(chc,doti,istwf_k,mpi_enreg,npw*nspinor,1,cwavef,ghc)
   write(message, '(a,i4,a,e14.6,a,e14.6,a)' ) 'Band nbr: ',iband,' 
Compute <C|H|C>=',chc,' Ha or ', chc*27.2107,' eV'
   call wrtout(ab_out,message,'COLL') !!! write the output in the 
***.out file only if prtvol == 10
   call wrtout(06,message,'COLL') !!! write the output in the log_file 
file only if prtvol == 10
end if
!!!!!!!!!!!!!!!!



Matthieu Verstraete wrote:

> Hello Goeffrey,
>
> you should try negative iscf (-2 -3) which actually call on the 
> getden. I am not sure that it will use your wfk, but it should. The 
> trouble then is that it will probably diagonalize the matrix it finds, 
> instead of feeding out the full H_mn. I tried to do this before with 
> different operators for H, but you may have to hack it so that it 
> outputs the matrix before calls to ZHPEV and so on (not hard, 
> especially for just 1 kpt). This could even be a new flag, which sets 
> nstep and iscf, and reads in DEN and WFK...
>
> Tell us how it goes,
>
> Matthieu

# Convergence with respect to the number of k points.

chkprim 0
prtvol 10 

#Definition of the unit cell
#Definition of the unit cell
acell 3*2.1203311987E+01          # This is equivalent to
#acell 3*21.3446668          # This is equivalent to GaAs lattice Cst
rprim  1.0  0.0  0.0   # FCC primitive vectors (to be scaled by acell)
       0.0  1.0  0.0
       0.0  0.0  1.0
xred
             -8.2239427067E-21  8.5221721693E-21  5.1564586525E-21
             -8.2239427067E-21  8.5221721693E-21  4.9999999999E-01
             -8.2239427067E-21  4.9999999999E-01  5.1564586525E-21
             -8.2239427067E-21  4.9999999999E-01  4.9999999999E-01
              4.9999999999E-01  8.5221721693E-21  5.1564586525E-21
              4.9999999999E-01  8.5221721693E-21  4.9999999999E-01
              4.9999999999E-01  4.9999999999E-01  5.1564586525E-21
              4.9999999999E-01  4.9999999999E-01  4.9999999999E-01
              4.5663886334E-04  2.4357916843E-01  2.4357916843E-01
             -4.5663886334E-04  2.4357916843E-01  7.5642083156E-01
             -4.5663886334E-04  7.5642083156E-01  2.4357916843E-01
              4.5663886334E-04  7.5642083156E-01  7.5642083156E-01
              4.9998342271E-01  2.4954775603E-01  2.4954775603E-01
              5.0001657727E-01  2.4954775603E-01  7.5045224395E-01
              5.0001657727E-01  7.5045224395E-01  2.4954775603E-01
              4.9998342271E-01  7.5045224395E-01  7.5045224395E-01
              2.4357916843E-01  4.5663886334E-04  2.4357916843E-01
              2.4357916843E-01 -4.5663886334E-04  7.5642083156E-01
              2.4954775603E-01  4.9998342271E-01  2.4954775603E-01
              2.4954775603E-01  5.0001657727E-01  7.5045224395E-01
              7.5642083156E-01 -4.5663886334E-04  2.4357916843E-01
              7.5642083156E-01  4.5663886334E-04  7.5642083156E-01
              7.5045224395E-01  5.0001657727E-01  2.4954775603E-01
              7.5045224395E-01  4.9998342271E-01  7.5045224395E-01
              2.4357916843E-01  2.4357916843E-01  4.5663886334E-04
              2.4954775603E-01  2.4954775603E-01  4.9998342271E-01
              2.4357916843E-01  7.5642083156E-01 -4.5663886334E-04
              2.4954775603E-01  7.5045224395E-01  5.0001657727E-01
              7.5642083156E-01  2.4357916843E-01 -4.5663886334E-04
              7.5045224395E-01  2.4954775603E-01  5.0001657727E-01
              7.5642083156E-01  7.5642083156E-01  4.5663886334E-04
              7.5045224395E-01  7.5045224395E-01  4.9998342271E-01
              3.7207521125E-01  1.2470729392E-01  3.7207521125E-01
              3.7401361666E-01  1.2495006418E-01  8.7504993580E-01
              3.7505083674E-01  6.2494916324E-01  3.7505083674E-01
              3.7207521125E-01  6.2792478874E-01  8.7529270606E-01
              8.7504993580E-01  1.2495006418E-01  3.7401361666E-01
              8.9458595135E-01  1.0541404863E-01  8.9458595135E-01
              8.7529270606E-01  6.2792478874E-01  3.7207521125E-01
              8.7504993580E-01  6.2598638333E-01  8.7504993580E-01
              1.0541404863E-01  1.0541404863E-01  1.0541404863E-01
              1.2495006418E-01  1.2495006418E-01  6.2598638333E-01
              1.2495006418E-01  6.2598638333E-01  1.2495006418E-01
              1.2470729392E-01  6.2792478874E-01  6.2792478874E-01
              6.2598638333E-01  1.2495006418E-01  1.2495006418E-01
              6.2792478874E-01  1.2470729392E-01  6.2792478874E-01
              6.2792478874E-01  6.2792478874E-01  1.2470729392E-01
              6.2494916324E-01  6.2494916324E-01  6.2494916324E-01
              1.2470729392E-01  3.7207521125E-01  3.7207521125E-01
              1.2495006418E-01  3.7401361666E-01  8.7504993580E-01
              1.2495006418E-01  8.7504993580E-01  3.7401361666E-01
              1.0541404863E-01  8.9458595135E-01  8.9458595135E-01
              6.2494916324E-01  3.7505083674E-01  3.7505083674E-01
              6.2792478874E-01  3.7207521125E-01  8.7529270606E-01
              6.2792478874E-01  8.7529270606E-01  3.7207521125E-01
              6.2598638333E-01  8.7504993580E-01  8.7504993580E-01
              3.7207521125E-01  3.7207521125E-01  1.2470729392E-01
              3.7505083674E-01  3.7505083674E-01  6.2494916324E-01
              3.7401361666E-01  8.7504993580E-01  1.2495006418E-01
              3.7207521125E-01  8.7529270606E-01  6.2792478874E-01
              8.7504993580E-01  3.7401361666E-01  1.2495006418E-01
              8.7529270606E-01  3.7207521125E-01  6.2792478874E-01
              8.9458595135E-01  8.9458595135E-01  1.0541404863E-01
              8.7504993580E-01  8.7504993580E-01  6.2598638333E-01

# Second: The BS...
getcell -1   # start from the previous acell values
getxred -1   # start from the previous xred values
#prtwf 0 # No wavefunction output is provided
iscf    -2   # no self-consistent calculations
getden  -1   # because start with the previous computed density
#getwfk  1
irdwfk 1 # start from the previous WF (wfi_WFK) which is here the GaAs at k=0 
W.F.
#kptopt  -2   # divide in two zones the BS
#nband   360  # 16*8=128 + 32 en plus pour avoir la conduction

#ndivk 1 1  #5 5  # 5 5 5 10, 12 and 17 divisions of the 3 segments, delimited

#kptbounds  0.02  0.02  0.02 # L point (0.5 0.5 0.5) divided by 25=> 0.02 
0.02 0.02 and 0.01 0.01 0.01
#           0.0  0.0  0.0 # Gamma point
#           0.04  0.0  0.0 # X point (1.0 0 0) => divided by 25 and take 3 
kpts 
# => 0.04 0.0 0.0 and 0.02 0.0 0.0 and 0.0 0.0 0.0
#         #1.0  1.0  1.0 # Gamma point in another cell

kptopt 0 # read directly nkpt, kpt, kptnrm and wtk  (corresponds to the usage 
before version 2.1) 
nband  256  # 8*32=256 ( 2atoms: 6 = 4 VB + 2 CB) en plus pour avoir la 
conduction
nkpt 1 #5 #nbre of kpts
kpt 
0.0 0.0 0.0
#0.02 0.02 0.02 ## list of kpts

tolwfr  1.0d2
enunit  1             # Will output the eigenenergies in eV 

#Definition of the atom types in the cristal
ntypat 3          # There are three types of atom
znucl 7 33 31     # The keyword "znucl" refers to the atomic number of the
                  # possible type(s) of atom. The pseudopotential(s)
                  # mentioned in the "files" file must correspond
                  # to the type(s) of atom. Here, the only type is Silicon.

#Definition of the atoms
natom 64           # There are 64 atoms
typat 1 2 2 2 2 2 2 2
      2 2 2 2 2 2 2 2
      2 2 2 2 2 2 2 2
      2 2 2 2 2 2 2 2
      3 3 3 3 3 3 3 3
      3 3 3 3 3 3 3 3
      3 3 3 3 3 3 3 3
      3 3 3 3 3 3 3 3

#Definition of the planewave basis set
#ecut 20.0         # Maximal kinetic energy cut-off, in Hartree
ecut 25.0
#ecut 27.5625
#ecut3 30.0
#ecut4 35.0
#ecut 32.0

#Definition of the SCF procedure
nstep 1          # Maximal number of SCF cycles
#toldff 5.0d-6 # A value ten times smaller than tolmxf is suggested (for 
example 5.0d-6 hartree/bohr)
#toldfe 1.0d-10    # Will stop when, twice in a row, the difference 
diemac 12.0       # Although this is not mandatory, it is worth to