TABLE OF CONTENTS
- ABINIT/m_matlu
- m_matlu/add_matlu
- m_matlu/checkdiag_matlu
- m_matlu/checkreal_matlu
- m_matlu/chg_repr_matlu
- m_matlu/conjg_matlu
- m_matlu/copy_matlu
- m_matlu/destroy_matlu
- m_matlu/diag_matlu
- m_matlu/diff_matlu
- m_matlu/fac_matlu
- m_matlu/gather_matlu
- m_matlu/identity_matlu
- m_matlu/init_matlu
- m_matlu/inverse_matlu
- m_matlu/ln_matlu
- m_matlu/matlu_type
- m_matlu/print_matlu
- m_matlu/printplot_matlu
- m_matlu/prod_matlu
- m_matlu/rotate_matlu
- m_matlu/shift_matlu
- m_matlu/slm2ylm_matlu
- m_matlu/sym_matlu
- m_matlu/trace_matlu
- m_matlu/zero_matlu
ABINIT/m_matlu [ Modules ]
NAME
m_matlu
FUNCTION
COPYRIGHT
Copyright (C) 2006-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
OUTPUT
NOTES
subroutines in this module must never call a subroutine of m_oper, m_green, m_self in order to avoid circular dependancies
SOURCE
24 #if defined HAVE_CONFIG_H 25 #include "config.h" 26 #endif 27 28 29 #include "abi_common.h" 30 31 MODULE m_matlu 32 33 use defs_basis 34 use m_errors 35 use m_abicore 36 37 use m_hide_lapack, only : xginv 38 use m_geometry, only : symredcart 39 40 implicit none 41 42 private 43 44 public :: init_matlu 45 public :: inverse_matlu 46 public :: destroy_matlu 47 public :: diff_matlu 48 public :: add_matlu 49 public :: print_matlu 50 public :: sym_matlu 51 public :: copy_matlu 52 public :: gather_matlu 53 public :: zero_matlu 54 public :: trace_matlu 55 public :: diag_matlu 56 public :: rotate_matlu 57 public :: shift_matlu 58 public :: checkdiag_matlu 59 public :: checkreal_matlu 60 public :: prod_matlu 61 public :: conjg_matlu 62 public :: ln_matlu 63 public :: slm2ylm_matlu 64 public :: fac_matlu 65 public :: printplot_matlu 66 public :: identity_matlu
m_matlu/add_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
add_matlu
FUNCTION
INPUTS
maltu1 <type(matlu_type)>= density matrix matlu1 in the local orbital basis and related variables maltu2 <type(matlu_type)>= density matrix matlu2 in the local orbital basis and related variables natom = number of atoms sign_matlu2= 1 add matlu1 and matlu2 -1 substract matlu2 to matlu1
OUTPUT
maltu3 <type(matlu_type)>= density matrix matlu3, sum/substract matlu1 and matlu2
SOURCE
978 subroutine add_matlu(matlu1,matlu2,matlu3,natom,sign_matlu2) 979 980 use defs_basis 981 use m_paw_dmft, only : paw_dmft_type 982 use m_crystal, only : crystal_t 983 implicit none 984 985 !Arguments ------------------------------------ 986 !type 987 integer,intent(in) :: natom,sign_matlu2 988 type(matlu_type), intent(in) :: matlu1(natom),matlu2(natom) 989 type(matlu_type), intent(inout) :: matlu3(natom) !vz_i 990 991 !local variables------------------------------- 992 integer :: iatom 993 ! ********************************************************************* 994 995 do iatom = 1 , natom 996 matlu3(iatom)%mat=matlu1(iatom)%mat+float(sign_matlu2)*matlu2(iatom)%mat 997 enddo ! natom 998 999 end subroutine add_matlu
m_matlu/checkdiag_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
checkdiag_matlu
FUNCTION
Check that matlu is diagonal in the orbital index with given precision
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
matlu(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: input quantity to rotate natom=number of atoms in cell. tol : precision
SIDE EFFECTS
NOTES
SOURCE
2250 subroutine checkdiag_matlu(matlu,natom,tol,nondiag) 2251 use defs_basis 2252 use defs_wvltypes 2253 use m_crystal, only : crystal_t 2254 implicit none 2255 2256 !Arguments ------------------------------------ 2257 !scalars 2258 real(dp),intent(in) :: tol 2259 integer, intent(in) :: natom 2260 logical, intent(out) :: nondiag 2261 !arrays 2262 type(matlu_type),intent(inout) :: matlu(natom) 2263 !Local variables------------------------------- 2264 !scalars 2265 integer :: iatom,im,im1,ispinor,ispinor1,isppol 2266 integer :: lpawu 2267 !arrays 2268 !************************************************************************ 2269 nondiag=.false. 2270 do iatom=1,natom 2271 lpawu=matlu(iatom)%lpawu 2272 if(lpawu.ne.-1) then 2273 do im=1,2*lpawu+1 2274 do im1=1,2*lpawu+1 2275 do isppol=1,matlu(1)%nsppol 2276 do ispinor=1,matlu(1)%nspinor 2277 ! if(im/=im1) write(std_out,*) "im,im1",im,im1,matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor) 2278 ! if(present(nondiag).eqv..false.) then 2279 ! if(im/=im1.and.(abs(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor))>tol)) then 2280 ! write(message,'(5i5)') im,im1,isppol,ispinor,ispinor 2281 ! call wrtout(std_out,message,'COLL') 2282 ! write(message,'(a,3e16.5)')" checkdiag_matlu: Warning ",matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor),tol 2283 ! call wrtout(std_out,message,'COLL') 2284 ! if(.not.present(opt)) ABI_ERROR("not present(opt)") 2285 ! if(matlu(1)%nspinor==1) ABI_ERROR("matlu%nspinor==1") 2286 ! endif 2287 ! endif 2288 do ispinor1=1,matlu(1)%nspinor 2289 ! if(present(nondiag)) then 2290 if((im/=im1.or.ispinor/=ispinor1)& 2291 & .and.abs(real(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)))>tol) then 2292 nondiag=.true. 2293 ! write(6,*) "NONDIAG", matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1) 2294 endif 2295 !if(ispinor/=ispinor1.and.(abs(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1))>tol)) then 2296 ! write(message,'(a,3e16.5)')" checkdiag_matlu :i Warning ",matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1),tol 2297 ! call wrtout(std_out,message,'COLL') 2298 ! write(message,'(5i5)') im,im1,isppol,ispinor,ispinor 2299 ! call wrtout(std_out,message,'COLL') 2300 ! if(matlu(1)%nspinor==1) ABI_ERROR("matlu%nspinor==1") 2301 !endif 2302 enddo 2303 enddo ! ispinor 2304 enddo ! isppol 2305 enddo ! im1 2306 enddo ! im 2307 endif ! lpawu 2308 enddo ! iatom 2309 2310 end subroutine checkdiag_matlu
m_matlu/checkreal_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
checkreal_matlu
FUNCTION
Check that matlu is real in the orbital index with given precision
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
matlu(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: input quantity to rotate natom=number of atoms in cell. tol : precision
SIDE EFFECTS
NOTES
SOURCE
2149 subroutine checkreal_matlu(matlu,natom,tol) 2150 use defs_basis 2151 use defs_wvltypes 2152 use m_crystal, only : crystal_t 2153 implicit none 2154 2155 !Arguments ------------------------------------ 2156 !scalars 2157 real(dp),intent(in) :: tol 2158 integer, intent(in) :: natom 2159 !arrays 2160 type(matlu_type),intent(inout) :: matlu(natom) 2161 !Local variables------------------------------- 2162 !scalars 2163 integer :: iatom,im,im1,ispinor,ispinor1,isppol 2164 integer :: lpawu 2165 character(len=500) :: message 2166 real(dp) :: maximag,maxoffdiag,maximagdiag 2167 !arrays 2168 !************************************************************************ 2169 maximag=zero 2170 maximagdiag=zero 2171 maxoffdiag=zero 2172 do iatom=1,natom 2173 lpawu=matlu(iatom)%lpawu 2174 if(lpawu.ne.-1) then 2175 do im=1,2*lpawu+1 2176 do im1=1,2*lpawu+1 2177 do isppol=1,matlu(1)%nsppol 2178 do ispinor=1,matlu(1)%nspinor 2179 do ispinor1=1,matlu(1)%nspinor 2180 if(abs(aimag(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)))>maximag) then 2181 maximag=abs(aimag(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1))) 2182 if(im==im1.and.ispinor==ispinor1) maximagdiag=abs(aimag(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1))) 2183 endif 2184 if((im/=im1.or.ispinor/=ispinor1).and.abs(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1))>maxoffdiag) then 2185 maxoffdiag=abs(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)) 2186 endif 2187 enddo 2188 enddo ! ispinor 2189 enddo ! isppol 2190 enddo ! im1 2191 enddo ! im 2192 endif ! lpawu 2193 enddo ! iatom 2194 if (maximagdiag>tol) then 2195 write(message,'(3x,2a,e12.4,a,e12.4,2a)') ch10,& 2196 & ' Diagonal part of the occupation matrix is complex: the imaginary part ',& 2197 & maximagdiag,' is larger than',tol,ch10 & 2198 & , "The calculation cannot handle it : check that your calculation is meaningfull" 2199 ABI_ERROR(message) 2200 endif 2201 if (maximag>tol) then 2202 write(message,'(3x,2a,e12.4,a,e12.4,2a)') ch10,& 2203 & ' Off diag occupation matrix is complex: the imaginary part ',maximag,' is larger than',tol,ch10& 2204 , "Check that your calculation is meaningfull" 2205 ABI_WARNING(message) 2206 else 2207 write(message,'(3x,2a,e12.4,a,e12.4,2a)') ch10,& 2208 & ' Occupation matrix is real: the imaginary part ',maximag,' is lower than',tol 2209 ABI_COMMENT(message) 2210 endif 2211 if (maxoffdiag>tol) then 2212 write(message,'(3x,2a,e12.4,a,e12.4,6a)') ch10,& 2213 & ' Occupation matrix is non diagonal : the maximum off-diag part ',maxoffdiag,' is larger than',tol,ch10& 2214 & , "The corresponding non diagonal elements will be neglected in the Weiss/Hybridization functions",ch10& 2215 & , "(Except if dmft_solv=8,9 where these elements are taken into accounts)",ch10& 2216 & , "This is an approximation" 2217 ABI_WARNING(message) 2218 else 2219 write(message,'(3x,2a,e12.4,a,e12.4,2a)') ch10,& 2220 & ' Occupation matrix is diagonal : the off-diag part ',maxoffdiag,' is lower than',tol 2221 ABI_COMMENT(message) 2222 endif 2223 2224 end subroutine checkreal_matlu
m_matlu/chg_repr_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
chg_repr_matlu
FUNCTION
Change representation of density matrix (useful for nspinor=2)
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
glocspsp(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: density matrix in the spin spin representation glocnm(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: density matrix in the magnetization representation natom=number of atoms in cell. option= 1 glocspsp is input, glocnm is computed option= -1 glocspsp is computed, glocnm is input prtopt= option to define level of printing
OUTPUT
glocspsp(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: density matrix in the spin spin representation glocnm(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: density matrix in the magnetization representation
SIDE EFFECTS
NOTES
SOURCE
1032 subroutine chg_repr_matlu(glocspsp,glocnm,natom,option,prtopt) 1033 use defs_basis 1034 use defs_wvltypes 1035 use m_crystal, only : crystal_t 1036 implicit none 1037 1038 !Arguments ------------------------------------ 1039 !scalars 1040 integer,intent(in) :: natom,option,prtopt 1041 !arrays 1042 type(matlu_type),intent(inout) :: glocspsp(natom) 1043 type(matlu_type),intent(inout) :: glocnm(natom) 1044 !Local variables------------------------------- 1045 !scalars 1046 integer :: iatom,isppol,lpawu,m1,m2,ndim,nsppol,mu 1047 complex(dpc) :: ci 1048 character(len=500) :: message 1049 1050 ! DBG_ENTER("COLL") 1051 1052 ci=j_dpc 1053 1054 !== Compute density matrix in density magnetization representation 1055 if (option==1) then 1056 nsppol=glocspsp(1)%nsppol 1057 do isppol=1,nsppol 1058 do iatom=1,natom 1059 if(glocspsp(iatom)%lpawu/=-1) then 1060 ndim=2*glocspsp(iatom)%lpawu+1 1061 do m1=1,ndim 1062 do m2=1,ndim 1063 glocnm(iatom)%mat(m1,m2,isppol,1,1)=glocspsp(iatom)%mat(m1,m2,isppol,1,1) & 1064 & +glocspsp(iatom)%mat(m1,m2,isppol,2,2) 1065 glocnm(iatom)%mat(m1,m2,isppol,4,1)=glocspsp(iatom)%mat(m1,m2,isppol,1,1) & 1066 & -glocspsp(iatom)%mat(m1,m2,isppol,2,2) 1067 glocnm(iatom)%mat(m1,m2,isppol,2,1)=glocspsp(iatom)%mat(m1,m2,isppol,1,2) & 1068 & +glocspsp(iatom)%mat(m1,m2,isppol,2,1) 1069 glocnm(iatom)%mat(m1,m2,isppol,3,1)= & 1070 & cmplx((aimag(glocspsp(iatom)%mat(m1,m2,isppol,2,1)) & 1071 & -aimag(glocspsp(iatom)%mat(m1,m2,isppol,1,2))), & 1072 & (-real(glocspsp(iatom)%mat(m1,m2,isppol,2,1))+ & 1073 & real(glocspsp(iatom)%mat(m1,m2,isppol,1,2))),kind=dp) 1074 enddo ! m2 1075 enddo ! m1 1076 if(abs(prtopt)>=3) then 1077 write(message,'(a)') " -- in n, m repr " 1078 call wrtout(std_out, message,'COLL') 1079 do mu=1,4 1080 do m1=1,ndim 1081 write(message,'(8x,(14(2f9.5,2x)))')(glocnm(iatom)%mat(m1,m2,isppol,mu,1),m2=1,ndim) 1082 call wrtout(std_out, message,'COLL') 1083 enddo ! m1 1084 write(message,'(a)') ch10 1085 call wrtout(std_out, message,'COLL') 1086 enddo ! mu 1087 endif ! prtopt >3 1088 endif ! lpawu/=-1 1089 enddo 1090 enddo 1091 1092 !== Compute back density matrix in upup dndn updn dnup representation 1093 else if (option==-1) then 1094 isppol=1 1095 do iatom=1,natom 1096 if(glocnm(iatom)%lpawu/=-1) then 1097 lpawu=glocnm(iatom)%lpawu 1098 ndim=2*glocnm(iatom)%lpawu+1 1099 do m1=1, 2*lpawu+1 1100 do m2=1, 2*lpawu+1 1101 glocspsp(iatom)%mat(m1,m2,isppol,1,1)=half*(glocnm(iatom)%mat(m1,m2,isppol,1,1)+glocnm(iatom)%mat(m1,m2,isppol,4,1)) 1102 glocspsp(iatom)%mat(m1,m2,isppol,2,2)=half*(glocnm(iatom)%mat(m1,m2,isppol,1,1)-glocnm(iatom)%mat(m1,m2,isppol,4,1)) 1103 glocspsp(iatom)%mat(m1,m2,isppol,1,2)=half*(glocnm(iatom)%mat(m1,m2,isppol,2,1)-ci*glocnm(iatom)%mat(m1,m2,isppol,3,1)) 1104 glocspsp(iatom)%mat(m1,m2,isppol,2,1)=half*(glocnm(iatom)%mat(m1,m2,isppol,2,1)+ci*glocnm(iatom)%mat(m1,m2,isppol,3,1)) 1105 end do ! m2 1106 end do ! m1 1107 if(abs(prtopt)>6) then 1108 write(message,'(a)') " -- in spin spin repr " 1109 call wrtout(std_out, message,'COLL') 1110 do mu=1,4 1111 do m1=1,ndim 1112 write(message,'(8x,14(2f9.5,2x))')(glocspsp(iatom)%mat(m1,m2,isppol,mu,1),m2=1,ndim) 1113 call wrtout(std_out, message,'COLL') 1114 enddo 1115 write(message,'(a)') ch10 1116 call wrtout(std_out, message,'COLL') 1117 enddo 1118 endif !prtopt>3 1119 endif ! lpawu/=-1 1120 end do ! iatom 1121 else 1122 message = "stop in chg_repr_matlu" 1123 ABI_ERROR(message) 1124 endif 1125 1126 1127 ! DBG_EXIT("COLL") 1128 1129 end subroutine chg_repr_matlu
m_matlu/conjg_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
conjg_matlu
FUNCTION
conjugate of input matlu
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
matlu1(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: input quantity
SIDE EFFECTS
NOTES
SOURCE
2405 subroutine conjg_matlu(matlu1,natom) 2406 use defs_basis 2407 use defs_wvltypes 2408 implicit none 2409 2410 !Arguments ------------------------------------ 2411 !scalars 2412 integer, intent(in) :: natom 2413 !arrays 2414 type(matlu_type), intent(inout) :: matlu1(natom) 2415 !Local variables------------------------------- 2416 !scalars 2417 integer :: iatom,im1,im2,ispinor2,ispinor1,isppol 2418 integer :: lpawu 2419 !arrays 2420 !************************************************************************ 2421 do iatom=1,natom 2422 lpawu=matlu1(iatom)%lpawu 2423 if(lpawu.ne.-1) then 2424 do isppol=1,matlu1(1)%nsppol 2425 do ispinor1=1,matlu1(1)%nspinor 2426 do ispinor2=1,matlu1(1)%nspinor 2427 do im1=1,2*lpawu+1 2428 do im2=1,2*lpawu+1 2429 matlu1(iatom)%mat(im1,im2,isppol,ispinor1,ispinor2)= & 2430 & conjg(matlu1(iatom)%mat(im1,im2,isppol,ispinor1,ispinor2)) 2431 enddo ! im2 2432 enddo ! im1 2433 enddo ! ispinor2 2434 enddo ! ispinor1 2435 enddo ! isppol 2436 endif ! lpawu 2437 enddo ! iatom 2438 2439 end subroutine conjg_matlu
m_matlu/copy_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
copy_matlu
FUNCTION
Copy matlu1 into matlu2
INPUTS
maltu1 <type(matlu_type)>= density matrix matlu1 in the local orbital basis and related variables natom = number of atoms
OUTPUT
maltu2 <type(matlu_type)>= density matrix matlu2 in the local orbital basis and related variables
SOURCE
290 subroutine copy_matlu(nmat1,nmat2,natom,opt_diag,opt_non_diag,opt_re) 291 292 use defs_basis 293 implicit none 294 295 !Arguments ------------------------------------ 296 !type 297 integer, intent(in) :: natom 298 type(matlu_type),intent(in) :: nmat1(natom) 299 type(matlu_type),intent(inout) :: nmat2(natom) !vz_i 300 integer, optional, intent(in) :: opt_diag,opt_non_diag,opt_re 301 302 !Local variables------------------------------- 303 integer :: iatom,isppol,im1,im2,ispinor,ispinor1,tndim 304 ! ********************************************************************* 305 306 do isppol=1,nmat1(1)%nsppol 307 do iatom=1,natom 308 if(nmat1(iatom)%lpawu.ne.-1) then 309 tndim=(2*nmat1(iatom)%lpawu+1) 310 do im1=1,tndim 311 do im2=1,tndim 312 do ispinor=1,nmat1(1)%nspinor 313 do ispinor1=1,nmat1(1)%nspinor 314 if(present(opt_diag)) then 315 nmat2(iatom)%mat(im1,im1,isppol,ispinor,ispinor)=nmat1(iatom)%mat(im1,im1,isppol,ispinor,ispinor) 316 else if(present(opt_non_diag)) then 317 if(im1/=im2.or.ispinor/=ispinor1) then 318 nmat2(iatom)%mat(im1,im2,isppol,ispinor,ispinor1)=nmat1(iatom)%mat(im1,im2,isppol,ispinor,ispinor1) 319 endif 320 else 321 nmat2(iatom)%mat(im1,im2,isppol,ispinor,ispinor1)=nmat1(iatom)%mat(im1,im2,isppol,ispinor,ispinor1) 322 if(present(opt_re)) nmat2(iatom)%mat(im1,im2,isppol,ispinor,ispinor1)=& 323 & cmplx(real(nmat1(iatom)%mat(im1,im2,isppol,ispinor,ispinor1)),0.d0,kind=dp) 324 endif 325 enddo 326 enddo 327 enddo 328 enddo 329 endif ! lpawu=-1 330 enddo ! iatom 331 enddo ! isppol 332 ! do iatom=1,natom 333 ! lpawu=nmat1(iatom)%lpawu 334 ! if(lpawu.ne.-1) then 335 ! nmat2(iatom)%mat=nmat1(iatom)%mat 336 ! endif 337 ! enddo 338 339 340 end subroutine copy_matlu
m_matlu/destroy_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
destroy_matlu
FUNCTION
deallocate matlu
INPUTS
maltu <type(matlu_type)>= density matrix in the local orbital basis and related variables natom = number of atoms
OUTPUT
SOURCE
249 subroutine destroy_matlu(matlu,natom) 250 251 use defs_basis 252 use m_crystal, only : crystal_t 253 implicit none 254 255 !Arguments ------------------------------------ 256 !scalars 257 integer, intent(in) :: natom 258 type(matlu_type),intent(inout) :: matlu(natom) 259 260 !Local variables------------------------------- 261 integer :: iatom 262 263 ! ********************************************************************* 264 265 do iatom=1,natom 266 if ( allocated(matlu(iatom)%mat) ) then 267 ABI_FREE(matlu(iatom)%mat) 268 end if 269 enddo 270 271 end subroutine destroy_matlu
m_matlu/diag_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
diag_matlu
FUNCTION
Diagonalize matlu matrix
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
matlu(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: input quantity to diagonalize natom=number of atoms prtopt= option to define level of printing nsppol_imp= if 1, one can diagonalize with the same matrix the Up and Dn matlu matrix. It is convenient because one can thus have the same interaction matrix for up and dn spins.
OUTPUT
matlu_diag(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: diagonalized density matrix eigvectmatlu(natom) <type(coeff2c_type)> = Eigenvectors corresponding to the diagonalization
SIDE EFFECTS
NOTES
SOURCE
1386 subroutine diag_matlu(matlu,matlu_diag,natom,prtopt,eigvectmatlu,nsppol_imp,checkstop,optreal,test) 1387 use defs_basis 1388 use defs_wvltypes 1389 use m_crystal, only : crystal_t 1390 use m_matrix, only : blockdiago_fordsyev,blockdiago_forzheev 1391 implicit none 1392 1393 !Arguments ------------------------------------ 1394 !scalars 1395 integer, intent(in) :: natom 1396 integer, intent(in) :: prtopt 1397 !arrays 1398 type(matlu_type),intent(inout) :: matlu(natom) 1399 type(matlu_type),intent(inout) :: matlu_diag(natom) !vz_i 1400 type(coeff2c_type),optional,intent(inout) :: eigvectmatlu(natom,matlu(1)%nsppol) !vz_i 1401 integer,optional,intent(in) :: nsppol_imp 1402 logical,optional,intent(in) :: checkstop 1403 integer,optional,intent(in) :: optreal 1404 integer,optional,intent(in) :: test 1405 !Local variables------------------------------- 1406 !scalars 1407 integer :: iatom,im1,im2,im3,imc,imc1,info,ispinor,ispinor1,isppol,lwork,tndim,lworkr 1408 integer :: nsppol,nsppolimp,nspinor 1409 logical :: checkstop_in,blockdiag 1410 character(len=500) :: message 1411 !arrays 1412 type(coeff2c_type),allocatable :: gathermatlu(:) 1413 real(dp),allocatable :: eig(:),rwork(:),work(:),valuer(:,:)!,valuer2(:,:) 1414 !real(dp),allocatable :: valuer3(:,:),valuer4(:,:) 1415 ! real(dp),allocatable :: eigvec(:,:) 1416 complex(dpc),allocatable :: zwork(:) 1417 logical :: donotdiag,print_temp_mat2 1418 complex(dpc),allocatable :: temp_mat(:,:) 1419 complex(dpc),allocatable :: temp_mat2(:,:) 1420 !debug complex(dpc),allocatable :: temp_mat3(:,:) 1421 !************************************************************************ 1422 1423 call zero_matlu(matlu_diag,natom) 1424 nsppol=matlu(1)%nsppol 1425 if(present(nsppol_imp)) then 1426 nsppolimp=nsppol_imp 1427 else 1428 nsppolimp=nsppol 1429 endif 1430 if(present(checkstop)) then 1431 checkstop_in=checkstop 1432 else 1433 checkstop_in=.true. 1434 endif 1435 if(present(test)) then 1436 blockdiag=(test==8) 1437 else 1438 blockdiag=.false. 1439 endif 1440 nspinor=matlu(1)%nspinor 1441 1442 donotdiag=.true. 1443 donotdiag=.false. 1444 ! =========================== 1445 ! Check is diagonalization is necessary and how 1446 ! =========================== 1447 do isppol=1,matlu(1)%nsppol 1448 do iatom=1,natom 1449 if(matlu(iatom)%lpawu.ne.-1) then 1450 tndim=(2*matlu(iatom)%lpawu+1) 1451 do im1=1,tndim 1452 do im2=1,tndim 1453 do ispinor=1,nspinor 1454 do ispinor1=1,nspinor 1455 if(abs(matlu(iatom)%mat(im1,im2,isppol,ispinor,ispinor1))>tol8.and.& 1456 & (im1/=im2.or.ispinor/=ispinor1)) then 1457 ! if matrix is diagonal: do not diagonalize 1458 donotdiag=.false. 1459 exit 1460 endif 1461 enddo 1462 enddo 1463 enddo 1464 enddo 1465 endif 1466 enddo 1467 enddo 1468 1469 if(donotdiag) then 1470 do isppol=1,matlu(1)%nsppol 1471 do iatom=1,natom 1472 if(matlu(iatom)%lpawu.ne.-1) then 1473 tndim=(2*matlu(iatom)%lpawu+1) 1474 eigvectmatlu(iatom,isppol)%value(:,:)=czero 1475 do im1=1,tndim 1476 eigvectmatlu(iatom,isppol)%value(im1,im1)=cone 1477 enddo 1478 endif 1479 enddo 1480 enddo 1481 call copy_matlu(matlu,matlu_diag,natom) 1482 write(message,'(a)') " Diagonalisation of matlu will not be performed" 1483 call wrtout(std_out,message,'COLL') 1484 return 1485 endif 1486 1487 1488 ! Below, gathermatlu is used to store eigenvectors (after diagonalization) 1489 ! For nsppol=2, and if nsppolimp=1, these eigenvectors computed for isppol=1, and applied to 1490 ! rotate matrix for isppol=2. It is the reason why the sum below is only 1491 ! from 1 to nsppolimp ! 1492 1493 1494 1495 do isppol=1,nsppolimp ! 1496 ! =========================== 1497 ! Define gathermatlu 1498 ! =========================== 1499 ABI_MALLOC(gathermatlu,(natom)) 1500 do iatom=1,natom 1501 if(matlu(iatom)%lpawu.ne.-1) then 1502 tndim=nspinor*(2*matlu(iatom)%lpawu+1) 1503 ABI_MALLOC(gathermatlu(iatom)%value,(tndim,tndim)) 1504 gathermatlu(iatom)%value=czero 1505 endif 1506 enddo 1507 if(nsppol==1.and.nspinor==2) then 1508 call gather_matlu(matlu,gathermatlu,natom,option=1,prtopt=1) 1509 else if((nsppol==2.or.nsppol==1).and.nspinor==1) then 1510 do iatom=1,natom 1511 if(matlu(iatom)%lpawu.ne.-1) then 1512 tndim=nspinor*(2*matlu(iatom)%lpawu+1) 1513 do im1=1,tndim 1514 do im2=1,tndim 1515 gathermatlu(iatom)%value(im1,im2)=matlu(iatom)%mat(im1,im2,isppol,1,1) 1516 enddo 1517 enddo 1518 endif 1519 enddo 1520 endif 1521 ! =========================== 1522 ! Diagonalize 1523 ! =========================== 1524 do iatom=1,natom 1525 if(matlu(iatom)%lpawu.ne.-1) then 1526 if(isppol<=nsppolimp) then 1527 tndim=nspinor*(2*matlu(iatom)%lpawu+1) 1528 !debug allocate(temp_mat2(tndim,tndim)) 1529 !debug temp_mat2=zero 1530 lwork=2*tndim-1 1531 ABI_MALLOC(rwork,(3*tndim-2)) 1532 rwork = zero 1533 1534 lworkr=tndim*(tndim+2)*2 1535 ABI_MALLOC(work,(lworkr)) 1536 work = zero 1537 ABI_MALLOC(valuer,(tndim,tndim)) 1538 ! ABI_MALLOC(valuer2,(tndim,tndim)) 1539 ! ABI_MALLOC(valuer3,(tndim,tndim)) 1540 ! ABI_MALLOC(valuer4,(tndim,tndim)) 1541 ABI_MALLOC(zwork,(lwork)) 1542 ! valuer2=zero 1543 ! valuer3=zero 1544 ! valuer4=zero 1545 zwork = czero 1546 ABI_MALLOC(eig,(tndim)) 1547 eig = zero 1548 info = 0 1549 if(prtopt>=4) then 1550 write(message,'(a,i4,a,i4)') " BEFORE DIAGONALIZATION for atom",iatom," and isppol",isppol 1551 call wrtout(std_out,message,'COLL') 1552 do im1=1,tndim 1553 write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')& 1554 & (gathermatlu(iatom)%value(im1,im2),im2=1,tndim) 1555 call wrtout(std_out,message,'COLL') 1556 end do 1557 endif 1558 !debug temp_mat2(:,:)=gathermatlu(iatom)%value(:,:) 1559 ! write(std_out,*)"diag" 1560 if(present(optreal).and.maxval(abs(aimag(gathermatlu(iatom)%value(:,:))))<tol6) then 1561 write(message,'(a,2x,a,e9.3,a)') ch10,"Imaginary part of Local Hamiltonian is lower than ",& 1562 & tol8, ": the real matrix is used" 1563 call wrtout(std_out,message,'COLL') 1564 valuer=real(gathermatlu(iatom)%value,kind=dp) 1565 ! write(message,'(a)') ch10 1566 ! call wrtout(std_out,message,'COLL') 1567 ! write(message,'(a,i4,a,i4)') "BEFORE valuer for atom",iatom," and isppol",isppol 1568 ! call wrtout(std_out,message,'COLL') 1569 ! do im1=1,tndim 1570 ! write(message,'(2(1x,18(1x,"(",f20.15,",",f20.15,")")))')& 1571 !& (valuer(im1,im2),im2=1,tndim) 1572 ! call wrtout(std_out,message,'COLL') 1573 ! end do 1574 ! do im1=1,tndim 1575 ! valuer(im1,im1)=real(im1,kind=dp)*0.00000000001_dp+valuer(im1,im1) 1576 ! enddo 1577 ! write(message,'(a)') ch10 1578 ! call wrtout(std_out,message,'COLL') 1579 ! write(message,'(a,i4,a,i4)') "BEFORE valuer for atom",iatom," and isppol",isppol 1580 ! call wrtout(std_out,message,'COLL') 1581 ! do im1=1,tndim 1582 ! write(message,'(2(1x,18(1x,f20.15,f20.15)))')& 1583 !& (valuer(im1,im2),im2=1,tndim) 1584 ! call wrtout(std_out,message,'COLL') 1585 ! end do 1586 !call dsyev('v','u',tndim,valuer,tndim,eig,work,lworkr,info) 1587 if (blockdiag) then 1588 call blockdiago_fordsyev(valuer,tndim,eig) 1589 else 1590 call dsyev('v','u',tndim,valuer,tndim,eig,work,lworkr,info) 1591 endif 1592 !! For reproductibility 1593 ! ! valuer2: eigenvector for the perturb matrix 1594 ! valuer2=real(gathermatlu(iatom)%value,kind=dp) 1595 ! do im1=1,tndim 1596 ! valuer2(im1,im1)=float(im1)*0.00000000001+valuer2(im1,im1) 1597 ! enddo 1598 ! call dsyev('v','u',tndim,valuer2,tndim,eig,work,lworkr,info) 1599 ! write(message,'(a)') ch10 1600 ! call wrtout(std_out,message,'COLL') 1601 ! write(message,'(a,i4,a,i4)') " valuer2 for atom",iatom," and isppol",isppol 1602 ! call wrtout(std_out,message,'COLL') 1603 ! do im1=1,tndim 1604 ! write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')& 1605 !& (valuer2(im1,im2),im2=1,tndim) 1606 ! call wrtout(std_out,message,'COLL') 1607 ! end do 1608 ! call dgemm('n','n',tndim,tndim,tndim,cone,valuer,tndim,& 1609 !& valuer2,tndim,czero,valuer3 ,tndim) 1610 ! call dgemm('c','n',tndim,tndim,tndim,cone,valuer2,tndim,& 1611 !& valuer3 ,tndim,czero,valuer4,tndim) 1612 ! ! valuer4: compute unpert matrix in the basis of the 1613 ! ! perturb basis 1614 ! write(message,'(a)') ch10 1615 ! call wrtout(std_out,message,'COLL') 1616 ! write(message,'(a,i4,a,i4)') "BEFORE valuer4 for atom",iatom," and isppol",isppol 1617 ! call wrtout(std_out,message,'COLL') 1618 ! do im1=1,tndim 1619 ! write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')& 1620 !& (valuer4(im1,im2),im2=1,tndim) 1621 ! call wrtout(std_out,message,'COLL') 1622 ! end do 1623 ! call dsyev('v','u',tndim,valuer4,tndim,eig,work,lworkr,info) 1624 ! ! valuer4: Diago valuer4 (nearly diag) 1625 ! write(message,'(a)') ch10 1626 ! call wrtout(std_out,message,'COLL') 1627 ! write(message,'(a,i4,a,i4)') "AFTER valuer4 for atom",iatom," and isppol",isppol 1628 ! call wrtout(std_out,message,'COLL') 1629 ! do im1=1,tndim 1630 ! write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')& 1631 !& (valuer4(im1,im2),im2=1,tndim) 1632 ! call wrtout(std_out,message,'COLL') 1633 ! end do 1634 ! call dgemm('n','n',tndim,tndim,tndim,cone,valuer2,tndim,& 1635 !& valuer4,tndim,czero,valuer ,tndim) 1636 !write(6,*) "INFO",info 1637 gathermatlu(iatom)%value=cmplx(valuer,0.d0,kind=dp) 1638 ! write(message,'(a,i4,a,i4)') "AFTER valuer for atom",iatom," and isppol",isppol 1639 ! call wrtout(std_out,message,'COLL') 1640 ! do im1=1,tndim 1641 ! write(message,'(2(1x,18(1x,"(",f20.15,",",f20.15,")")))')& 1642 !& (valuer(im1,im2),im2=1,tndim) 1643 ! call wrtout(std_out,message,'COLL') 1644 ! end do 1645 else 1646 if(present(optreal).and.maxval(abs(aimag(gathermatlu(iatom)%value(:,:))))>tol8) then 1647 write(message,'(a)') " Local hamiltonian in correlated basis is complex" 1648 ABI_COMMENT(message) 1649 endif 1650 call zheev('v','u',tndim,gathermatlu(iatom)%value,tndim,eig,zwork,lwork,rwork,info) 1651 !call blockdiago_forzheev(gathermatlu(iatom)%value,tndim,eig) 1652 endif 1653 if(prtopt>=3) then 1654 write(message,'(a)') ch10 1655 call wrtout(std_out,message,'COLL') 1656 write(message,'(a,i4,a,i4)') " EIGENVECTORS for atom",iatom," and isppol",isppol 1657 call wrtout(std_out,message,'COLL') 1658 do im1=1,tndim 1659 write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')& 1660 & (gathermatlu(iatom)%value(im1,im2),im2=1,tndim) 1661 call wrtout(std_out,message,'COLL') 1662 end do 1663 ! do im1=1,tndim 1664 ! xcheck=czero 1665 ! do im3=1,tndim 1666 ! do im2=1,tndim 1667 ! xcheck=xcheck+gathermatlu(iatom)%value(im1,im2)*conjg(gathermatlu(iatom)%value(im2,im3)) 1668 ! end do 1669 ! end do 1670 ! write(6,*) "check",im3,im1,xcheck 1671 ! end do 1672 endif 1673 ! write(std_out,*) "eig",eig 1674 ! =========================== 1675 ! Put eigenvalue in matlu_diag 1676 ! =========================== 1677 imc=0 1678 do ispinor=1,nspinor 1679 do im1=1,2*matlu(iatom)%lpawu+1 1680 imc=imc+1 1681 matlu_diag(iatom)%mat(im1,im1,isppol,ispinor,ispinor)=eig(imc) 1682 !debug temp_mat2(imc,imc)=eig(imc) 1683 enddo 1684 enddo 1685 if(prtopt>=2) then 1686 ! write(message,'(a,12(1x,18(1x,"(",f9.3,",",f9.3,")")))')& 1687 !& ch10,(eig(im1),im1=1,tndim) 1688 ! call wrtout(std_out,message,'COLL') 1689 !call wrtout(std_out,message,'COLL') 1690 !write(std_out,*) "EIG", eig 1691 endif 1692 ABI_FREE(zwork) 1693 ABI_FREE(rwork) 1694 ABI_FREE(work) 1695 ABI_FREE(valuer) 1696 ! ABI_FREE(valuer2) 1697 ! ABI_FREE(valuer3) 1698 ! ABI_FREE(valuer4) 1699 ABI_FREE(eig) 1700 ! endif 1701 ! enddo 1702 ! =========================== 1703 ! Keep eigenvectors gathermatlu 1704 ! =========================== 1705 if (present(eigvectmatlu)) then 1706 tndim=nspinor*(2*matlu(iatom)%lpawu+1) 1707 eigvectmatlu(iatom,isppol)%value(:,:)=gathermatlu(iatom)%value(:,:) 1708 ! write(std_out,*) "eigvect in diag_matlu" 1709 ! do im1=1,tndim 1710 ! write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')& 1711 !& (gathermatlu(iatom)%value(im1,im2),im2=1,tndim) 1712 ! call wrtout(std_out,message,'COLL') 1713 ! end do 1714 endif 1715 1716 1717 ! ================================================================== 1718 endif 1719 if(nsppolimp==1.and.matlu(1)%nsppol==2) then 1720 ! If necessary rotate levels for this other spin, assuming the same 1721 ! rotation matrix: it have to be checked afterwards that the matrix is 1722 ! diagonal 1723 ! =================================================================== 1724 ABI_MALLOC(temp_mat,(tndim,tndim)) 1725 ABI_MALLOC(temp_mat2,(tndim,tndim)) 1726 temp_mat(:,:)=czero 1727 ! input matrix: gathermatlu 1728 ! rotation matrix: eigvectmatlu 1729 ! intermediate matrix: temp_mat 1730 ! result matrix: temp_mat2 1731 do im1=1,tndim 1732 do im2=1,tndim 1733 gathermatlu(iatom)%value(im1,im2)=matlu(iatom)%mat(im1,im2,2,1,1) 1734 enddo 1735 enddo 1736 if(prtopt>=3) then 1737 write(message,'(a,i4,a,i4)') " GATHERMATLU for atom",iatom," inside if nsppolimp==1" 1738 call wrtout(std_out,message,'COLL') 1739 do im1=1,tndim 1740 write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')& 1741 & (gathermatlu(iatom)%value(im1,im2),im2=1,tndim) 1742 call wrtout(std_out,message,'COLL') 1743 end do 1744 endif 1745 call zgemm('n','n',tndim,tndim,tndim,cone,gathermatlu(iatom)%value,tndim,& 1746 & eigvectmatlu(iatom,1)%value,tndim,czero,temp_mat ,tndim) 1747 call zgemm('c','n',tndim,tndim,tndim,cone,eigvectmatlu(iatom,1)%value,tndim,& 1748 & temp_mat ,tndim,czero,temp_mat2,tndim) 1749 eigvectmatlu(iatom,2)%value=eigvectmatlu(iatom,1)%value 1750 imc=0 1751 print_temp_mat2=.false. 1752 do ispinor=1,nspinor 1753 do im1=1,2*matlu(iatom)%lpawu+1 1754 imc=imc+1 1755 imc1=0 1756 do ispinor1=1,nspinor 1757 do im2=1,2*matlu(iatom)%lpawu+1 1758 imc1=imc1+1 1759 matlu_diag(iatom)%mat(im1,im2,2,ispinor,ispinor1)=temp_mat2(imc,imc1) 1760 if (imc/=imc1.and.(abs(temp_mat2(imc,imc1))>tol5)) then 1761 write(message,'(3a,i4,2f16.4)') ch10,'diag_matlu= Matrix for spin number 2 obtained with', & 1762 & ' eigenvectors from diagonalization for spin nb 1 is non diagonal for atom:',iatom,& 1763 & abs(temp_mat2(imc,imc1)),tol5 1764 call wrtout(std_out,message,'COLL') 1765 if(.not.checkstop_in.and.(abs(temp_mat2(imc,imc1))>0.1)) print_temp_mat2=.true. 1766 if(checkstop_in) print_temp_mat2=.true. 1767 endif 1768 ! write(std_out,*) temp_mat2(imc,imc1) 1769 !debug temp_mat2(imc,imc)=eig(imc) 1770 enddo 1771 enddo 1772 enddo 1773 enddo 1774 if(print_temp_mat2.and.prtopt>=3) then 1775 write(message,'(a)') " temp_mat2" 1776 call wrtout(std_out,message,'COLL') 1777 do im1=1,tndim 1778 write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')& 1779 & (temp_mat2(im1,im3),im3=1,tndim) 1780 call wrtout(std_out,message,'COLL') 1781 end do 1782 if(iatom==2) then 1783 ABI_ERROR("iatom==2") 1784 end if 1785 endif 1786 ABI_FREE(temp_mat) 1787 ABI_FREE(temp_mat2) 1788 endif 1789 1790 1791 endif ! end lpawu/=-1 1792 1793 !! for check only 1794 !debug if(matlu(iatom)%lpawu.ne.-1) then 1795 !debug allocate(temp_mat(tndim,tndim)) 1796 !debug allocate(temp_mat3(tndim,tndim)) 1797 !debug do im1=1,tndim 1798 !debug do im2=1,tndim 1799 !debug! rot_mat(iatom,isppol)%value(im1,im2)=rot_mat_orig(iatom,isppol)%value(im1,im2) 1800 !debug temp_mat3(im1,im2)=conjg(gathermatlu(iatom)%value(im2,im1)) 1801 !debug enddo 1802 !debug enddo 1803 !debug temp_mat(:,:)=czero 1804 !debug! input matrix: temp_mat2 1805 !debug! rotation matrix: gathermatlu 1806 !debug! intermediate matrix: temp_mat 1807 !debug! result matrix: temp_mat2 1808 !debug call zgemm('n','c',tndim,tndim,tndim,cone,temp_mat2 ,tndim,& 1809 !debug& temp_mat3,tndim,czero,temp_mat ,tndim) 1810 !debug call zgemm('n','n',tndim,tndim,tndim,cone,temp_mat3,tndim,& 1811 !debug& temp_mat ,tndim,czero,temp_mat2,tndim) 1812 !debug! call zgemm('n','c',tndim,tndim,tndim,cone,temp_mat2 ,tndim,& 1813 !debug!& gathermatlu(iatom)%value,tndim,czero,temp_mat ,tndim) 1814 !debug! call zgemm('n','n',tndim,tndim,tndim,cone,gathermatlu(iatom)%value,tndim,& 1815 !debug!& temp_mat ,tndim,czero,temp_mat2,tndim) 1816 !debug write(std_out,*) "result" 1817 !debug do im1=1,tndim 1818 !debug write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')& 1819 !debug& (temp_mat2(im1,im2),im2=1,tndim) 1820 !debug call wrtout(std_out,message,'COLL') 1821 !debug end do 1822 !debug deallocate(temp_mat) 1823 !debug deallocate(temp_mat3) 1824 !debug endif ! lpawu 1825 1826 enddo ! iatom 1827 ! End loop over atoms 1828 ! =========================== 1829 do iatom=1,natom 1830 if(matlu(iatom)%lpawu.ne.-1) then 1831 !debug deallocate(temp_mat2) 1832 ABI_FREE(gathermatlu(iatom)%value) 1833 endif 1834 enddo 1835 ABI_FREE(gathermatlu) 1836 enddo ! isppol 1837 1838 end subroutine diag_matlu
m_matlu/diff_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
diff_matlu
FUNCTION
INPUTS
char1 = character describing matlu1 char2 = character describing matlu2 matlu1(natom) <type(matlu_type)>= density matrix 1 in the local orbital basis and related variables matlu2(natom) <type(matlu_type)>= density matrix 2 in the local orbital basis and related variables natom = number of atoms option =1 if diff> toldiff , stop 0 print diff and toldiff else do not test and do not print toldiff = maximum value for the difference between matlu1 and matlu2
OUTPUT
SOURCE
870 subroutine diff_matlu(char1,char2,matlu1,matlu2,natom,option,toldiff,ierr,zero_or_one) 871 872 use defs_basis 873 use m_paw_dmft, only : paw_dmft_type 874 use m_crystal, only : crystal_t 875 use m_io_tools, only : flush_unit 876 implicit none 877 878 !Arguments ------------------------------------ 879 !type 880 integer,intent(in) :: natom,option 881 type(matlu_type), intent(in) :: matlu1(natom),matlu2(natom) 882 character(len=*), intent(in) :: char1,char2 883 real(dp),intent(in) :: toldiff 884 integer,intent(out), optional :: ierr 885 integer,intent(in), optional :: zero_or_one 886 887 !local variables------------------------------- 888 integer :: iatom,idiff,ispinor,ispinor1,isppol,m1,m,lpawu,nspinor,nsppol 889 real(dp) :: matludiff 890 character(len=500) :: message 891 ! ********************************************************************* 892 nsppol=matlu1(1)%nsppol 893 nspinor=matlu1(1)%nspinor 894 895 matludiff=zero 896 idiff=0 897 do iatom = 1 , natom 898 lpawu=matlu1(iatom)%lpawu 899 if(lpawu/=-1) then 900 do isppol = 1 , nsppol 901 do ispinor = 1 , nspinor 902 do ispinor1 = 1, nspinor 903 do m1 = 1 , 2*lpawu+1 904 do m = 1 , 2*lpawu+1 905 idiff=idiff+1 906 matludiff=matludiff+ & 907 & sqrt( real(matlu1(iatom)%mat(m1,m,isppol,ispinor,ispinor1) & 908 & - matlu2(iatom)%mat(m1,m,isppol,ispinor,ispinor1))**2 & 909 & +aimag(matlu1(iatom)%mat(m1,m,isppol,ispinor,ispinor1) & 910 & - matlu2(iatom)%mat(m1,m,isppol,ispinor,ispinor1))**2 ) 911 ! write(std_out,*) m,m1,matlu1(iatom)%mat(m1,m,isppol,ispinor,ispinor1),matlu2(iatom)%mat(m1,m,isppol,ispinor,ispinor1),matludiff 912 enddo 913 enddo 914 end do ! ispinor1 915 end do ! ispinor 916 enddo ! isppol 917 endif ! lpawu/=1 918 enddo ! natom 919 if(.not.present(zero_or_one)) matludiff=matludiff/float(idiff) 920 921 if(option==1.or.option==0) then 922 if( matludiff < toldiff ) then 923 write(message,'(5a,6x,3a,4x,e12.4,a,e12.4)') ch10,& 924 & ' ** Differences between ',trim(char1),' and ',ch10,trim(char2),' are small enough:',& 925 & ch10,matludiff,' is lower than',toldiff 926 call wrtout(std_out,message,'COLL') 927 if(present(ierr)) ierr=0 928 else 929 write(message,'(5a,3x,3a,3x,e12.4,a,e12.4)') ch10,& 930 & 'Differences between ',trim(char1),' and ',ch10,trim(char2),' is too large:',& 931 & ch10,matludiff,' is larger than',toldiff 932 ABI_WARNING(message) 933 ! write(message,'(8a,4x,e12.4,a,e12.4)') ch10," Matrix for ",trim(char1) 934 write(message,'(a,3x,a)') ch10,trim(char1) 935 call wrtout(std_out,message,'COLL') 936 call print_matlu(matlu1,natom,prtopt=1,opt_diag=-1) 937 write(message,'(a,3x,a)') ch10,trim(char2) 938 call wrtout(std_out,message,'COLL') 939 call print_matlu(matlu2,natom,prtopt=1,opt_diag=-1) 940 if (present(zero_or_one).and.(mod(matludiff,1.d0)< toldiff)) then 941 write(message,'(a,3x,a)') ch10," The norm is not identity for this k-point but& 942 & is compatible with a high symmetry point" 943 call wrtout(std_out,message,'COLL') 944 else if(present(zero_or_one)) then 945 write(message,'(a,3x,a)') ch10," The norm is not identity for this k-point but& 946 & might be compatible with a high symmetry point: it should be checked" 947 call wrtout(std_out,message,'COLL') 948 else 949 if(option==1) then 950 call flush_unit(std_out) 951 ABI_ERROR("option==1, aborting now!") 952 end if 953 end if 954 if(present(ierr)) ierr=-1 955 endif 956 endif 957 958 end subroutine diff_matlu
m_matlu/fac_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
fac_matlu
FUNCTION
shift matlu matrix
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
matlu(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: input quantity to rotate natom=number of atoms in cell. shift= shift of the diagonal part.
OUTPUT
matlu(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: shifted matrix
SIDE EFFECTS
NOTES
SOURCE
2666 subroutine fac_matlu(matlu,natom,fac) 2667 use defs_basis 2668 implicit none 2669 2670 !Arguments ------------------------------------ 2671 !scalars 2672 integer, intent(in) :: natom 2673 !arrays 2674 type(matlu_type),intent(inout) :: matlu(natom) 2675 complex(dpc),intent(in) :: fac 2676 !Local variables------------------------------- 2677 !scalars 2678 integer :: iatom,im,im1,ispinor,ispinor1,isppol 2679 integer :: lpawu 2680 ! character(len=500) :: message 2681 !arrays 2682 !************************************************************************ 2683 do iatom=1,natom 2684 lpawu=matlu(iatom)%lpawu 2685 if(lpawu.ne.-1) then 2686 do im=1,2*lpawu+1 2687 do im1=1,2*lpawu+1 2688 do isppol=1,matlu(1)%nsppol 2689 do ispinor=1,matlu(1)%nspinor 2690 do ispinor1=1,matlu(1)%nspinor 2691 matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)=& 2692 & matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)*fac 2693 enddo 2694 enddo 2695 enddo 2696 enddo 2697 enddo 2698 endif 2699 enddo 2700 2701 end subroutine fac_matlu
m_matlu/gather_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
gather_matlu
FUNCTION
Create new array from matlu
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
gloc(natom) <type(matlu_type)> = density matrix in the spin spin representation gatherloc(natom) <type(coeff2c_type)> = density matrix where spin and angular momentum are gathered in the same index natom=number of atoms in cell. option= 1 go from gloc to gathergloc option= -1 go from gathergloc to gloc prtopt= option to define level of printing
OUTPUT
gloc(natom) <type(matlu_type)> = density matrix in the spin spin representation gatherloc(natom) <type(coeff2c_type)> = density matrix where spin and angular momentum are gathered in the same index
SIDE EFFECTS
SOURCE
1271 subroutine gather_matlu(gloc,gathergloc,natom,option,prtopt) 1272 use defs_basis 1273 use defs_wvltypes 1274 use m_crystal, only : crystal_t 1275 implicit none 1276 1277 ! type matlus_type 1278 ! SEQUENCE 1279 ! complex(dpc), pointer :: mat(:,:) 1280 ! end type matlus_type 1281 1282 !Arguments ------------------------------------ 1283 !scalars 1284 integer,intent(in) :: natom,option,prtopt 1285 type(coeff2c_type), intent(inout) :: gathergloc(natom) 1286 type(matlu_type),intent(inout) :: gloc(natom) 1287 !Local variables------------------------------- 1288 !scalars 1289 integer :: iatom,im1,im2,ispinor,ispinor1,isppol,isppol1 1290 integer :: jc1,jc2,ml1,ml2,ndim,nspinor,nsppol,tndim 1291 character(len=500) :: message 1292 1293 ! DBG_ENTER("COLL") 1294 nsppol=gloc(1)%nsppol 1295 nspinor=gloc(1)%nspinor 1296 1297 do iatom=1,natom 1298 if(gloc(iatom)%lpawu.ne.-1) then 1299 !==------------------------------------- 1300 1301 ndim=2*gloc(iatom)%lpawu+1 1302 tndim=nsppol*nspinor*ndim 1303 1304 !== Put norm into array "gathergloc" 1305 jc1=0 1306 do isppol=1,nsppol 1307 do ispinor=1,nspinor 1308 do ml1=1,ndim 1309 jc1=jc1+1 1310 jc2=0 1311 do isppol1=1,nsppol 1312 do ispinor1=1,nspinor 1313 do ml2=1,ndim 1314 jc2=jc2+1 1315 if(option==1) then 1316 if(isppol==isppol1) then 1317 gathergloc(iatom)%value(jc1,jc2)=gloc(iatom)%mat(ml1,ml2,isppol,ispinor,ispinor1) 1318 endif 1319 else if(option==-1) then 1320 if(isppol==isppol1) then 1321 gloc(iatom)%mat(ml1,ml2,isppol,ispinor,ispinor1)=gathergloc(iatom)%value(jc1,jc2) 1322 endif 1323 endif 1324 enddo 1325 enddo ! ispinor1 1326 enddo ! isppol1 1327 enddo 1328 enddo !ispinor 1329 enddo ! isppol 1330 endif 1331 enddo ! iatom 1332 if(option==1.and.prtopt==3) then 1333 do iatom=1,natom 1334 if(gloc(iatom)%lpawu.ne.-1) then 1335 tndim=nsppol*nspinor*(2*gloc(iatom)%lpawu+1) 1336 write(message,'(2a,i5)') ch10,' (gathermatlu:) For atom', iatom 1337 call wrtout(std_out,message,'COLL') 1338 do im1=1,tndim 1339 write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')& 1340 & (gathergloc(iatom)%value(im1,im2),im2=1,tndim) 1341 call wrtout(std_out,message,'COLL') 1342 end do 1343 endif 1344 enddo ! iatom 1345 else if(option==-1.and.prtopt==3) then 1346 call print_matlu(gloc,natom,prtopt) 1347 endif 1348 1349 1350 1351 ! DBG_EXIT("COLL") 1352 1353 end subroutine gather_matlu
m_matlu/identity_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
identity_matlu
FUNCTION
Make the matlu the identity
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
matlu(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: input quantity to rotate natom=number of atoms in cell.
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
2802 subroutine identity_matlu(matlu,natom) 2803 use defs_basis 2804 implicit none 2805 2806 !Arguments ------------------------------------ 2807 !scalars 2808 integer, intent(in) :: natom 2809 !arrays 2810 type(matlu_type),intent(inout) :: matlu(natom) 2811 !Local variables------------------------------- 2812 !scalars 2813 integer :: iatom,im,ispinor,isppol 2814 integer :: ndim 2815 ! character(len=500) :: message 2816 !arrays 2817 !************************************************************************ 2818 ! not yet tested and used 2819 do iatom=1,natom 2820 if(matlu(iatom)%lpawu.ne.-1) then 2821 ndim=2*matlu(iatom)%lpawu+1 2822 do isppol=1,matlu(iatom)%nsppol 2823 do im=1,ndim 2824 do ispinor=1,matlu(iatom)%nspinor 2825 matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)= cone 2826 enddo ! ispinor 2827 enddo ! im 2828 enddo 2829 endif ! lpawu 2830 enddo ! iatom 2831 end subroutine identity_matlu
m_matlu/init_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
init_matlu
FUNCTION
Allocate variables used in type matlu_type.
INPUTS
natom = number of atoms nspinor = number of spinorial components nsppol = number of polarisation components lpawu_natom(natom) = value of lpawu for all atoms maltu <type(matlu_type)>= density matrix in the local orbital basis and related variables OUTPUTS maltu <type(matlu_type)>= density matrix in the local orbital basis and related variables
SOURCE
136 subroutine init_matlu(natom,nspinor,nsppol,lpawu_natom,matlu) 137 138 use defs_basis 139 use m_crystal, only : crystal_t 140 implicit none 141 142 !Arguments ------------------------------------ 143 !scalars 144 integer :: natom,nspinor,nsppol 145 !type 146 integer, intent(in) :: lpawu_natom(natom) 147 type(matlu_type), intent(inout) :: matlu(natom) 148 !Local variables ------------------------------------ 149 integer :: iatom,lpawu 150 151 !************************************************************************ 152 153 ! matlu%mband = mband 154 ! matlu%dmftbandf = dmftbandf 155 ! matlu%dmftbandi = dmftbandi 156 ! matlu%nkpt = nkpt 157 ! matlu%mbandc = 0 158 matlu%nsppol = nsppol 159 matlu%nspinor = nspinor 160 do iatom=1,natom 161 lpawu=lpawu_natom(iatom) 162 matlu(iatom)%lpawu=lpawu 163 if(lpawu.ne.-1) then 164 ABI_MALLOC(matlu(iatom)%mat,(2*lpawu+1,2*lpawu+1,nsppol,nspinor,nspinor)) 165 matlu(iatom)%mat=czero 166 else 167 ABI_MALLOC(matlu(iatom)%mat,(0,0,nsppol,nspinor,nspinor)) 168 endif 169 enddo 170 171 end subroutine init_matlu
m_matlu/inverse_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
inverse_matlu
FUNCTION
Inverse local quantity.
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
matlu(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: input quantity to inverse natom=number of atoms in cell. prtopt= option to define level of printing
OUTPUT
matlu(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: inverse of input matrix
SIDE EFFECTS
NOTES
SOURCE
798 subroutine inverse_matlu(matlu,natom,prtopt) 799 use defs_basis 800 use defs_wvltypes 801 use m_crystal, only : crystal_t 802 implicit none 803 804 !Arguments ------------------------------------ 805 !scalars 806 integer, intent(in) :: natom 807 integer, intent(in) :: prtopt 808 !arrays 809 type(matlu_type),intent(inout) :: matlu(natom) 810 !Local variables------------------------------- 811 integer :: iatom,tndim 812 integer :: nsppol,nspinor 813 !scalars 814 type(coeff2c_type),allocatable :: gathermatlu(:) 815 !************************************************************************ 816 817 818 nspinor=matlu(1)%nspinor 819 nsppol=matlu(1)%nsppol 820 if(prtopt>0) then 821 endif 822 ABI_MALLOC(gathermatlu,(natom)) 823 do iatom=1,natom 824 if(matlu(iatom)%lpawu.ne.-1) then 825 tndim=nsppol*nspinor*(2*matlu(iatom)%lpawu+1) 826 ABI_MALLOC(gathermatlu(iatom)%value,(tndim,tndim)) 827 gathermatlu(iatom)%value=czero 828 endif 829 enddo 830 831 call gather_matlu(matlu,gathermatlu,natom,option=1,prtopt=1) 832 do iatom=1,natom 833 if(matlu(iatom)%lpawu.ne.-1) then 834 tndim=nsppol*nspinor*(2*matlu(iatom)%lpawu+1) 835 !call matcginv_dpc(gathermatlu(iatom)%value,tndim,tndim) 836 call xginv(gathermatlu(iatom)%value,tndim) 837 endif 838 enddo 839 call gather_matlu(matlu,gathermatlu,natom,option=-1,prtopt=1) 840 841 do iatom=1,natom 842 if(matlu(iatom)%lpawu.ne.-1) then 843 ABI_FREE(gathermatlu(iatom)%value) 844 endif 845 enddo 846 ABI_FREE(gathermatlu) 847 end subroutine inverse_matlu
m_matlu/ln_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
ln_matlu
FUNCTION
Compute the logarithm of matlu (only if diagonal for the moment)
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
matlu1(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: input quantity
SIDE EFFECTS
NOTES
SOURCE
2463 subroutine ln_matlu(matlu1,natom) 2464 use defs_basis 2465 use defs_wvltypes 2466 implicit none 2467 2468 !Arguments ------------------------------------ 2469 !scalars 2470 integer, intent(in) :: natom 2471 !arrays 2472 type(matlu_type), intent(inout) :: matlu1(natom) 2473 !Local variables------------------------------- 2474 !scalars 2475 integer :: iatom,im,ispinor,isppol 2476 integer :: lpawu 2477 character(len=500) :: message 2478 !arrays 2479 !************************************************************************ 2480 !call checkdiag_matlu(matlu1,natom,tol8) 2481 do iatom=1,natom 2482 lpawu=matlu1(iatom)%lpawu 2483 if(lpawu.ne.-1) then 2484 do isppol=1,matlu1(1)%nsppol 2485 do ispinor=1,matlu1(1)%nspinor 2486 do im=1,2*lpawu+1 2487 if( real(matlu1(iatom)%mat(im,im,isppol,ispinor,ispinor))<zero) then 2488 write(message,'(2a,2es13.5,a)') ch10," ln_matlu: PROBLEM " & 2489 & , matlu1(iatom)%mat(im,im,isppol,ispinor,ispinor) 2490 ABI_ERROR(message) 2491 endif 2492 matlu1(iatom)%mat(im,im,isppol,ispinor,ispinor)= & 2493 & log(matlu1(iatom)%mat(im,im,isppol,ispinor,ispinor)) 2494 enddo ! im 2495 enddo ! ispinor 2496 enddo ! isppol 2497 endif ! lpawu 2498 enddo ! iatom 2499 2500 end subroutine ln_matlu
m_matlu/matlu_type [ Types ]
NAME
matlu_type
FUNCTION
This structured datatype contains an matrix for the correlated subspace
SOURCE
79 type, public :: matlu_type ! for each atom 80 81 integer :: lpawu 82 ! value of the angular momentum for correlated electrons 83 84 ! integer :: natom 85 ! number of atoms (given for each atom, not useful..could be changed) 86 ! 87 ! integer :: mband 88 ! ! Number of bands 89 ! 90 ! integer :: mbandc 91 ! ! Total number of bands in the Kohn-Sham Basis for PAW+DMFT 92 ! 93 ! integer :: natpawu ! Number of correlated atoms 94 ! 95 ! integer :: nkpt 96 ! ! Number of k-point in the IBZ. 97 character(len=12) :: whichmatlu 98 ! describe the type of local matrix computed (greenDFT, etc..) 99 ! 100 integer :: nspinor 101 ! Number of spinorial component 102 ! 103 integer :: nsppol 104 ! Number of polarization 105 106 complex(dpc), allocatable :: mat(:,:,:,:,:) 107 ! local quantity 108 109 end type matlu_type 110 111 !---------------------------------------------------------------------- 112 113 114 CONTAINS !========================================================================================
m_matlu/print_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
print_matlu
FUNCTION
INPUTS
maltu <type(matlu_type)>= density matrix in the local orbital basis and related variables natom= number of atoms prtopt= option for printing opt_diag= 0 print non diagonal matrix (real or complex according to nspinor) -1 print non diagonal complex matrix >=1 print diagonal matrix (real or complex according to nspinor) opt_ab_out= 0 print matrix on std_out /=0 print matrix on ab_out
OUTPUT
SOURCE
363 subroutine print_matlu(matlu,natom,prtopt,opt_diag,opt_ab_out,opt_exp,argout,compl) 364 365 use defs_basis 366 use m_crystal, only : crystal_t 367 implicit none 368 369 !Arguments ------------------------------------ 370 !type 371 integer, intent(in):: natom 372 type(matlu_type),intent(in) :: matlu(natom) 373 integer, intent(in) :: prtopt 374 integer, optional, intent(in) :: opt_diag,opt_ab_out,opt_exp,argout,compl 375 !Local variables------------------------------- 376 integer :: iatom,ispinor,ispinor1,isppol,m1,m,lpawu,nspinor,nsppol,optdiag,optab_out,arg_out 377 character(len=500) :: message 378 character(len=4) :: mode_paral 379 type(matlu_type), allocatable :: matlu_tmp(:) 380 character(len=9),parameter :: dspinm(2,2)= RESHAPE((/"n ","mx ","my ","mz "/),(/2,2/)) 381 logical :: testcmplx 382 real(dp) :: noccspin 383 ! ********************************************************************* 384 mode_paral='COLL' 385 if(present(opt_diag)) then 386 optdiag=opt_diag 387 else 388 optdiag=0 389 endif 390 if(present(opt_ab_out)) then 391 optab_out=opt_ab_out 392 else 393 optab_out=0 394 endif 395 if(optab_out==0) then 396 arg_out=std_out 397 else 398 arg_out=ab_out 399 endif 400 if(present(argout)) then 401 arg_out=argout 402 mode_paral='PERS' 403 endif 404 nspinor=matlu(1)%nspinor 405 nsppol=matlu(1)%nsppol 406 testcmplx=(nspinor==2) 407 if(present(compl)) testcmplx=(nspinor==2).or.(compl==1) 408 409 do iatom = 1 , natom 410 lpawu=matlu(iatom)%lpawu 411 if(lpawu/=-1) then 412 write(message,'(2a,i4)') ch10,' -------> For Correlated Atom', iatom 413 call wrtout(arg_out, message,mode_paral) 414 do isppol = 1 , nsppol 415 if(present(opt_ab_out).and.nsppol==2) then 416 noccspin=zero 417 do m1=1, 2*lpawu +1 418 noccspin=noccspin+REAL(matlu(iatom)%mat(m1,m1,isppol,1,1)) 419 enddo 420 !write(message,fmt='(7x,a,i3,a,f10.5)') ". Occ. for lpawu and for spin",isppol," =",noccspin 421 !call wrtout(arg_out, message,mode_paral) 422 endif 423 enddo 424 425 do isppol = 1 , nsppol 426 if(nspinor == 1) then 427 write(message,'(a,10x,a,i3,i3)') ch10,'-- polarization spin component',isppol 428 call wrtout(arg_out, message,mode_paral) 429 endif 430 do ispinor = 1 , nspinor 431 do ispinor1 = 1, nspinor 432 if(nspinor == 2) then 433 write(message,'(a,10x,a,i3,i3)') ch10,'-- spin components',ispinor,ispinor1 434 call wrtout(arg_out, message,mode_paral) 435 endif 436 if(optdiag<=0) then 437 do m1=1, 2*lpawu+1 438 if(optdiag==0) then 439 if(nspinor==1.and.abs(prtopt)>0) then 440 if(present(opt_exp)) then 441 write(message,'(5x,20e24.14)') (REAL(matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1) 442 ! call wrtout(arg_out, message,mode_paral) 443 ! write(message,'(5x,20e20.14)') (REAL(sqrt(matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1))),m=1,2*lpawu+1) 444 ! call wrtout(arg_out, message,mode_paral) 445 ! write(message,'(5x,20e20.14)') (REAL(1.d0/sqrt(matlu(iatom)%mat(m,m,isppol,ispinor,ispinor1))),m=1,2*lpawu+1) 446 else 447 write(message,'(5x,20f10.5)') (REAL(matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1) 448 endif 449 endif 450 if(testcmplx.and.abs(prtopt)>0) then 451 if(present(opt_exp)) then 452 if(opt_exp==2) then 453 write(message,'(5x,14(2e18.10,1x))') ((matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1) 454 else 455 write(message,'(5x,14(2e14.4,2x))') ((matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1) 456 endif 457 else 458 write(message,'(5x,14(2f9.5,2x))')((matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1) 459 endif 460 !& write(message,'(5x,14(2f15.11,2x))')((matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1) 461 endif 462 else if(optdiag==-1) then 463 write(message,'(5x,14(2f10.5,2x))')((matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1) 464 endif 465 call wrtout(arg_out, message,mode_paral) 466 enddo 467 elseif (optdiag>=1) then 468 if(nspinor==1.and.abs(prtopt)>0) & 469 & write(message,'(5x,20f10.5)') (REAL(matlu(iatom)%mat(m,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1) 470 if(testcmplx.and.abs(prtopt)>0) & 471 & write(message,'(5x,14(2f9.5,2x))')((matlu(iatom)%mat(m,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1) 472 ! write(std_out,'(5x,14(2f9.5,2x))')((matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1) 473 call wrtout(arg_out, message,mode_paral) 474 endif 475 end do ! ispinor1 476 end do ! ispinor 477 if(nspinor==2.and.prtopt>=5) then 478 ABI_MALLOC(matlu_tmp,(0:natom)) 479 call init_matlu(natom,nspinor,nsppol,matlu(1:natom)%lpawu,matlu_tmp) 480 matlu_tmp(iatom)%mat(m1,m,isppol,1,1)= matlu(iatom)%mat(m1,m,isppol,1,1)+matlu(iatom)%mat(m1,m,isppol,2,2) 481 matlu_tmp(iatom)%mat(m1,m,isppol,2,2)= matlu(iatom)%mat(m1,m,isppol,1,1)+matlu(iatom)%mat(m1,m,isppol,2,2) 482 matlu_tmp(iatom)%mat(m1,m,isppol,1,2)= matlu(iatom)%mat(m1,m,isppol,1,2)+matlu(iatom)%mat(m1,m,isppol,2,1) 483 matlu_tmp(iatom)%mat(m1,m,isppol,2,1)= & 484 & (matlu(iatom)%mat(m1,m,isppol,1,2)+matlu(iatom)%mat(m1,m,isppol,2,1))*cmplx(0_dp,1_dp,kind=dp) 485 do ispinor = 1 , nspinor 486 do ispinor1 = 1, nspinor 487 write(message,'(a,10x,a,a)') ch10,'-- spin components',dspinm(ispinor,ispinor1) 488 call wrtout(arg_out, message,mode_paral) 489 write(message,'(5x,14(2f9.5,2x))')((matlu_tmp(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1) 490 call wrtout(arg_out, message,mode_paral) 491 end do ! ispinor1 492 end do ! ispinor 493 call destroy_matlu(matlu_tmp,natom) 494 ABI_FREE(matlu_tmp) 495 endif 496 enddo ! isppol 497 ! if(nsppol==1.and.nspinor==1) then 498 ! write(message,'(a,10x,a,i3,a)') ch10,'-- polarization spin component',isppol+1,' is identical' 499 ! call wrtout(arg_out, message,mode_paral) 500 ! endif 501 endif ! lpawu/=1 502 enddo ! natom 503 504 505 end subroutine print_matlu
m_matlu/printplot_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
printplot_matlu
FUNCTION
shift matlu matrix
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
matlu(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: input quantity to rotate natom=number of atoms in cell. shift= shift of the diagonal part.
OUTPUT
matlu(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: shifted matrix
SIDE EFFECTS
NOTES
SOURCE
2730 subroutine printplot_matlu(matlu,natom,freq,char1,units,imre) 2731 use defs_basis 2732 use m_fstrings, only : int2char4 2733 implicit none 2734 2735 !Arguments ------------------------------------ 2736 !scalars 2737 integer, intent(in) :: natom,units 2738 integer, optional, intent(in) :: imre 2739 real(dp), intent(in) :: freq 2740 !arrays 2741 type(matlu_type),intent(inout) :: matlu(natom) 2742 character(len=*), intent(in) :: char1 2743 !Local variables------------------------------- 2744 !scalars 2745 integer :: iatom,im,im1,ispinor,ispinor1,isppol 2746 integer :: lpawu,unitnb 2747 character(len=4) :: tag_at 2748 character(len=fnlen) :: tmpfil,tmpfilre,tmpfilim 2749 ! character(len=500) :: message 2750 !arrays 2751 !************************************************************************ 2752 ! not yet tested and used 2753 do iatom=1,natom 2754 lpawu=matlu(iatom)%lpawu 2755 if(lpawu.ne.-1) then 2756 unitnb=units+iatom 2757 call int2char4(iatom,tag_at) 2758 if(present(imre)) then 2759 tmpfilre = trim(char1)//tag_at//"re" 2760 tmpfilim = trim(char1)//tag_at//"im" 2761 open (unit=unitnb+10,file=trim(tmpfilre),status='unknown',form='formatted') 2762 open (unit=unitnb+20,file=trim(tmpfilim),status='unknown',form='formatted') 2763 write(unitnb+10,'(400e26.16)') freq,(((((real(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)),& 2764 im=1,2*lpawu+1),ispinor=1,matlu(1)%nspinor),im1=1,2*lpawu+1),ispinor1=1,matlu(1)%nspinor),isppol=1,matlu(1)%nsppol) 2765 write(unitnb+20,'(400e26.16)') freq,(((((aimag(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)),& 2766 im=1,2*lpawu+1),ispinor=1,matlu(1)%nspinor),im1=1,2*lpawu+1),ispinor1=1,matlu(1)%nspinor),isppol=1,matlu(1)%nsppol) 2767 else 2768 tmpfil = trim(char1)//tag_at 2769 open (unit=unitnb,file=trim(tmpfil),status='unknown',form='formatted') 2770 write(unitnb,'(400e26.16)') freq,(((((matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1),& 2771 im=1,2*lpawu+1),ispinor=1,matlu(1)%nspinor),im1=1,2*lpawu+1),ispinor1=1,matlu(1)%nspinor),isppol=1,matlu(1)%nsppol) 2772 endif 2773 endif 2774 enddo 2775 end subroutine printplot_matlu
m_matlu/prod_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
prod_matlu
FUNCTION
Do the matrix product of two matlus
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
matlu1(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: input quantity matlu2(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: input quantity
OUTPUT
matlu3(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: output quantity
SIDE EFFECTS
NOTES
SOURCE
2338 subroutine prod_matlu(matlu1,matlu2,matlu3,natom) 2339 use defs_basis 2340 use defs_wvltypes 2341 use m_crystal, only : crystal_t 2342 implicit none 2343 2344 !Arguments ------------------------------------ 2345 !scalars 2346 integer, intent(in) :: natom 2347 !arrays 2348 type(matlu_type), intent(in) :: matlu1(natom),matlu2(natom) 2349 type(matlu_type), intent(inout) :: matlu3(natom) 2350 !Local variables------------------------------- 2351 !scalars 2352 integer :: iatom,im1,im2,im3,ispinor1,ispinor2,ispinor3,isppol 2353 integer :: lpawu 2354 !arrays 2355 !************************************************************************ 2356 call zero_matlu(matlu3,natom) 2357 do iatom=1,natom 2358 lpawu=matlu1(iatom)%lpawu 2359 if(lpawu.ne.-1) then 2360 do isppol=1,matlu1(1)%nsppol 2361 do ispinor1=1,matlu1(1)%nspinor 2362 do ispinor2=1,matlu1(1)%nspinor 2363 do ispinor3=1,matlu1(1)%nspinor 2364 do im1=1,2*lpawu+1 2365 do im2=1,2*lpawu+1 2366 do im3=1,2*lpawu+1 2367 matlu3(iatom)%mat(im1,im2,isppol,ispinor1,ispinor2)= & 2368 & matlu3(iatom)%mat(im1,im2,isppol,ispinor1,ispinor2)+ & 2369 & matlu1(iatom)%mat(im1,im3,isppol,ispinor1,ispinor3)*& 2370 & matlu2(iatom)%mat(im3,im2,isppol,ispinor3,ispinor2) 2371 enddo ! im3 2372 enddo ! im2 2373 enddo ! im1 2374 enddo ! ispinor3 2375 enddo ! ispinor2 2376 enddo ! ispinor1 2377 enddo ! isppol 2378 endif ! lpawu 2379 enddo ! iatom 2380 2381 end subroutine prod_matlu
m_matlu/rotate_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
rotate_matlu
FUNCTION
Rotate matlu matrix
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
matlu(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: input quantity to rotate rot_mat(natom) <type(coeff2c_type)> = Rotation matrix (usually from diag_matlu) natom=number of atoms in cell. prtopt= option to define level of printing
OUTPUT
matlu(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: Rotated matrix
SIDE EFFECTS
NOTES
SOURCE
1868 subroutine rotate_matlu(matlu,rot_mat,natom,prtopt,inverse) 1869 use defs_basis 1870 use defs_wvltypes 1871 use m_crystal, only : crystal_t 1872 implicit none 1873 1874 !Arguments ------------------------------------ 1875 !scalars 1876 integer, intent(in) :: natom 1877 integer, intent(in) :: prtopt 1878 integer, intent(in) :: inverse 1879 !arrays 1880 type(matlu_type),intent(inout) :: matlu(natom) 1881 type(coeff2c_type),optional,intent(inout) :: rot_mat(natom,matlu(1)%nsppol) 1882 !Local variables------------------------------- 1883 !scalars 1884 integer :: iatom,im1,im2,isppol 1885 integer :: nsppol,nspinor,tndim 1886 !arrays 1887 type(coeff2c_type),allocatable :: gathermatlu(:) 1888 ! type(coeff2c_type),allocatable :: rot_mat_orig(:,:) 1889 type(coeff2c_type),allocatable :: rot_mat_orig(:) 1890 complex(dpc),allocatable :: temp_mat(:,:) 1891 !************************************************************************ 1892 if(prtopt==1) then 1893 endif 1894 nsppol=matlu(1)%nsppol 1895 nspinor=matlu(1)%nspinor 1896 ! ABI_MALLOC(rot_mat_orig,(natom,matlu(1)%nsppol)) 1897 1898 do isppol=1,nsppol 1899 1900 ! =========================== 1901 ! Define gathermatlu and rot_mat_orig and allocate 1902 ! =========================== 1903 ABI_MALLOC(rot_mat_orig,(natom)) 1904 ABI_MALLOC(gathermatlu,(natom)) 1905 do iatom=1,natom 1906 if(matlu(iatom)%lpawu.ne.-1) then 1907 tndim=nspinor*(2*matlu(iatom)%lpawu+1) 1908 ABI_MALLOC(gathermatlu(iatom)%value,(tndim,tndim)) 1909 gathermatlu(iatom)%value=czero 1910 ! ABI_MALLOC(rot_mat_orig(iatom,isppol)%value,(tndim,tndim)) 1911 ! rot_mat_orig(iatom,isppol)%value(:,:)=rot_mat(iatom,isppol)%value(:,:) 1912 ABI_MALLOC(rot_mat_orig(iatom)%value,(tndim,tndim)) 1913 rot_mat_orig(iatom)%value(:,:)=rot_mat(iatom,isppol)%value(:,:) 1914 endif 1915 enddo 1916 if(nsppol==1.and.nspinor==2) then 1917 call gather_matlu(matlu,gathermatlu,natom,option=1,prtopt=1) 1918 else if((nsppol==2.or.nsppol==1).and.nspinor==1) then 1919 do iatom=1,natom 1920 if(matlu(iatom)%lpawu.ne.-1) then 1921 tndim=nspinor*(2*matlu(iatom)%lpawu+1) 1922 do im1=1,tndim 1923 do im2=1,tndim 1924 gathermatlu(iatom)%value(im1,im2)=matlu(iatom)%mat(im1,im2,isppol,1,1) 1925 enddo 1926 enddo 1927 endif 1928 enddo 1929 endif 1930 ! write(std_out,*) "gathermatlu in rotate matlu" 1931 ! do im1=1,tndim 1932 ! write(message,'(12(1x,18(1x,"(",e17.10,",",e17.10,")")))')& 1933 ! (gathermatlu(1)%value(im1,im2),im2=1,tndim) 1934 ! call wrtout(std_out,message,'COLL') 1935 ! end do 1936 1937 ! =========================== 1938 ! If necessary, invert rot_mat 1939 ! =========================== 1940 if(inverse==1) then 1941 do iatom=1,natom 1942 if(matlu(iatom)%lpawu.ne.-1) then 1943 tndim=nspinor*(2*matlu(iatom)%lpawu+1) 1944 do im1=1,tndim 1945 do im2=1,tndim 1946 ! rot_mat(iatom,isppol)%value(im1,im2)=conjg(rot_mat_orig(iatom,isppol)%value(im2,im1)) 1947 rot_mat(iatom,isppol)%value(im1,im2)=conjg(rot_mat_orig(iatom)%value(im2,im1)) 1948 enddo 1949 enddo 1950 endif ! lpawu 1951 enddo ! iatom 1952 endif 1953 ! write(std_out,*) "rot_mat_orig " 1954 ! do im1=1,tndim 1955 ! write(message,'(12(1x,18(1x,"(",e18.10,",",e18.10,")")))')& 1956 ! & (rot_mat_orig(1)%value(im1,im2),im2=1,tndim) 1957 ! call wrtout(std_out,message,'COLL') 1958 ! end do 1959 ! write(std_out,*) "rot_mat " 1960 ! do im1=1,tndim 1961 ! write(message,'(12(1x,18(1x,"(",e18.10,",",e18.10,")")))')& 1962 ! & (rot_mat(1,1)%value(im1,im2),im2=1,tndim) 1963 ! call wrtout(std_out,message,'COLL') 1964 ! end do 1965 1966 ! =========================== 1967 ! Rotate 1968 ! =========================== 1969 ABI_MALLOC(temp_mat,(tndim,tndim)) 1970 do iatom=1,natom 1971 if(matlu(iatom)%lpawu.ne.-1) then 1972 tndim=nspinor*(2*matlu(iatom)%lpawu+1) 1973 temp_mat(:,:)=czero 1974 ! input matrix: gathermatlu 1975 ! rotation matrix: rot_mat 1976 ! intermediate matrix: temp_mat 1977 ! result matrix: gathermatlu 1978 ! temp_mat = gathermatlu * conjg(rot_mat) 1979 call zgemm('n','c',tndim,tndim,tndim,cone,gathermatlu(iatom)%value ,tndim,& 1980 & rot_mat(iatom,isppol)%value,tndim,czero,temp_mat ,tndim) 1981 ! gathermatlu = rot_mat * temp_mat = rot_mat * gathermatlu * conjg(rot_mat) 1982 call zgemm('n','n',tndim,tndim,tndim,cone,rot_mat(iatom,isppol)%value,tndim,& 1983 & temp_mat ,tndim,czero,gathermatlu(iatom)%value,tndim) 1984 endif ! lpawu 1985 enddo ! iatom 1986 !do iatom=1,natom 1987 ! if(matlu(iatom)%lpawu.ne.-1) then 1988 ! write(std_out,*) "temp_mat in rotate_matlu 2" 1989 ! do im1=1,tndim 1990 ! write(message,'(12(1x,18(1x,"(",f17.10,",",f17.10,")")))')& 1991 !& (temp_mat(im1,im2),im2=1,tndim) 1992 ! call wrtout(std_out,message,'COLL') 1993 ! end do 1994 ! endif 1995 !enddo 1996 !do iatom=1,natom 1997 ! if(matlu(iatom)%lpawu.ne.-1) then 1998 ! write(std_out,*) "gathermatlu in rotate_matlu 2" 1999 ! do im1=1,tndim 2000 ! write(message,'(12(1x,18(1x,"(",f17.10,",",f17.10,")")))')& 2001 !& (gathermatlu(iatom)%value(im1,im2),im2=1,tndim) 2002 ! call wrtout(std_out,message,'COLL') 2003 ! end do 2004 ! endif 2005 !enddo 2006 ABI_FREE(temp_mat) 2007 !ABI_ERROR("Aborting now") 2008 2009 ! Choose inverse rotation: reconstruct correct rot_mat from rot_mat_orig 2010 ! ======================================================================== 2011 if(inverse==1) then 2012 do iatom=1,natom 2013 if(matlu(iatom)%lpawu.ne.-1) then 2014 tndim=nspinor*(2*matlu(iatom)%lpawu+1) 2015 do im1=1,tndim 2016 do im2=1,tndim 2017 ! rot_mat(iatom,isppol)%value(im1,im2)=rot_mat_orig(iatom,isppol)%value(im1,im2) 2018 rot_mat(iatom,isppol)%value(im1,im2)=rot_mat_orig(iatom)%value(im1,im2) 2019 enddo 2020 enddo 2021 endif ! lpawu 2022 enddo ! iatom 2023 endif 2024 2025 ! =========================== 2026 ! Put data into matlu(iatom) 2027 ! =========================== 2028 if(nsppol==1.and.nspinor==2) then 2029 call gather_matlu(matlu,gathermatlu,natom,option=-1,prtopt=1) 2030 else if((nsppol==2.or.nsppol==1).and.nspinor==1) then 2031 do iatom=1,natom 2032 if(matlu(iatom)%lpawu.ne.-1) then 2033 tndim=nspinor*(2*matlu(iatom)%lpawu+1) 2034 do im1=1,tndim 2035 do im2=1,tndim 2036 matlu(iatom)%mat(im1,im2,isppol,1,1)= gathermatlu(iatom)%value(im1,im2) 2037 enddo 2038 enddo 2039 endif 2040 enddo 2041 endif ! test nsppol/nspinor 2042 ! =========================== 2043 ! Deallocations 2044 ! =========================== 2045 do iatom=1,natom 2046 if(matlu(iatom)%lpawu.ne.-1) then 2047 ABI_FREE(gathermatlu(iatom)%value) 2048 ! ABI_FREE(rot_mat_orig(iatom,isppol)%value) 2049 ABI_FREE(rot_mat_orig(iatom)%value) 2050 endif 2051 enddo 2052 ABI_FREE(gathermatlu) 2053 ABI_FREE(rot_mat_orig) 2054 enddo ! isppol 2055 2056 end subroutine rotate_matlu
m_matlu/shift_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
shift_matlu
FUNCTION
shift matlu matrix
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
matlu(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: input quantity to rotate natom=number of atoms in cell. shift= shift of the diagonal part.
OUTPUT
matlu(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: shifted matrix
SIDE EFFECTS
NOTES
SOURCE
2085 subroutine shift_matlu(matlu,natom,shift,signe) 2086 use defs_basis 2087 use defs_wvltypes 2088 use m_crystal, only : crystal_t 2089 implicit none 2090 2091 !Arguments ------------------------------------ 2092 !scalars 2093 integer, intent(in) :: natom 2094 !arrays 2095 type(matlu_type),intent(inout) :: matlu(natom) 2096 complex(dpc),intent(in) :: shift(natom) 2097 integer, optional,intent(in) :: signe 2098 !Local variables------------------------------- 2099 !scalars 2100 integer :: iatom,im,ispinor,isppol 2101 integer :: lpawu,signe_used 2102 ! character(len=500) :: message 2103 !arrays 2104 !************************************************************************ 2105 signe_used=1 2106 if(present(signe)) then 2107 if(signe==-1) signe_used=-1 2108 endif 2109 do iatom=1,natom 2110 lpawu=matlu(iatom)%lpawu 2111 if(lpawu.ne.-1) then 2112 do im=1,2*lpawu+1 2113 do isppol=1,matlu(1)%nsppol 2114 do ispinor=1,matlu(1)%nspinor 2115 matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)=& 2116 & matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)-signe_used*shift(iatom) 2117 enddo 2118 enddo 2119 enddo 2120 endif 2121 enddo 2122 2123 end subroutine shift_matlu
m_matlu/slm2ylm_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
slm2ylm_matlu
FUNCTION
Transform mat from Slm to Ylm basis or vice versa
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
matlu1(natom)%(nsppol,nspinor,nspinor,ndim,ndim) :: input quantity natom :: number of atoms option=1 go from Slm to Ylm basis option=2 go from Ylm to Slm basis
SIDE EFFECTS
NOTES
SOURCE
2526 subroutine slm2ylm_matlu(matlu,natom,option,optprt) 2527 use defs_basis 2528 use defs_wvltypes 2529 implicit none 2530 2531 !Arguments ------------------------------------ 2532 !scalars 2533 integer, intent(in) :: natom,option,optprt 2534 !arrays 2535 type(matlu_type), intent(inout) :: matlu(natom) 2536 !Local variables------------------------------- 2537 !scalars 2538 integer :: iatom,im,ispinor,isppol,ispinor2 2539 integer :: lpawu,ll,mm,jm,ii,jj,im1,im2 2540 character(len=500) :: message 2541 real(dp) :: onem 2542 complex(dpc),allocatable :: slm2ylm(:,:) 2543 complex(dpc),allocatable :: mat_inp_c(:,:) 2544 complex(dpc),allocatable :: mat_out_c(:,:) 2545 complex(dpc) :: tmp2 2546 real(dp),parameter :: invsqrt2=one/sqrt2 2547 !arrays 2548 !************************************************************************ 2549 2550 do iatom=1,natom 2551 lpawu=matlu(iatom)%lpawu 2552 if(lpawu.ne.-1) then 2553 ll=lpawu 2554 ABI_MALLOC(slm2ylm,(2*ll+1,2*ll+1)) 2555 slm2ylm=czero 2556 do im=1,2*ll+1 2557 mm=im-ll-1;jm=-mm+ll+1 2558 onem=dble((-1)**mm) 2559 if (mm> 0) then 2560 slm2ylm(im,im)= cmplx(onem*invsqrt2,zero,kind=dp) 2561 slm2ylm(jm,im)= cmplx(invsqrt2, zero,kind=dp) 2562 end if 2563 if (mm==0) then 2564 slm2ylm(im,im)=cone 2565 end if 2566 if (mm< 0) then 2567 slm2ylm(im,im)= cmplx(zero, invsqrt2,kind=dp) 2568 slm2ylm(jm,im)=-cmplx(zero,onem*invsqrt2,kind=dp) 2569 end if 2570 end do 2571 if(optprt>2) then 2572 write(message,'(2a)') ch10,"SLM2YLM matrix" 2573 call wrtout(std_out,message,'COLL') 2574 do im1=1,ll*2+1 2575 write(message,'(12(1x,9(1x,"(",f9.5,",",f9.5,")")))')& 2576 & (slm2ylm(im1,im2),im2=1,ll*2+1) 2577 call wrtout(std_out,message,'COLL') 2578 end do 2579 endif 2580 do isppol=1,matlu(1)%nsppol 2581 do ispinor=1,matlu(1)%nspinor 2582 do ispinor2=1,matlu(1)%nspinor 2583 ABI_MALLOC(mat_out_c,(2*ll+1,2*ll+1)) 2584 ABI_MALLOC(mat_inp_c,(2*ll+1,2*ll+1)) 2585 mat_inp_c(:,:) = matlu(iatom)%mat(:,:,isppol,ispinor,ispinor2) 2586 mat_out_c=czero 2587 2588 if(optprt>2) then 2589 write(message,'(2a, i2, a, i2, a, i2)') ch10,"SLM input matrix, isppol=", isppol, ", ispinor=", ispinor,& 2590 & ", ispinor2=", ispinor2 2591 call wrtout(std_out,message,'COLL') 2592 do im1=1,ll*2+1 2593 write(message,'(12(1x,9(1x,"(",f9.5,",",f9.5,")")))')& 2594 & (mat_inp_c(im1,im2),im2=1,ll*2+1) 2595 call wrtout(std_out,message,'COLL') 2596 end do 2597 endif 2598 2599 do jm=1,2*ll+1 2600 do im=1,2*ll+1 2601 tmp2=czero 2602 do ii=1,2*ll+1 2603 do jj=1,2*ll+1 2604 if(option==1) then 2605 tmp2=tmp2+mat_inp_c(ii,jj)*(slm2ylm(im,ii))*CONJG(slm2ylm(jm,jj)) 2606 else if(option==2) then 2607 tmp2=tmp2+mat_inp_c(ii,jj)*CONJG(slm2ylm(ii,im))*(slm2ylm(jj,jm)) 2608 end if 2609 end do 2610 end do 2611 mat_out_c(im,jm)=tmp2 2612 end do 2613 end do 2614 2615 if(optprt>2) then 2616 write(message,'(2a, i2, a, i2, a, i2)') ch10,"YLM output matrix, isppol=", isppol, ", ispinor=", ispinor,& 2617 & ", ispinor2=", ispinor2 2618 call wrtout(std_out,message,'COLL') 2619 do im1=1,ll*2+1 2620 write(message,'(12(1x,9(1x,"(",f9.5,",",f9.5,")")))')& 2621 & (mat_out_c(im1,im2),im2=1,ll*2+1) 2622 call wrtout(std_out,message,'COLL') 2623 end do 2624 endif 2625 2626 matlu(iatom)%mat(:,:,isppol,ispinor,ispinor2)=mat_out_c(:,:) 2627 ABI_FREE(mat_out_c) 2628 ABI_FREE(mat_inp_c) 2629 enddo ! im 2630 enddo ! ispinor 2631 enddo ! isppol 2632 ABI_FREE(slm2ylm) 2633 endif ! lpawu 2634 enddo ! iatom 2635 2636 2637 end subroutine slm2ylm_matlu
m_matlu/sym_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
sym_matlu
FUNCTION
Symetrise local quantity.
COPYRIGHT
Copyright (C) 2005-2022 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
cryst_struc <type(crystal_t)>=crystal structure data gloc(natom) <type(matlu_type)>= density matrix in the local orbital basis and related variables pawang <type(pawang)>=paw angular mesh and related data paw_dmft <type(paw_dmft_type)>= paw+dmft related data
OUTPUT
gloc(natom) <type(matlu_type)>= density matrix symetrized in the local orbital basis and related variables
SIDE EFFECTS
NOTES
SOURCE
535 subroutine sym_matlu(cryst_struc,gloc,pawang,paw_dmft) 536 537 use defs_basis 538 ! use defs_wvltypes 539 use m_pawang, only : pawang_type 540 use m_crystal, only : crystal_t 541 use m_paw_dmft, only: paw_dmft_type 542 543 implicit none 544 545 !Arguments ------------------------------------ 546 !scalars 547 type(crystal_t),intent(in) :: cryst_struc 548 type(pawang_type),intent(in) :: pawang 549 type(paw_dmft_type), intent(in) :: paw_dmft 550 !arrays 551 type(matlu_type),intent(inout) :: gloc(cryst_struc%natom) 552 !scalars 553 !Local variables------------------------------- 554 !scalars 555 integer :: at_indx,iatom,irot,ispinor,ispinor1,isppol,lpawu,m1,m2,m3,m4,mu 556 integer :: natom,ndim,nsppol,nspinor,nu,t2g,m1s,m2s,m3s,m4s,lpawu_zarot,x2my2d 557 complex(dpc) :: sumrho,summag(3),rotmag(3),ci 558 real(dp) :: zarot2 559 !arrays 560 ! complex(dpc),allocatable :: glocnm(:,:,:,:,:) 561 type(matlu_type),allocatable :: glocnm(:) 562 ! complex(dpc),allocatable :: glocnms(:,:,:,:,:) 563 type(matlu_type),allocatable :: glocnms(:) 564 type(matlu_type),allocatable :: glocsym(:) 565 real(dp),allocatable :: symrec_cart(:,:,:) 566 integer :: mt2g(3),mx2my2d 567 mt2g(1)=1 568 mt2g(2)=2 569 mt2g(3)=4 570 mx2my2d=5 571 t2g=paw_dmft%dmftqmc_t2g 572 x2my2d=paw_dmft%dmftqmc_x2my2d 573 574 ! DBG_ENTER("COLL") 575 576 ci=cone 577 nspinor=gloc(1)%nspinor 578 nsppol=gloc(1)%nsppol 579 natom=cryst_struc%natom 580 581 ABI_MALLOC(glocnm,(natom)) 582 ABI_MALLOC(glocnms,(natom)) 583 ABI_MALLOC(glocsym,(natom)) 584 call init_matlu(natom,nspinor,nsppol,gloc(1:natom)%lpawu,glocnm) 585 call init_matlu(natom,nspinor,nsppol,gloc(1:natom)%lpawu,glocnms) 586 call init_matlu(natom,nspinor,nsppol,gloc(1:natom)%lpawu,glocsym) 587 588 589 !========= Case nspinor ==1 ======================== 590 591 if (nspinor==1) then 592 ispinor=1 593 ispinor1=1 594 do iatom=1,cryst_struc%natom 595 do isppol=1,nsppol 596 if(gloc(iatom)%lpawu/=-1) then 597 lpawu=gloc(iatom)%lpawu 598 do m1=1, 2*lpawu+1 599 do m2=1, 2*lpawu+1 600 do irot=1,cryst_struc%nsym 601 at_indx=cryst_struc%indsym(4,irot,iatom) 602 do m3=1, 2*lpawu+1 603 do m4=1, 2*lpawu+1 604 if(t2g==1) then 605 m1s=mt2g(m1) 606 m2s=mt2g(m2) 607 m3s=mt2g(m3) 608 m4s=mt2g(m4) 609 lpawu_zarot=2 610 else if (x2my2d==1) then 611 m1s=mx2my2d 612 m2s=mx2my2d 613 m3s=mx2my2d 614 m4s=mx2my2d 615 lpawu_zarot=2 616 else 617 m1s=m1 618 m2s=m2 619 m3s=m3 620 m4s=m4 621 lpawu_zarot=lpawu 622 endif 623 zarot2=pawang%zarot(m3s,m1s,lpawu_zarot+1,irot)*pawang%zarot(m4s,m2s,lpawu_zarot+1,irot) 624 glocsym(iatom)%mat(m1,m2,isppol,ispinor,ispinor1)=& 625 & glocsym(iatom)%mat(m1,m2,isppol,ispinor,ispinor1)& 626 & +gloc(at_indx)%mat(m3,m4,isppol,ispinor,ispinor1)*zarot2 627 end do ! m3 628 end do ! m4 629 end do ! irot 630 glocsym(iatom)%mat(m1,m2,isppol,ispinor,ispinor1)=& 631 & glocsym(iatom)%mat(m1,m2,isppol,ispinor,ispinor1)/real(cryst_struc%nsym,kind=dp) 632 end do ! m2 633 end do ! m1 634 endif ! lpawu/=-1 635 end do ! isppol 636 end do ! iatom 637 !== Put glocsym into gloc 638 do iatom=1,cryst_struc%natom 639 if(gloc(iatom)%lpawu/=-1) then 640 gloc(iatom)%mat=glocsym(iatom)%mat 641 ! gloc(iatom)%mat(:,:,1,:,:)=(glocsym(iatom)%mat(:,:,1,:,:) & 642 !& + glocsym(iatom)%mat(:,:,2,:,:))/two 643 ! gloc(iatom)%mat(:,:,2,:,:)= gloc(iatom)%mat(:,:,1,:,:) 644 ! write(std_out,*) "WARNING: SYM non mag" 645 ! write(ab_out,*) "WARNING: SYM non mag" 646 endif 647 end do ! iatom 648 649 !========= Case nspinor ==2 ======================== 650 651 else if (nspinor==2) then 652 653 !== Allocate temporary arrays 654 do iatom=1,cryst_struc%natom 655 if(gloc(iatom)%lpawu/=-1) then 656 ndim=2*gloc(iatom)%lpawu+1 657 ABI_FREE(glocnm(iatom)%mat) 658 ABI_FREE(glocnms(iatom)%mat) 659 ABI_FREE(glocsym(iatom)%mat) 660 ABI_MALLOC(glocnm(iatom)%mat,(ndim,ndim,nsppol,4,1)) 661 ABI_MALLOC(glocnms(iatom)%mat,(ndim,ndim,nsppol,4,1)) 662 ABI_MALLOC(glocsym(iatom)%mat,(ndim,ndim,nsppol,2,2)) 663 endif 664 enddo 665 ABI_MALLOC(symrec_cart,(3,3,cryst_struc%nsym)) 666 667 !== Compute symrec_cart 668 do irot=1,cryst_struc%nsym 669 call symredcart(cryst_struc%gprimd,cryst_struc%rprimd,symrec_cart(:,:,irot),cryst_struc%symrec(:,:,irot)) 670 end do 671 672 !== Compute density matrix in density and magnetization representation 673 call chg_repr_matlu(gloc,glocnm,cryst_struc%natom,option=1,prtopt=1) 674 675 !== Do the sum over symetrized density matrix (in n,m repr) 676 isppol=1 677 do iatom=1,cryst_struc%natom 678 if(gloc(iatom)%lpawu/=-1) then 679 lpawu=gloc(iatom)%lpawu 680 ndim=2*gloc(iatom)%lpawu+1 681 do m1=1, 2*lpawu+1 682 do m2=1, 2*lpawu+1 683 sumrho=czero 684 rotmag=czero 685 do irot=1,cryst_struc%nsym 686 summag=czero 687 at_indx=cryst_struc%indsym(4,irot,iatom) 688 do m3=1, 2*lpawu+1 689 do m4=1, 2*lpawu+1 690 if(t2g==1) then 691 m1s=mt2g(m1) 692 m2s=mt2g(m2) 693 m3s=mt2g(m3) 694 m4s=mt2g(m4) 695 lpawu_zarot=2 696 else if (x2my2d==1) then 697 m1s=mx2my2d 698 m2s=mx2my2d 699 m3s=mx2my2d 700 m4s=mx2my2d 701 lpawu_zarot=2 702 else 703 m1s=m1 704 m2s=m2 705 m3s=m3 706 m4s=m4 707 lpawu_zarot=lpawu 708 endif 709 zarot2=pawang%zarot(m3s,m2s,lpawu_zarot+1,irot)*pawang%zarot(m4s,m1s,lpawu_zarot+1,irot) 710 sumrho=sumrho + glocnm(at_indx)%mat(m4,m3,isppol,1,1) * zarot2 711 do mu=1,3 712 summag(mu)=summag(mu) + glocnm(at_indx)%mat(m4,m3,isppol,mu+1,1) * zarot2 713 enddo 714 end do ! m3 715 end do !m4 716 717 ! == special case of magnetization 718 do nu=1,3 719 do mu=1,3 720 rotmag(mu)=rotmag(mu)+symrec_cart(mu,nu,irot)*summag(nu) 721 end do 722 end do 723 ! write(std_out,'(a,3i4,2x,3(2f10.5,2x))') "rotmag",irot,m1,m2,(rotmag(mu),mu=1,3) 724 end do ! irot 725 726 ! == Normalizes sum 727 sumrho=sumrho/real(cryst_struc%nsym,kind=dp) 728 ! sumrho=glocnm(isppol,1,iatom,m1,m2) ! test without sym 729 glocnms(iatom)%mat(m1,m2,isppol,1,1)=sumrho 730 do mu=1,3 731 rotmag(mu)=rotmag(mu)/real(cryst_struc%nsym,kind=dp) 732 ! rotmag(mu)=glocnm(isppol,mu+1,iatom,m1,m2) ! test without sym 733 glocnms(iatom)%mat(m1,m2,isppol,mu+1,1)=rotmag(mu) 734 enddo 735 end do ! m2 736 end do ! m1 737 endif ! lpawu/=-1 738 end do ! iatom 739 740 !== Compute back density matrix in upup dndn updn dnup representation 741 call chg_repr_matlu(glocsym,glocnms,cryst_struc%natom,option=-1,prtopt=1) 742 743 !== Put glocsym into gloc 744 do iatom=1,cryst_struc%natom 745 if(gloc(iatom)%lpawu/=-1) then 746 gloc(iatom)%mat=glocsym(iatom)%mat 747 ! gloc(iatom)%mat(:,:,1,:,:)=(glocsym(iatom)%mat(:,:,1,:,:) & 748 !& + glocsym(iatom)%mat(:,:,2,:,:))/two 749 ! gloc(iatom)%mat(:,:,2,:,:)= gloc(iatom)%mat(:,:,1,:,:) 750 ! write(std_out,*) "WARNING: SYM non mag" 751 ! write(ab_out,*) "WARNING: SYM non mag" 752 endif 753 end do ! iatom 754 755 ABI_FREE(symrec_cart) 756 endif 757 758 call destroy_matlu(glocnm,cryst_struc%natom) 759 call destroy_matlu(glocnms,cryst_struc%natom) 760 call destroy_matlu(glocsym,cryst_struc%natom) 761 ABI_FREE(glocnm) 762 ABI_FREE(glocnms) 763 ABI_FREE(glocsym) 764 !==============end of nspinor==2 case =========== 765 766 767 ! DBG_EXIT("COLL") 768 769 end subroutine sym_matlu
m_matlu/trace_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
trace_matlu
FUNCTION
Compute the trace of the matlu matrix
INPUTS
maltu(natom) <type(matlu_type)>= density matrix in the local orbital basis and related variables natom = number of atoms
OUTPUT
trace_loc(natom,nsppol+1)= trace for each atoms and each polarization, trace_loc(iatom,nsppol+1) is the full trace over polarization also.
SOURCE
1150 subroutine trace_matlu(matlu,natom,trace_loc,itau) 1151 1152 use defs_basis 1153 implicit none 1154 1155 !Arguments ------------------------------------ 1156 !type 1157 integer, intent(in) :: natom 1158 type(matlu_type), intent(in) :: matlu(natom) 1159 real(dp),intent(inout),target,optional :: trace_loc(natom,matlu(1)%nsppol+1) 1160 integer, intent(in),optional :: itau 1161 1162 !local variables------------------------------- 1163 integer :: iatom,isppol,ispinor,m,lpawu 1164 integer :: nsppol,nspinor 1165 real(dp), pointer :: traceloc(:,:)=>null() 1166 character(len=500) :: message 1167 ! ********************************************************************* 1168 nsppol=matlu(1)%nsppol 1169 nspinor=matlu(1)%nspinor 1170 if(present(trace_loc)) then 1171 traceloc=>trace_loc 1172 else 1173 ABI_MALLOC(traceloc,(natom,matlu(1)%nsppol+1)) 1174 endif 1175 1176 traceloc=zero 1177 do iatom = 1 , natom 1178 lpawu=matlu(iatom)%lpawu 1179 if(lpawu/=-1) then 1180 write(message,'(2a,i4)') ch10,' -------> For Correlated Atom', iatom 1181 if(.not.present(itau)) then 1182 call wrtout(std_out, message,'COLL') 1183 end if 1184 if(present(itau)) then 1185 if (itau>0) then 1186 call wrtout(std_out, message,'COLL') 1187 end if 1188 endif 1189 do isppol = 1 , nsppol 1190 do ispinor = 1 , nspinor 1191 do m = 1 , 2*lpawu+1 1192 traceloc(iatom,isppol)=traceloc(iatom,isppol)+& 1193 & matlu(iatom)%mat(m,m,isppol,ispinor,ispinor) 1194 enddo 1195 enddo 1196 traceloc(iatom,nsppol+1)=traceloc(iatom,nsppol+1)+traceloc(iatom,isppol) 1197 enddo 1198 if(nsppol==1.and.nspinor==1) traceloc(iatom,nsppol+1)=traceloc(iatom,nsppol+1)*two 1199 if(.not.present(itau)) then 1200 write(message,'(8x,a,f12.6)') 'Nb of Corr. elec. from G(w) is:'& 1201 & ,traceloc(iatom,nsppol+1) 1202 call wrtout(std_out, message,'COLL') 1203 endif 1204 if(present(itau)) then 1205 if(itau==1) then 1206 write(message,'(8x,a,f12.6)') 'Nb of Corr. elec. from G(tau) is:'& 1207 & ,traceloc(iatom,nsppol+1) 1208 call wrtout(std_out, message,'COLL') 1209 else if(itau==-1) then 1210 write(message,'(8x,a,f12.6)') 'Nb: Sum of the values of G0(tau=0-) is:'& 1211 & ,traceloc(iatom,nsppol+1) 1212 call wrtout(std_out, message,'COLL') 1213 else if(itau==4) then 1214 write(message,'(8x,a,f12.6)') 'Trace of matlu matrix is:'& 1215 & ,traceloc(iatom,nsppol+1) 1216 call wrtout(std_out, message,'COLL') 1217 endif 1218 endif 1219 endif 1220 enddo 1221 do iatom = 1 , natom 1222 lpawu=matlu(iatom)%lpawu 1223 if(lpawu/=-1) then 1224 !! MAG 1225 if(nsppol>1) then 1226 ! if(nsppol>1.and.present(itau)) then 1227 ! if(itau==1) then 1228 write(message,'(8x,a,f12.6)') 'DMFT Cor. El. Mag: ',traceloc(iatom,2)-traceloc(iatom,1) 1229 call wrtout(std_out, message,'COLL') 1230 ! endif 1231 endif 1232 1233 endif 1234 enddo 1235 if(.not.present(trace_loc)) then 1236 ABI_FREE(traceloc) 1237 traceloc => null() 1238 endif 1239 1240 end subroutine trace_matlu
m_matlu/zero_matlu [ Functions ]
[ Top ] [ m_matlu ] [ Functions ]
NAME
zero_matlu
FUNCTION
zero_matlu
INPUTS
maltu <type(matlu_type)>= density matrix in the local orbital basis and related variables natom = number of atoms
OUTPUT
maltu <type(matlu_type)>= density matrix in the local orbital basis and related variables
SOURCE
190 subroutine zero_matlu(matlu,natom,onlynondiag) 191 192 use defs_basis 193 implicit none 194 195 !Arguments ------------------------------------ 196 !scalars 197 integer, intent(in) :: natom 198 type(matlu_type),intent(inout) :: matlu(natom) 199 integer, optional, intent(in) :: onlynondiag 200 !Local variables------------------------------- 201 integer :: iatom,im,im1,ispinor,ispinor1,isppol,ndim 202 203 !********************************************************************* 204 205 if(present(onlynondiag)) then 206 do iatom=1,natom 207 if(matlu(iatom)%lpawu.ne.-1) then 208 do ispinor=1,matlu(iatom)%nspinor 209 ndim=(2*matlu(iatom)%lpawu+1) 210 do im=1,ndim 211 do im1=1,ndim 212 do ispinor1=1,matlu(iatom)%nspinor 213 if(im/=im1.or.ispinor/=ispinor1) then 214 do isppol=1,matlu(iatom)%nsppol 215 matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)=czero 216 enddo 217 end if 218 end do 219 end do 220 end do 221 end do 222 endif 223 enddo 224 else 225 do iatom=1,natom 226 matlu(iatom)%mat=czero 227 enddo 228 endif 229 230 231 end subroutine zero_matlu