TABLE OF CONTENTS


ABINIT/m_matlu [ Modules ]

[ Top ] [ Modules ]

NAME

  m_matlu

FUNCTION

COPYRIGHT

 Copyright (C) 2006-2018 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

PARENTS

CHILDREN

NOTES

  subroutines in this module must never call
   a subroutine of m_oper, m_green, m_self
   in order to avoid circular dependancies

SOURCE

29 #if defined HAVE_CONFIG_H
30 #include "config.h"
31 #endif
32 
33 
34 #include "abi_common.h"
35 
36 MODULE m_matlu
37 
38  use defs_basis
39  use m_errors
40  use m_abicore
41 
42  use m_hide_lapack,  only : xginv
43  use m_geometry, only : symredcart
44 
45  implicit none
46 
47  private
48 
49  public :: init_matlu
50  public :: inverse_matlu
51  public :: destroy_matlu
52  public :: diff_matlu
53  public :: add_matlu
54  public :: print_matlu
55  public :: sym_matlu
56  public :: copy_matlu
57  public :: gather_matlu
58  public :: zero_matlu
59  public :: trace_matlu
60  public :: diag_matlu
61  public :: rotate_matlu
62  public :: shift_matlu
63  public :: checkdiag_matlu
64  public :: checkreal_matlu
65  public :: prod_matlu
66  public :: conjg_matlu
67  public :: ln_matlu
68  public :: slm2ylm_matlu
69  public :: fac_matlu
70  public :: printplot_matlu
71  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

PARENTS

      dyson,hybridization_asymptotic_coefficient,m_green
      psichi_renormalization,qmc_prep_ctqmc

CHILDREN

SOURCE

1005 subroutine add_matlu(matlu1,matlu2,matlu3,natom,sign_matlu2)
1006 
1007  use defs_basis
1008  use m_paw_dmft, only : paw_dmft_type
1009  use m_crystal, only : crystal_t
1010 
1011 !This section has been created automatically by the script Abilint (TD).
1012 !Do not modify the following lines by hand.
1013 #undef ABI_FUNC
1014 #define ABI_FUNC 'add_matlu'
1015 !End of the abilint section
1016 
1017  implicit none
1018 
1019 !Arguments ------------------------------------
1020 !type
1021  integer,intent(in) :: natom,sign_matlu2
1022  type(matlu_type), intent(in) :: matlu1(natom),matlu2(natom)
1023  type(matlu_type), intent(inout) :: matlu3(natom) !vz_i
1024 
1025 !local variables-------------------------------
1026  integer :: iatom
1027 ! *********************************************************************
1028 
1029  do iatom = 1 , natom
1030    matlu3(iatom)%mat=matlu1(iatom)%mat+float(sign_matlu2)*matlu2(iatom)%mat
1031  enddo ! natom
1032 
1033 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-2018 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

PARENTS

      compute_levels

CHILDREN

SOURCE

2356  subroutine checkdiag_matlu(matlu,natom,tol,nondiag)
2357  use defs_basis
2358  use defs_wvltypes
2359  use m_crystal, only : crystal_t
2360 
2361 !This section has been created automatically by the script Abilint (TD).
2362 !Do not modify the following lines by hand.
2363 #undef ABI_FUNC
2364 #define ABI_FUNC 'checkdiag_matlu'
2365 !End of the abilint section
2366 
2367  implicit none
2368 
2369 !Arguments ------------------------------------
2370 !scalars
2371  real(dp),intent(in) :: tol
2372  integer, intent(in) :: natom
2373  logical, intent(out) :: nondiag
2374 !arrays
2375  type(matlu_type),intent(inout) :: matlu(natom)
2376 !Local variables-------------------------------
2377 !scalars
2378  integer :: iatom,im,im1,ispinor,ispinor1,isppol
2379  integer :: lpawu
2380 !arrays
2381 !************************************************************************
2382  nondiag=.false.
2383  do iatom=1,natom
2384    lpawu=matlu(iatom)%lpawu
2385    if(lpawu.ne.-1) then
2386      do im=1,2*lpawu+1
2387        do im1=1,2*lpawu+1
2388          do isppol=1,matlu(1)%nsppol
2389            do ispinor=1,matlu(1)%nspinor
2390 !              if(im/=im1) write(std_out,*) "im,im1",im,im1,matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor)
2391  !            if(present(nondiag).eqv..false.) then
2392  !              if(im/=im1.and.(abs(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor))>tol))  then
2393  !                write(message,'(5i5)') im,im1,isppol,ispinor,ispinor
2394  !                call wrtout(std_out,message,'COLL')
2395  !                write(message,'(a,3e16.5)')" checkdiag_matlu: Warning ",matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor),tol
2396  !                call wrtout(std_out,message,'COLL')
2397  !                if(.not.present(opt)) MSG_ERROR("not present(opt)")
2398  !                if(matlu(1)%nspinor==1) MSG_ERROR("matlu%nspinor==1")
2399  !              endif
2400 !             endif
2401              do ispinor1=1,matlu(1)%nspinor
2402  !              if(present(nondiag)) then
2403                  if((im/=im1.or.ispinor/=ispinor1)&
2404 &                         .and.abs(real(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)))>tol) then
2405                    nondiag=.true.
2406                   ! write(6,*) "NONDIAG", matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)
2407                  endif
2408                !if(ispinor/=ispinor1.and.(abs(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1))>tol))  then
2409                !  write(message,'(a,3e16.5)')" checkdiag_matlu :i Warning ",matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1),tol
2410                !  call wrtout(std_out,message,'COLL')
2411                !  write(message,'(5i5)') im,im1,isppol,ispinor,ispinor
2412                !  call wrtout(std_out,message,'COLL')
2413                !  if(matlu(1)%nspinor==1) MSG_ERROR("matlu%nspinor==1")
2414                !endif
2415              enddo
2416            enddo ! ispinor
2417          enddo ! isppol
2418        enddo ! im1
2419      enddo ! im
2420    endif ! lpawu
2421  enddo ! iatom
2422 
2423  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-2018 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

PARENTS

      qmc_prep_ctqmc

CHILDREN

SOURCE

2243  subroutine checkreal_matlu(matlu,natom,tol)
2244  use defs_basis
2245  use defs_wvltypes
2246  use m_crystal, only : crystal_t
2247 
2248 !This section has been created automatically by the script Abilint (TD).
2249 !Do not modify the following lines by hand.
2250 #undef ABI_FUNC
2251 #define ABI_FUNC 'checkreal_matlu'
2252 !End of the abilint section
2253 
2254  implicit none
2255 
2256 !Arguments ------------------------------------
2257 !scalars
2258  real(dp),intent(in) :: tol
2259  integer, intent(in) :: natom
2260 !arrays
2261  type(matlu_type),intent(inout) :: matlu(natom)
2262 !Local variables-------------------------------
2263 !scalars
2264  integer :: iatom,im,im1,ispinor,ispinor1,isppol
2265  integer :: lpawu
2266  character(len=500) :: message
2267  real(dp) :: maximag,maxoffdiag,maximagdiag
2268 !arrays
2269 !************************************************************************
2270  maximag=zero
2271  maximagdiag=zero
2272  maxoffdiag=zero
2273  do iatom=1,natom
2274    lpawu=matlu(iatom)%lpawu
2275    if(lpawu.ne.-1) then
2276      do im=1,2*lpawu+1
2277        do im1=1,2*lpawu+1
2278          do isppol=1,matlu(1)%nsppol
2279            do ispinor=1,matlu(1)%nspinor
2280              do ispinor1=1,matlu(1)%nspinor
2281                if(abs(aimag(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)))>maximag) then
2282                  maximag=abs(aimag(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)))
2283                  if(im==im1.and.ispinor==ispinor1) maximagdiag=abs(aimag(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)))
2284                endif
2285                if((im/=im1.or.ispinor/=ispinor1).and.abs(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1))>maxoffdiag) then
2286                  maxoffdiag=abs(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1))
2287                endif
2288              enddo
2289            enddo ! ispinor
2290          enddo ! isppol
2291        enddo ! im1
2292      enddo ! im
2293    endif ! lpawu
2294  enddo ! iatom
2295  if (maximagdiag>tol) then
2296    write(message,'(3x,2a,e12.4,a,e12.4,2a)') ch10,&
2297 &   ' Diagonal part of the occupation matrix is complex: the imaginary part ',&
2298 &     maximagdiag,' is larger than',tol,ch10  &
2299 &    , "The calculation cannot handle it : check that your calculation is meaningfull"
2300    MSG_ERROR(message)
2301  endif
2302  if (maximag>tol) then
2303    write(message,'(3x,2a,e12.4,a,e12.4,2a)') ch10,&
2304 &   ' Off diag occupation matrix is complex: the imaginary part ',maximag,' is larger than',tol,ch10&
2305     , "Check that your calculation is meaningfull"
2306    MSG_WARNING(message)
2307  else
2308    write(message,'(3x,2a,e12.4,a,e12.4,2a)') ch10,&
2309 &   ' Occupation matrix is real: the imaginary part ',maximag,' is lower than',tol
2310    MSG_COMMENT(message)
2311  endif
2312  if (maxoffdiag>tol) then
2313    write(message,'(3x,2a,e12.4,a,e12.4,6a)') ch10,&
2314 &   ' Occupation matrix is non diagonal : the maximum off-diag part ',maxoffdiag,' is larger than',tol,ch10&
2315 &    , "The corresponding non diagonal elements will be neglected in the Weiss/Hybridization functions",ch10&
2316 &    , "(Except if dmft_solv=8 where these elements are taken into accounts)",ch10&
2317 &    , "This is an approximation"
2318    MSG_WARNING(message)
2319  else
2320    write(message,'(3x,2a,e12.4,a,e12.4,2a)') ch10,&
2321 &   ' Occupation matrix is diagonal : the off-diag part ',maxoffdiag,' is lower than',tol
2322    MSG_COMMENT(message)
2323  endif
2324 
2325  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-2018 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

PARENTS

      m_matlu

CHILDREN

SOURCE

1071  subroutine chg_repr_matlu(glocspsp,glocnm,natom,option,prtopt)
1072  use defs_basis
1073  use defs_wvltypes
1074  use m_crystal, only : crystal_t
1075 
1076 !This section has been created automatically by the script Abilint (TD).
1077 !Do not modify the following lines by hand.
1078 #undef ABI_FUNC
1079 #define ABI_FUNC 'chg_repr_matlu'
1080 !End of the abilint section
1081 
1082  implicit none
1083 
1084 !Arguments ------------------------------------
1085 !scalars
1086  integer,intent(in) :: natom,option,prtopt
1087 !arrays
1088  type(matlu_type),intent(inout) :: glocspsp(natom)
1089  type(matlu_type),intent(inout) :: glocnm(natom)
1090 !Local variables-------------------------------
1091 !scalars
1092  integer :: iatom,isppol,lpawu,m1,m2,ndim,nsppol,mu
1093  complex(dpc) :: ci
1094  character(len=500) :: message
1095 
1096 ! DBG_ENTER("COLL")
1097 
1098  ci=j_dpc
1099 
1100 !==  Compute density matrix in density magnetization representation
1101  if (option==1) then
1102   nsppol=glocspsp(1)%nsppol
1103   do isppol=1,nsppol
1104    do iatom=1,natom
1105     if(glocspsp(iatom)%lpawu/=-1) then
1106      ndim=2*glocspsp(iatom)%lpawu+1
1107      do m1=1,ndim
1108       do m2=1,ndim
1109        glocnm(iatom)%mat(m1,m2,isppol,1,1)=glocspsp(iatom)%mat(m1,m2,isppol,1,1) &
1110 &                                         +glocspsp(iatom)%mat(m1,m2,isppol,2,2)
1111        glocnm(iatom)%mat(m1,m2,isppol,4,1)=glocspsp(iatom)%mat(m1,m2,isppol,1,1) &
1112 &                                         -glocspsp(iatom)%mat(m1,m2,isppol,2,2)
1113        glocnm(iatom)%mat(m1,m2,isppol,2,1)=glocspsp(iatom)%mat(m1,m2,isppol,1,2) &
1114 &                                         +glocspsp(iatom)%mat(m1,m2,isppol,2,1)
1115        glocnm(iatom)%mat(m1,m2,isppol,3,1)= &
1116 &                               cmplx((aimag(glocspsp(iatom)%mat(m1,m2,isppol,2,1))   &
1117 &                                     -aimag(glocspsp(iatom)%mat(m1,m2,isppol,1,2))), &
1118 &                                    (-real(glocspsp(iatom)%mat(m1,m2,isppol,2,1))+  &
1119 &                                      real(glocspsp(iatom)%mat(m1,m2,isppol,1,2))),kind=dp)
1120       enddo  ! m2
1121      enddo ! m1
1122      if(abs(prtopt)>=3) then
1123       write(message,'(a)') "        -- in n, m repr "
1124       call wrtout(std_out,  message,'COLL')
1125       do mu=1,4
1126        do m1=1,ndim
1127         write(message,'(8x,(14(2f9.5,2x)))')(glocnm(iatom)%mat(m1,m2,isppol,mu,1),m2=1,ndim)
1128         call wrtout(std_out,  message,'COLL')
1129        enddo ! m1
1130        write(message,'(a)') ch10
1131        call wrtout(std_out,  message,'COLL')
1132       enddo ! mu
1133      endif ! prtopt >3
1134     endif ! lpawu/=-1
1135    enddo
1136   enddo
1137 
1138 !==  Compute back density matrix in upup dndn updn dnup representation
1139  else  if (option==-1) then
1140   isppol=1
1141   do iatom=1,natom
1142    if(glocnm(iatom)%lpawu/=-1) then
1143     lpawu=glocnm(iatom)%lpawu
1144     ndim=2*glocnm(iatom)%lpawu+1
1145     do m1=1, 2*lpawu+1
1146      do m2=1, 2*lpawu+1
1147       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))
1148       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))
1149       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))
1150       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))
1151      end do  ! m2
1152     end do ! m1
1153     if(abs(prtopt)>6) then
1154      write(message,'(a)') "        -- in spin spin repr "
1155      call wrtout(std_out,  message,'COLL')
1156      do mu=1,4
1157       do m1=1,ndim
1158        write(message,'(8x,14(2f9.5,2x))')(glocspsp(iatom)%mat(m1,m2,isppol,mu,1),m2=1,ndim)
1159        call wrtout(std_out,  message,'COLL')
1160       enddo
1161       write(message,'(a)') ch10
1162       call wrtout(std_out,  message,'COLL')
1163      enddo
1164     endif !prtopt>3
1165    endif ! lpawu/=-1
1166   end do ! iatom
1167  else
1168   message = "stop in chg_repr_matlu"
1169   MSG_ERROR(message)
1170  endif
1171 
1172 
1173 ! DBG_EXIT("COLL")
1174 
1175  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-2018 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

PARENTS

CHILDREN

SOURCE

2534  subroutine conjg_matlu(matlu1,natom)
2535  use defs_basis
2536  use defs_wvltypes
2537 
2538 !This section has been created automatically by the script Abilint (TD).
2539 !Do not modify the following lines by hand.
2540 #undef ABI_FUNC
2541 #define ABI_FUNC 'conjg_matlu'
2542 !End of the abilint section
2543 
2544  implicit none
2545 
2546 !Arguments ------------------------------------
2547 !scalars
2548  integer, intent(in) :: natom
2549 !arrays
2550  type(matlu_type), intent(inout) :: matlu1(natom)
2551 !Local variables-------------------------------
2552 !scalars
2553  integer :: iatom,im1,im2,ispinor2,ispinor1,isppol
2554  integer :: lpawu
2555 !arrays
2556 !************************************************************************
2557  do iatom=1,natom
2558    lpawu=matlu1(iatom)%lpawu
2559    if(lpawu.ne.-1) then
2560      do isppol=1,matlu1(1)%nsppol
2561        do ispinor1=1,matlu1(1)%nspinor
2562          do ispinor2=1,matlu1(1)%nspinor
2563            do im1=1,2*lpawu+1
2564              do im2=1,2*lpawu+1
2565                matlu1(iatom)%mat(im1,im2,isppol,ispinor1,ispinor2)= &
2566 &               conjg(matlu1(iatom)%mat(im1,im2,isppol,ispinor1,ispinor2))
2567              enddo ! im2
2568            enddo ! im1
2569          enddo ! ispinor2
2570        enddo ! ispinor1
2571      enddo ! isppol
2572    endif ! lpawu
2573  enddo ! iatom
2574 
2575  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

PARENTS

      m_datafordmft,hubbard_one,m_green,m_matlu,m_oper,m_self,qmc_prep_ctqmc
      spectral_function

CHILDREN

SOURCE

315 subroutine copy_matlu(nmat1,nmat2,natom,opt_diag,opt_non_diag,opt_re)
316 
317  use defs_basis
318 
319 !This section has been created automatically by the script Abilint (TD).
320 !Do not modify the following lines by hand.
321 #undef ABI_FUNC
322 #define ABI_FUNC 'copy_matlu'
323 !End of the abilint section
324 
325  implicit none
326 
327 !Arguments ------------------------------------
328 !type
329  integer, intent(in) :: natom
330  type(matlu_type),intent(in) :: nmat1(natom)
331  type(matlu_type),intent(inout) :: nmat2(natom) !vz_i
332  integer, optional, intent(in) :: opt_diag,opt_non_diag,opt_re
333 
334 !Local variables-------------------------------
335  integer :: iatom,isppol,im1,im2,ispinor,ispinor1,tndim
336 ! *********************************************************************
337 
338  do isppol=1,nmat1(1)%nsppol
339    do iatom=1,natom
340      if(nmat1(iatom)%lpawu.ne.-1) then
341        tndim=(2*nmat1(iatom)%lpawu+1)
342        do im1=1,tndim
343          do im2=1,tndim
344            do ispinor=1,nmat1(1)%nspinor
345              do ispinor1=1,nmat1(1)%nspinor
346                if(present(opt_diag)) then
347                  nmat2(iatom)%mat(im1,im1,isppol,ispinor,ispinor)=nmat1(iatom)%mat(im1,im1,isppol,ispinor,ispinor)
348                else if(present(opt_non_diag)) then
349                  if(im1/=im2.or.ispinor/=ispinor1) then
350                    nmat2(iatom)%mat(im1,im2,isppol,ispinor,ispinor1)=nmat1(iatom)%mat(im1,im2,isppol,ispinor,ispinor1)
351                  endif
352                else
353                  nmat2(iatom)%mat(im1,im2,isppol,ispinor,ispinor1)=nmat1(iatom)%mat(im1,im2,isppol,ispinor,ispinor1)
354                  if(present(opt_re)) nmat2(iatom)%mat(im1,im2,isppol,ispinor,ispinor1)=&
355 &                         cmplx(real(nmat1(iatom)%mat(im1,im2,isppol,ispinor,ispinor1)),0.d0,kind=dp)
356                endif
357              enddo
358            enddo
359          enddo
360        enddo
361      endif ! lpawu=-1
362    enddo ! iatom
363  enddo ! isppol
364 !   do iatom=1,natom
365 !    lpawu=nmat1(iatom)%lpawu
366 !    if(lpawu.ne.-1) then
367 !     nmat2(iatom)%mat=nmat1(iatom)%mat
368 !    endif
369 !   enddo
370 
371 
372 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

PARENTS

      m_datafordmft,hubbard_one,m_green,m_matlu,m_oper,qmc_prep_ctqmc

CHILDREN

SOURCE

261 subroutine destroy_matlu(matlu,natom)
262 
263  use defs_basis
264  use m_crystal, only : crystal_t
265 
266 !This section has been created automatically by the script Abilint (TD).
267 !Do not modify the following lines by hand.
268 #undef ABI_FUNC
269 #define ABI_FUNC 'destroy_matlu'
270 !End of the abilint section
271 
272  implicit none
273 
274 !Arguments ------------------------------------
275 !scalars
276  integer, intent(in) :: natom
277  type(matlu_type),intent(inout) :: matlu(natom)
278 
279 !Local variables-------------------------------
280  integer :: iatom
281 
282 ! *********************************************************************
283 
284  do iatom=1,natom
285   if ( allocated(matlu(iatom)%mat) )   then
286     ABI_DEALLOCATE(matlu(iatom)%mat)
287   end if
288  enddo
289 
290 end subroutine destroy_matlu

m_matlu/diag_matlu [ Functions ]

[ Top ] [ m_matlu ] [ Functions ]

NAME

 diag_matlu

FUNCTION

 Diagonalize matlu matrix

COPYRIGHT

 Copyright (C) 2005-2018 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

PARENTS

      hubbard_one,qmc_prep_ctqmc

CHILDREN

SOURCE

1447  subroutine diag_matlu(matlu,matlu_diag,natom,prtopt,eigvectmatlu,nsppol_imp,checkstop,optreal,test)
1448  use defs_basis
1449  use defs_wvltypes
1450  use m_crystal, only : crystal_t
1451  use m_matrix,         only : blockdiago_fordsyev,blockdiago_forzheev
1452 
1453 !This section has been created automatically by the script Abilint (TD).
1454 !Do not modify the following lines by hand.
1455 #undef ABI_FUNC
1456 #define ABI_FUNC 'diag_matlu'
1457 !End of the abilint section
1458 
1459  implicit none
1460 
1461 !Arguments ------------------------------------
1462 !scalars
1463  integer, intent(in) :: natom
1464  integer, intent(in) :: prtopt
1465 !arrays
1466  type(matlu_type),intent(inout) :: matlu(natom)
1467  type(matlu_type),intent(inout) :: matlu_diag(natom) !vz_i
1468  type(coeff2c_type),optional,intent(inout) :: eigvectmatlu(natom,matlu(1)%nsppol) !vz_i
1469  integer,optional,intent(in) :: nsppol_imp
1470  logical,optional,intent(in) :: checkstop
1471  integer,optional,intent(in) :: optreal
1472  integer,optional,intent(in) :: test
1473 !Local variables-------------------------------
1474 !scalars
1475  integer :: iatom,im1,im2,im3,imc,imc1,info,ispinor,ispinor1,isppol,lwork,tndim,lworkr
1476  integer :: nsppol,nsppolimp,nspinor
1477  logical :: checkstop_in,blockdiag
1478  character(len=500) :: message
1479 !arrays
1480  type(coeff2c_type),allocatable :: gathermatlu(:)
1481  real(dp),allocatable :: eig(:),rwork(:),work(:),valuer(:,:)!,valuer2(:,:)
1482  !real(dp),allocatable :: valuer3(:,:),valuer4(:,:)
1483 ! real(dp),allocatable :: eigvec(:,:)
1484  complex(dpc),allocatable :: zwork(:)
1485  logical :: donotdiag,print_temp_mat2
1486  complex(dpc),allocatable :: temp_mat(:,:)
1487  complex(dpc),allocatable :: temp_mat2(:,:)
1488 !debug complex(dpc),allocatable :: temp_mat3(:,:)
1489 !************************************************************************
1490 
1491  call zero_matlu(matlu_diag,natom)
1492  nsppol=matlu(1)%nsppol
1493  if(present(nsppol_imp)) then
1494    nsppolimp=nsppol_imp
1495  else
1496    nsppolimp=nsppol
1497  endif
1498  if(present(checkstop)) then
1499   checkstop_in=checkstop
1500  else
1501   checkstop_in=.true.
1502  endif
1503  if(present(test)) then
1504   blockdiag=(test==8)
1505  else
1506   blockdiag=.false.
1507  endif
1508  nspinor=matlu(1)%nspinor
1509 
1510  donotdiag=.true.
1511  donotdiag=.false.
1512 ! ===========================
1513 ! Check is diagonalization is necessary and how
1514 ! ===========================
1515  do isppol=1,matlu(1)%nsppol
1516    do iatom=1,natom
1517       if(matlu(iatom)%lpawu.ne.-1) then
1518        tndim=(2*matlu(iatom)%lpawu+1)
1519         do im1=1,tndim
1520           do im2=1,tndim
1521             do ispinor=1,nspinor
1522               do ispinor1=1,nspinor
1523                 if(abs(matlu(iatom)%mat(im1,im2,isppol,ispinor,ispinor1))>tol8.and.&
1524 &                 (im1/=im2.or.ispinor/=ispinor1)) then
1525 !                 if matrix is diagonal: do not diagonalize
1526                   donotdiag=.false.
1527                   exit
1528                 endif
1529               enddo
1530             enddo
1531           enddo
1532         enddo
1533       endif
1534    enddo
1535  enddo
1536 
1537  if(donotdiag) then
1538    do isppol=1,matlu(1)%nsppol
1539      do iatom=1,natom
1540         if(matlu(iatom)%lpawu.ne.-1) then
1541           tndim=(2*matlu(iatom)%lpawu+1)
1542           eigvectmatlu(iatom,isppol)%value(:,:)=czero
1543           do im1=1,tndim
1544             eigvectmatlu(iatom,isppol)%value(im1,im1)=cone
1545           enddo
1546         endif
1547      enddo
1548    enddo
1549    call copy_matlu(matlu,matlu_diag,natom)
1550    write(message,'(a)')  "   Diagonalisation of matlu will not be performed"
1551    call wrtout(std_out,message,'COLL')
1552    return
1553  endif
1554 
1555 
1556 ! Below, gathermatlu is used to store eigenvectors (after diagonalization)
1557 ! For nsppol=2, and if nsppolimp=1, these eigenvectors computed for isppol=1, and applied to
1558 ! rotate matrix for isppol=2. It is the reason why the sum below is only
1559 ! from 1 to nsppolimp !
1560 
1561 
1562 
1563  do isppol=1,nsppolimp !
1564 ! ===========================
1565 ! Define gathermatlu
1566 ! ===========================
1567    ABI_DATATYPE_ALLOCATE(gathermatlu,(natom))
1568    do iatom=1,natom
1569      if(matlu(iatom)%lpawu.ne.-1) then
1570        tndim=nspinor*(2*matlu(iatom)%lpawu+1)
1571        ABI_ALLOCATE(gathermatlu(iatom)%value,(tndim,tndim))
1572        gathermatlu(iatom)%value=czero
1573      endif
1574    enddo
1575    if(nsppol==1.and.nspinor==2) then
1576      call gather_matlu(matlu,gathermatlu,natom,option=1,prtopt=1)
1577    else if((nsppol==2.or.nsppol==1).and.nspinor==1) then
1578      do iatom=1,natom
1579        if(matlu(iatom)%lpawu.ne.-1) then
1580          do im1=1,tndim
1581            do im2=1,tndim
1582              gathermatlu(iatom)%value(im1,im2)=matlu(iatom)%mat(im1,im2,isppol,1,1)
1583            enddo
1584          enddo
1585        endif
1586      enddo
1587    endif
1588 ! ===========================
1589 ! Diagonalize
1590 ! ===========================
1591    do iatom=1,natom
1592      if(matlu(iatom)%lpawu.ne.-1) then
1593        if(isppol<=nsppolimp) then
1594          tndim=nspinor*(2*matlu(iatom)%lpawu+1)
1595 !debug       allocate(temp_mat2(tndim,tndim))
1596 !debug       temp_mat2=zero
1597          lwork=2*tndim-1
1598          ABI_ALLOCATE(rwork,(3*tndim-2))
1599          rwork = zero
1600 
1601          lworkr=tndim*(tndim+2)*2
1602          ABI_ALLOCATE(work,(lworkr))
1603          work = zero
1604          ABI_ALLOCATE(valuer,(tndim,tndim))
1605 !         ABI_ALLOCATE(valuer2,(tndim,tndim))
1606 !         ABI_ALLOCATE(valuer3,(tndim,tndim))
1607 !         ABI_ALLOCATE(valuer4,(tndim,tndim))
1608          ABI_ALLOCATE(zwork,(lwork))
1609 !         valuer2=zero
1610 !         valuer3=zero
1611 !         valuer4=zero
1612          zwork = czero
1613          ABI_ALLOCATE(eig,(tndim))
1614          eig = zero
1615          info = 0
1616          if(prtopt>=4) then
1617            write(message,'(a,i4,a,i4)')  "       BEFORE DIAGONALIZATION for atom",iatom,"  and isppol",isppol
1618            call wrtout(std_out,message,'COLL')
1619            do im1=1,tndim
1620              write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')&
1621 &             (gathermatlu(iatom)%value(im1,im2),im2=1,tndim)
1622              call wrtout(std_out,message,'COLL')
1623            end do
1624          endif
1625 !debug       temp_mat2(:,:)=gathermatlu(iatom)%value(:,:)
1626 !           write(std_out,*)"diag"
1627          if(present(optreal).and.maxval(abs(aimag(gathermatlu(iatom)%value(:,:))))<tol6) then
1628            write(message,'(a,2x,a,e9.3,a)') ch10,"Imaginary part of Local Hamiltonian is lower than ",&
1629 &                   tol8, ": the real matrix is used"
1630            call wrtout(std_out,message,'COLL')
1631            valuer=real(gathermatlu(iatom)%value,kind=dp)
1632 !           write(message,'(a)') ch10
1633 !           call wrtout(std_out,message,'COLL')
1634 !           write(message,'(a,i4,a,i4)')  "BEFORE valuer for atom",iatom,"  and isppol",isppol
1635 !           call wrtout(std_out,message,'COLL')
1636 !           do im1=1,tndim
1637 !             write(message,'(2(1x,18(1x,"(",f20.15,",",f20.15,")")))')&
1638 !&             (valuer(im1,im2),im2=1,tndim)
1639 !             call wrtout(std_out,message,'COLL')
1640 !           end do
1641 !           do im1=1,tndim
1642 !             valuer(im1,im1)=real(im1,kind=dp)*0.00000000001_dp+valuer(im1,im1)
1643 !           enddo
1644 !           write(message,'(a)') ch10
1645 !           call wrtout(std_out,message,'COLL')
1646 !           write(message,'(a,i4,a,i4)')  "BEFORE valuer for atom",iatom,"  and isppol",isppol
1647 !           call wrtout(std_out,message,'COLL')
1648 !           do im1=1,tndim
1649 !             write(message,'(2(1x,18(1x,f20.15,f20.15)))')&
1650 !&             (valuer(im1,im2),im2=1,tndim)
1651 !             call wrtout(std_out,message,'COLL')
1652 !           end do
1653            !call dsyev('v','u',tndim,valuer,tndim,eig,work,lworkr,info)
1654            if (blockdiag) then
1655              call blockdiago_fordsyev(valuer,tndim,eig)
1656            else
1657              call dsyev('v','u',tndim,valuer,tndim,eig,work,lworkr,info)
1658            endif
1659 !!       For reproductibility
1660 !           ! valuer2: eigenvector for the perturb matrix
1661 !           valuer2=real(gathermatlu(iatom)%value,kind=dp)
1662 !           do im1=1,tndim
1663 !             valuer2(im1,im1)=float(im1)*0.00000000001+valuer2(im1,im1)
1664 !           enddo
1665 !           call dsyev('v','u',tndim,valuer2,tndim,eig,work,lworkr,info)
1666 !           write(message,'(a)') ch10
1667 !           call wrtout(std_out,message,'COLL')
1668 !           write(message,'(a,i4,a,i4)')  "       valuer2 for atom",iatom,"  and isppol",isppol
1669 !           call wrtout(std_out,message,'COLL')
1670 !           do im1=1,tndim
1671 !             write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')&
1672 !&             (valuer2(im1,im2),im2=1,tndim)
1673 !             call wrtout(std_out,message,'COLL')
1674 !           end do
1675 !           call dgemm('n','n',tndim,tndim,tndim,cone,valuer,tndim,&
1676 !&            valuer2,tndim,czero,valuer3                ,tndim)
1677 !           call dgemm('c','n',tndim,tndim,tndim,cone,valuer2,tndim,&
1678 !&            valuer3                   ,tndim,czero,valuer4,tndim)
1679 !           ! valuer4: compute unpert matrix in the basis of the
1680 !           ! perturb basis
1681 !           write(message,'(a)') ch10
1682 !           call wrtout(std_out,message,'COLL')
1683 !           write(message,'(a,i4,a,i4)')  "BEFORE valuer4 for atom",iatom,"  and isppol",isppol
1684 !           call wrtout(std_out,message,'COLL')
1685 !           do im1=1,tndim
1686 !             write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')&
1687 !&             (valuer4(im1,im2),im2=1,tndim)
1688 !             call wrtout(std_out,message,'COLL')
1689 !           end do
1690 !           call dsyev('v','u',tndim,valuer4,tndim,eig,work,lworkr,info)
1691 !           ! valuer4: Diago valuer4 (nearly diag)
1692 !           write(message,'(a)') ch10
1693 !           call wrtout(std_out,message,'COLL')
1694 !           write(message,'(a,i4,a,i4)')  "AFTER  valuer4 for atom",iatom,"  and isppol",isppol
1695 !           call wrtout(std_out,message,'COLL')
1696 !           do im1=1,tndim
1697 !             write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')&
1698 !&             (valuer4(im1,im2),im2=1,tndim)
1699 !             call wrtout(std_out,message,'COLL')
1700 !           end do
1701 !           call dgemm('n','n',tndim,tndim,tndim,cone,valuer2,tndim,&
1702 !&            valuer4,tndim,czero,valuer                ,tndim)
1703            write(6,*) "INFO",info
1704            gathermatlu(iatom)%value=cmplx(valuer,0.d0,kind=dp)
1705 !           write(message,'(a,i4,a,i4)')  "AFTER valuer for atom",iatom,"  and isppol",isppol
1706 !           call wrtout(std_out,message,'COLL')
1707 !           do im1=1,tndim
1708 !             write(message,'(2(1x,18(1x,"(",f20.15,",",f20.15,")")))')&
1709 !&             (valuer(im1,im2),im2=1,tndim)
1710 !             call wrtout(std_out,message,'COLL')
1711 !           end do
1712          else
1713            if(present(optreal).and.maxval(abs(aimag(gathermatlu(iatom)%value(:,:))))>tol8) then
1714              write(message,'(a)') " Local hamiltonian in correlated basis is complex"
1715              MSG_COMMENT(message)
1716            endif
1717            call zheev('v','u',tndim,gathermatlu(iatom)%value,tndim,eig,zwork,lwork,rwork,info)
1718            !call blockdiago_forzheev(gathermatlu(iatom)%value,tndim,eig)
1719          endif
1720          if(prtopt>=3) then
1721            write(message,'(a)') ch10
1722            call wrtout(std_out,message,'COLL')
1723            write(message,'(a,i4,a,i4)')  "       EIGENVECTORS for atom",iatom,"  and isppol",isppol
1724            call wrtout(std_out,message,'COLL')
1725            do im1=1,tndim
1726              write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')&
1727 &             (gathermatlu(iatom)%value(im1,im2),im2=1,tndim)
1728              call wrtout(std_out,message,'COLL')
1729            end do
1730           ! do im1=1,tndim
1731           !   xcheck=czero
1732           !   do im3=1,tndim
1733           !     do im2=1,tndim
1734           !       xcheck=xcheck+gathermatlu(iatom)%value(im1,im2)*conjg(gathermatlu(iatom)%value(im2,im3))
1735           !     end do
1736           !   end do
1737           !   write(6,*) "check",im3,im1,xcheck
1738           ! end do
1739          endif
1740 !       write(std_out,*) "eig",eig
1741 ! ===========================
1742 ! Put eigenvalue in matlu_diag
1743 ! ===========================
1744          imc=0
1745          do ispinor=1,nspinor
1746            do im1=1,2*matlu(iatom)%lpawu+1
1747              imc=imc+1
1748              matlu_diag(iatom)%mat(im1,im1,isppol,ispinor,ispinor)=eig(imc)
1749 !debug             temp_mat2(imc,imc)=eig(imc)
1750            enddo
1751          enddo
1752          if(prtopt>=2) then
1753 !            write(message,'(a,12(1x,18(1x,"(",f9.3,",",f9.3,")")))')&
1754 !&            ch10,(eig(im1),im1=1,tndim)
1755 !             call wrtout(std_out,message,'COLL')
1756            !call wrtout(std_out,message,'COLL')
1757             !write(std_out,*) "EIG", eig
1758          endif
1759          ABI_DEALLOCATE(zwork)
1760          ABI_DEALLOCATE(rwork)
1761          ABI_DEALLOCATE(work)
1762          ABI_DEALLOCATE(valuer)
1763 !         ABI_DEALLOCATE(valuer2)
1764 !         ABI_DEALLOCATE(valuer3)
1765 !         ABI_DEALLOCATE(valuer4)
1766          ABI_DEALLOCATE(eig)
1767 !     endif
1768 !   enddo
1769 ! ===========================
1770 ! Keep eigenvectors gathermatlu
1771 ! ===========================
1772          if (present(eigvectmatlu)) then
1773            tndim=nspinor*(2*matlu(iatom)%lpawu+1)
1774            eigvectmatlu(iatom,isppol)%value(:,:)=gathermatlu(iatom)%value(:,:)
1775 !           write(std_out,*) "eigvect in diag_matlu"
1776 !           do im1=1,tndim
1777 !             write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')&
1778 !&             (gathermatlu(iatom)%value(im1,im2),im2=1,tndim)
1779 !             call wrtout(std_out,message,'COLL')
1780 !           end do
1781          endif
1782 
1783 
1784 ! ==================================================================
1785        endif
1786        if(nsppolimp==1.and.matlu(1)%nsppol==2) then
1787 ! If necessary rotate levels for this other spin, assuming the same
1788 ! rotation matrix: it have to be checked afterwards that the matrix is
1789 ! diagonal
1790 ! ===================================================================
1791          ABI_ALLOCATE(temp_mat,(tndim,tndim))
1792          ABI_ALLOCATE(temp_mat2,(tndim,tndim))
1793          temp_mat(:,:)=czero
1794 !        input matrix: gathermatlu
1795 !        rotation matrix: eigvectmatlu
1796 !        intermediate matrix: temp_mat
1797 !        result matrix: temp_mat2
1798          do im1=1,tndim
1799            do im2=1,tndim
1800              gathermatlu(iatom)%value(im1,im2)=matlu(iatom)%mat(im1,im2,2,1,1)
1801            enddo
1802          enddo
1803          if(prtopt>=3) then
1804            write(message,'(a,i4,a,i4)')  "       GATHERMATLU for atom",iatom," inside if nsppolimp==1"
1805            call wrtout(std_out,message,'COLL')
1806            do im1=1,tndim
1807              write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')&
1808 &             (gathermatlu(iatom)%value(im1,im2),im2=1,tndim)
1809              call wrtout(std_out,message,'COLL')
1810            end do
1811          endif
1812          call zgemm('n','n',tndim,tndim,tndim,cone,gathermatlu(iatom)%value,tndim,&
1813 &          eigvectmatlu(iatom,1)%value,tndim,czero,temp_mat                ,tndim)
1814          call zgemm('c','n',tndim,tndim,tndim,cone,eigvectmatlu(iatom,1)%value,tndim,&
1815 &          temp_mat                   ,tndim,czero,temp_mat2,tndim)
1816          eigvectmatlu(iatom,2)%value=eigvectmatlu(iatom,1)%value
1817          imc=0
1818          print_temp_mat2=.false.
1819          do ispinor=1,nspinor
1820            do im1=1,2*matlu(iatom)%lpawu+1
1821              imc=imc+1
1822              imc1=0
1823              do ispinor1=1,nspinor
1824                do im2=1,2*matlu(iatom)%lpawu+1
1825                  imc1=imc1+1
1826                  matlu_diag(iatom)%mat(im1,im2,2,ispinor,ispinor1)=temp_mat2(imc,imc1)
1827                  if (imc/=imc1.and.(abs(temp_mat2(imc,imc1))>tol5)) then
1828                    write(message,'(3a,i4,2f16.4)') ch10,'diag_matlu= Matrix for spin number 2 obtained with', &
1829 &                   ' eigenvectors from diagonalization for spin nb 1 is non diagonal for atom:',iatom,&
1830 &                    abs(temp_mat2(imc,imc1)),tol5
1831                    call wrtout(std_out,message,'COLL')
1832                    if(.not.checkstop_in.and.(abs(temp_mat2(imc,imc1))>0.1)) print_temp_mat2=.true.
1833                    if(checkstop_in) print_temp_mat2=.true.
1834                  endif
1835 !                 write(std_out,*) temp_mat2(imc,imc1)
1836 !debug             temp_mat2(imc,imc)=eig(imc)
1837                enddo
1838              enddo
1839            enddo
1840          enddo
1841          if(print_temp_mat2.and.prtopt>=3) then
1842            write(message,'(a)')  "       temp_mat2"
1843            call wrtout(std_out,message,'COLL')
1844            do im1=1,tndim
1845              write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')&
1846 &             (temp_mat2(im1,im3),im3=1,tndim)
1847              call wrtout(std_out,message,'COLL')
1848            end do
1849            if(iatom==2) then
1850              MSG_ERROR("iatom==2")
1851            end if
1852          endif
1853          ABI_DEALLOCATE(temp_mat)
1854          ABI_DEALLOCATE(temp_mat2)
1855        endif
1856 
1857 
1858      endif ! end lpawu/=-1
1859 
1860 !!  for check only
1861 !debug     if(matlu(iatom)%lpawu.ne.-1) then
1862 !debug       allocate(temp_mat(tndim,tndim))
1863 !debug       allocate(temp_mat3(tndim,tndim))
1864 !debug           do im1=1,tndim
1865 !debug             do im2=1,tndim
1866 !debug!               rot_mat(iatom,isppol)%value(im1,im2)=rot_mat_orig(iatom,isppol)%value(im1,im2)
1867 !debug               temp_mat3(im1,im2)=conjg(gathermatlu(iatom)%value(im2,im1))
1868 !debug             enddo
1869 !debug           enddo
1870 !debug       temp_mat(:,:)=czero
1871 !debug!      input matrix: temp_mat2
1872 !debug!      rotation matrix: gathermatlu
1873 !debug!      intermediate matrix: temp_mat
1874 !debug!      result matrix: temp_mat2
1875 !debug       call zgemm('n','c',tndim,tndim,tndim,cone,temp_mat2   ,tndim,&
1876 !debug&        temp_mat3,tndim,czero,temp_mat                ,tndim)
1877 !debug       call zgemm('n','n',tndim,tndim,tndim,cone,temp_mat3,tndim,&
1878 !debug&        temp_mat                   ,tndim,czero,temp_mat2,tndim)
1879 !debug!       call zgemm('n','c',tndim,tndim,tndim,cone,temp_mat2   ,tndim,&
1880 !debug!&        gathermatlu(iatom)%value,tndim,czero,temp_mat                ,tndim)
1881 !debug!       call zgemm('n','n',tndim,tndim,tndim,cone,gathermatlu(iatom)%value,tndim,&
1882 !debug!&        temp_mat                   ,tndim,czero,temp_mat2,tndim)
1883 !debug         write(std_out,*) "result"
1884 !debug         do im1=1,tndim
1885 !debug           write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')&
1886 !debug&           (temp_mat2(im1,im2),im2=1,tndim)
1887 !debug           call wrtout(std_out,message,'COLL')
1888 !debug         end do
1889 !debug       deallocate(temp_mat)
1890 !debug       deallocate(temp_mat3)
1891 !debug     endif ! lpawu
1892 
1893    enddo  ! iatom
1894 ! End loop over atoms
1895 ! ===========================
1896    do iatom=1,natom
1897      if(matlu(iatom)%lpawu.ne.-1) then
1898 !debug       deallocate(temp_mat2)
1899        ABI_DEALLOCATE(gathermatlu(iatom)%value)
1900      endif
1901    enddo
1902    ABI_DATATYPE_DEALLOCATE(gathermatlu)
1903  enddo ! isppol
1904 
1905  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

PARENTS

      m_datafordmft,m_green,m_oper,qmc_prep_ctqmc

CHILDREN

SOURCE

896 subroutine diff_matlu(char1,char2,matlu1,matlu2,natom,option,toldiff,ierr)
897 
898  use defs_basis
899  use m_paw_dmft, only : paw_dmft_type
900  use m_crystal, only : crystal_t
901  use m_io_tools,           only : flush_unit
902 
903 !This section has been created automatically by the script Abilint (TD).
904 !Do not modify the following lines by hand.
905 #undef ABI_FUNC
906 #define ABI_FUNC 'diff_matlu'
907 !End of the abilint section
908 
909  implicit none
910 
911 !Arguments ------------------------------------
912 !type
913  integer,intent(in) :: natom,option
914  type(matlu_type), intent(in) :: matlu1(natom),matlu2(natom)
915  character(len=*), intent(in) :: char1,char2
916  real(dp),intent(in) :: toldiff
917  integer,intent(out), optional :: ierr
918 
919 !local variables-------------------------------
920  integer :: iatom,idiff,ispinor,ispinor1,isppol,m1,m,lpawu,nspinor,nsppol
921  real(dp) :: matludiff
922  character(len=500) :: message
923 ! *********************************************************************
924  nsppol=matlu1(1)%nsppol
925  nspinor=matlu1(1)%nspinor
926 
927  matludiff=zero
928  idiff=0
929  do iatom = 1 , natom
930   lpawu=matlu1(iatom)%lpawu
931   if(lpawu/=-1) then
932    do isppol = 1 , nsppol
933     do ispinor = 1 , nspinor
934      do ispinor1 = 1, nspinor
935       do m1 = 1 , 2*lpawu+1
936        do m = 1 ,  2*lpawu+1
937         idiff=idiff+1
938         matludiff=matludiff+ &
939 &        sqrt( real(matlu1(iatom)%mat(m1,m,isppol,ispinor,ispinor1)      &
940 &            -      matlu2(iatom)%mat(m1,m,isppol,ispinor,ispinor1))**2  &
941 &            +aimag(matlu1(iatom)%mat(m1,m,isppol,ispinor,ispinor1)      &
942 &            -      matlu2(iatom)%mat(m1,m,isppol,ispinor,ispinor1))**2  )
943 !       write(std_out,*) m,m1,matlu1(iatom)%mat(m1,m,isppol,ispinor,ispinor1),matlu2(iatom)%mat(m1,m,isppol,ispinor,ispinor1),matludiff
944        enddo
945       enddo
946      end do ! ispinor1
947     end do ! ispinor
948    enddo ! isppol
949   endif ! lpawu/=1
950  enddo ! natom
951  matludiff=matludiff/float(idiff)
952  if(option==1.or.option==0) then
953   if( matludiff < toldiff ) then
954    write(message,'(5a,6x,3a,4x,e12.4,a,e12.4)') ch10,&
955 &   '   ** Differences between ',trim(char1),' and ',ch10,trim(char2),' are small enough:',&
956 &   ch10,matludiff,' is lower than',toldiff
957    call wrtout(std_out,message,'COLL')
958    if(present(ierr)) ierr=0
959   else
960    write(message,'(5a,3x,3a,3x,e12.4,a,e12.4)') ch10,&
961 &   'Differences between ',trim(char1),' and ',ch10,trim(char2),' is too large:',&
962 &   ch10,matludiff,' is larger than',toldiff
963    MSG_WARNING(message)
964 !   write(message,'(8a,4x,e12.4,a,e12.4)') ch10,"  Matrix for ",trim(char1)
965    write(message,'(a,3x,a)') ch10,trim(char1)
966    call wrtout(std_out,message,'COLL')
967    call print_matlu(matlu1,natom,prtopt=1,opt_diag=-1)
968    write(message,'(a,3x,a)') ch10,trim(char2)
969    call wrtout(std_out,message,'COLL')
970    call print_matlu(matlu2,natom,prtopt=1,opt_diag=-1)
971    if(option==1) then
972      call flush_unit(std_out)
973      MSG_ERROR("option==1, aborting now!")
974    end if
975    if(present(ierr)) ierr=-1
976   endif
977  endif
978 
979 end subroutine diff_matlu

m_matlu/fac_matlu [ Functions ]

[ Top ] [ m_matlu ] [ Functions ]

NAME

 fac_matlu

FUNCTION

 shift matlu matrix

COPYRIGHT

 Copyright (C) 2005-2018 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

PARENTS

      qmc_prep_ctqmc

CHILDREN

SOURCE

2828  subroutine fac_matlu(matlu,natom,fac)
2829  use defs_basis
2830 
2831 !This section has been created automatically by the script Abilint (TD).
2832 !Do not modify the following lines by hand.
2833 #undef ABI_FUNC
2834 #define ABI_FUNC 'fac_matlu'
2835 !End of the abilint section
2836 
2837  implicit none
2838 
2839 !Arguments ------------------------------------
2840 !scalars
2841  integer, intent(in) :: natom
2842 !arrays
2843  type(matlu_type),intent(inout) :: matlu(natom)
2844  complex(dpc),intent(in)      :: fac
2845 !Local variables-------------------------------
2846 !scalars
2847  integer :: iatom,im,im1,ispinor,ispinor1,isppol
2848  integer :: lpawu
2849 ! character(len=500) :: message
2850 !arrays
2851 !************************************************************************
2852  do iatom=1,natom
2853    lpawu=matlu(iatom)%lpawu
2854    if(lpawu.ne.-1) then
2855      do im=1,2*lpawu+1
2856        do im1=1,2*lpawu+1
2857          do isppol=1,matlu(1)%nsppol
2858            do ispinor=1,matlu(1)%nspinor
2859              do ispinor1=1,matlu(1)%nspinor
2860                matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)=&
2861 &                matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)*fac
2862              enddo
2863            enddo
2864          enddo
2865        enddo
2866      enddo
2867    endif
2868  enddo
2869 
2870  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-2018 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

PARENTS

      m_matlu,psichi_renormalization

CHILDREN

SOURCE

1320  subroutine gather_matlu(gloc,gathergloc,natom,option,prtopt)
1321  use defs_basis
1322  use defs_wvltypes
1323  use m_crystal, only : crystal_t
1324 
1325 !This section has been created automatically by the script Abilint (TD).
1326 !Do not modify the following lines by hand.
1327 #undef ABI_FUNC
1328 #define ABI_FUNC 'gather_matlu'
1329 !End of the abilint section
1330 
1331  implicit none
1332 
1333 ! type  matlus_type
1334 !  SEQUENCE
1335 !  complex(dpc), pointer :: mat(:,:)
1336 ! end type matlus_type
1337 
1338 !Arguments ------------------------------------
1339 !scalars
1340  integer,intent(in) :: natom,option,prtopt
1341  type(coeff2c_type), intent(inout) :: gathergloc(natom)
1342  type(matlu_type),intent(inout) :: gloc(natom)
1343 !Local variables-------------------------------
1344 !scalars
1345  integer :: iatom,im1,im2,ispinor,ispinor1,isppol,isppol1
1346  integer :: jc1,jc2,ml1,ml2,ndim,nspinor,nsppol,tndim
1347  character(len=500) :: message
1348 
1349 ! DBG_ENTER("COLL")
1350  nsppol=gloc(1)%nsppol
1351  nspinor=gloc(1)%nspinor
1352 
1353  do iatom=1,natom
1354    if(gloc(iatom)%lpawu.ne.-1) then
1355 !==-------------------------------------
1356 
1357      ndim=2*gloc(iatom)%lpawu+1
1358      tndim=nsppol*nspinor*ndim
1359 
1360 !== Put norm into array "gathergloc"
1361      jc1=0
1362      do isppol=1,nsppol
1363        do ispinor=1,nspinor
1364          do ml1=1,ndim
1365            jc1=jc1+1
1366            jc2=0
1367            do isppol1=1,nsppol
1368              do ispinor1=1,nspinor
1369                do ml2=1,ndim
1370                  jc2=jc2+1
1371                  if(option==1) then
1372                    if(isppol==isppol1) then
1373                      gathergloc(iatom)%value(jc1,jc2)=gloc(iatom)%mat(ml1,ml2,isppol,ispinor,ispinor1)
1374                    endif
1375                  else if(option==-1) then
1376                    if(isppol==isppol1) then
1377                      gloc(iatom)%mat(ml1,ml2,isppol,ispinor,ispinor1)=gathergloc(iatom)%value(jc1,jc2)
1378                    endif
1379                  endif
1380                enddo
1381              enddo ! ispinor1
1382            enddo ! isppol1
1383          enddo
1384        enddo !ispinor
1385      enddo ! isppol
1386    endif
1387  enddo ! iatom
1388  if(option==1.and.prtopt==3) then
1389    do iatom=1,natom
1390      if(gloc(iatom)%lpawu.ne.-1) then
1391        tndim=nsppol*nspinor*(2*gloc(iatom)%lpawu+1)
1392        write(message,'(2a,i5)') ch10,' (gathermatlu:) For atom', iatom
1393        call wrtout(std_out,message,'COLL')
1394        do im1=1,tndim
1395          write(message,'(12(1x,18(1x,"(",f9.3,",",f9.3,")")))')&
1396 &         (gathergloc(iatom)%value(im1,im2),im2=1,tndim)
1397          call wrtout(std_out,message,'COLL')
1398        end do
1399      endif
1400    enddo ! iatom
1401  else if(option==-1.and.prtopt==3) then
1402    call print_matlu(gloc,natom,prtopt)
1403  endif
1404 
1405 
1406 
1407 ! DBG_EXIT("COLL")
1408 
1409  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-2018 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

PARENTS

      qmc_prep_ctqmc

CHILDREN

SOURCE

2988  subroutine identity_matlu(matlu,natom)
2989  use defs_basis
2990 
2991 !This section has been created automatically by the script Abilint (TD).
2992 !Do not modify the following lines by hand.
2993 #undef ABI_FUNC
2994 #define ABI_FUNC 'identity_matlu'
2995 !End of the abilint section
2996 
2997  implicit none
2998 
2999 !Arguments ------------------------------------
3000 !scalars
3001  integer, intent(in) :: natom
3002 !arrays
3003  type(matlu_type),intent(inout) :: matlu(natom)
3004 !Local variables-------------------------------
3005 !scalars
3006  integer :: iatom,im,ispinor,isppol
3007  integer :: ndim
3008 ! character(len=500) :: message
3009 !arrays
3010 !************************************************************************
3011  ! not yet tested and used
3012  do iatom=1,natom
3013    if(matlu(iatom)%lpawu.ne.-1) then
3014      ndim=2*matlu(iatom)%lpawu+1
3015      do isppol=1,matlu(iatom)%nsppol
3016        do im=1,ndim
3017          do ispinor=1,matlu(iatom)%nspinor
3018            matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)= cone
3019          enddo ! ispinor
3020        enddo ! im
3021      enddo
3022    endif ! lpawu
3023  enddo ! iatom
3024  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

PARENTS

      m_datafordmft,hubbard_one,m_green,m_matlu,m_oper,qmc_prep_ctqmc

CHILDREN

SOURCE

146 subroutine init_matlu(natom,nspinor,nsppol,lpawu_natom,matlu)
147 
148  use defs_basis
149  use m_crystal, only : crystal_t
150 
151 !This section has been created automatically by the script Abilint (TD).
152 !Do not modify the following lines by hand.
153 #undef ABI_FUNC
154 #define ABI_FUNC 'init_matlu'
155 !End of the abilint section
156 
157  implicit none
158 
159 !Arguments ------------------------------------
160 !scalars
161  integer :: natom,nspinor,nsppol
162 !type
163  integer, intent(in) :: lpawu_natom(natom)
164  type(matlu_type), intent(inout) :: matlu(natom)
165 !Local variables ------------------------------------
166  integer :: iatom,lpawu
167 
168 !************************************************************************
169 
170 ! matlu%mband       = mband
171 ! matlu%dmftbandf   = dmftbandf
172 ! matlu%dmftbandi   = dmftbandi
173 ! matlu%nkpt        = nkpt
174 ! matlu%mbandc  = 0
175  matlu%nsppol      = nsppol
176  matlu%nspinor     = nspinor
177  do iatom=1,natom
178   lpawu=lpawu_natom(iatom)
179   matlu(iatom)%lpawu=lpawu
180   if(lpawu.ne.-1) then
181    ABI_ALLOCATE(matlu(iatom)%mat,(2*lpawu+1,2*lpawu+1,nsppol,nspinor,nspinor))
182    matlu(iatom)%mat=czero
183   else
184    ABI_ALLOCATE(matlu(iatom)%mat,(0,0,nsppol,nspinor,nspinor))
185   endif
186  enddo
187 
188 end subroutine init_matlu

m_matlu/inverse_matlu [ Functions ]

[ Top ] [ m_matlu ] [ Functions ]

NAME

 inverse_matlu

FUNCTION

 Inverse local quantity.

COPYRIGHT

 Copyright (C) 2005-2018 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

PARENTS

      m_oper

CHILDREN

SOURCE

812  subroutine inverse_matlu(matlu,natom,prtopt)
813  use defs_basis
814  use defs_wvltypes
815  use m_crystal, only : crystal_t
816 
817 !This section has been created automatically by the script Abilint (TD).
818 !Do not modify the following lines by hand.
819 #undef ABI_FUNC
820 #define ABI_FUNC 'inverse_matlu'
821 !End of the abilint section
822 
823  implicit none
824 
825 !Arguments ------------------------------------
826 !scalars
827  integer, intent(in) :: natom
828  integer, intent(in) :: prtopt
829 !arrays
830  type(matlu_type),intent(inout) :: matlu(natom)
831 !Local variables-------------------------------
832  integer :: iatom,tndim
833  integer :: nsppol,nspinor
834 !scalars
835  type(coeff2c_type),allocatable :: gathermatlu(:)
836  !************************************************************************
837 
838 
839  nspinor=matlu(1)%nspinor
840  nsppol=matlu(1)%nsppol
841  if(prtopt>0) then
842  endif
843  ABI_DATATYPE_ALLOCATE(gathermatlu,(natom))
844  do iatom=1,natom
845    if(matlu(iatom)%lpawu.ne.-1) then
846      tndim=nsppol*nspinor*(2*matlu(iatom)%lpawu+1)
847      ABI_ALLOCATE(gathermatlu(iatom)%value,(tndim,tndim))
848      gathermatlu(iatom)%value=czero
849    endif
850  enddo
851 
852  call gather_matlu(matlu,gathermatlu,natom,option=1,prtopt=1)
853  do iatom=1,natom
854    if(matlu(iatom)%lpawu.ne.-1) then
855      tndim=nsppol*nspinor*(2*matlu(iatom)%lpawu+1)
856      !call matcginv_dpc(gathermatlu(iatom)%value,tndim,tndim)
857      call xginv(gathermatlu(iatom)%value,tndim)
858    endif
859  enddo
860  call gather_matlu(matlu,gathermatlu,natom,option=-1,prtopt=1)
861 
862  do iatom=1,natom
863    if(matlu(iatom)%lpawu.ne.-1) then
864      ABI_DEALLOCATE(gathermatlu(iatom)%value)
865    endif
866  enddo
867  ABI_DATATYPE_DEALLOCATE(gathermatlu)
868  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-2018 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

PARENTS

CHILDREN

SOURCE

2603  subroutine ln_matlu(matlu1,natom)
2604  use defs_basis
2605  use defs_wvltypes
2606 
2607 !This section has been created automatically by the script Abilint (TD).
2608 !Do not modify the following lines by hand.
2609 #undef ABI_FUNC
2610 #define ABI_FUNC 'ln_matlu'
2611 !End of the abilint section
2612 
2613  implicit none
2614 
2615 !Arguments ------------------------------------
2616 !scalars
2617  integer, intent(in) :: natom
2618 !arrays
2619  type(matlu_type), intent(inout) :: matlu1(natom)
2620 !Local variables-------------------------------
2621 !scalars
2622  integer :: iatom,im,ispinor,isppol
2623  integer :: lpawu
2624  character(len=500) :: message
2625 !arrays
2626 !************************************************************************
2627  !call checkdiag_matlu(matlu1,natom,tol8)
2628  do iatom=1,natom
2629    lpawu=matlu1(iatom)%lpawu
2630    if(lpawu.ne.-1) then
2631      do isppol=1,matlu1(1)%nsppol
2632        do ispinor=1,matlu1(1)%nspinor
2633          do im=1,2*lpawu+1
2634            if( real(matlu1(iatom)%mat(im,im,isppol,ispinor,ispinor))<zero) then
2635              write(message,'(2a,2es13.5,a)') ch10," ln_matlu: PROBLEM " &
2636 &             , matlu1(iatom)%mat(im,im,isppol,ispinor,ispinor)
2637              MSG_ERROR(message)
2638            endif
2639            matlu1(iatom)%mat(im,im,isppol,ispinor,ispinor)= &
2640 &           log(matlu1(iatom)%mat(im,im,isppol,ispinor,ispinor))
2641          enddo ! im
2642        enddo ! ispinor
2643      enddo ! isppol
2644    endif ! lpawu
2645  enddo ! iatom
2646 
2647  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

 84  type, public :: matlu_type ! for each atom
 85 
 86   integer :: lpawu
 87   ! value of the angular momentum for correlated electrons
 88 
 89 !  integer :: natom
 90    ! number of atoms (given for each atom, not useful..could be changed)
 91 !
 92 !  integer :: mband
 93 !  ! Number of bands
 94 !
 95 !  integer :: mbandc
 96 !  ! Total number of bands in the Kohn-Sham Basis for PAW+DMFT
 97 !
 98 !  integer :: natpawu         ! Number of correlated atoms
 99 !
100 !  integer :: nkpt
101 !  ! Number of k-point in the IBZ.
102   character(len=12) :: whichmatlu
103   ! describe the type of local matrix computed (greenLDA, etc..)
104 !
105   integer :: nspinor
106   ! Number of spinorial component
107 !
108   integer :: nsppol
109   ! Number of polarization
110 
111   complex(dpc), allocatable :: mat(:,:,:,:,:)
112   ! local quantity
113 
114  end type matlu_type
115 
116 !----------------------------------------------------------------------
117 
118 
119 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

PARENTS

      compute_levels,m_datafordmft,hubbard_one,impurity_solve,m_green,m_matlu
      m_oper,m_self,psichi_renormalization,qmc_prep_ctqmc

CHILDREN

SOURCE

401 subroutine print_matlu(matlu,natom,prtopt,opt_diag,opt_ab_out,opt_exp,argout,compl)
402 
403  use defs_basis
404  use m_crystal, only : crystal_t
405 
406 !This section has been created automatically by the script Abilint (TD).
407 !Do not modify the following lines by hand.
408 #undef ABI_FUNC
409 #define ABI_FUNC 'print_matlu'
410 !End of the abilint section
411 
412  implicit none
413 
414 !Arguments ------------------------------------
415 !type
416  integer, intent(in):: natom
417  type(matlu_type),intent(in) :: matlu(natom)
418  integer, intent(in) :: prtopt
419  integer, optional, intent(in) :: opt_diag,opt_ab_out,opt_exp,argout,compl
420 !Local variables-------------------------------
421  integer :: iatom,ispinor,ispinor1,isppol,m1,m,lpawu,nspinor,nsppol,optdiag,optab_out,arg_out
422  character(len=500) :: message
423  character(len=4) :: mode_paral
424  type(matlu_type), allocatable :: matlu_tmp(:)
425  character(len=9),parameter :: dspinm(2,2)= RESHAPE((/"n        ","mx       ","my       ","mz       "/),(/2,2/))
426  logical :: testcmplx
427  real(dp) :: noccspin
428 ! *********************************************************************
429  mode_paral='COLL'
430  if(present(opt_diag)) then
431    optdiag=opt_diag
432  else
433    optdiag=0
434  endif
435  if(present(opt_ab_out)) then
436    optab_out=opt_ab_out
437  else
438    optab_out=0
439  endif
440  if(optab_out==0) then
441    arg_out=std_out
442  else
443    arg_out=ab_out
444  endif
445  if(present(argout)) then
446   arg_out=argout
447   mode_paral='PERS'
448  endif
449  nspinor=matlu(1)%nspinor
450  nsppol=matlu(1)%nsppol
451  testcmplx=(nspinor==2)
452  if(present(compl)) testcmplx=(nspinor==2).or.(compl==1)
453 
454  do iatom = 1 , natom
455    lpawu=matlu(iatom)%lpawu
456    if(lpawu/=-1) then
457      write(message,'(2a,i4)')  ch10,'   -------> For Correlated Atom', iatom
458      call wrtout(arg_out,  message,mode_paral)
459      do isppol = 1 , nsppol
460        if(present(opt_ab_out).and.nsppol==2) then
461          noccspin=zero
462          do m1=1, 2*lpawu +1
463            noccspin=noccspin+REAL(matlu(iatom)%mat(m1,m1,isppol,1,1))
464          enddo
465          !write(message,fmt='(7x,a,i3,a,f10.5)') ". Occ. for lpawu and for spin",isppol," =",noccspin
466          !call wrtout(arg_out, message,mode_paral)
467        endif
468      enddo
469 
470      do isppol = 1 , nsppol
471        if(nspinor == 1) then
472          write(message,'(a,10x,a,i3,i3)')  ch10,'-- polarization spin component',isppol
473          call wrtout(arg_out,  message,mode_paral)
474        endif
475        do ispinor = 1 , nspinor
476          do ispinor1 = 1, nspinor
477            if(nspinor == 2) then
478              write(message,'(a,10x,a,i3,i3)')  ch10,'-- spin components',ispinor,ispinor1
479              call wrtout(arg_out,  message,mode_paral)
480            endif
481            if(optdiag<=0) then
482              do m1=1, 2*lpawu+1
483                if(optdiag==0) then
484                  if(nspinor==1.and.abs(prtopt)>0)  then
485                    if(present(opt_exp)) then
486                      write(message,'(5x,20e24.14)') (REAL(matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1)
487 !                     call wrtout(arg_out,  message,mode_paral)
488 !                     write(message,'(5x,20e20.14)') (REAL(sqrt(matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1))),m=1,2*lpawu+1)
489 !                     call wrtout(arg_out,  message,mode_paral)
490 !                     write(message,'(5x,20e20.14)') (REAL(1.d0/sqrt(matlu(iatom)%mat(m,m,isppol,ispinor,ispinor1))),m=1,2*lpawu+1)
491                    else
492                      write(message,'(5x,20f10.5)') (REAL(matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1)
493                    endif
494                  endif
495                  if(testcmplx.and.abs(prtopt)>0) then
496                    if(present(opt_exp)) then
497                      if(opt_exp==2) then
498                        write(message,'(5x,14(2e18.10,1x))') ((matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1)
499                      else
500                        write(message,'(5x,14(2e14.4,2x))') ((matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1)
501                      endif
502                    else
503                      write(message,'(5x,14(2f9.5,2x))')((matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1)
504                    endif
505 !&                  write(message,'(5x,14(2f15.11,2x))')((matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1)
506                  endif
507                else if(optdiag==-1) then
508                  write(message,'(5x,14(2f10.5,2x))')((matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1)
509                endif
510                call wrtout(arg_out,  message,mode_paral)
511              enddo
512            elseif (optdiag>=1) then
513              if(nspinor==1.and.abs(prtopt)>0) &
514 &             write(message,'(5x,20f10.5)') (REAL(matlu(iatom)%mat(m,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1)
515              if(testcmplx.and.abs(prtopt)>0) &
516 &             write(message,'(5x,14(2f9.5,2x))')((matlu(iatom)%mat(m,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1)
517 !            write(std_out,'(5x,14(2f9.5,2x))')((matlu(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1)
518              call wrtout(arg_out,  message,mode_paral)
519            endif
520          end do ! ispinor1
521        end do ! ispinor
522        if(nspinor==2.and.prtopt>=5) then
523          ABI_DATATYPE_ALLOCATE(matlu_tmp,(0:natom))
524          call init_matlu(natom,nspinor,nsppol,matlu(1:natom)%lpawu,matlu_tmp)
525          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)
526          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)
527          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)
528          matlu_tmp(iatom)%mat(m1,m,isppol,2,1)= &
529 &           (matlu(iatom)%mat(m1,m,isppol,1,2)+matlu(iatom)%mat(m1,m,isppol,2,1))*cmplx(0_dp,1_dp,kind=dp)
530          do ispinor = 1 , nspinor
531            do ispinor1 = 1, nspinor
532              write(message,'(a,10x,a,a)')  ch10,'-- spin components',dspinm(ispinor,ispinor1)
533              call wrtout(arg_out,  message,mode_paral)
534              write(message,'(5x,14(2f9.5,2x))')((matlu_tmp(iatom)%mat(m1,m,isppol,ispinor,ispinor1)),m=1,2*lpawu+1)
535              call wrtout(arg_out,  message,mode_paral)
536            end do ! ispinor1
537          end do ! ispinor
538          call destroy_matlu(matlu_tmp,natom)
539          ABI_DATATYPE_DEALLOCATE(matlu_tmp)
540        endif
541      enddo ! isppol
542 !     if(nsppol==1.and.nspinor==1) then
543 !       write(message,'(a,10x,a,i3,a)')  ch10,'-- polarization spin component',isppol+1,' is identical'
544 !       call wrtout(arg_out,  message,mode_paral)
545 !     endif
546    endif ! lpawu/=1
547  enddo ! natom
548 
549 
550 end subroutine print_matlu

m_matlu/printplot_matlu [ Functions ]

[ Top ] [ m_matlu ] [ Functions ]

NAME

 printplot_matlu

FUNCTION

 shift matlu matrix

COPYRIGHT

 Copyright (C) 2005-2018 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

PARENTS

      qmc_prep_ctqmc

CHILDREN

SOURCE

2904  subroutine printplot_matlu(matlu,natom,freq,char1,units,imre)
2905  use defs_basis
2906  use m_fstrings,       only : int2char4
2907 
2908 !This section has been created automatically by the script Abilint (TD).
2909 !Do not modify the following lines by hand.
2910 #undef ABI_FUNC
2911 #define ABI_FUNC 'printplot_matlu'
2912 !End of the abilint section
2913 
2914  implicit none
2915 
2916 !Arguments ------------------------------------
2917 !scalars
2918  integer, intent(in) :: natom,units
2919  integer, optional, intent(in) :: imre
2920  real(dp), intent(in) :: freq
2921 !arrays
2922  type(matlu_type),intent(inout) :: matlu(natom)
2923  character(len=*), intent(in) :: char1
2924 !Local variables-------------------------------
2925 !scalars
2926  integer :: iatom,im,im1,ispinor,ispinor1,isppol
2927  integer :: lpawu,unitnb
2928  character(len=4) :: tag_at
2929  character(len=fnlen) :: tmpfil,tmpfilre,tmpfilim
2930 ! character(len=500) :: message
2931 !arrays
2932 !************************************************************************
2933  ! not yet tested and used
2934  do iatom=1,natom
2935    lpawu=matlu(iatom)%lpawu
2936    if(lpawu.ne.-1) then
2937      unitnb=units+iatom
2938      call int2char4(iatom,tag_at)
2939      if(present(imre)) then
2940        tmpfilre = trim(char1)//tag_at//"re"
2941        tmpfilim = trim(char1)//tag_at//"im"
2942        open (unit=unitnb+10,file=trim(tmpfilre),status='unknown',form='formatted')
2943        open (unit=unitnb+20,file=trim(tmpfilim),status='unknown',form='formatted')
2944        write(unitnb+10,'(400e26.16)') freq,(((((real(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)),&
2945        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)
2946        write(unitnb+20,'(400e26.16)') freq,(((((aimag(matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)),&
2947        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)
2948      else
2949        tmpfil = trim(char1)//tag_at
2950        open (unit=unitnb,file=trim(tmpfil),status='unknown',form='formatted')
2951        write(unitnb,'(400e26.16)') freq,(((((matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1),&
2952        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)
2953      endif
2954    endif
2955  enddo
2956  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-2018 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

PARENTS

      m_oper,qmc_prep_ctqmc

CHILDREN

SOURCE

2456  subroutine prod_matlu(matlu1,matlu2,matlu3,natom)
2457  use defs_basis
2458  use defs_wvltypes
2459  use m_crystal, only : crystal_t
2460 
2461 !This section has been created automatically by the script Abilint (TD).
2462 !Do not modify the following lines by hand.
2463 #undef ABI_FUNC
2464 #define ABI_FUNC 'prod_matlu'
2465 !End of the abilint section
2466 
2467  implicit none
2468 
2469 !Arguments ------------------------------------
2470 !scalars
2471  integer, intent(in) :: natom
2472 !arrays
2473  type(matlu_type), intent(in) :: matlu1(natom),matlu2(natom)
2474  type(matlu_type), intent(inout) :: matlu3(natom)
2475 !Local variables-------------------------------
2476 !scalars
2477  integer :: iatom,im1,im2,im3,ispinor1,ispinor2,ispinor3,isppol
2478  integer :: lpawu
2479 !arrays
2480 !************************************************************************
2481  call zero_matlu(matlu3,natom)
2482  do iatom=1,natom
2483    lpawu=matlu1(iatom)%lpawu
2484    if(lpawu.ne.-1) then
2485      do isppol=1,matlu1(1)%nsppol
2486        do ispinor1=1,matlu1(1)%nspinor
2487          do ispinor2=1,matlu1(1)%nspinor
2488            do ispinor3=1,matlu1(1)%nspinor
2489              do im1=1,2*lpawu+1
2490                do im2=1,2*lpawu+1
2491                  do im3=1,2*lpawu+1
2492                    matlu3(iatom)%mat(im1,im2,isppol,ispinor1,ispinor2)= &
2493 &                    matlu3(iatom)%mat(im1,im2,isppol,ispinor1,ispinor2)+ &
2494 &                    matlu1(iatom)%mat(im1,im3,isppol,ispinor1,ispinor3)*&
2495 &                    matlu2(iatom)%mat(im3,im2,isppol,ispinor3,ispinor2)
2496                  enddo ! im3
2497                enddo ! im2
2498              enddo ! im1
2499            enddo ! ispinor3
2500          enddo ! ispinor2
2501        enddo ! ispinor1
2502      enddo ! isppol
2503    endif ! lpawu
2504  enddo ! iatom
2505 
2506  end subroutine prod_matlu

m_matlu/rotate_matlu [ Functions ]

[ Top ] [ m_matlu ] [ Functions ]

NAME

 rotate_matlu

FUNCTION

 Rotate matlu matrix

COPYRIGHT

 Copyright (C) 2005-2018 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

PARENTS

      hubbard_one,qmc_prep_ctqmc

CHILDREN

SOURCE

1940  subroutine rotate_matlu(matlu,rot_mat,natom,prtopt,inverse)
1941  use defs_basis
1942  use defs_wvltypes
1943  use m_crystal, only : crystal_t
1944 
1945 !This section has been created automatically by the script Abilint (TD).
1946 !Do not modify the following lines by hand.
1947 #undef ABI_FUNC
1948 #define ABI_FUNC 'rotate_matlu'
1949 !End of the abilint section
1950 
1951  implicit none
1952 
1953 !Arguments ------------------------------------
1954 !scalars
1955  integer, intent(in) :: natom
1956  integer, intent(in) :: prtopt
1957  integer, intent(in) :: inverse
1958 !arrays
1959  type(matlu_type),intent(inout) :: matlu(natom)
1960  type(coeff2c_type),optional,intent(inout) :: rot_mat(natom,matlu(1)%nsppol)
1961 !Local variables-------------------------------
1962 !scalars
1963  integer :: iatom,im1,im2,isppol
1964  integer :: nsppol,nspinor,tndim
1965 !arrays
1966  type(coeff2c_type),allocatable :: gathermatlu(:)
1967 ! type(coeff2c_type),allocatable :: rot_mat_orig(:,:)
1968  type(coeff2c_type),allocatable :: rot_mat_orig(:)
1969  complex(dpc),allocatable :: temp_mat(:,:)
1970 !************************************************************************
1971  if(prtopt==1) then
1972  endif
1973  nsppol=matlu(1)%nsppol
1974  nspinor=matlu(1)%nspinor
1975 ! ABI_DATATYPE_ALLOCATE(rot_mat_orig,(natom,matlu(1)%nsppol))
1976 
1977  do isppol=1,nsppol
1978 
1979 ! ===========================
1980 ! Define gathermatlu and rot_mat_orig and allocate
1981 ! ===========================
1982    ABI_DATATYPE_ALLOCATE(rot_mat_orig,(natom))
1983    ABI_DATATYPE_ALLOCATE(gathermatlu,(natom))
1984    do iatom=1,natom
1985      if(matlu(iatom)%lpawu.ne.-1) then
1986        tndim=nspinor*(2*matlu(iatom)%lpawu+1)
1987        ABI_ALLOCATE(gathermatlu(iatom)%value,(tndim,tndim))
1988        gathermatlu(iatom)%value=czero
1989 !       ABI_ALLOCATE(rot_mat_orig(iatom,isppol)%value,(tndim,tndim))
1990 !       rot_mat_orig(iatom,isppol)%value(:,:)=rot_mat(iatom,isppol)%value(:,:)
1991        ABI_ALLOCATE(rot_mat_orig(iatom)%value,(tndim,tndim))
1992        rot_mat_orig(iatom)%value(:,:)=rot_mat(iatom,isppol)%value(:,:)
1993      endif
1994    enddo
1995    if(nsppol==1.and.nspinor==2) then
1996      call gather_matlu(matlu,gathermatlu,natom,option=1,prtopt=1)
1997    else if((nsppol==2.or.nsppol==1).and.nspinor==1) then
1998      do iatom=1,natom
1999        if(matlu(iatom)%lpawu.ne.-1) then
2000          do im1=1,tndim
2001            do im2=1,tndim
2002              gathermatlu(iatom)%value(im1,im2)=matlu(iatom)%mat(im1,im2,isppol,1,1)
2003            enddo
2004          enddo
2005        endif
2006      enddo
2007    endif
2008         ! write(std_out,*) "gathermatlu in rotate matlu"
2009         ! do im1=1,tndim
2010         !   write(message,'(12(1x,18(1x,"(",e17.10,",",e17.10,")")))')&
2011         !    (gathermatlu(1)%value(im1,im2),im2=1,tndim)
2012         !   call wrtout(std_out,message,'COLL')
2013         ! end do
2014 
2015 ! ===========================
2016 ! If necessary, invert rot_mat
2017 ! ===========================
2018    if(inverse==1) then
2019      do iatom=1,natom
2020        if(matlu(iatom)%lpawu.ne.-1) then
2021          tndim=nspinor*(2*matlu(iatom)%lpawu+1)
2022            do im1=1,tndim
2023              do im2=1,tndim
2024 !               rot_mat(iatom,isppol)%value(im1,im2)=conjg(rot_mat_orig(iatom,isppol)%value(im2,im1))
2025                rot_mat(iatom,isppol)%value(im1,im2)=conjg(rot_mat_orig(iatom)%value(im2,im1))
2026              enddo
2027            enddo
2028        endif ! lpawu
2029      enddo ! iatom
2030    endif
2031         ! write(std_out,*) "rot_mat_orig "
2032         ! do im1=1,tndim
2033         !   write(message,'(12(1x,18(1x,"(",e18.10,",",e18.10,")")))')&
2034         ! &   (rot_mat_orig(1)%value(im1,im2),im2=1,tndim)
2035         !   call wrtout(std_out,message,'COLL')
2036         ! end do
2037         ! write(std_out,*) "rot_mat "
2038         ! do im1=1,tndim
2039         !   write(message,'(12(1x,18(1x,"(",e18.10,",",e18.10,")")))')&
2040         ! &   (rot_mat(1,1)%value(im1,im2),im2=1,tndim)
2041         !   call wrtout(std_out,message,'COLL')
2042         ! end do
2043 
2044 ! ===========================
2045 ! Rotate
2046 ! ===========================
2047    ABI_ALLOCATE(temp_mat,(tndim,tndim))
2048    do iatom=1,natom
2049      if(matlu(iatom)%lpawu.ne.-1) then
2050        tndim=nspinor*(2*matlu(iatom)%lpawu+1)
2051        temp_mat(:,:)=czero
2052 !      input matrix: gathermatlu
2053 !      rotation matrix: rot_mat
2054 !      intermediate matrix: temp_mat
2055 !      result matrix: gathermatlu
2056        ! temp_mat = gathermatlu * conjg(rot_mat)
2057        call zgemm('n','c',tndim,tndim,tndim,cone,gathermatlu(iatom)%value   ,tndim,&
2058 &        rot_mat(iatom,isppol)%value,tndim,czero,temp_mat                ,tndim)
2059        ! gathermatlu = rot_mat * temp_mat = rot_mat * gathermatlu * conjg(rot_mat)
2060        call zgemm('n','n',tndim,tndim,tndim,cone,rot_mat(iatom,isppol)%value,tndim,&
2061 &        temp_mat                   ,tndim,czero,gathermatlu(iatom)%value,tndim)
2062      endif ! lpawu
2063    enddo ! iatom
2064    !do iatom=1,natom
2065    !  if(matlu(iatom)%lpawu.ne.-1) then
2066    !    write(std_out,*) "temp_mat in rotate_matlu 2"
2067    !    do im1=1,tndim
2068    !      write(message,'(12(1x,18(1x,"(",f17.10,",",f17.10,")")))')&
2069    !&       (temp_mat(im1,im2),im2=1,tndim)
2070    !      call wrtout(std_out,message,'COLL')
2071    !    end do
2072    !  endif
2073    !enddo
2074    !do iatom=1,natom
2075    !  if(matlu(iatom)%lpawu.ne.-1) then
2076    !    write(std_out,*) "gathermatlu in rotate_matlu 2"
2077    !    do im1=1,tndim
2078    !      write(message,'(12(1x,18(1x,"(",f17.10,",",f17.10,")")))')&
2079    !&       (gathermatlu(iatom)%value(im1,im2),im2=1,tndim)
2080    !      call wrtout(std_out,message,'COLL')
2081    !    end do
2082    !  endif
2083    !enddo
2084    ABI_DEALLOCATE(temp_mat)
2085      !MSG_ERROR("Aborting now")
2086 
2087 ! Choose inverse rotation: reconstruct correct rot_mat from rot_mat_orig
2088 ! ========================================================================
2089    if(inverse==1) then
2090      do iatom=1,natom
2091        if(matlu(iatom)%lpawu.ne.-1) then
2092          tndim=nspinor*(2*matlu(iatom)%lpawu+1)
2093            do im1=1,tndim
2094              do im2=1,tndim
2095 !               rot_mat(iatom,isppol)%value(im1,im2)=rot_mat_orig(iatom,isppol)%value(im1,im2)
2096                rot_mat(iatom,isppol)%value(im1,im2)=rot_mat_orig(iatom)%value(im1,im2)
2097              enddo
2098            enddo
2099        endif ! lpawu
2100      enddo ! iatom
2101    endif
2102 
2103 ! ===========================
2104 ! Put data into matlu(iatom)
2105 ! ===========================
2106    if(nsppol==1.and.nspinor==2) then
2107      call gather_matlu(matlu,gathermatlu,natom,option=-1,prtopt=1)
2108    else if((nsppol==2.or.nsppol==1).and.nspinor==1) then
2109      do iatom=1,natom
2110        if(matlu(iatom)%lpawu.ne.-1) then
2111          do im1=1,tndim
2112            do im2=1,tndim
2113              matlu(iatom)%mat(im1,im2,isppol,1,1)= gathermatlu(iatom)%value(im1,im2)
2114            enddo
2115          enddo
2116        endif
2117      enddo
2118    endif ! test nsppol/nspinor
2119 ! ===========================
2120 ! Deallocations
2121 ! ===========================
2122    do iatom=1,natom
2123      if(matlu(iatom)%lpawu.ne.-1) then
2124        ABI_DEALLOCATE(gathermatlu(iatom)%value)
2125 !       ABI_DEALLOCATE(rot_mat_orig(iatom,isppol)%value)
2126        ABI_DEALLOCATE(rot_mat_orig(iatom)%value)
2127      endif
2128    enddo
2129    ABI_DATATYPE_DEALLOCATE(gathermatlu)
2130    ABI_DATATYPE_DEALLOCATE(rot_mat_orig)
2131  enddo ! isppol
2132 
2133  end subroutine rotate_matlu

m_matlu/shift_matlu [ Functions ]

[ Top ] [ m_matlu ] [ Functions ]

NAME

 shift_matlu

FUNCTION

 shift matlu matrix

COPYRIGHT

 Copyright (C) 2005-2018 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

PARENTS

      m_green,m_self,qmc_prep_ctqmc

CHILDREN

SOURCE

2167  subroutine shift_matlu(matlu,natom,shift,signe)
2168  use defs_basis
2169  use defs_wvltypes
2170  use m_crystal, only : crystal_t
2171 
2172 !This section has been created automatically by the script Abilint (TD).
2173 !Do not modify the following lines by hand.
2174 #undef ABI_FUNC
2175 #define ABI_FUNC 'shift_matlu'
2176 !End of the abilint section
2177 
2178  implicit none
2179 
2180 !Arguments ------------------------------------
2181 !scalars
2182  integer, intent(in) :: natom
2183 !arrays
2184  type(matlu_type),intent(inout) :: matlu(natom)
2185  complex(dpc),intent(in) :: shift(natom)
2186  integer, optional,intent(in) :: signe
2187 !Local variables-------------------------------
2188 !scalars
2189  integer :: iatom,im,ispinor,isppol
2190  integer :: lpawu,signe_used
2191 ! character(len=500) :: message
2192 !arrays
2193 !************************************************************************
2194  signe_used=1
2195  if(present(signe)) then
2196    if(signe==-1) signe_used=-1
2197  endif
2198  do iatom=1,natom
2199    lpawu=matlu(iatom)%lpawu
2200    if(lpawu.ne.-1) then
2201      do im=1,2*lpawu+1
2202        do isppol=1,matlu(1)%nsppol
2203          do ispinor=1,matlu(1)%nspinor
2204            matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)=&
2205 &           matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)-signe_used*shift(iatom)
2206          enddo
2207        enddo
2208      enddo
2209    endif
2210  enddo
2211 
2212  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-2018 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

PARENTS

      qmc_prep_ctqmc

CHILDREN

SOURCE

2678  subroutine slm2ylm_matlu(matlu,natom,option,optprt)
2679  use defs_basis
2680  use defs_wvltypes
2681 
2682 !This section has been created automatically by the script Abilint (TD).
2683 !Do not modify the following lines by hand.
2684 #undef ABI_FUNC
2685 #define ABI_FUNC 'slm2ylm_matlu'
2686 !End of the abilint section
2687 
2688  implicit none
2689 
2690 !Arguments ------------------------------------
2691 !scalars
2692  integer, intent(in) :: natom,option,optprt
2693 !arrays
2694  type(matlu_type), intent(inout) :: matlu(natom)
2695 !Local variables-------------------------------
2696 !scalars
2697  integer :: iatom,im,ispinor,isppol,ispinor2
2698  integer :: lpawu,ll,mm,jm,ii,jj,im1,im2
2699  character(len=500) :: message
2700  real(dp) :: onem
2701  complex(dpc),allocatable :: slm2ylm(:,:)
2702  complex(dpc),allocatable :: mat_inp_c(:,:)
2703  complex(dpc),allocatable :: mat_out_c(:,:)
2704  complex(dpc) :: tmp2
2705  real(dp),parameter :: invsqrt2=one/sqrt2
2706 !arrays
2707 !************************************************************************
2708 
2709  do iatom=1,natom
2710    lpawu=matlu(iatom)%lpawu
2711    if(lpawu.ne.-1) then
2712      ll=lpawu
2713      ABI_ALLOCATE(slm2ylm,(2*ll+1,2*ll+1))
2714      slm2ylm=czero
2715      do im=1,2*ll+1
2716        mm=im-ll-1;jm=-mm+ll+1
2717        onem=dble((-1)**mm)
2718        if (mm> 0) then
2719          slm2ylm(im,im)= cmplx(onem*invsqrt2,zero,kind=dp)
2720          slm2ylm(jm,im)= cmplx(invsqrt2,     zero,kind=dp)
2721        end if
2722        if (mm==0) then
2723          slm2ylm(im,im)=cone
2724        end if
2725        if (mm< 0) then
2726          slm2ylm(im,im)= cmplx(zero,     invsqrt2,kind=dp)
2727          slm2ylm(jm,im)=-cmplx(zero,onem*invsqrt2,kind=dp)
2728        end if
2729      end do
2730      if(optprt>2) then
2731        write(message,'(2a)') ch10,"SLM2YLM matrix"
2732        call wrtout(std_out,message,'COLL')
2733        do im1=1,ll*2+1
2734          write(message,'(12(1x,9(1x,"(",f9.5,",",f9.5,")")))')&
2735 &         (slm2ylm(im1,im2),im2=1,ll*2+1)
2736          call wrtout(std_out,message,'COLL')
2737        end do
2738      endif
2739      do isppol=1,matlu(1)%nsppol
2740        do ispinor=1,matlu(1)%nspinor
2741          do ispinor2=1,matlu(1)%nspinor
2742            ABI_ALLOCATE(mat_out_c,(2*ll+1,2*ll+1))
2743            ABI_ALLOCATE(mat_inp_c,(2*ll+1,2*ll+1))
2744            mat_inp_c(:,:) = matlu(iatom)%mat(:,:,isppol,ispinor,ispinor2)
2745            mat_out_c=czero
2746 
2747            if(optprt>2) then
2748              write(message,'(2a)') ch10,"SLM input matrix"
2749              call wrtout(std_out,message,'COLL')
2750              do im1=1,ll*2+1
2751                write(message,'(12(1x,9(1x,"(",f9.5,",",f9.5,")")))')&
2752 &               (mat_inp_c(im1,im2),im2=1,ll*2+1)
2753                call wrtout(std_out,message,'COLL')
2754              end do
2755            endif
2756 
2757            do jm=1,2*ll+1
2758              do im=1,2*ll+1
2759                tmp2=czero
2760                do ii=1,2*ll+1
2761                  do jj=1,2*ll+1
2762                    if(option==1) then
2763                      tmp2=tmp2+mat_inp_c(ii,jj)*(slm2ylm(im,ii))*CONJG(slm2ylm(jm,jj))
2764                    else if(option==2) then
2765                      tmp2=tmp2+mat_inp_c(ii,jj)*CONJG(slm2ylm(ii,im))*(slm2ylm(jj,jm))
2766                    end if
2767                  end do
2768                end do
2769                mat_out_c(im,jm)=tmp2
2770              end do
2771            end do
2772 
2773            if(optprt>2) then
2774              write(message,'(2a)') ch10,"YLM output matrix"
2775              call wrtout(std_out,message,'COLL')
2776              do im1=1,ll*2+1
2777                write(message,'(12(1x,9(1x,"(",f9.5,",",f9.5,")")))')&
2778       &         (mat_out_c(im1,im2),im2=1,ll*2+1)
2779                call wrtout(std_out,message,'COLL')
2780              end do
2781            endif
2782 
2783            matlu(iatom)%mat(:,:,isppol,ispinor,ispinor2)=mat_out_c(:,:)
2784            ABI_DEALLOCATE(mat_out_c)
2785            ABI_DEALLOCATE(mat_inp_c)
2786          enddo ! im
2787        enddo ! ispinor
2788      enddo ! isppol
2789    endif ! lpawu
2790  enddo ! iatom
2791 
2792  ABI_DEALLOCATE(slm2ylm)
2793 
2794  end subroutine slm2ylm_matlu

m_matlu/sym_matlu [ Functions ]

[ Top ] [ m_matlu ] [ Functions ]

NAME

 sym_matlu

FUNCTION

 Symetrise local quantity.

COPYRIGHT

 Copyright (C) 2005-2018 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

OUTPUT

  gloc(natom) <type(matlu_type)>= density matrix symetrized in the local orbital basis and related variables

SIDE EFFECTS

NOTES

PARENTS

      compute_levels,m_datafordmft,hybridization_asymptotic_coefficient,m_green
      psichi_renormalization,qmc_prep_ctqmc

CHILDREN

SOURCE

586  subroutine sym_matlu(cryst_struc,gloc,pawang)
587 
588  use defs_basis
589 ! use defs_wvltypes
590  use m_pawang, only : pawang_type
591  use m_crystal, only : crystal_t
592 
593 !This section has been created automatically by the script Abilint (TD).
594 !Do not modify the following lines by hand.
595 #undef ABI_FUNC
596 #define ABI_FUNC 'sym_matlu'
597 !End of the abilint section
598 
599  implicit none
600 
601 !Arguments ------------------------------------
602 !scalars
603  type(crystal_t),intent(in) :: cryst_struc
604  type(pawang_type),intent(in) :: pawang
605 !arrays
606  type(matlu_type),intent(inout) :: gloc(cryst_struc%natom)
607 !Local variables-------------------------------
608 !scalars
609  integer :: at_indx,iatom,irot,ispinor,ispinor1,isppol,lpawu,m1,m2,m3,m4,mu
610  integer :: natom,ndim,nsppol,nspinor,nu
611  complex(dpc) :: sumrho,summag(3),rotmag(3),ci
612  real(dp) :: zarot2
613 !arrays
614 ! complex(dpc),allocatable :: glocnm(:,:,:,:,:)
615  type(matlu_type),allocatable :: glocnm(:)
616 ! complex(dpc),allocatable :: glocnms(:,:,:,:,:)
617  type(matlu_type),allocatable :: glocnms(:)
618  type(matlu_type),allocatable :: glocsym(:)
619  real(dp),allocatable :: symrec_cart(:,:,:)
620 
621 ! DBG_ENTER("COLL")
622 
623  ci=cone
624  nspinor=gloc(1)%nspinor
625  nsppol=gloc(1)%nsppol
626  natom=cryst_struc%natom
627 
628  ABI_DATATYPE_ALLOCATE(glocnm,(natom))
629  ABI_DATATYPE_ALLOCATE(glocnms,(natom))
630  ABI_DATATYPE_ALLOCATE(glocsym,(natom))
631  call init_matlu(natom,nspinor,nsppol,gloc(1:natom)%lpawu,glocnm)
632  call init_matlu(natom,nspinor,nsppol,gloc(1:natom)%lpawu,glocnms)
633  call init_matlu(natom,nspinor,nsppol,gloc(1:natom)%lpawu,glocsym)
634 
635 
636 !=========  Case nspinor ==1 ========================
637 
638  if (nspinor==1) then
639   ispinor=1
640   ispinor1=1
641   do iatom=1,cryst_struc%natom
642    do isppol=1,nsppol
643     if(gloc(iatom)%lpawu/=-1) then
644      lpawu=gloc(iatom)%lpawu
645      do m1=1, 2*lpawu+1
646       do m2=1, 2*lpawu+1
647        do irot=1,cryst_struc%nsym
648         at_indx=cryst_struc%indsym(4,irot,iatom)
649         do m3=1, 2*lpawu+1
650          do m4=1, 2*lpawu+1
651           zarot2=pawang%zarot(m3,m1,lpawu+1,irot)*pawang%zarot(m4,m2,lpawu+1,irot)
652           glocsym(iatom)%mat(m1,m2,isppol,ispinor,ispinor1)=&
653 &          glocsym(iatom)%mat(m1,m2,isppol,ispinor,ispinor1)&
654 &          +gloc(at_indx)%mat(m3,m4,isppol,ispinor,ispinor1)*zarot2
655          end do  ! m3
656         end do  ! m4
657        end do  ! irot
658        glocsym(iatom)%mat(m1,m2,isppol,ispinor,ispinor1)=&
659 &       glocsym(iatom)%mat(m1,m2,isppol,ispinor,ispinor1)/real(cryst_struc%nsym,kind=dp)
660       end do ! m2
661      end do ! m1
662     endif ! lpawu/=-1
663    end do ! isppol
664   end do ! iatom
665 !==  Put glocsym into gloc
666   do iatom=1,cryst_struc%natom
667     if(gloc(iatom)%lpawu/=-1) then
668       gloc(iatom)%mat=glocsym(iatom)%mat
669 !      gloc(iatom)%mat(:,:,1,:,:)=(glocsym(iatom)%mat(:,:,1,:,:) &
670 !&      + glocsym(iatom)%mat(:,:,2,:,:))/two
671 !      gloc(iatom)%mat(:,:,2,:,:)= gloc(iatom)%mat(:,:,1,:,:)
672 !      write(std_out,*) "WARNING: SYM non mag"
673 !      write(ab_out,*) "WARNING: SYM non mag"
674     endif
675   end do ! iatom
676 
677 !=========  Case nspinor ==2 ========================
678 
679  else if (nspinor==2) then
680 
681 !== Allocate temporary arrays
682   do iatom=1,cryst_struc%natom
683    if(gloc(iatom)%lpawu/=-1) then
684     ndim=2*gloc(iatom)%lpawu+1
685     ABI_DEALLOCATE(glocnm(iatom)%mat)
686     ABI_DEALLOCATE(glocnms(iatom)%mat)
687     ABI_DEALLOCATE(glocsym(iatom)%mat)
688     ABI_ALLOCATE(glocnm(iatom)%mat,(ndim,ndim,nsppol,4,1))
689     ABI_ALLOCATE(glocnms(iatom)%mat,(ndim,ndim,nsppol,4,1))
690     ABI_ALLOCATE(glocsym(iatom)%mat,(ndim,ndim,nsppol,2,2))
691    endif
692   enddo
693   ABI_ALLOCATE(symrec_cart,(3,3,cryst_struc%nsym))
694 
695 !==  Compute symrec_cart
696   do irot=1,cryst_struc%nsym
697    call symredcart(cryst_struc%gprimd,cryst_struc%rprimd,symrec_cart(:,:,irot),cryst_struc%symrec(:,:,irot))
698   end do
699 
700 !==  Compute density matrix in density and magnetization representation
701   call chg_repr_matlu(gloc,glocnm,cryst_struc%natom,option=1,prtopt=1)
702 
703 !==  Do the sum over symetrized density matrix (in n,m repr)
704   isppol=1
705   do iatom=1,cryst_struc%natom
706    if(gloc(iatom)%lpawu/=-1) then
707     lpawu=gloc(iatom)%lpawu
708     ndim=2*gloc(iatom)%lpawu+1
709     do m1=1, 2*lpawu+1
710      do m2=1, 2*lpawu+1
711       sumrho=czero
712       rotmag=czero
713       do irot=1,cryst_struc%nsym
714        summag=czero
715        at_indx=cryst_struc%indsym(4,irot,iatom)
716        do m3=1, 2*lpawu+1
717         do m4=1, 2*lpawu+1
718          zarot2=pawang%zarot(m3,m2,lpawu+1,irot)*pawang%zarot(m4,m1,lpawu+1,irot)
719          sumrho=sumrho +  glocnm(at_indx)%mat(m4,m3,isppol,1,1)  * zarot2
720          do mu=1,3
721           summag(mu)=summag(mu) + glocnm(at_indx)%mat(m4,m3,isppol,mu+1,1) * zarot2
722          enddo
723         end do ! m3
724        end do !m4
725 
726 !       ==  special case of magnetization
727        do nu=1,3
728         do mu=1,3
729          rotmag(mu)=rotmag(mu)+symrec_cart(mu,nu,irot)*summag(nu)
730         end do
731        end do
732 !      write(std_out,'(a,3i4,2x,3(2f10.5,2x))') "rotmag",irot,m1,m2,(rotmag(mu),mu=1,3)
733       end do ! irot
734 
735 !       ==  Normalizes sum
736       sumrho=sumrho/real(cryst_struc%nsym,kind=dp)
737 !        sumrho=glocnm(isppol,1,iatom,m1,m2) ! test without sym
738       glocnms(iatom)%mat(m1,m2,isppol,1,1)=sumrho
739       do mu=1,3
740        rotmag(mu)=rotmag(mu)/real(cryst_struc%nsym,kind=dp)
741 !          rotmag(mu)=glocnm(isppol,mu+1,iatom,m1,m2) ! test without sym
742        glocnms(iatom)%mat(m1,m2,isppol,mu+1,1)=rotmag(mu)
743       enddo
744      end do  ! m2
745     end do ! m1
746    endif ! lpawu/=-1
747   end do ! iatom
748 
749 !==  Compute back density matrix in upup dndn updn dnup representation
750   call chg_repr_matlu(glocsym,glocnms,cryst_struc%natom,option=-1,prtopt=1)
751 
752 !==  Put glocsym into gloc
753   do iatom=1,cryst_struc%natom
754     if(gloc(iatom)%lpawu/=-1) then
755       gloc(iatom)%mat=glocsym(iatom)%mat
756 !      gloc(iatom)%mat(:,:,1,:,:)=(glocsym(iatom)%mat(:,:,1,:,:) &
757 !&      + glocsym(iatom)%mat(:,:,2,:,:))/two
758 !      gloc(iatom)%mat(:,:,2,:,:)= gloc(iatom)%mat(:,:,1,:,:)
759 !      write(std_out,*) "WARNING: SYM non mag"
760 !      write(ab_out,*) "WARNING: SYM non mag"
761     endif
762   end do ! iatom
763 
764   ABI_DEALLOCATE(symrec_cart)
765  endif
766 
767  call destroy_matlu(glocnm,cryst_struc%natom)
768  call destroy_matlu(glocnms,cryst_struc%natom)
769  call destroy_matlu(glocsym,cryst_struc%natom)
770  ABI_DATATYPE_DEALLOCATE(glocnm)
771  ABI_DATATYPE_DEALLOCATE(glocnms)
772  ABI_DATATYPE_DEALLOCATE(glocsym)
773 !==============end of nspinor==2 case ===========
774 
775 
776 ! DBG_EXIT("COLL")
777 
778  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.

PARENTS

      m_green,m_oper

CHILDREN

SOURCE

1201  subroutine trace_matlu(matlu,natom,trace_loc,itau)
1202 
1203  use defs_basis
1204 
1205 !This section has been created automatically by the script Abilint (TD).
1206 !Do not modify the following lines by hand.
1207 #undef ABI_FUNC
1208 #define ABI_FUNC 'trace_matlu'
1209 !End of the abilint section
1210 
1211  implicit none
1212 
1213 !Arguments ------------------------------------
1214 !type
1215  integer, intent(in) :: natom
1216  type(matlu_type), intent(in) :: matlu(natom)
1217  real(dp),intent(inout),target,optional :: trace_loc(natom,matlu(1)%nsppol+1)
1218  integer, intent(in),optional :: itau
1219 
1220 !local variables-------------------------------
1221  integer :: iatom,isppol,ispinor,m,lpawu
1222  integer :: nsppol,nspinor
1223  real(dp), pointer :: traceloc(:,:)=>null()
1224  character(len=500) :: message
1225 ! *********************************************************************
1226  nsppol=matlu(1)%nsppol
1227  nspinor=matlu(1)%nspinor
1228  if(present(trace_loc)) then
1229    traceloc=>trace_loc
1230  else
1231    ABI_ALLOCATE(traceloc,(natom,matlu(1)%nsppol+1))
1232  endif
1233 
1234  traceloc=zero
1235  do iatom = 1 , natom
1236    lpawu=matlu(iatom)%lpawu
1237    if(lpawu/=-1) then
1238       write(message,'(2a,i4)')  ch10,'   -------> For Correlated Atom', iatom
1239      if(.not.present(itau)) then
1240        call wrtout(std_out,  message,'COLL')
1241      end if
1242      if(present(itau)) then
1243        if (itau>0) then
1244          call wrtout(std_out,  message,'COLL')
1245        end if
1246      endif
1247      do isppol = 1 , nsppol
1248        do ispinor = 1 , nspinor
1249          do m = 1 ,  2*lpawu+1
1250            traceloc(iatom,isppol)=traceloc(iatom,isppol)+&
1251 &           matlu(iatom)%mat(m,m,isppol,ispinor,ispinor)
1252          enddo
1253        enddo
1254       traceloc(iatom,nsppol+1)=traceloc(iatom,nsppol+1)+traceloc(iatom,isppol)
1255      enddo
1256      if(nsppol==1.and.nspinor==1)  traceloc(iatom,nsppol+1)=traceloc(iatom,nsppol+1)*two
1257      if(.not.present(itau)) then
1258        write(message,'(8x,a,f12.6)')   'Nb of Corr. elec. from G(w) is:'&
1259 &       ,traceloc(iatom,nsppol+1)
1260        call wrtout(std_out,  message,'COLL')
1261      endif
1262      if(present(itau)) then
1263        if(itau==1) then
1264          write(message,'(8x,a,f12.6)')   'Nb of Corr. elec. from G(tau) is:'&
1265 &         ,traceloc(iatom,nsppol+1)
1266          call wrtout(std_out,  message,'COLL')
1267        else if(itau==-1) then
1268          write(message,'(8x,a,f12.6)')   'Nb: Sum of the values of G0(tau=0-) is:'&
1269 &       ,traceloc(iatom,nsppol+1)
1270          call wrtout(std_out,  message,'COLL')
1271        else if(itau==4) then
1272          write(message,'(8x,a,f12.6)')   'Trace of matlu matrix is:'&
1273 &         ,traceloc(iatom,nsppol+1)
1274          call wrtout(std_out,  message,'COLL')
1275        endif
1276      endif
1277    endif
1278  enddo
1279  if(.not.present(trace_loc)) then
1280   ABI_DEALLOCATE(traceloc)
1281   traceloc => null()
1282  endif
1283 
1284  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

PARENTS

      m_green,m_matlu

CHILDREN

SOURCE

212 subroutine zero_matlu(matlu,natom)
213 
214  use defs_basis
215 
216 !This section has been created automatically by the script Abilint (TD).
217 !Do not modify the following lines by hand.
218 #undef ABI_FUNC
219 #define ABI_FUNC 'zero_matlu'
220 !End of the abilint section
221 
222  implicit none
223 
224 !Arguments ------------------------------------
225 !scalars
226  integer, intent(in) :: natom
227  type(matlu_type),intent(inout) :: matlu(natom)
228 !Local variables-------------------------------
229  integer :: iatom
230 
231 !*********************************************************************
232 
233  do iatom=1,natom
234   matlu(iatom)%mat=czero
235  enddo
236 
237 
238 end subroutine zero_matlu