TABLE OF CONTENTS


ABINIT/m_matlu [ Modules ]

[ Top ] [ 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 ]

[ Top ] [ m_matlu ] [ 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