TABLE OF CONTENTS


ABINIT/m_melemts [ Modules ]

[ Top ] [ Modules ]

NAME

  m_melemts

FUNCTION

  This module defines an object used as database to store matrix
  elements of several potentials and operators between two Bloch states.
  These values are used in the GW part of abinit to evaluate QP energies
  using the perturbative approach.

COPYRIGHT

 Copyright (C) 2008-2018 ABINIT group (MG)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

NOTES

  * This module is supposed to be used only in the GW part to facilitate
    further developments. The object might change in the future. Thus
    contact Matteo Giantomassi if you wish to use this piece of code
    for your developments. In particular we might decide to switch
    to ragged arrays Mels(nkibz,nsppol*nspinor**2)%data

  * Routines tagged with "@type_name" are tightly connected to the definition of the data type.
    Tightly connected means that the proper functioning of the implementation relies on the
    assumption that the tagged procedure is consistent with the type declaration.
    Every time a developer changes the structure "type_name" adding new entries, he/she has to make sure
    that all the tightly connected routines are changed accordingly to accommodate the modification of the data type.
    Typical examples of tightly connected routines are creation, destruction or reset methods.

TODO

  This module can be moved to a higher level directory.

PARENTS

CHILDREN

SOURCE

41 #if defined HAVE_CONFIG_H
42 #include "config.h"
43 #endif
44 
45 #include "abi_common.h"
46 
47 MODULE m_melemts
48 
49  use defs_basis
50  use m_errors
51  use m_xmpi
52  use m_abicore
53 
54  use m_fstrings,       only : tolower, sjoin
55  use m_numeric_tools,  only : print_arr
56 
57  implicit none
58 
59  private

m_melemts/melements_free [ Functions ]

[ Top ] [ m_melemts ] [ Functions ]

NAME

 melements_free

FUNCTION

  Free all dynamic memory of the database

INPUTS

  Mels<melements_t>=The database to be freed

OUTPUT

  See side effects

PARENTS

      sigma

CHILDREN

      my_select_melements

SOURCE

325 subroutine melements_free(Mels)
326 
327 
328 !This section has been created automatically by the script Abilint (TD).
329 !Do not modify the following lines by hand.
330 #undef ABI_FUNC
331 #define ABI_FUNC 'melements_free'
332 !End of the abilint section
333 
334  implicit none
335 
336 !Arguments ------------------------------------
337  type(melements_t),intent(inout) :: Mels
338 
339 ! *************************************************************************
340 
341  DBG_ENTER("COLL")
342 
343  ! @melements_t
344 
345 !integer arrays
346  if (allocated(Mels%bands_idx)) then
347    ABI_FREE(Mels%bands_idx)
348  end if
349  if (allocated(Mels%iscalc))   then
350    ABI_FREE(Mels%iscalc)
351  end if
352 
353 !real arrays
354  if (allocated(Mels%kibz)) then
355    ABI_FREE(Mels%kibz)
356  end if
357 
358 !complex arrays
359  if (allocated(Mels%hbare))   then
360    ABI_FREE(Mels%hbare)
361  end if
362  if (allocated(Mels%sxcore)) then
363    ABI_FREE(Mels%sxcore)
364  end if
365  if (allocated(Mels%vhartree))   then
366    ABI_FREE(Mels%vhartree)
367  end if
368  if (allocated(Mels%vlexx)) then
369    ABI_FREE(Mels%vlexx)
370  end if
371  if (allocated(Mels%vu)) then
372    ABI_FREE(Mels%vu)
373  end if
374  if (allocated(Mels%vxc)) then
375    ABI_FREE(Mels%vxc)
376  end if
377  if (allocated(Mels%vxcval)) then
378    ABI_FREE(Mels%vxcval)
379  end if
380  if (allocated(Mels%vxcval_hybrid)) then
381    ABI_FREE(Mels%vxcval_hybrid)
382  end if
383 
384  ! * Reset all has_* flags.
385  call melflags_reset(Mels%flags)
386 
387  DBG_EXIT("COLL")
388 
389 end subroutine melements_free

m_melemts/melements_herm [ Functions ]

[ Top ] [ m_melemts ] [ Functions ]

NAME

 melements_herm

FUNCTION

  Reconstruc the lower triangle of all calculated arrays.
  Assuming Hermitian operator. Works both for collinear and  non-collinear cases.

INPUTS

  Mels=The database
  [aname]=The name of the array to be symmetrized, by default
    all calculated arrays are filled.

SIDE EFFECTS

  All arrays whose flag is 2, are filled assuming an Hermitian operator.

PARENTS

      calc_vhxc_me

CHILDREN

      my_select_melements

SOURCE

633 subroutine melements_herm(Mels,aname)
634 
635 
636 !This section has been created automatically by the script Abilint (TD).
637 !Do not modify the following lines by hand.
638 #undef ABI_FUNC
639 #define ABI_FUNC 'melements_herm'
640 !End of the abilint section
641 
642  implicit none
643 
644 !Arguments ------------------------------------
645 !scalars
646  character(len=*),optional,intent(in) :: aname
647  type(melements_t),intent(inout) :: Mels
648 
649 !Local variables-------------------------------
650  integer :: is,ik,ib,jb,iab,iab_tr,iname
651  integer,pointer :: flag_p
652  character(len=NAMELEN) :: key
653 !arrays
654  integer,parameter :: trsp_idx(2:4) = [2,4,3]
655  complex(dpc),ABI_CONTIGUOUS pointer :: arr_p(:,:,:,:)
656 
657 ! *************************************************************************
658 
659  ! === Symmetrize matrix elements ===
660  ! * In the collinear case, generate the lower triangle by just doing a complex conjugate.
661  ! * In the noncollinear case do also a transposition since A_{12}^{ab} = A_{21}^{ba}^*
662  !   2-->2, 3-->4, 4-->3
663  !
664  do iname=1,NNAMES
665    key = ANAMES(iname)
666    if (PRESENT(aname)) then
667      if (key /= aname) CYCLE
668    end if
669 
670    call my_select_melements(Mels,key,flag_p,arr_p)
671 
672    if (flag_p>0) then
673      do ik=1,Mels%nkibz
674        do is=1,Mels%nsppol
675 
676          do jb=Mels%bmin,Mels%bmax
677            do ib=Mels%bmin,jb ! Upper triangle
678 
679              if (ib/=jb) then
680                arr_p(jb,ib,ik,is)=CONJG(arr_p(ib,jb,ik,is))
681                if (Mels%nspinor==2) then
682                  do iab=2,4
683                    iab_tr=trsp_idx(iab)
684                    arr_p(jb,ib,ik,iab)=CONJG(arr_p(ib,jb,ik,iab_tr))
685                  end do
686                end if
687              else ! For ib==jb force real-valued
688                arr_p(jb,ib,ik,is)=half*(arr_p(jb,ib,ik,is)+CONJG(arr_p(jb,ib,ik,is)))
689                if (Mels%nspinor==2) arr_p(jb,ib,ik,2)=half*(arr_p(ib,jb,ik,2)+CONJG(arr_p(ib,jb,ik,2)))
690              end if
691 
692            end do !ib
693          end do !jb
694 
695        end do !is
696      end do !ik
697    end if
698 
699  end do !inames
700 
701 end subroutine melements_herm

m_melemts/melements_init [ Functions ]

[ Top ] [ m_melemts ] [ Functions ]

NAME

 melements_init

FUNCTION

  Initialize the database, allocate arrays according to
  Mflags_in, zeroing the content of the allocated arrays.

INPUTS

  nspinor=Number of spinor components
  nsppol=Number of independent spin polarizations.
  nspden=Number of spin-density components
  b1,b2=min and Max band index over spin and k-points.

OUTPUT

  Mels=The initialized database with dimensions and allocated memory.

PARENTS

      calc_vhxc_me

CHILDREN

      my_select_melements

SOURCE

497 subroutine melements_init(Mels,Mflags_in,nsppol,nspden,nspinor,nkibz,kibz,bands_idx)
498 
499 
500 !This section has been created automatically by the script Abilint (TD).
501 !Do not modify the following lines by hand.
502 #undef ABI_FUNC
503 #define ABI_FUNC 'melements_init'
504 !End of the abilint section
505 
506  implicit none
507 
508 !Arguments ------------------------------------
509 !scalars
510  integer,intent(in) :: nspinor,nspden,nsppol,nkibz
511  type(melflags_t),intent(in) :: Mflags_in
512  type(melements_t),intent(out) :: Mels
513 !arrays
514  integer,intent(in) :: bands_idx(2,nkibz,nsppol)
515  real(dp),intent(in) :: kibz(3,nkibz)
516 
517 !Local variables-------------------------------
518  integer :: ikibz,isppol,bmin,bmax,b1,b2
519 
520 ! *************************************************************************
521 
522  DBG_ENTER("COLL")
523 
524  !@melements_t
525 
526  ! * Copy flags.
527  call melflags_copy(Mflags_in, Mels%flags)
528 
529  ! * Copy dimensions.
530  Mels%nkibz   = nkibz
531  Mels%nsppol  = nsppol
532  Mels%nspinor = nspinor
533  Mels%nspden  = nspden
534 
535  ABI_MALLOC(Mels%bands_idx,(2,nkibz,nsppol))
536  Mels%bands_idx=bands_idx
537 
538  ABI_MALLOC(Mels%iscalc,(nkibz,nsppol))
539  Mels%iscalc=0
540 
541  bmin = HUGE(1); bmax =-HUGE(1)
542  do isppol=1,Mels%nsppol
543    do ikibz=1,Mels%nkibz
544      if (ANY(Mels%bands_idx(:,ikibz,isppol)/=0)) then
545        b1 = Mels%bands_idx(1,ikibz,isppol)
546        b2 = Mels%bands_idx(2,ikibz,isppol)
547        Mels%iscalc(ikibz,isppol)=1
548        bmin = MIN(bmin,b1)
549        bmax = MAX(bmax,b2)
550        ABI_CHECK(b2>=b1 .and. b1>0,"Wrong b1, b2")
551      end if
552    end do
553  end do
554 
555  if (bmin==HUGE(1).or.bmax==-HUGE(1)) then
556    MSG_BUG("Wrong bands_idx")
557  end if
558 
559  Mels%bmin = bmin
560  Mels%bmax = bmax
561 
562  b1 = Mels%bmin
563  b2 = Mels%bmax
564 
565 ! real arrays
566  ABI_MALLOC(Mels%kibz,(3,nkibz))
567  Mels%kibz = kibz
568 
569 ! complex arrays
570  if (Mels%flags%has_hbare == 1) then
571    ABI_CALLOC(Mels%hbare,(b1:b2,b1:b2,nkibz,nsppol*nspinor**2))
572  end if
573 
574  if (Mels%flags%has_sxcore == 1) then
575    ABI_CALLOC(Mels%sxcore,(b1:b2,b1:b2,nkibz,nsppol*nspinor**2))
576  end if
577 
578  if (Mels%flags%has_vhartree == 1) then
579    ABI_CALLOC(Mels%vhartree,(b1:b2,b1:b2,nkibz,nsppol*nspinor**2))
580  end if
581 
582  if (Mels%flags%has_lexexch == 1) then
583    ABI_CALLOC(Mels%vlexx,(b1:b2,b1:b2,nkibz,nsppol*nspinor**2))
584  end if
585 
586  if (Mels%flags%has_vu == 1) then
587    ABI_CALLOC(Mels%vu,(b1:b2,b1:b2,nkibz,nsppol*nspinor**2))
588  end if
589 
590  if (Mels%flags%has_vxc == 1) then
591    ABI_CALLOC(Mels%vxc,(b1:b2,b1:b2,nkibz,nsppol*nspinor**2))
592  end if
593 
594  if (Mels%flags%has_vxcval == 1) then
595    ABI_CALLOC(Mels%vxcval,(b1:b2,b1:b2,nkibz,nsppol*nspinor**2))
596  end if
597 
598  if (Mels%flags%has_vxcval_hybrid == 1) then
599    ABI_CALLOC(Mels%vxcval_hybrid,(b1:b2,b1:b2,nkibz,nsppol*nspinor**2))
600  end if
601 
602  DBG_EXIT("COLL")
603 
604 end subroutine melements_init

m_melemts/melements_mpisum [ Functions ]

[ Top ] [ m_melemts ] [ Functions ]

NAME

 melements_mpisum

FUNCTION

  Perform a collective SUM within the MPI communicator comm
  of the matrix elements stored in the database.

INPUTS

  Mels=The database
  [aname]=The name of a particular array to be summed, by default
    all allocated arrays are considered.

SIDE EFFECTS

  All arrays whose flag==1 are summed within the MPI communicator comm.
  In output the corresponding flas is set to 2.

PARENTS

      calc_vhxc_me

CHILDREN

      my_select_melements

SOURCE

731 subroutine melements_mpisum(Mels,comm,aname)
732 
733 
734 !This section has been created automatically by the script Abilint (TD).
735 !Do not modify the following lines by hand.
736 #undef ABI_FUNC
737 #define ABI_FUNC 'melements_mpisum'
738 !End of the abilint section
739 
740  implicit none
741 
742 !Arguments ------------------------------------
743 !scalars
744  integer,intent(in) :: comm
745  character(len=*),optional,intent(in) :: aname
746  type(melements_t),intent(inout) :: Mels
747 
748 !Local variables-------------------------------
749  integer :: iname,ierr
750  integer,pointer :: flag_p
751  character(len=NAMELEN) :: key
752  character(len=500) :: msg
753 !arrays
754  complex(dpc),ABI_CONTIGUOUS pointer :: arr_p(:,:,:,:)
755 
756 ! *************************************************************************
757 
758  do iname=1,NNAMES
759    key = ANAMES(iname)
760    if (PRESENT(aname)) then
761     if (key /= aname) CYCLE
762    end if
763 
764    call my_select_melements(Mels,key,flag_p,arr_p)
765 
766    if (flag_p==1) then
767      call xmpi_sum(arr_p,comm,ierr)
768      if (ierr/=0) then
769        write(msg,'(a,i4,2a)')" xmpi_sum reported ierr= ",ierr," for key= ",TRIM(key)
770        MSG_ERROR(msg)
771      end if
772      flag_p=2 ! Tag this array as calculated
773    end if
774  end do
775 
776 end subroutine melements_mpisum

m_melemts/melements_print [ Functions ]

[ Top ] [ m_melemts ] [ Functions ]

NAME

 melements_print

FUNCTION

  Printout of the content of all calculated array.
  Optionally, it is possible to print the content of a single entry of the database.

INPUTS

  Mels=The database
  [unit]=the unit number for output, defaults to std_out
  [prtvol]=verbosity level, defaults to 0
  [mode_paral]=either "COLL" or "PERS", default to "COLL"
  [header]=title for info

OUTPUT

  Only writing

PARENTS

      sigma

CHILDREN

      my_select_melements

SOURCE

807 subroutine melements_print(Mels,names_list,header,unit,prtvol,mode_paral)
808 
809 
810 !This section has been created automatically by the script Abilint (TD).
811 !Do not modify the following lines by hand.
812 #undef ABI_FUNC
813 #define ABI_FUNC 'melements_print'
814 !End of the abilint section
815 
816  implicit none
817 
818 !Arguments ------------------------------------
819 !scalars
820  type(melements_t),intent(in) :: Mels
821  integer,optional,intent(in) :: prtvol,unit
822  character(len=*),optional,intent(in) :: names_list(:)
823  character(len=*),optional,intent(in) :: header
824  character(len=4),optional,intent(in) :: mode_paral
825 
826 !Local variables-------------------------------
827  integer :: my_unt,my_prtvol,max_r,max_c,ii
828  integer :: isppol,ikibz,iab,ib,b1,b2,my_nkeys,ikey
829  integer,pointer :: flag_p
830  character(len=4) :: my_mode
831  character(len=NAMELEN) :: key
832  character(len=500) :: msg,str,fmt
833 !arrays
834  integer,allocatable :: tab(:)
835  character(len=NAMELEN),allocatable :: my_keys(:)
836  complex(dpc),allocatable :: mat(:,:)
837 
838  type rarr_dpc4
839    complex(dpc),ABI_CONTIGUOUS pointer :: arr_p(:,:,:,:)
840  end type rarr_dpc4
841  type(rarr_dpc4),allocatable :: data_p(:)
842 
843 ! *************************************************************************
844 
845  !@melements_t
846  my_unt   =std_out; if (PRESENT(unit      )) my_unt   =unit
847  my_prtvol=0      ; if (PRESENT(prtvol    )) my_prtvol=prtvol
848  my_mode  ='COLL' ; if (PRESENT(mode_paral)) my_mode  =mode_paral
849 
850  if (Mels%nspinor == 2) MSG_WARNING("nspinor=2 not coded")
851 
852  if (PRESENT(names_list)) then
853    my_nkeys=SIZE(names_list)
854    ABI_MALLOC(my_keys,(my_nkeys))
855    my_keys = names_list
856  else
857    my_nkeys=NNAMES
858    ABI_MALLOC(my_keys,(NNAMES))
859    my_keys = ANAMES
860  end if
861 
862  ABI_DT_MALLOC(data_p,(my_nkeys))
863  ABI_MALLOC(tab,(my_nkeys))
864  tab=0
865 
866  my_nkeys=0; str = "  ib"; ii=4
867  do ikey=1,SIZE(my_keys)
868    key = my_keys(ikey)
869    call my_select_melements(Mels,key,flag_p,data_p(ikey)%arr_p)
870    if (flag_p==2) then
871      my_nkeys=my_nkeys+1
872      tab(my_nkeys)=ikey
873      str(ii+1:)=" "//TRIM(tolower(key))
874      ii = ii+MAX(1+LEN_TRIM(key),10)
875      ABI_CHECK(ii<490,"I'm gonna SIGFAULT!")
876    end if
877  end do
878 
879  write(msg,'(2a)')ch10,' === Matrix Elements stored in Mels% [eV] === '
880  if (PRESENT(header)) write(msg,'(4a)')ch10,' === '//TRIM(ADJUSTL(header))//' [eV] === '
881  call wrtout(my_unt,msg,my_mode)
882  if (Mels%nspinor == 2) then
883    call wrtout(my_unt, "Sum_ab M_ab, M_11, M_22, Re(M_12), IM(Re_12)" ,my_mode)
884  end if
885 
886  if (my_nkeys==0) GOTO 10
887  write(fmt,'(a,i4,a)')'(1x,i3,',my_nkeys,'(1x,f9.5))' ! width of 10 chars
888 
889  do isppol=1,Mels%nsppol
890    do ikibz=1,Mels%nkibz
891     if (Mels%iscalc(ikibz,isppol)/=1) CYCLE
892 
893     write(msg,'(a,3es16.8,a,i2,a)')" kpt= (",Mels%kibz(:,ikibz),") spin=",isppol,":"
894     call wrtout(my_unt,msg,my_mode)
895 
896     b1 = Mels%bands_idx(1,ikibz,isppol)
897     b2 = Mels%bands_idx(2,ikibz,isppol)
898 
899     if (Mels%flags%only_diago==1.or.my_prtvol==0) then
900       ! Print only the diagonal.
901       write(msg,'(a)')str
902       call wrtout(my_unt,msg,my_mode)
903       do ib=b1,b2
904         if (Mels%nspinor == 1) then
905           write(msg,fmt)ib,(REAL(data_p(tab(ikey))%arr_p(ib,ib,ikibz,1))*Ha_eV, ikey=1,my_nkeys)
906           call wrtout(my_unt,msg,my_mode)
907         else
908           ! Write sum_ab, then diagonal elements, finally Re_12, Im_12
909           write(msg,fmt)ib,(real(sum(data_p(tab(ikey))%arr_p(ib,ib,ikibz,:)))*Ha_eV, ikey=1,my_nkeys)
910           call wrtout(my_unt,msg,my_mode)
911           if (my_prtvol > 0) then
912             write(msg,fmt)ib,(real(data_p(tab(ikey))%arr_p(ib,ib,ikibz,1))*Ha_eV, ikey=1,my_nkeys)
913             call wrtout(my_unt,msg,my_mode)
914             write(msg,fmt)ib,(real(data_p(tab(ikey))%arr_p(ib,ib,ikibz,2))*Ha_eV, ikey=1,my_nkeys)
915             call wrtout(my_unt,msg,my_mode)
916             write(msg,fmt)ib,(real(data_p(tab(ikey))%arr_p(ib,ib,ikibz,3))*Ha_eV, ikey=1,my_nkeys)
917             call wrtout(my_unt,msg,my_mode)
918             write(msg,fmt)ib,(aimag(data_p(tab(ikey))%arr_p(ib,ib,ikibz,3))*Ha_eV, ikey=1,my_nkeys)
919             call wrtout(my_unt,msg,my_mode)
920           end if
921         end if
922       end do
923 
924     else
925       ! Print full matrix.
926       max_r = b2-b1+1
927       max_c = MIN(b2-b1+1, 9)
928       ABI_MALLOC(mat,(b1:b2,b1:b2))
929       do ikey=1,my_nkeys
930         write(msg,'(3a)')" **** Off-diagonal elements of ",TRIM(my_keys(tab(ikey)))," **** "
931         call wrtout(my_unt,msg,my_mode)
932         do iab=1,Mels%nspinor**2
933           mat = data_p(tab(ikey))%arr_p(b1:b2,b1:b2,ikibz,iab)*Ha_eV
934           call print_arr(mat,max_r,max_c,my_unt,my_mode)
935         end do
936         write(msg,'(a)')ch10
937         call wrtout(my_unt,msg,my_mode)
938       end do
939       ABI_FREE(mat)
940     end if
941 
942    end do !ikibz
943  end do ! isppol
944 
945 10 continue
946  ABI_FREE(my_keys)
947  ABI_DT_FREE(data_p)
948  ABI_FREE(tab)
949 
950 end subroutine melements_print

m_melemts/melements_t [ Types ]

[ Top ] [ m_melemts ] [ Types ]

NAME

FUNCTION

  Structure defining a database to store the matrix elements of operators
  needed for GW calculations.

SOURCE

106  type,public :: melements_t
107 
108   integer :: nkibz
109   ! Number of k-points.
110 
111   integer :: nsppol
112   ! Number of independent spin-polarizations.
113 
114   integer :: nspinor
115   ! 1 for collinear, 2 for noncollinear.
116 
117   integer :: nspden
118   ! Number of independent spin-density components.
119 
120   integer :: bmin,bmax
121   ! min and Max band index over k-points and spin.
122   ! Used to dimension the arrays below.
123 
124   integer, allocatable :: bands_idx(:,:,:)
125   ! bands_idx(2,nkibz,nsppol)
126   ! min and Max band index for each k-point and spin.
127 
128   integer, allocatable :: iscalc(:,:)
129   ! stat_k(nkibz,nsppol)
130   ! 1 if this k-point and spin has been calculated, 0 otherwise.
131 
132   real(dp), allocatable :: kibz(:,:)
133   ! kibz(3,nkibz)
134   ! The list of k-points in reduced coordinates.
135 
136   complex(dpc), allocatable :: hbare(:,:,:,:)
137   ! hbare(b1:b2,b1:b2,nkibz,nsppol*nspinor**2)
138   ! Matrix elements of the bare Hamiltonian.
139 
140   complex(dpc), allocatable :: sxcore(:,:,:,:)
141   ! sxcore(b1:b2,b1:b2,nkibz,nsppol*nspinor**2)
142   ! Matrix elements of the Fock operator generated by core electrons.
143 
144   complex(dpc), allocatable :: vhartree(:,:,:,:)
145   ! vhartree(b1:b2,b1:b2,nkibz,nsppol*nspinor**2)
146   ! Matrix elements of the Hartree potential.
147 
148   complex(dpc), allocatable :: vlexx(:,:,:,:)
149   ! vlexx(b1:b2,b1:b2,nkibz,nsppol*nspinor**2)
150   ! Matrix elements of the local exact exchange potential.
151 
152   complex(dpc), allocatable :: vu(:,:,:,:)
153   ! vu(b1:b2,b1:b2,nkibz,nsppol*nspinor**2)
154   ! Matrix elements of the U Hamiltonian.
155 
156   complex(dpc), allocatable :: vxc(:,:,:,:)
157   ! vxc(b1:b2,b1:b2,nkibz,nsppol*nspinor**2)
158   ! Matrix elements of XC potential, including core if present.
159 
160   complex(dpc), allocatable :: vxcval(:,:,:,:)
161   ! vxcval(b1:b2,b1:b2,nkibz,nsppol*nspinor**2)
162   ! Matrix elements of the XC potential, valence-only contribution.
163 
164   complex(dpc), allocatable :: vxcval_hybrid(:,:,:,:)
165   ! vxcval_hybrid(b1:b2,b1:b2,nkibz,nsppol*nspinor**2)
166   ! Matrix elements of the XC potential for hybrid calculations, valence-only contribution.
167 
168   type(melflags_t) :: flags
169 
170  end type melements_t
171 
172  public :: melements_init      ! Initialize the object
173  public :: melements_free      ! Free memory
174  public :: melements_herm      ! Construct the lower triangle from the upper triangle
175  public :: melements_mpisum    ! Perform a collective SUM within the MPI communicator comm
176  public :: melements_print     ! Print matrix elements
177  public :: melements_zero      ! Set matrix elements connecting states with different irrep to zero.
178  !public :: mels_get_exene_core

m_melemts/melements_zero [ Functions ]

[ Top ] [ m_melemts ] [ Functions ]

NAME

 melements_zero

FUNCTION

  Set matrix elements connecting states with different irreducible representation to zero.

INPUTS

   irrep_tab=Array used to select the entries that have to be set to zero.
     irrep_tab(ib,ik,is)=gives the index of the irreducible representation associated to state (ib,ik,is).
  [aname]=The name of the array to be symmetrized, by default
    all calculated arrays are filled.

SIDE EFFECTS

  Mels= All arrays elements connecting states belonging to different irreps are set to zero.

PARENTS

      sigma

CHILDREN

      my_select_melements

SOURCE

 979 subroutine melements_zero(Mels,irrep_tab,aname)
 980 
 981 
 982 !This section has been created automatically by the script Abilint (TD).
 983 !Do not modify the following lines by hand.
 984 #undef ABI_FUNC
 985 #define ABI_FUNC 'melements_zero'
 986 !End of the abilint section
 987 
 988  implicit none
 989 
 990 !Arguments ------------------------------------
 991 !scalars
 992  character(len=*),optional,intent(in) :: aname
 993  type(melements_t),intent(inout) :: Mels
 994 !arrays
 995  integer,intent(in) :: irrep_tab(Mels%bmin:Mels%bmax,Mels%nkibz,Mels%nsppol)
 996 
 997 !Local variables-------------------------------
 998  integer :: is,ik,ib,jb,iname,irrep_j,irrep_i
 999  integer,pointer :: flag_p
1000  character(len=NAMELEN) :: key
1001 !arrays
1002  complex(dpc),ABI_CONTIGUOUS pointer :: arr_p(:,:,:,:)
1003 
1004 ! *************************************************************************
1005 
1006  do iname=1,NNAMES
1007    key = ANAMES(iname)
1008    if (PRESENT(aname)) then
1009      if (key /= aname) CYCLE
1010    end if
1011 
1012    call my_select_melements(Mels,key,flag_p,arr_p)
1013 
1014    if (flag_p>0) then
1015      do is=1,Mels%nsppol
1016        do ik=1,Mels%nkibz
1017 
1018          do jb=Mels%bmin,Mels%bmax
1019            irrep_j = irrep_tab(jb,ik,is)
1020            do ib=Mels%bmin,Mels%bmax
1021              irrep_i = irrep_tab(ib,ik,is)
1022              !
1023              ! Set this matrix element to zero if the irreps are known and they differ.
1024              if (irrep_i/=irrep_j .and. ALL((/irrep_i,irrep_j/) /=0) ) then
1025                !write(std_out,*)"setting to zero ",ib,jb,ik,is
1026                if (Mels%nspinor==2) then
1027                  arr_p(ib,jb,ik,is)=czero
1028                else
1029                  arr_p(ib,jb,ik,:)=czero
1030                end if
1031              end if
1032 
1033            end do !ib
1034          end do !jb
1035 
1036        end do !is
1037      end do !ik
1038    end if
1039 
1040  end do !inames
1041 
1042 end subroutine melements_zero

m_melemts/melflags_copy [ Functions ]

[ Top ] [ m_melemts ] [ Functions ]

NAME

 melflags_copy

FUNCTION

  Copy an object storing the flags.

INPUTS

  Mflags_in=The flags to be copied.

OUTPUT

  Mflags_out=The new set of flags.

PARENTS

      m_melemts

CHILDREN

      my_select_melements

SOURCE

269 subroutine melflags_copy(Mflags_in, Mflags_out)
270 
271 
272 !This section has been created automatically by the script Abilint (TD).
273 !Do not modify the following lines by hand.
274 #undef ABI_FUNC
275 #define ABI_FUNC 'melflags_copy'
276 !End of the abilint section
277 
278  implicit none
279 
280 !Arguments ------------------------------------
281  type(melflags_t),intent(in)    :: Mflags_in
282  type(melflags_t),intent(inout) :: Mflags_out
283 
284 ! *************************************************************************
285 
286  call melflags_reset(Mflags_out)
287 
288  ! @melflags_t
289  Mflags_out%has_hbare           = Mflags_in%has_hbare
290  Mflags_out%has_sxcore          = Mflags_in%has_sxcore
291  Mflags_out%has_vhartree        = Mflags_in%has_vhartree
292  Mflags_out%has_vu              = Mflags_in%has_vu
293  Mflags_out%has_vxc             = Mflags_in%has_vxc
294  Mflags_out%has_vxcval          = Mflags_in%has_vxcval
295  Mflags_out%has_vxcval_hybrid   = Mflags_in%has_vxcval_hybrid
296  Mflags_out%has_lexexch         = Mflags_in%has_lexexch
297  Mflags_out%only_diago          = Mflags_in%only_diago
298 
299 end subroutine melflags_copy

m_melemts/melflags_reset [ Functions ]

[ Top ] [ m_melemts ] [ Functions ]

NAME

 melflags_reset

FUNCTION

  Set all flags in melflags_t to 0.

INPUTS

PARENTS

      m_melemts,sigma

CHILDREN

      my_select_melements

SOURCE

216 subroutine melflags_reset(Mflags)
217 
218 
219 !This section has been created automatically by the script Abilint (TD).
220 !Do not modify the following lines by hand.
221 #undef ABI_FUNC
222 #define ABI_FUNC 'melflags_reset'
223 !End of the abilint section
224 
225  implicit none
226 
227 !Arguments ------------------------------------
228  type(melflags_t),intent(inout) :: Mflags
229 
230 ! *************************************************************************
231 
232  ! @melflags_t
233  Mflags%has_hbare           = 0
234  Mflags%has_sxcore          = 0
235  Mflags%has_vhartree        = 0
236  Mflags%has_vu              = 0
237  Mflags%has_vxc             = 0
238  Mflags%has_vxcval          = 0
239  Mflags%has_vxcval_hybrid   = 0
240  Mflags%has_lexexch         = 0
241  Mflags%only_diago          = 0
242 
243 end subroutine melflags_reset

m_melemts/melflags_t [ Types ]

[ Top ] [ m_melemts ] [ Types ]

NAME

FUNCTION

  Container for the flags defining the status of the corresponding
  pointer defined in the type melements_t. Possible values are:
    * 0 if the correspondending array is not allocated.
    * 1 if allocated but not yet calculated.
    * 2 if allocated and calculated

SOURCE

76  type,public :: melflags_t
77 
78   integer :: has_hbare=0
79   integer :: has_lexexch=0
80   integer :: has_sxcore=0
81   integer :: has_vhartree=0
82   integer :: has_vu=0
83   integer :: has_vxc=0
84   integer :: has_vxcval=0
85   integer :: has_vxcval_hybrid=0
86   integer :: only_diago=0
87   ! 1 if only diagonal elements are calculated
88 
89  end type melflags_t
90 
91  public :: melflags_reset  ! Reset the value of the flags.
92  public :: melflags_copy   ! Copy the object

m_melemts/my_select_melements [ Functions ]

[ Top ] [ m_melemts ] [ Functions ]

NAME

  my_select_melements

FUNCTION

  Helper function returning a pointer to the array "aname" as well as the status of the array

INPUTS

  Mels<melements_t>=The database.
  aname=String with the name of the array.

OUTPUT

  flag_p=Pointer to the integer defining the status of the array, see melflags_t.
  arr_p=The pointer to the array.

PARENTS

      m_melemts

CHILDREN

      my_select_melements

SOURCE

417 subroutine my_select_melements(Mels,aname,flag_p,arr_p)
418 
419 
420 !This section has been created automatically by the script Abilint (TD).
421 !Do not modify the following lines by hand.
422 #undef ABI_FUNC
423 #define ABI_FUNC 'my_select_melements'
424 !End of the abilint section
425 
426  implicit none
427 
428 !Arguments ------------------------------------
429 !scalars
430  integer,pointer :: flag_p
431  character(len=*),intent(in) :: aname
432  type(melements_t),target,intent(in) :: Mels
433 !arrays
434  complex(dpc),ABI_CONTIGUOUS pointer :: arr_p(:,:,:,:)
435 
436 ! *************************************************************************
437 
438  SELECT CASE (tolower(aname))
439  CASE ("hbare")
440    flag_p => Mels%flags%has_hbare
441    arr_p  => Mels%hbare
442  CASE ("sxcore")
443    flag_p => Mels%flags%has_sxcore
444    arr_p  => Mels%sxcore
445  CASE ("vhartree")
446    flag_p => Mels%flags%has_vhartree
447    arr_p  => Mels%vhartree
448  CASE ("vlexx")
449    flag_p => Mels%flags%has_lexexch
450    arr_p  => Mels%vlexx
451  CASE ("vu")
452    flag_p => Mels%flags%has_vu
453    arr_p  => Mels%vu
454  CASE ("vxc")
455    flag_p => Mels%flags%has_vxc
456    arr_p  => Mels%vxc
457  CASE ("vxcval")
458    flag_p => Mels%flags%has_vxcval
459    arr_p  => Mels%vxcval
460  CASE ("vxcval_hybrid")
461    flag_p => Mels%flags%has_vxcval_hybrid
462    arr_p  => Mels%vxcval_hybrid
463  CASE DEFAULT
464    MSG_ERROR(sjoin("Wrong aname: ", aname))
465  END SELECT
466 
467 end subroutine my_select_melements

m_sigma/mels_get_exene_core [ Functions ]

[ Top ] [ m_sigma ] [ Functions ]

NAME

  mels_get_exene_core

FUNCTION

  Compute exchange energy.

INPUTS

  mels<melements_t>=Matrix elements.
  kmesh<kmesh_t>=BZ sampling.
  bands<band_t>=Bands with occupation factors

PARENTS

SOURCE

1063 !pure function mels_get_exene_core(mels,kmesh,bands) result(ex_energy)
1064 !
1065 !!This section has been created automatically by the script Abilint (TD).
1066 !!Do not modify the following lines by hand.
1067 !#undef ABI_FUNC
1068 !#define ABI_FUNC 'mels_get_exene_core'
1069 !!End of the abilint section
1070 !
1071 ! implicit none
1072 !
1073 !!Arguments ------------------------------------
1074 !!scalars
1075 ! real(dp) :: ex_energy
1076 ! type(melements_t),intent(in) :: mels
1077 ! type(kmesh_t),intent(in) :: kmesh
1078 ! type(ebands_t),intent(in) :: bands
1079 !
1080 !!Local variables-------------------------------
1081 !!scalars
1082 ! integer :: ik,ib,spin
1083 ! real(dp) :: wtk,occ_bks
1084 !
1085 !! *************************************************************************
1086 !
1087 ! ex_energy = zero
1088 !
1089 ! do spin=1,mels%nsppol
1090 !   do ik=1,mels%nkibz
1091 !     wtk = kmesh%wt(ik)
1092 !     do ib=mels%bmin,mels%bmax
1093 !       occ_bks = bands%occ(ib,ik,spin)
1094 !       if (mels%nspinor==1) then
1095 !         ex_energy = ex_energy + half * occ_bks * wtk * mels%sxcore(ib,ib,ik,spin)
1096 !       else
1097 !         ex_energy = ex_energy + half * occ_bks wtk *SUM(mels%sxcore(ib,ib,ik,:))
1098 !       end if
1099 !     end do
1100 !   end do
1101 ! end do
1102 !
1103 !end function mels_get_exene_core
1104 !!!***
1105 
1106 END MODULE m_melemts