TABLE OF CONTENTS
ABINIT/pawxenergy [ Functions ]
NAME
pawxenergy
FUNCTION
Compute contributions to energy for PAW+ local exact exchange calculations
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (BA,FJ) 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/Infos/contributors.
INPUTS
pawprtvol=control print volume and debugging output for PAW pawrhoij <type(pawrhoij_type)>= paw rhoij occupancies and related data pawtab <type(pawtab_type)>=paw tabulated starting data: %lexexch=l used for local exact-exchange %vex(2*lexexch+1*4)=screened coulomb matrix
SIDE EFFECTS
eexex=energy is updated with the contribution of the cuyrrent atom
PARENTS
pawdenpot
CHILDREN
wrtout
SOURCE
34 #if defined HAVE_CONFIG_H 35 #include "config.h" 36 #endif 37 38 #include "abi_common.h" 39 40 subroutine pawxenergy(eexex,pawprtvol,pawrhoij,pawtab) 41 42 use m_profiling_abi 43 use defs_basis 44 use m_errors 45 46 use m_pawtab, only : pawtab_type 47 use m_pawrhoij, only : pawrhoij_type 48 49 !This section has been created automatically by the script Abilint (TD). 50 !Do not modify the following lines by hand. 51 #undef ABI_FUNC 52 #define ABI_FUNC 'pawxenergy' 53 use interfaces_14_hidewrite 54 !End of the abilint section 55 56 implicit none 57 58 !Arguments --------------------------------------------- 59 !scalars 60 integer,intent(in) :: pawprtvol 61 real(dp),intent(inout) :: eexex 62 type(pawrhoij_type),intent(in) :: pawrhoij 63 type(pawtab_type),intent(in) :: pawtab 64 65 !Local variables --------------------------------------- 66 !scalars 67 integer :: irhoij,irhoij1,ispden,jrhoij,jrhoij1,klmn,klmn1,lexexch,ll,m11,m21,m31,m41,n1 68 integer :: n2,n3,n4,nk,nn1,nn2 69 real(dp) :: eexextemp 70 character(len=500) :: message 71 !arrays 72 integer :: indn(3,3) 73 real(dp) :: factnk(6) 74 75 ! ***************************************************** 76 77 DBG_ENTER("COLL") 78 79 lexexch=pawtab%lexexch 80 if (pawtab%nproju==1) nk=1 81 if (pawtab%nproju==2) nk=6 82 factnk(1)=one;factnk(2)=one;factnk(3)=one 83 factnk(4)=two;factnk(5)=two;factnk(6)=two 84 indn(1,1)=1;indn(1,2)=4;indn(1,3)=5 85 indn(2,1)=4;indn(2,2)=2;indn(2,3)=6 86 indn(3,1)=5;indn(3,2)=6;indn(3,3)=3 87 88 !====================================================== 89 !Compute local exact exchange Energy 90 !----------------------------------------------------- 91 eexextemp=zero 92 93 do ispden=1,pawrhoij%nspden 94 jrhoij=1 95 do irhoij=1,pawrhoij%nrhoijsel 96 klmn=pawrhoij%rhoijselect(irhoij) 97 if(pawtab%indklmn(3,klmn)==0.and.pawtab%indklmn(4,klmn)==2*lexexch) then 98 m11=pawtab%klmntomn(1,klmn);m21=pawtab%klmntomn(2,klmn) 99 n1=pawtab%klmntomn(3,klmn);n2=pawtab%klmntomn(4,klmn) 100 nn1=(n1*n2)/2+1 101 jrhoij1=1 102 do irhoij1=1,pawrhoij%nrhoijsel 103 klmn1=pawrhoij%rhoijselect(irhoij1) 104 if(pawtab%indklmn(3,klmn1)==0.and.pawtab%indklmn(4,klmn1)==2*lexexch) then 105 m31=pawtab%klmntomn(1,klmn1);m41=pawtab%klmntomn(2,klmn1) 106 n3=pawtab%klmntomn(3,klmn1);n4=pawtab%klmntomn(4,klmn1) 107 nn2=(n3*n4)/2+1 108 do ll=1,lexexch+1 109 ! eexextemp=eexextemp-pawtab%vex(m11,m31,m41,m21,ll)*factnk(indn(nn1,nn2))*pawtab%fk(indn(nn1,nn2),ll)& 110 ! & *pawrhoij%rhoijp(jrhoij,ispden)*pawrhoij%rhoijp(jrhoij1,ispden) 111 eexextemp=eexextemp-pawtab%vex(m11,m31,m41,m21,ll)*pawtab%dltij(klmn)*pawtab%fk(indn(nn1,nn2),ll)& 112 & *pawtab%dltij(klmn1)*pawrhoij%rhoijp(jrhoij,ispden)*pawrhoij%rhoijp(jrhoij1,ispden) 113 end do 114 end if 115 jrhoij1=jrhoij1+pawrhoij%cplex 116 end do !irhoij1 117 end if 118 jrhoij=jrhoij+pawrhoij%cplex 119 end do !irhoij 120 end do ! ispden 121 eexextemp=eexextemp/two 122 eexex=eexex+eexextemp*pawtab%exchmix 123 124 if (abs(pawprtvol)>=2) then 125 write(message, '(a)' )" Contributions to the direct expression of energy:" 126 call wrtout(std_out,message,'COLL') 127 write(message,fmt='(a,e20.10,a)') " HF exchange energy =",eexextemp,ch10 128 call wrtout(std_out,message,'COLL') 129 end if 130 131 DBG_EXIT("COLL") 132 133 end subroutine pawxenergy