TABLE OF CONTENTS
ABINIT/scfeig [ Functions ]
NAME
scfeig
FUNCTION
Compute the largest eigenvalue and eigenvector of the SCF cycle. A brute force algorithm is presently used.
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (XG,GMR) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
istep= number of the step in the SCF cycle nfft=(effective) number of FFT grid points (for this processor) nspden=number of spin-density components
OUTPUT
(see side effects)
SIDE EFFECTS
vtrial0(nfft,nspden)= contains vtrial at istep == 1 vtrial(nfft,nspden)= at input, it is the trial potential that gave vresid . at output, it is an updated trial potential vrespc(nfft,nspden)=the input preconditioned residual potential work(nfft,nspden,2)=work space
PARENTS
m_ab7_mixing
CHILDREN
wrtout
SOURCE
41 #if defined HAVE_CONFIG_H 42 #include "config.h" 43 #endif 44 45 #include "abi_common.h" 46 47 48 subroutine scfeig(istep,nfft,nspden,vrespc,vtrial,vtrial0,work,errid,errmess) 49 50 use defs_basis 51 use m_profiling_abi 52 53 !This section has been created automatically by the script Abilint (TD). 54 !Do not modify the following lines by hand. 55 #undef ABI_FUNC 56 #define ABI_FUNC 'scfeig' 57 use interfaces_14_hidewrite 58 !End of the abilint section 59 60 implicit none 61 62 !Arguments ------------------------------------ 63 !scalars 64 integer,intent(in) :: istep,nfft,nspden 65 integer,intent(out) :: errid 66 character(len = 500), intent(out) :: errmess 67 !arrays 68 real(dp),intent(inout) :: vtrial0(nfft,nspden),work(nfft,nspden,2) 69 real(dp),intent(inout) :: vrespc(nfft,nspden) 70 real(dp), intent(inout) :: vtrial(nfft,nspden) 71 72 !Local variables------------------------------- 73 !scalars 74 integer :: ifft,isp 75 real(dp) :: eigen_scf,factor,fix_resid,resid_new,resid_old 76 character(len=500) :: message 77 78 ! ************************************************************************* 79 80 errid = AB7_NO_ERROR 81 82 if(nspden==4)then 83 errid = AB7_ERROR_MIXING_ARG 84 write(errmess, *) ' scfeig : does not work yet for nspden=4' 85 return 86 end if 87 88 !Set a fixed residual square for normalization of eigenvectors 89 fix_resid=1.0d-4 90 91 !A few initialisations for the first istep 92 if(istep==1)then 93 94 write(message, '(a,es12.4,a,a,a,a,a,a,a)' )& 95 & ' scfeig : fixed PC_residual square =',fix_resid,ch10,& 96 & ' Note that fixed resid should always be much larger',ch10,& 97 & ' than initial PC resid square, still sufficiently',ch10,& 98 & ' small to reduce anharmonic effects ',ch10 99 call wrtout(std_out,message,'COLL') 100 101 ! Compute the preconditioned residual 102 resid_old=0.0_dp 103 do isp=1,nspden 104 do ifft=1,nfft 105 resid_old=resid_old+vrespc(ifft,isp)**2 106 end do 107 end do 108 write(message, '(a,es12.4)' )& 109 & ' scfeig : initial PC_residual square =',resid_old 110 call wrtout(std_out,message,'COLL') 111 if(resid_old>1.0d-8)then 112 errid = AB7_ERROR_MIXING_ARG 113 write(errmess,'(a,a,a,a,a,a,a,a,a,a)') ch10,& 114 & ' scfeig : ERROR -',ch10,& 115 & ' This value is not good enough to allow',ch10,& 116 & ' the computation of the eigenvectors of the SCF cycle.',ch10,& 117 & ' It should be better than 1.0d-8 .',ch10,& 118 & ' Action : improve the accuracy of your starting wavefunctions.' 119 return 120 end if 121 122 ! Also transfer vtrial in vtrial_old 123 vtrial0(:,:)=vtrial(:,:) 124 125 ! In order to start the search for eigenvectors, 126 ! use the tiny residual vector, renormalized 127 factor=sqrt(fix_resid/resid_old) 128 work(:,:,1)=vrespc(:,:)*factor 129 vtrial(:,:)=vtrial0(:,:)+work(:,:,1) 130 131 ! If istep is not equal to 1 132 else if(istep>=2)then 133 ! 134 ! Compute the corresponding operator expectation value 135 ! And put the residual vector minus the difference 136 ! between vtrial and vtrial_old 137 ! (this is actually the action of the operator !) in vect(*,2) 138 eigen_scf=0.0_dp 139 do isp=1,nspden 140 do ifft=1,nfft 141 eigen_scf=eigen_scf+& 142 & work(ifft,isp,1) * vrespc(ifft,isp) 143 end do 144 end do 145 146 do isp=1,nspden 147 do ifft=1,nfft 148 vrespc(ifft,isp)=vrespc(ifft,isp)& 149 & +vtrial(ifft,isp)-vtrial0(ifft,isp) 150 work(ifft,isp,2)=vrespc(ifft,isp) 151 end do 152 end do 153 eigen_scf=eigen_scf/fix_resid 154 write(message, '(a,es12.4,a)' ) & 155 & ' scfeig : Operator expectation value ',eigen_scf,' (extremal eigenvalue * diemix)' 156 call wrtout(std_out,message,'COLL') 157 call wrtout(ab_out,message,'COLL') 158 ! 159 ! Compute residual of vect(*,2) 160 resid_new=zero 161 do isp=1,min(nspden,2) 162 do ifft=1,nfft 163 resid_new=resid_new+ work(ifft,isp,2) ** 2 164 end do 165 end do 166 if (nspden==4) then 167 do ifft=1,nfft 168 resid_new=resid_new+two*(work(ifft,3,2)**2+work(ifft,4,2)**2) 169 end do 170 end if 171 factor=sqrt(fix_resid/resid_new) 172 if(eigen_scf<zero) then 173 factor=-factor ! the new vector MAY be oposite to the old one 174 ! if(factor<-one) factor=-factor ! the new vector is not opposed to the old one 175 end if 176 write(message, '(a,es12.4)' ) & 177 & ' scfeig : Inverse of renormalization factor ',one/factor 178 call wrtout(std_out,message,'COLL') 179 call wrtout(ab_out,message,'COLL') 180 write(message, '(a,es12.4)' ) & 181 & ' scfeig : Convergence criterion value (->0 at convergency) ',one/factor-eigen_scf-one 182 call wrtout(std_out,message,'COLL') 183 call wrtout(ab_out,message,'COLL') 184 185 work(:,:,1)=work(:,:,2)*factor 186 vtrial(:,:)=vtrial0(:,:)+work(:,:,1) 187 ! End the different istep cases 188 end if 189 190 end subroutine scfeig