TABLE OF CONTENTS


ABINIT/m_hexc [ Modules ]

[ Top ] [ Modules ]

NAME

 m_hexc

FUNCTION

 module for excitonic hamiltonian for Haydock

COPYRIGHT

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

PARENTS

      m_haydock

CHILDREN

SOURCE

22 #if defined HAVE_CONFIG_H
23 #include "config.h"
24 #endif
25 
26 #include "abi_common.h"
27 
28 MODULE m_hexc
29 
30  use defs_basis
31  use m_abicore
32  use m_bs_defs
33  use m_xmpi
34  use m_errors
35  use m_nctk
36  use m_haydock_io
37  use m_linalg_interfaces
38 #ifdef HAVE_NETCDF
39  use netcdf
40 #endif
41 
42  use m_time,              only : timab
43  use m_fstrings,          only : indent, strcat, sjoin, itoa
44  use defs_datatypes,      only : ebands_t, pseudopotential_type
45  use m_hide_blas,         only : xdotc, xgemv
46  use m_numeric_tools,     only : print_arr, symmetrize, hermitianize, wrap2_pmhalf
47  use m_crystal,           only : crystal_t
48  use m_bz_mesh,           only : kmesh_t, findqg0, get_bz_item
49  use m_double_grid,       only : double_grid_t, get_kpt_from_indices_coarse, compute_corresp
50  use m_wfd,               only : wfd_t,wfd_sym_ur,wfd_get_ur, wfd_change_ngfft
51  use m_bse_io,            only : exc_read_rcblock, exc_write_optme, exc_ham_ncwrite
52  use m_pawtab,            only : pawtab_type
53  use m_vcoul,             only : vcoul_t
54  use m_bseinterp,            only : interpolator_t, interpolator_init, interpolator_normalize, &
55 &                    interpolator_free, int_alloc_work, int_free_work
56 
57  implicit none
58 
59  private

m_haydock/hexc_t [ Types ]

[ Top ] [ m_haydock ] [ Types ]

NAME

 hexc_t

FUNCTION

  Store the excitonic hamiltonian and other related information

SOURCE

 71  type,public :: hexc_t
 72 
 73  !scalars
 74     integer :: comm
 75     ! MPI communicator
 76 
 77     integer :: hsize_coarse
 78     ! Size of the coarse hamiltonian
 79 
 80     integer :: hsize
 81     ! Size of the hamiltonian and the kets
 82     ! (= hsize_coarse without interpolation, = hsize_dense with interpolation)
 83 
 84     integer :: nbz
 85     ! Number of kpoints for the full problem
 86     ! (= nbz_coarse without interpolation, = nbz_dense with interpolation)
 87 
 88     integer :: my_t1
 89     ! Lower limit of MPI paral
 90 
 91     integer :: my_t2
 92     ! Upper limit of MPI paral
 93 
 94     integer :: my_nt
 95     ! Number of transitions treat by node
 96     ! = my_t2 - my_t1 + 1
 97 
 98     integer :: nbnd_coarse
 99     ! Product of number of bands conduction X valence
100 
101     ! Pointers to data that are already in memory
102     type(excparam),pointer :: bsp => null()
103     ! parameters for BS
104 
105     type(excfiles),pointer :: bs_files => null()
106     ! files for BSE
107 
108     type(crystal_t),pointer :: crystal => null()
109     ! crystal info
110 
111     type(kmesh_t),pointer :: kmesh_coarse => null()
112     ! kmesh of the coarse mesh
113 
114     type(kmesh_t),pointer :: kmesh => null()
115     ! kmesh of the full problem
116 
117     type(ebands_t),pointer :: ks_bst => null()
118     type(ebands_t),pointer :: qp_bst => null()
119     ! band structures of the full problem
120 
121     type(wfd_t),pointer :: wfd_coarse => null()
122     ! Wfd of the coarse problem
123 
124     type(wfd_t),pointer :: wfd => null()
125     ! wfd of the full problem
126 
127  !arrays
128     complex(dpc),allocatable :: hreso(:,:)
129     ! Resonant part of the hamiltonian
130 
131     complex(dpc),allocatable :: hcoup(:,:)
132     ! Coupling part of the hamiltonian
133 
134     complex(dpc),allocatable :: diag_coarse(:)
135     ! Diagonal part of the hamiltonian with transition energies
136 
137  end type hexc_t

m_hexc/hexc_build_hinterp [ Functions ]

[ Top ] [ m_hexc ] [ Functions ]

NAME

 hexc_build_hinterp

FUNCTION

 Pre-compute interpolated hamiltonian and store it in memory

INPUTS

OUTPUT

SIDE EFFECTS

 hexc, hexc_i
   Pre-compute info to save CPU time when computing matmul

PARENTS

      m_haydock

CHILDREN

      timab

SOURCE

623 subroutine hexc_build_hinterp(hexc,hexc_i)
624 
625 
626 !This section has been created automatically by the script Abilint (TD).
627 !Do not modify the following lines by hand.
628 #undef ABI_FUNC
629 #define ABI_FUNC 'hexc_build_hinterp'
630 !End of the abilint section
631 
632  implicit none
633 
634 !Arguments ---------------------------
635  type(hexc_t),intent(inout) :: hexc
636  type(hexc_interp_t),intent(inout) :: hexc_i
637 
638 !Local variables ---------------------
639  integer :: ierr,ncerr,ncid
640  character(len=500) :: msg
641 
642 !*****************************************************************************
643 
644  write(msg,"(a,f8.1,a)")"Memory needed for hinterp = ",one*(hexc_i%hsize_dense**2)*2*dpc*b2Mb," Mb"
645  call wrtout(std_out,msg,"COLL")
646 
647  ABI_STAT_MALLOC(hexc_i%hinterp,(hexc_i%hsize_dense,hexc_i%hsize_dense), ierr)
648  ABI_CHECK(ierr==0, 'Out of memory in hinterp')
649 
650  call hexc_compute_hinterp(hexc%BSp, hexc%hsize_coarse, hexc_i%hsize_dense, hexc_i%all_hmat, &
651 &  hexc_i%interpolator%double_grid,hexc%nbnd_coarse, hexc_i%interpolator, &
652 &  hexc_i%kdense2div, hexc_i%all_acoeffs,hexc_i%all_bcoeffs, hexc_i%all_ccoeffs, &
653 &  hexc_i%Kmesh_dense, hexc_i%Vcp_dense, hexc%crystal%gmet, hexc_i%hinterp, hexc_i%m3_width)
654 
655  if( allocated(hexc_i%all_acoeffs) ) then
656    ABI_FREE(hexc_i%all_acoeffs)
657  end if
658  if( allocated(hexc_i%all_bcoeffs) ) then
659    ABI_FREE(hexc_i%all_bcoeffs)
660  end if
661  if( allocated(hexc_i%all_ccoeffs) ) then
662    ABI_FREE(hexc_i%all_ccoeffs)
663  end if
664 
665 
666  if (hexc%BSp%prt_ncham) then
667 #ifdef HAVE_NETCDF
668    MSG_COMMENT("Printing HEXC_I.nc file")
669    ncerr = nctk_open_create(ncid, trim(hexc%BS_files%out_basename)//"_HEXC_I.nc", xmpi_comm_self)
670    NCF_CHECK_MSG(ncerr, "Creating HEXC_I file")
671    call exc_ham_ncwrite(ncid, hexc_i%Kmesh_dense, hexc%BSp, hexc_i%hsize_dense, hexc%BSp%nreh_interp, &
672 &       hexc%BSp%vcks2t_interp, hexc_i%hinterp, hexc_i%diag_dense)
673    NCF_CHECK(nf90_close(ncid))
674 #else
675    ABI_UNUSED(ncid)
676 #endif
677  end if
678 
679 end subroutine hexc_build_hinterp

m_hexc/hexc_compute_hinterp [ Functions ]

[ Top ] [ m_hexc ] [ Functions ]

NAME

 hexc_compute_hinterp

FUNCTION

 Compute interpolated matrix elements for methods 2 and 3

INPUTS

 BSp<type(excparam)=The parameter for the Bethe-Salpeter run.
 hsize_coarse=Size of the coarse Hamiltonian
 hsize_dense=Size of the dense Hamiltonian
 hmat(hsize_coarse,hsize_coarse,8)=Excitonic matrix
 grid<double_grid_t> = Correspondence between coarse and dense k-mesh.
 nbnd_coarse = Total number of bands
 interpolator<interpolator_t> = Interpolator
 kdense2div = Mapping kdense2div
 acoeffs, bcoeffs, ccoeffs = decomposition "W = a/q^2 + b/q + c"
 Kmesh_dense<type(kmesh_t)>=The list of k-points in the BZ, IBZ and symmetry tables.
 Vcp_dense<vcoul_t>=Coulomb interation in G-space on the dense Q-mesh
 gmet(3,3)=Metric tensor in G-space

OUTPUT

   hinterp = Interpolated hamiltonian

PARENTS

      m_hexc

CHILDREN

      timab

SOURCE

 849 subroutine hexc_compute_hinterp(BSp,hsize_coarse,hsize_dense,hmat,grid,nbnd_coarse,&
 850 &  interpolator,kdense2div,acoeffs,bcoeffs,ccoeffs,Kmesh_dense,Vcp_dense,gmet,hinterp,&
 851 &  m3_width)
 852 
 853 
 854 !This section has been created automatically by the script Abilint (TD).
 855 !Do not modify the following lines by hand.
 856 #undef ABI_FUNC
 857 #define ABI_FUNC 'hexc_compute_hinterp'
 858 !End of the abilint section
 859 
 860  implicit none
 861 
 862 !Arguments ------------------------------------
 863 !scalars
 864  integer,intent(in) :: hsize_coarse,hsize_dense,nbnd_coarse !,ntrans
 865  real(dp),intent(in) :: m3_width
 866  type(excparam),intent(in) :: BSp
 867  type(double_grid_t),intent(in) :: grid
 868  type(vcoul_t),intent(in) :: Vcp_dense
 869  type(kmesh_t),intent(in) :: Kmesh_dense
 870  type(interpolator_t),target,intent(inout) :: interpolator
 871 !arrays
 872  integer,intent(in) :: kdense2div(grid%nbz_dense)
 873  real(dp),intent(in) :: gmet(3,3)
 874  complex(dpc),intent(in) :: hmat(hsize_coarse,hsize_coarse)
 875  complex(dpc),intent(in) :: acoeffs(hsize_coarse,hsize_coarse)
 876  complex(dpc),intent(in) :: bcoeffs(hsize_coarse,hsize_coarse)
 877  complex(dpc),intent(in) :: ccoeffs(hsize_coarse,hsize_coarse)
 878  complex(dpc),intent(out) :: hinterp(hsize_dense,hsize_dense)
 879 
 880 !Local variables ------------------------------
 881 !scalars
 882  integer,parameter :: spin1=1, spin2=1
 883  integer :: ic,iv,iv1,ic1,ik_dense,ik_coarse,it_coarse,it_dense,idiv,ibnd_coarse,ibnd_coarse1,ineighbour
 884  integer :: icp,ivp,ikp_dense,ikp_coarse,itp_coarse,itp_dense,idivp,ibndp_coarse,ibndp_coarse1,ineighbourp,itp_coarse1
 885  integer :: itc,it_dense1,indwithnb, corresp_ind
 886  integer :: limitnbz, inb,ierr
 887  real(dp) :: factor,vc_sqrt_qbz,qnorm
 888  complex(dpc) :: term
 889  logical :: newway, use_herm
 890 !arrays
 891  real(dp) :: kmkp(3),q2(3),shift(3),qinred(3),tsec(2)
 892  complex(dpc),allocatable :: Cmat(:,:,:) !Temp matrices for optimized version
 893  complex(dpc),allocatable :: tmp_Cmat(:)
 894  complex(dpc),allocatable :: work_coeffs(:,:)
 895  integer,allocatable :: band2it(:)
 896  complex(dpc),ABI_CONTIGUOUS pointer :: btemp(:),ctemp(:)
 897 
 898 !************************************************************************
 899 
 900  call timab(696,1,tsec)
 901 
 902  newway = .True.
 903  use_herm = .True.
 904 
 905  if (any(BSp%interp_mode == [2,3])) then
 906    if(Vcp_dense%mode /= 'CRYSTAL' .and. Vcp_dense%mode /= 'AUXILIARY_FUNCTION') then
 907      MSG_BUG('Vcp_dense%mode not implemented yet !')
 908    end if
 909  end if
 910 
 911  if(BSp%nsppol > 1) then
 912    MSG_BUG("nsppol > 1 not yet implemented")
 913  end if
 914 
 915  factor = one/grid%ndiv
 916 
 917  hinterp = czero; term = czero
 918 
 919  ABI_STAT_MALLOC(Cmat,(nbnd_coarse,nbnd_coarse,interpolator%nvert), ierr)
 920  ABI_CHECK(ierr==0, "out of memory in Cmat")
 921  Cmat = czero
 922 
 923  ABI_MALLOC(band2it,(nbnd_coarse))
 924  call int_alloc_work(interpolator,nbnd_coarse*interpolator%nvert)
 925 
 926  btemp => interpolator%btemp
 927  ctemp => interpolator%ctemp
 928 
 929  if(newway) then
 930    ABI_MALLOC(tmp_Cmat,(nbnd_coarse))
 931    ABI_MALLOC(work_coeffs,(nbnd_coarse,interpolator%nvert))
 932  end if
 933 
 934  do ik_dense = 1,grid%nbz_dense
 935    write(std_out,*) "Kdense = ",ik_dense,"/",grid%nbz_dense
 936    ik_coarse = grid%dense_to_coarse(ik_dense)
 937    if(use_herm) then
 938      limitnbz = ik_dense
 939    else
 940      limitnbz = grid%nbz_dense
 941    end if
 942 
 943    do ikp_dense = 1,limitnbz
 944      ikp_coarse = grid%dense_to_coarse(ikp_dense)
 945 
 946      do iv1 = BSp%lomo_spin(spin2),BSp%homo_spin(spin2)
 947        do ic1 = BSp%lumo_spin(spin2),BSp%humo_spin(spin2)
 948          itp_coarse1 = BSp%vcks2t(iv1,ic1,ikp_coarse,spin2)
 949          ibndp_coarse1 = (iv1-BSp%lomo_spin(spin2))*BSp%maxnbndc+(ic1-BSp%lumo_spin(spin2)+1)
 950 
 951          band2it(ibndp_coarse1) = itp_coarse1
 952        end do
 953      end do
 954 
 955 
 956      if (any(BSp%interp_mode == [2,3])) then
 957        ! Check if we are along the diagonal
 958        kmkp = Kmesh_dense%bz(:,ik_dense) - Kmesh_dense%bz(:,ikp_dense)
 959 
 960        call wrap2_pmhalf(kmkp(:),q2(:),shift(:))
 961        qinred = MATMUL(grid%kptrlatt_coarse,q2)
 962 
 963        ! We are outside the diagonal
 964        if (BSp%interp_mode==3 .and. ANY((ABS(qinred)-tol7) > m3_width)) cycle
 965 
 966        qnorm = two_pi*SQRT(DOT_PRODUCT(q2,MATMUL(gmet,q2)))
 967 
 968        if(ALL(ABS(q2(:)) < 1.e-3)) then
 969          vc_sqrt_qbz = SQRT(Vcp_dense%i_sz)
 970        else
 971          vc_sqrt_qbz = SQRT(four_pi/qnorm**2)
 972        end if
 973 
 974        !!DEBUG CHK !
 975        !!COMPUTE Qpoint
 976        !call findqg0(iq_bz,g0,kmkp,Qmesh_dense%nbz,Qmesh_dense%bz,BSp%mG0)
 977 
 978        !! * Get iq_ibz, and symmetries from iq_bz
 979        !call get_BZ_item(Qmesh_dense,iq_bz,qbz,iq_ibz,isym_q,itim_q)
 980 
 981        !if(iq_ibz > 1 .and. ABS(vc_sqrt_qbz - Vcp_dense%vc_sqrt(1,iq_ibz)) > 1.e-3) then
 982        !   write(*,*) "vc_sqrt_qbz = ",vc_sqrt_qbz
 983        !   write(*,*) "Vcp_dense%vc_sqrt(1,iq_ibz) = ",Vcp_dense%vc_sqrt(1,iq_ibz)
 984        !   MSG_ERROR("vcp are not the same !")
 985        !else if(iq_ibz == 1 .and. ABS(vc_sqrt_qbz - SQRT(Vcp_dense%i_sz)) > 1.e-3) then
 986        !   write(*,*) "vc_sqrt_qbz = ",vc_sqrt_qbz
 987        !   write(*,*) "SQRT(Vcp_dense%i_sz) = ",SQRT(Vcp_dense%i_sz)
 988        !   MSG_ERROR("vcp are not the same !")
 989        !end if
 990        !!END DEBUG CHK !
 991      end if
 992 
 993      if(newway) then
 994 
 995        Cmat = czero
 996 
 997        work_coeffs = czero
 998 
 999        do ineighbour = 1,interpolator%nvert
1000 
1001          ! Loop over the (c, v) part of the left transition
1002          do iv = BSp%lomo_spin(spin1),BSp%homo_spin(spin1)
1003            do ic = BSp%lumo_spin(spin1),BSp%humo_spin(spin1)
1004 
1005              it_dense = BSp%vcks2t_interp(iv,ic,ik_dense,spin1)
1006              it_coarse = BSp%vcks2t(iv,ic,ik_coarse,spin1)
1007              ibnd_coarse = (iv-BSp%lomo_spin(spin1))*BSp%maxnbndc+(ic-BSp%lumo_spin(spin1)+1)
1008 
1009              itc = interpolator%corresp(it_coarse,ineighbour,spin1)
1010 
1011              if (any(BSp%interp_mode == [1,3,4])) then
1012 
1013                !work_coeffs(:,:) = hmat(itc,band2it(:),:)
1014                do inb = 1,interpolator%nvert
1015                  work_coeffs(:,inb) = hmat(itc,interpolator%corresp(band2it(:),inb,spin2))
1016                end do
1017 
1018                call hexc_compute_subhinterp(BSp,grid,nbnd_coarse,&
1019 &        interpolator,kdense2div,&
1020 &        work_coeffs,tmp_Cmat,ikp_dense,&
1021 &        interpolator%overlaps(:,:,:,ikp_dense,spin2))
1022 
1023                if(any(BSp%interp_mode == [1,4])) then
1024                  Cmat(ibnd_coarse,:,ineighbour) = tmp_Cmat
1025                else if (BSp%interp_mode == 3) then
1026                  Cmat(ibnd_coarse,:,ineighbour) = -tmp_Cmat
1027                end if
1028              end if
1029 
1030 
1031              if (any(BSp%interp_mode == [2,3])) then
1032                !work_coeffs(:,:) = acoeffs(itc,band2it(:),:)
1033                do inb = 1,interpolator%nvert
1034                  work_coeffs(:,inb) = acoeffs(itc,interpolator%corresp(band2it(:),inb,spin2))
1035                end do
1036 
1037                call hexc_compute_subhinterp(BSp,grid,nbnd_coarse,&
1038 &        interpolator,kdense2div,&
1039 &        work_coeffs,tmp_Cmat,ikp_dense,&
1040 &        interpolator%overlaps(:,:,:,ikp_dense,spin2))
1041 
1042                tmp_Cmat = tmp_Cmat * (vc_sqrt_qbz**2)
1043                Cmat(ibnd_coarse,:,ineighbour) = Cmat(ibnd_coarse,:,ineighbour) + tmp_Cmat
1044              end if
1045 
1046 
1047              if (any(BSp%interp_mode == [2,3])) then
1048                !work_coeffs(:,:) = bcoeffs(itc,band2it(:),:)
1049                do inb = 1,interpolator%nvert
1050                  work_coeffs(:,inb) = bcoeffs(itc,interpolator%corresp(band2it(:),inb,spin2))
1051                end do
1052 
1053                call hexc_compute_subhinterp(BSp,grid,nbnd_coarse,&
1054 &        interpolator,kdense2div,&
1055 &        work_coeffs,tmp_Cmat,ikp_dense,&
1056 &        interpolator%overlaps(:,:,:,ikp_dense,spin2))
1057 
1058                tmp_Cmat = tmp_Cmat * (vc_sqrt_qbz)
1059                Cmat(ibnd_coarse,:,ineighbour) = Cmat(ibnd_coarse,:,ineighbour) + tmp_Cmat
1060              end if
1061 
1062 
1063              if (any(BSp%interp_mode == [2,3])) then
1064                !work_coeffs(:,:) = ccoeffs(itc,band2it(:),:)
1065                do inb = 1,interpolator%nvert
1066                  work_coeffs(:,inb) = ccoeffs(itc,interpolator%corresp(band2it(:),inb,spin2))
1067                end do
1068 
1069                call hexc_compute_subhinterp(BSp,grid,nbnd_coarse,&
1070 &        interpolator,kdense2div,&
1071 &        work_coeffs,tmp_Cmat,ikp_dense,&
1072 &        interpolator%overlaps(:,:,:,ikp_dense,spin2))
1073 
1074                Cmat(ibnd_coarse,:,ineighbour) = Cmat(ibnd_coarse,:,ineighbour) + tmp_Cmat
1075              end if
1076            end do ! ic
1077          end do ! iv
1078        end do ! ineighbour
1079 
1080      else
1081        ! Loop over the (c, v) part of the left transition
1082        do iv = BSp%lomo_spin(spin1),BSp%homo_spin(spin1)
1083          do ic = BSp%lumo_spin(spin1),BSp%humo_spin(spin1)
1084 
1085            it_dense = BSp%vcks2t_interp(iv,ic,ik_dense,spin1)
1086            it_coarse = BSp%vcks2t(iv,ic,ik_coarse,spin1)
1087            ibnd_coarse = (iv-BSp%lomo_spin(spin1))*BSp%maxnbndc+(ic-BSp%lumo_spin(spin1)+1)
1088 
1089            ! Loop over the (c', v') part of the right transition
1090            do ivp = BSp%lomo_spin(spin2),BSp%homo_spin(spin2)
1091              do icp = BSp%lumo_spin(spin2),BSp%humo_spin(spin2)
1092 
1093                itp_dense = BSp%vcks2t_interp(ivp,icp,ikp_dense,spin2)
1094                itp_coarse = BSp%vcks2t(ivp,icp,ikp_coarse,spin2)
1095                ibndp_coarse = (ivp-Bsp%lomo_spin(spin2))*BSp%maxnbndc+(icp-BSp%lumo_spin(spin2)+1)
1096                ! Now we now it_dense, and itp_dense
1097 
1098                idivp = kdense2div(ikp_dense)
1099 
1100                btemp = czero; ctemp = czero
1101 
1102                ! MG TODO: This way of looping is not optimal
1103                do ineighbour = 1,interpolator%nvert
1104                  itc = interpolator%corresp(it_coarse,ineighbour,spin1)
1105 
1106                  do ineighbourp = 1,interpolator%nvert
1107 
1108                    do iv1 = BSp%lomo_spin(spin2),BSp%homo_spin(spin2)
1109                      do ic1 = BSp%lumo_spin(spin2),BSp%humo_spin(spin2)
1110 
1111                        ibndp_coarse1 = (iv1-BSp%lomo_spin(spin2))*BSp%maxnbndc+(ic1-BSp%lumo_spin(spin2)+1)
1112                        indwithnb = (ineighbourp-1)*nbnd_coarse+ibndp_coarse1
1113                        itp_coarse1 = BSp%vcks2t(iv1,ic1,ikp_coarse,spin2)
1114                        corresp_ind = interpolator%corresp(itp_coarse1,ineighbourp,spin2)
1115 
1116                        select case (BSp%interp_mode)
1117                        case (1,4)
1118                          interpolator%btemp(indwithnb) = hmat(itc,corresp_ind)
1119                        case (2)
1120                          interpolator%btemp(indwithnb) = acoeffs(itc,corresp_ind)*(vc_sqrt_qbz**2) &
1121 &                                         + bcoeffs(itc,corresp_ind)*(vc_sqrt_qbz) &
1122 &                                         + ccoeffs(itc,corresp_ind)
1123                        case (3)
1124                          ! Diff between divergence and hmat
1125                          interpolator%btemp(indwithnb) = acoeffs(itc,corresp_ind)*(vc_sqrt_qbz**2) &
1126 &                                         + bcoeffs(itc,corresp_ind)*(vc_sqrt_qbz) &
1127 &                                         + ccoeffs(itc,corresp_ind) &
1128 &                                         - hmat(itc,corresp_ind)
1129                        case default
1130                          MSG_ERROR("Wrong Bsp%interp_mode")
1131                        end select
1132 
1133                        ctemp(indwithnb) = &
1134 &                        interpolator%overlaps(iv1,ivp,ineighbourp,ikp_dense,spin2) &
1135 &                        * GWPC_CONJG(interpolator%overlaps(ic1,icp,ineighbourp,ikp_dense,spin2)) &
1136 &                        *interpolator%interp_factors(ineighbourp,idivp)
1137                      end do ! ic1
1138                    end do !iv1
1139 
1140                  end do !ineighbourp
1141                  Cmat(ibnd_coarse,ibndp_coarse,ineighbour) = xdotc(interpolator%nvert*nbnd_coarse,&
1142 &                                   ctemp,1,btemp,1)
1143                end do !ineighbour
1144 
1145              end do !icp
1146            end do !ivp
1147 
1148          end do !ic
1149        end do !iv
1150 
1151      end if
1152 
1153      do iv = BSp%lomo_spin(spin1),BSp%homo_spin(spin1)
1154        do ic = BSp%lumo_spin(spin1),BSp%humo_spin(spin1)
1155          it_dense = BSp%vcks2t_interp(iv,ic,ik_dense,spin1)
1156          it_coarse = BSp%vcks2t(iv,ic,ik_coarse,spin1)
1157          ibnd_coarse = (iv-BSp%lomo_spin(spin1))*BSp%maxnbndc+(ic-BSp%lumo_spin(spin1)+1)
1158 
1159          idiv = kdense2div(ik_dense)
1160 
1161          do ivp = BSp%lomo_spin(spin2),BSp%homo_spin(spin2)
1162            do icp = BSp%lumo_spin(spin2),BSp%humo_spin(spin2)
1163              itp_dense = BSp%vcks2t_interp(ivp,icp,ikp_dense,spin2)
1164              itp_coarse = BSp%vcks2t(ivp,icp,ikp_coarse,spin2)
1165              ibndp_coarse = (ivp-Bsp%lomo_spin(spin2))*BSp%maxnbndc+(icp-BSp%lumo_spin(spin2)+1)
1166 
1167              ! Hermicity
1168              if(use_herm .and. it_dense < itp_dense) then
1169                continue
1170              end if
1171 
1172              !btemp = czero; ctemp = czero
1173 
1174              do ineighbour = 1,interpolator%nvert
1175                do iv1 = BSp%lomo_spin(spin1),BSp%homo_spin(spin1)
1176                  do ic1 = BSp%lumo_spin(spin1),BSp%humo_spin(spin1)
1177                    ibnd_coarse1 = (iv1-BSp%lomo_spin(spin1))*BSp%maxnbndc+(ic1-BSp%lumo_spin(spin1)+1)
1178                    it_dense1 = BSp%vcks2t_interp(iv1,ic1,ik_dense,spin1)
1179                    indwithnb = (ineighbour-1)*nbnd_coarse+ibnd_coarse1
1180 
1181                    btemp(indwithnb) = Cmat(ibnd_coarse1,ibndp_coarse,ineighbour)
1182 
1183                    ctemp(indwithnb) = GWPC_CONJG(interpolator%overlaps(iv1,iv,ineighbour,ik_dense,spin1)) &
1184 &                                    *interpolator%overlaps(ic1,ic,ineighbour,ik_dense,spin1) &
1185 &                                    *interpolator%interp_factors(ineighbour,idiv)
1186                  end do !ic1
1187                end do !iv1
1188              end do !ineighbour
1189 
1190              ! Save interpolated value.
1191              hinterp(it_dense,itp_dense) = xdotc(interpolator%nvert*nbnd_coarse,ctemp,1,btemp,1)
1192              !DOT_PRODUCT(ctemp,btemp)
1193 
1194            end do !icp
1195          end do !ivp
1196 
1197        end do !ic
1198      end do !iv
1199 
1200    end do !ikp
1201  end do !ik
1202 
1203  ! Enforce hermiticity
1204  if(use_herm) then
1205    do itp_dense = 1,BSp%nreh_interp(spin2)
1206      do it_dense = itp_dense,BSp%nreh_interp(spin1)
1207        if(it_dense == itp_dense) then
1208          hinterp(itp_dense,it_dense) = DBLE(hinterp(it_dense,itp_dense))
1209        else
1210          hinterp(itp_dense,it_dense) = CONJG(hinterp(it_dense,itp_dense))
1211        end if
1212      end do
1213    end do
1214  end if
1215 
1216  ABI_FREE(Cmat)
1217  ABI_FREE(band2it)
1218 
1219  if(newway) then
1220    ABI_FREE(tmp_Cmat)
1221    ABI_FREE(work_coeffs)
1222  end if
1223 
1224  nullify(ctemp)
1225  nullify(btemp)
1226  call int_free_work(interpolator)
1227 
1228  hinterp = hinterp*factor
1229 
1230  call timab(696,2,tsec)
1231 
1232 end subroutine hexc_compute_hinterp

m_hexc/hexc_compute_subhinterp [ Functions ]

[ Top ] [ m_hexc ] [ Functions ]

NAME

 hexc_compute_subhinterp

FUNCTION

 Compute the interpolation for work_coeffs in dense mesh

INPUTS

 BSp<excparam>=Parameters for BS run
 grid<double_grid_t>=Double grid info
 nbnd_coarse=Number of bands (= nbndv * nbndc)
 interpolator<interpolator_t>=Interpolation info
 kdense2div=Mapping between dense point and coarse point
 work_coeffs=Coefficients to be interpolated
 ikp_dense=Current kpoint
 overlaps=Wavefunction overlaps

OUTPUT

 Cmat(nbnd_coarse) = Interpolated coefficients

PARENTS

      m_hexc

CHILDREN

      timab

SOURCE

712 subroutine hexc_compute_subhinterp(BSp,grid,nbnd_coarse,&
713 &  interpolator,kdense2div,work_coeffs,Cmat,ikp_dense,overlaps)
714 
715 
716 !This section has been created automatically by the script Abilint (TD).
717 !Do not modify the following lines by hand.
718 #undef ABI_FUNC
719 #define ABI_FUNC 'hexc_compute_subhinterp'
720 !End of the abilint section
721 
722  implicit none
723 
724 !Arguments ------------------------------------
725 !scalars
726  integer,intent(in) :: nbnd_coarse
727  integer,intent(in) :: ikp_dense
728  type(excparam),intent(in) :: BSp
729  type(double_grid_t),intent(in) :: grid
730  type(interpolator_t),target,intent(inout) :: interpolator
731 !arrays
732  integer,intent(in) :: kdense2div(grid%nbz_dense)
733  complex(gwpc),intent(in) :: overlaps(interpolator%mband_coarse,interpolator%mband_dense,interpolator%nvert)
734  complex(dpc),intent(in) :: work_coeffs(nbnd_coarse,interpolator%nvert)
735  complex(dpc),intent(out) :: Cmat(nbnd_coarse)
736 
737 !Local variables ------------------------------
738 !scalars
739  integer,parameter :: spin1=1,spin2=1
740  integer :: iv1,ic1
741  integer :: icp,ivp,idivp,ibndp_coarse,ibndp_coarse1,ineighbourp
742  integer :: indwithnb
743  integer :: lumo2,lomo2,humo2,homo2
744  complex(dpc) :: tmp_val, tmp2, tmp4
745 !arrays
746  complex(dpc),ABI_CONTIGUOUS pointer :: btemp(:),ctemp(:)
747 
748 !*********************************************************************
749 
750  btemp => interpolator%btemp
751  ctemp => interpolator%ctemp
752 
753  btemp = czero
754  ctemp = czero
755 
756  lumo2 = BSp%lumo_spin(spin2)
757  lomo2 = BSp%lomo_spin(spin2)
758  humo2 = BSp%humo_spin(spin2)
759  homo2 = BSp%homo_spin(spin2)
760 
761  Cmat = czero
762 
763  idivp = kdense2div(ikp_dense)
764 
765  do ineighbourp = 1,interpolator%nvert
766 
767    btemp(((ineighbourp-1)*nbnd_coarse+1):(ineighbourp*nbnd_coarse)) = &
768 &       interpolator%interp_factors(ineighbourp,idivp)*work_coeffs(:,ineighbourp)
769 
770  end do !ineighbourp
771 
772  ! Loop over the (c', v') part of the right transition
773  do ivp = lomo2,homo2
774    do icp = lumo2,humo2
775 
776      ibndp_coarse = (ivp-lomo2)*BSp%maxnbndc+(icp-lumo2+1)
777      ! Now we now it_dense, and itp_dense
778 
779      do ineighbourp = 1,interpolator%nvert
780 
781        do iv1 = lomo2, homo2
782          ! BSp%lumo_spin(spin2),BSp%humo_spin(spin2)
783 
784          tmp4 = overlaps(iv1,ivp,ineighbourp)
785 
786          do ic1 = lumo2,humo2
787            tmp2 = GWPC_CONJG(overlaps(ic1,icp,ineighbourp))
788            ! BSp%lomo_spin(spin2),BSp%homo_spin(spin2)
789 
790            ibndp_coarse1 = (iv1-lomo2)*BSp%maxnbndc+(ic1-lumo2+1)
791            indwithnb = (ineighbourp-1)*nbnd_coarse+ibndp_coarse1
792 
793            ctemp(indwithnb) = &
794 &             tmp4 &
795 &            *tmp2
796 
797          end do ! iv1
798        end do !ic1
799 
800      end do !ineighbourp
801 
802      tmp_val = xdotc(interpolator%nvert*nbnd_coarse,ctemp,1,btemp,1)
803      !tmp_val = DOT_PRODUCT(ctemp,btemp)
804 
805      Cmat(ibndp_coarse) = tmp_val
806    end do !ivp
807  end do !icp
808 
809  nullify(btemp)
810  nullify(ctemp)
811 
812 end subroutine hexc_compute_subhinterp

m_hexc/hexc_free [ Functions ]

[ Top ] [ m_hexc ] [ Functions ]

NAME

 hexc_free

FUNCTION

 Destroy the interpolator object in memory

PARENTS

      m_haydock

CHILDREN

      timab

SOURCE

1252 subroutine hexc_free(hexc)
1253 
1254 
1255 !This section has been created automatically by the script Abilint (TD).
1256 !Do not modify the following lines by hand.
1257 #undef ABI_FUNC
1258 #define ABI_FUNC 'hexc_free'
1259 !End of the abilint section
1260 
1261  implicit none
1262 
1263 !Arguments ---------------------------
1264  type(hexc_t),intent(inout) :: hexc
1265 
1266 !Local variables ---------------------
1267 
1268 !*****************************************************************************
1269 
1270  if( associated(hexc%bsp) ) then
1271    nullify(hexc%bsp)
1272  end if
1273 
1274  if( associated(hexc%crystal) ) then
1275    nullify(hexc%crystal)
1276  end if
1277 
1278  if( associated(hexc%kmesh_coarse) ) then
1279    nullify(hexc%kmesh_coarse)
1280  end if
1281 
1282  if( allocated(hexc%hreso) ) then
1283    ABI_FREE(hexc%hreso)
1284  end if
1285 
1286  if( allocated(hexc%hcoup) ) then
1287    ABI_FREE(hexc%hcoup)
1288  end if
1289 
1290  if( allocated(hexc%diag_coarse) ) then
1291    ABI_FREE(hexc%diag_coarse)
1292  end if
1293 
1294 end subroutine hexc_free

m_hexc/hexc_init [ Functions ]

[ Top ] [ m_hexc ] [ Functions ]

NAME

 hexc_init

FUNCTION

 Construct the hexc object

INPUTS

 BSp<excparam>=Parameters of BS
 BS_files<excparam>=Files for BS
 Cryst<crystal_t>=Info on the crystalline structure
 Kmesh_coarse<kmesh_t>=Kmesh info
 Wfd_coarse<wfd_t>=Wavefunction descriptor
 KS_BSt<ebands_t>=Kohn-Sham band structure
 QP_BSt<ebands_t>=Quasi-Particle band structure
 comm=communicator

OUTPUT

 hexc<hexc_t>=Excitonic Hamiltonian

PARENTS

      m_haydock

CHILDREN

      timab

SOURCE

248 subroutine hexc_init(hexc, BSp, BS_files, Cryst, Kmesh_coarse, Wfd_coarse, KS_BSt, QP_BSt, comm)
249 
250 
251 !This section has been created automatically by the script Abilint (TD).
252 !Do not modify the following lines by hand.
253 #undef ABI_FUNC
254 #define ABI_FUNC 'hexc_init'
255 !End of the abilint section
256 
257  implicit none
258 
259 !Arguments ---------------------------
260 !scalars
261  integer,intent(in) :: comm
262  type(hexc_t),intent(inout) :: hexc
263  type(excparam),intent(in),target :: BSp
264  type(excfiles),intent(in),target :: BS_files
265  type(crystal_t),intent(in),target :: Cryst
266  type(kmesh_t),intent(in),target :: Kmesh_coarse
267  type(wfd_t),intent(in),target :: Wfd_coarse
268  type(ebands_t),intent(in),target :: KS_BSt, QP_BSt
269 !arrays
270 
271 !Local variables ---------------------
272 !scalars
273  integer :: ierr,ncid,ncerr
274  integer :: max_r, max_c
275  integer :: hsize
276  integer :: spin, spad, itt ! For diagonal !
277  logical :: is_resonant, diago_is_real, use_mpio=.FALSE.
278  character(len=fnlen) :: hreso_fname, hcoup_fname
279  character(len=500) :: msg
280 !arrays
281  complex(dpc),allocatable :: test(:,:)
282 
283 !*****************************************************************************
284 
285  hexc%bsp => BSp
286  hexc%bs_files => BS_files
287  hexc%crystal => Cryst
288  hexc%kmesh_coarse => Kmesh_coarse
289 
290 
291  hexc%comm = comm
292  hsize = SUM(BSp%nreh)
293 
294  hexc%ks_bst => KS_BSt
295  hexc%qp_bst => QP_BSt
296  hexc%wfd => Wfd_coarse
297  hexc%wfd_coarse => Wfd_coarse
298  hexc%kmesh => Kmesh_coarse
299  hexc%hsize_coarse = hsize
300  hexc%hsize = hsize
301  hexc%nbz = Kmesh_coarse%nbz
302 
303  hexc%nbnd_coarse = BSp%maxnbndv*BSp%maxnbndc
304 
305  ! Divide the columns of the Hamiltonian among the nodes.
306  call xmpi_split_work(hsize,comm,hexc%my_t1,hexc%my_t2,msg,ierr)
307  if (ierr/=0) then
308    MSG_WARNING(msg)
309  end if
310 
311  hexc%my_nt = hexc%my_t2 - hexc%my_t1 + 1
312  ABI_CHECK(hexc%my_nt>0,"found processor with 0 rows")
313 
314  ABI_STAT_MALLOC(hexc%hreso,(hsize,hexc%my_t1:hexc%my_t2), ierr)
315  ABI_CHECK(ierr==0, "out of memory in hreso")
316 
317  ! Read the resonant block from file.
318  if (BS_files%in_hreso /= BSE_NOFILE) then
319    hreso_fname = BS_files%in_hreso
320  else
321    hreso_fname = BS_files%out_hreso
322  end if
323 
324  is_resonant=.TRUE.; diago_is_real=(.not.BSp%have_complex_ene)
325  call exc_read_rcblock(hreso_fname,Bsp,is_resonant,diago_is_real,BSp%nsppol,BSp%nreh,hsize,&
326 &   hexc%my_t1,hexc%my_t2,hexc%hreso,use_mpio,comm)
327 
328  !BEGIN DEBUG
329  if (use_mpio) then
330    MSG_WARNING("Testing MPI-IO routines")
331    ABI_STAT_MALLOC(test,(hsize,hexc%my_t1:hexc%my_t2), ierr)
332    ABI_CHECK(ierr==0, "out of memory in hreso")
333    diago_is_real=(.not.BSp%have_complex_ene)
334    call exc_read_rcblock(hreso_fname,Bsp,is_resonant,diago_is_real,Bsp%nsppol,Bsp%nreh,hsize,&
335 &     hexc%my_t1,hexc%my_t2,test,.FALSE.,comm)
336    test = test-hexc%hreso
337    write(std_out,*)"DEBUG: Diff MPI-IO - Fortran ",MAXVAL(ABS(test))
338    max_r=20; max_c=10
339    write(std_out,*)" **** Testing resonant block **** "
340    call print_arr(test,max_r=max_r,max_c=max_c,unit=std_out)
341    if (BSp%nsppol==2) then
342      write(std_out,*)" **** D down down ****"
343      call print_arr(test(hsize/2+1:,hsize/2+1:),max_r=max_r,max_c=max_c,unit=std_out)
344      write(std_out,*)" **** V up down ****"
345      call print_arr(test(1:hsize/2,hsize/2+1:),max_r=max_r,max_c=max_c,unit=std_out)
346      write(std_out,*)" **** V down up ****"
347      call print_arr(test(hsize/2+1:,1:hsize/2),max_r=max_r,max_c=max_c,unit=std_out)
348    end if
349    ABI_FREE(test)
350  end if
351  !END DEBUG
352 
353  !
354  ! Read coupling block.
355  if (BSp%use_coupling>0) then
356    ABI_CHECK(.not. Bsp%use_interp,"interpolation with coupling not coded!")
357    if (BS_files%in_hcoup /= BSE_NOFILE) then
358      hcoup_fname = BS_files%in_hcoup
359    else
360      hcoup_fname = BS_files%out_hcoup
361    end if
362 
363    ABI_STAT_MALLOC(hexc%hcoup,(hsize,hexc%my_t1:hexc%my_t2), ierr)
364    ABI_CHECK(ierr==0, "out of memory in hcoup")
365    is_resonant=.FALSE.; diago_is_real=.FALSE.
366    call exc_read_rcblock(hcoup_fname,Bsp,is_resonant,diago_is_real,BSp%nsppol,BSp%nreh,hsize,&
367 &     hexc%my_t1,hexc%my_t2,hexc%hcoup,use_mpio,comm)
368    !call symmetrize(hcoup,"ALL")
369 
370    if (use_mpio) then
371      MSG_WARNING("Testing MPI-IO routines")
372      ABI_STAT_MALLOC(test,(hsize,hexc%my_t1:hexc%my_t2), ierr)
373      ABI_CHECK(ierr==0, "out of memory in text")
374      diago_is_real=.FALSE.
375      call exc_read_rcblock(hcoup_fname,Bsp,is_resonant,diago_is_real,BSp%nsppol,Bsp%nreh,hsize,&
376 &       hexc%my_t1,hexc%my_t2,test,.FALSE.,comm)
377      test = test-hexc%hcoup
378      write(std_out,*)"DEBUG: Diff MPI-IO - Fortran ",MAXVAL(ABS(test))
379      max_r=20; max_c=10
380      write(std_out,*)" **** Testing coupling block **** "
381      call print_arr(test,max_r=max_r,max_c=max_c,unit=std_out)
382      if (BSp%nsppol==2) then
383        write(std_out,*)" **** D down down ****"
384        call print_arr(test(hsize/2+1:,hsize/2+1:),max_r=max_r,max_c=max_c,unit=std_out)
385        write(std_out,*)" **** V up down ****"
386        call print_arr(test(1:hsize/2,hsize/2+1:),max_r=max_r,max_c=max_c,unit=std_out)
387        write(std_out,*)" **** V down up ****"
388        call print_arr(test(hsize/2+1:,1:hsize/2),max_r=max_r,max_c=max_c,unit=std_out)
389      end if
390      ABI_FREE(test)
391    end if
392  end if
393 
394  if(BSp%prt_ncham .or. BSp%use_interp) then
395    ! I want to store the diagonal part for future use (printing or interpolation) !
396    ABI_MALLOC(hexc%diag_coarse,(hexc%hsize_coarse))
397    spad=0
398    do spin=1,BSp%nsppol
399      if(spin==2) spad=BSp%nreh(1)
400      do itt=1,BSp%nreh(spin) ! 1 is for spin 1
401        hexc%diag_coarse(spad+itt) = Bsp%Trans(itt,spin)%en
402      end do
403    end do
404 
405    if (BSp%prt_ncham) then
406 #ifdef HAVE_NETCDF
407      ncerr = nctk_open_create(ncid, trim(hexc%BS_files%out_basename)//"_HEXC.nc", xmpi_comm_self)
408      NCF_CHECK_MSG(ncerr, "Creating HEXC file")
409      call exc_ham_ncwrite(ncid, hexc%Kmesh_coarse, hexc%BSp, hexc%hsize_coarse, hexc%BSp%nreh, &
410 &         hexc%BSp%vcks2t,hexc%hreso,hexc%diag_coarse)
411      NCF_CHECK(nf90_close(ncid))
412 #else
413      ABI_UNUSED(ncid)
414 #endif
415    end if
416  end if
417 
418 end subroutine hexc_init

m_hexc/hexc_interp_free [ Functions ]

[ Top ] [ m_hexc ] [ Functions ]

NAME

 hexc_interp_free

FUNCTION

  Destroy the interpolator object in memory

PARENTS

      m_haydock

CHILDREN

      timab

SOURCE

1314 subroutine hexc_interp_free(hexc_i)
1315 
1316 
1317 !This section has been created automatically by the script Abilint (TD).
1318 !Do not modify the following lines by hand.
1319 #undef ABI_FUNC
1320 #define ABI_FUNC 'hexc_interp_free'
1321 !End of the abilint section
1322 
1323  implicit none
1324 
1325 !Arguments ---------------------------
1326  type(hexc_interp_t),intent(inout) :: hexc_i
1327 
1328 !Local variables ---------------------
1329 
1330 !*****************************************************************************
1331 
1332  if( allocated(hexc_i%kdense2div) ) then
1333    ABI_FREE(hexc_i%kdense2div)
1334  end if
1335 
1336  if( allocated(hexc_i%div2kdense) ) then
1337    ABI_FREE(hexc_i%div2kdense)
1338  end if
1339 
1340  if( allocated(hexc_i%diag_dense) ) then
1341    ABI_FREE(hexc_i%diag_dense)
1342  end if
1343 
1344  if( allocated(hexc_i%hinterp) ) then
1345    ABI_FREE(hexc_i%hinterp)
1346  end if
1347 
1348  if( allocated(hexc_i%all_hmat) ) then
1349    ABI_FREE(hexc_i%all_hmat)
1350  end if
1351 
1352  if( allocated(hexc_i%all_acoeffs) ) then
1353    ABI_FREE(hexc_i%all_acoeffs)
1354  end if
1355 
1356  if( allocated(hexc_i%all_bcoeffs) ) then
1357    ABI_FREE(hexc_i%all_bcoeffs)
1358  end if
1359 
1360  if( allocated(hexc_i%all_ccoeffs) ) then
1361    ABI_FREE(hexc_i%all_ccoeffs)
1362  end if
1363 
1364  if( associated(hexc_i%kmesh_dense) ) then
1365    nullify(hexc_i%kmesh_dense)
1366  end if
1367 
1368  if( associated(hexc_i%vcp_dense) ) then
1369    nullify(hexc_i%vcp_dense)
1370  end if
1371 
1372  call interpolator_free(hexc_i%interpolator)
1373 
1374 end subroutine hexc_interp_free

m_hexc/hexc_interp_init [ Functions ]

[ Top ] [ m_hexc ] [ Functions ]

NAME

 hexc_interp_init

FUNCTION

 Construct the hexc_interp object

INPUTS

 hexc<hexc_t>=Excitonic hamiltonian
 Kmesh_dense<kmesh_t>=Kmesh info
 Vcp_dense<vcoul_t>=Dense mesh info about coulomb
 double_grid<double_grid_t>=Link between dense and coarse mesh
 Wfd_dense<wfd_t>=Wavefunction descriptor
 KS_BSt_dense<ebands_t>=Kohn-Sham band structure
 QP_BSt_dense<ebands_t>=Quasi-Particle band structure
 Psps <type(pseudopotential_type)>=variables related to pseudopotentials.
 Pawtab(Cryst%ntypat*usepaw)<pawtab_type>=PAW tabulated starting data.
 comm=communicator

OUTPUT

 hexc_i<hexc_interp_t>=Interpolated excitonic hamiltonian

SIDE EFFECTS

 hexc
   Will be modified so that the size of the problem is full interpolated hamiltonian
 Wfd, Wfd_dense
   The memory might be modified by computing wavefunctions

PARENTS

      m_haydock

CHILDREN

      timab

SOURCE

460 subroutine hexc_interp_init(hexc_i, hexc, m3_width, method, Kmesh_dense, Vcp_dense, &
461 &    double_grid, Wfd_dense, KS_BSt_dense, QP_BSt_dense, Psps, Pawtab)
462 
463 
464 !This section has been created automatically by the script Abilint (TD).
465 !Do not modify the following lines by hand.
466 #undef ABI_FUNC
467 #define ABI_FUNC 'hexc_interp_init'
468 !End of the abilint section
469 
470  implicit none
471 
472 !Arguments ---------------------------
473 !scalars
474  integer,intent(in) :: method
475  real(dp),intent(in) :: m3_width
476  type(hexc_t),intent(inout) :: hexc
477  type(hexc_interp_t),intent(inout) :: hexc_i
478  type(double_grid_t),intent(in),target :: double_grid
479  type(wfd_t),intent(inout),target :: Wfd_dense !, Wfd
480  type(kmesh_t),intent(in),target :: Kmesh_dense
481  type(pseudopotential_type),intent(in) :: Psps
482  type(vcoul_t),intent(in),target :: Vcp_dense
483  type(ebands_t),intent(in),target :: KS_BSt_dense, QP_BSt_dense
484 !arrays
485  type(pawtab_type),intent(in) :: Pawtab(hexc%crystal%ntypat*hexc%Wfd_coarse%usepaw)
486 
487 !Local variables ---------------------
488 !scalars
489  integer,parameter :: spin1 = 1
490  integer :: nsppol,ierr,ii,itt,nproc, my_rank,hsize
491  real(dp),parameter :: threshold = 0.1_dp
492  type(excparam) :: BSp
493  logical :: is_resonant, diago_is_real, use_mpio
494 !arrays
495  character(len=fnlen) :: tmpfname, hreso_fname
496 
497 !*****************************************************************************
498 
499  BSp = hexc%bsp
500  ABI_CHECK(BSp%nsppol == 1,"nsppol > 1 not implemented yet")
501 
502  hsize = hexc%hsize_coarse
503 
504  nproc  = xmpi_comm_size(hexc%comm); my_rank= xmpi_comm_rank(hexc%comm)
505  nsppol = hexc%Bsp%nsppol
506 
507  ABI_CHECK(nproc == 1,"Parallelization not available in interpolation")
508 
509  hexc_i%m3_width = m3_width
510 
511  hexc_i%kmesh_dense => Kmesh_dense
512  hexc_i%vcp_dense => Vcp_dense
513  hexc_i%hsize_dense = SUM(BSp%nreh_interp)
514  hexc%hsize = hexc_i%hsize_dense
515  hexc%nbz = Kmesh_dense%nbz
516 
517  hexc%ks_bst => KS_BSt_dense
518  hexc%qp_bst => QP_BSt_dense
519  hexc%wfd => Wfd_dense
520  hexc%kmesh => Kmesh_dense
521 
522  ! No parallelization !
523  hexc%my_t1 = 1
524  hexc%my_t2 = hexc_i%hsize_dense
525 
526  ! Initialize the interpolator
527  call interpolator_init(hexc_i%interpolator, double_grid, Wfd_dense, hexc%Wfd_coarse, Kmesh_dense, &
528 &    hexc%Kmesh_coarse, hexc%BSp, hexc%crystal, Psps, Pawtab, method)
529 
530  if(BSp%sum_overlaps) then
531    call interpolator_normalize(hexc_i%interpolator)
532  end if
533 
534  ABI_MALLOC(hexc_i%kdense2div,(double_grid%nbz_dense))
535  ABI_MALLOC(hexc_i%div2kdense,(double_grid%nbz_coarse,double_grid%ndiv))
536 
537  call compute_corresp(double_grid,hexc_i%div2kdense,hexc_i%kdense2div)
538 
539  if (any(BSp%interp_mode == [2,3])) then
540    ! Read a, b, c coefficient matrices from file.
541    ! For the time being, we read the full matrix in a temporary array, and
542    ! then we store the data in a form suitable for the interpolation.
543    is_resonant=.TRUE.; diago_is_real=(.not.BSp%have_complex_ene); use_mpio=.FALSE.
544 
545    if (hexc%bs_files%in_hreso /= BSE_NOFILE) then
546      hreso_fname = hexc%bs_files%in_hreso
547    else
548      hreso_fname = hexc%bs_files%out_hreso
549    end if
550 
551    tmpfname = hreso_fname; ii = LEN_TRIM(hreso_fname)
552 
553    ! TODO: Write new IO routines to read MPI-distributed data in a format suitable for the interpolation
554    tmpfname(ii-2:ii+1) = 'ABSR'
555    ABI_STAT_MALLOC(hexc_i%all_acoeffs,(hsize,hsize), ierr)
556    ABI_CHECK(ierr==0, "out of memory in all_acoeffs")
557    call exc_read_rcblock(tmpfname,Bsp,is_resonant,diago_is_real,nsppol,BSp%nreh,hsize,1,hsize,&
558 &     hexc_i%all_acoeffs,use_mpio,hexc%comm)
559 
560    tmpfname(ii-2:ii+1) = 'BBSR'
561    ABI_STAT_MALLOC(hexc_i%all_bcoeffs,(hsize,hsize), ierr)
562    ABI_CHECK(ierr==0, "out of memory in all_bcoeffs")
563    call exc_read_rcblock(tmpfname,Bsp,is_resonant,diago_is_real,nsppol,BSp%nreh,hsize,1,hsize,&
564 &     hexc_i%all_bcoeffs,use_mpio,hexc%comm)
565 
566    tmpfname(ii-2:ii+1) = 'CBSR'
567    ABI_STAT_MALLOC(hexc_i%all_ccoeffs,(hsize,hsize), ierr)
568    ABI_CHECK(ierr==0, "out of memory in all_ccoeffs")
569    call exc_read_rcblock(tmpfname,Bsp,is_resonant,diago_is_real,nsppol,BSp%nreh,hsize,1,hsize,&
570 &     hexc_i%all_ccoeffs,use_mpio,hexc%comm)
571  end if
572 
573  ! Compute overlaps & compute all hmat
574  ABI_STAT_MALLOC(hexc_i%all_hmat,(hsize,hsize), ierr)
575  ABI_CHECK(ierr==0, "out of memory in all_hmat")
576 
577  hexc_i%all_hmat(:,:) = hexc%hreso(:,:)
578 
579  do itt=1,hsize
580    hexc_i%all_hmat(itt,itt) = hexc_i%all_hmat(itt,itt) - hexc%diag_coarse(itt)
581  end do
582 
583  ! I don't need the diag_coarse any more
584  ABI_FREE(hexc%diag_coarse)
585 
586  ! Compute diagonal part of the dense Ham
587  ABI_MALLOC(hexc_i%diag_dense,(hexc_i%hsize_dense))
588  do itt=1,BSp%nreh_interp(spin1) ! 1 is for spin 1
589    hexc_i%diag_dense(itt) = Bsp%Trans_interp(itt,spin1)%en
590  end do
591 
592 end subroutine hexc_interp_init

m_hexc/hexc_interp_matmul [ Functions ]

[ Top ] [ m_hexc ] [ Functions ]

NAME

 hexc_interp_matmul

FUNCTION

 Compute matrix-vector product Hmat * phi by interpolating coarse Hmat

INPUTS

  BSp<type(excparam)>=Parameters defining the BS calculation
  hsize_coarse = Size of the coarse hamiltonian
  hsize_dense = Size of the dense hamiltonian
  hmat(hsize_coarse,hsize_coarse,8) = coarse hamiltonian
  phi(hsize_dense) = ket on which apply the matrix
  grid <double_grid_t> = Correspondence between coarse and dense k-mesh.
  nbnd_coarse = Total number of bands
  interpolator<interpolator_t> = Interpolator
  div2kdense = Mapping from coarse and division -> dense mesh
  kdense2div = Mapping from dense mesh -> division

OUTPUT

  hphi(hsize_dense) = Interp(hmat)*phi

PARENTS

      m_hexc

CHILDREN

      timab

SOURCE

1409 subroutine hexc_interp_matmul(BSp,hsize_coarse,hsize_dense,hmat,phi,hphi,grid,&
1410 &   nbnd_coarse,interpolator,div2kdense,kdense2div)
1411 
1412 
1413 !This section has been created automatically by the script Abilint (TD).
1414 !Do not modify the following lines by hand.
1415 #undef ABI_FUNC
1416 #define ABI_FUNC 'hexc_interp_matmul'
1417 !End of the abilint section
1418 
1419  implicit none
1420 
1421 !Arguments ------------------------------------
1422 !scalars
1423  integer,intent(in) :: hsize_coarse,hsize_dense,nbnd_coarse !,ntrans
1424  type(excparam),intent(in) :: BSp
1425  type(double_grid_t),intent(in) :: grid
1426  type(interpolator_t),intent(in) :: interpolator
1427 !arrays
1428  integer,intent(in) :: div2kdense(grid%nbz_coarse,grid%ndiv), kdense2div(grid%nbz_dense)
1429  complex(dpc),intent(in) :: phi(hsize_dense)
1430  complex(dpc),intent(in) :: hmat(hsize_coarse,hsize_coarse)
1431  complex(dpc),intent(inout) :: hphi(hsize_dense)
1432 
1433 !Local variables ------------------------------
1434 !scalars
1435  integer :: itt,ik_dense,ik_coarse,it_coarse
1436  integer :: ic,iv,iv1,ic1, ibnd_coarse
1437  integer :: ibnd_coarse1
1438  integer :: ineighbour,idense,ikpt
1439  integer :: my_k1,my_k2,ind_with_nb,is, is1
1440  real(dp) :: factor
1441  complex(dpc) :: tmp
1442  logical,parameter :: use_blas=.True.
1443 !arrays
1444  integer :: allindices(nbnd_coarse)
1445  complex(dpc) :: allp(hsize_coarse,interpolator%nvert), test(hsize_coarse)
1446  complex(dpc) :: ophi(grid%nbz_dense,interpolator%nvert,nbnd_coarse)
1447  complex(dpc),allocatable :: b(:), c(:),A(:,:)
1448  complex(dpc),allocatable :: tmp_array(:), tmp_array2(:,:)
1449 
1450 !************************************************************************
1451 
1452  factor = one/grid%ndiv
1453 
1454  !hphi = czero
1455 
1456  ! Outer index : k point in the dense zone
1457  ! Sum over vc
1458  ! Index of result : k point in the dense zone, v2,c2,neighbour
1459 
1460  ! Parallelization on nbz in the coarse mesh !
1461  my_k1 = 1
1462  my_k2 = grid%nbz_coarse
1463 
1464  ABI_MALLOC(A,(interpolator%nvert*nbnd_coarse,nbnd_coarse))
1465  ABI_MALLOC(b,(nbnd_coarse))
1466  ABI_MALLOC(c,(interpolator%nvert*nbnd_coarse))
1467 
1468  c = czero; ophi = czero
1469 
1470  do ik_dense = 1,grid%nbz_dense
1471    ! if( ik_dense is not in my set of k-points)
1472    !   ! continue
1473    !
1474    do is1 = 1, BSp%nsppol
1475      do iv1 = BSp%lomo_spin(is1),Bsp%homo_spin(is1)
1476        do ic1 = BSp%lumo_spin(is1),Bsp%humo_spin(is1)
1477          ibnd_coarse = (iv1-BSp%lomo_spin(is1))*BSp%maxnbndc+(ic1-BSp%lumo_spin(is1)+1)
1478          itt = BSp%vcks2t_interp(iv1,ic1,ik_dense,is1)
1479          allindices(ibnd_coarse) = itt
1480        end do !ic1
1481      end do !iv1
1482    end do !is1
1483 
1484    b(:) = phi(allindices(:))
1485 
1486    do is = 1, BSp%nsppol
1487      do iv = BSp%lomo_spin(is),Bsp%homo_spin(is)
1488        do ic = BSp%lumo_spin(is),Bsp%humo_spin(is)
1489          ibnd_coarse = (iv-BSp%lomo_spin(is))*BSp%maxnbndc+(ic-BSp%lumo_spin(is)+1)
1490          idense = Bsp%vcks2t_interp(iv,ic,ik_dense,is)
1491 
1492          do ineighbour = 1,interpolator%nvert
1493            ind_with_nb = (ineighbour-1)*(nbnd_coarse)+ibnd_coarse
1494 
1495            !A(ind_with_nb,:) = overlaps(allindices(:),ibnd_coarse,ineighbour)
1496 
1497            ! Should be optimized !!!
1498            do iv1 = BSp%lomo_spin(is),Bsp%homo_spin(is)
1499              do ic1 = BSp%lumo_spin(is),Bsp%humo_spin(is)
1500                ibnd_coarse1 = (iv1-BSp%lomo_spin(is))*BSp%maxnbndc+(ic1-BSp%lumo_spin(is)+1)
1501                A(ind_with_nb,ibnd_coarse1) = GWPC_CONJG(interpolator%overlaps(iv,iv1,ineighbour,ik_dense,is)) &
1502 &                                          *interpolator%overlaps(ic,ic1,ineighbour,ik_dense,is)
1503              end do !ic1
1504            end do !iv1
1505          end do !ineighbour
1506        end do !ic
1507      end do !iv
1508    end do !is
1509 
1510    if(use_blas) then
1511      call xgemv('N',interpolator%nvert*nbnd_coarse,nbnd_coarse,cone,A,interpolator%nvert*nbnd_coarse,b,1,czero,c,1)
1512    else
1513      c = MATMUL(A,b)
1514    end if
1515 
1516    do is = 1, BSp%nsppol
1517      do iv = BSp%lomo_spin(is),BSp%homo_spin(is)
1518        do ic = BSp%lumo_spin(is),BSp%humo_spin(is)
1519          ibnd_coarse = (iv-BSp%lomo_spin(is))*BSp%maxnbndc+(ic-BSp%lumo_spin(is)+1)
1520          do ineighbour = 1,interpolator%nvert
1521            ind_with_nb = (ineighbour-1)*(nbnd_coarse)+ibnd_coarse
1522            ophi(ik_dense,ineighbour,ibnd_coarse) = c(ind_with_nb)
1523          end do !ineighbour
1524        end do !ic
1525      end do !iv
1526    end do !is
1527 
1528  end do !ik_dense
1529 
1530  ABI_FREE(A)
1531  ABI_FREE(b)
1532  ABI_FREE(c)
1533 
1534  !call xmpi_sum_(ophi,comm,ierr)
1535 
1536  ! Outer index : k,v,c in the coarse zone, ineighbour
1537  ! Sum over all k-dense relative to one coarse point
1538  ! Index of result : k,v,c in the coarse zone, ineighbour
1539 
1540  ABI_MALLOC(b,(grid%ndiv))
1541  ABI_MALLOC(c,(grid%ndiv))
1542 
1543  allp = czero
1544 
1545  do is = 1, BSp%nsppol
1546    do ineighbour = 1,interpolator%nvert
1547 
1548      do it_coarse = 1, BSp%nreh(is)
1549        ibnd_coarse = (Bsp%trans(it_coarse,is)%v-BSp%lomo_spin(is))*BSp%maxnbndc+&
1550 &            (BSp%Trans(it_coarse,is)%c-BSp%lumo_spin(is)+1)
1551        ik_coarse = BSp%trans(it_coarse,is)%k
1552        !b(:) = interp_factors(it_coarse,ineighbour,:)
1553        b(:) = interpolator%interp_factors(ineighbour,:)
1554        !c(:) = ophi(indices(it_coarse,:),ineighbour,ibnd_coarse)
1555        c(:) = ophi(div2kdense(ik_coarse,:),ineighbour,ibnd_coarse)
1556        tmp = DOT_PRODUCT(b,c)
1557        allp(it_coarse,ineighbour) = tmp
1558      end do
1559 
1560    end do
1561  end do
1562 
1563  !call xmpi_sum_(allp,comm,ierr)
1564 
1565  ABI_FREE(b)
1566  ABI_FREE(c)
1567 
1568  ABI_MALLOC(tmp_array,(hsize_coarse))
1569  ABI_MALLOC(tmp_array2,(hsize_coarse,hsize_coarse))
1570  tmp_array(:) = czero
1571  tmp_array2(:,:) = czero
1572 
1573  test = czero
1574 
1575  ! Second step : Multiplication by hmat
1576  do ineighbour = 1,interpolator%nvert
1577    if(use_blas) then
1578      !call xgemv('N',hsize_coarse,hsize_coarse,cone,factor*(hmat(:,:,ineighbour)),hsize_coarse,allp(:,ineighbour),1,czero,tmp_array,1)
1579      !tmp_array2 = hmat(:,:,ineighbour)
1580      tmp_array2 = hmat(:,interpolator%corresp(:,ineighbour,1)) ! 1 is for spin
1581      tmp_array2 = factor*tmp_array2
1582      call xgemv('N',hsize_coarse,hsize_coarse,cone,tmp_array2,hsize_coarse,allp(:,ineighbour),1,czero,tmp_array,1)
1583      test = test + tmp_array
1584    else
1585      test = test+MATMUL(factor*(hmat(:,interpolator%corresp(:,ineighbour,1))),allp(:,ineighbour))
1586    end if
1587  end do
1588 
1589  ABI_FREE(tmp_array)
1590  ABI_FREE(tmp_array2)
1591 
1592  ! Outer index : ineighbour
1593  ! Sum over all v c
1594  ! Index of result : ineighbour, k_dense, v,c
1595  ABI_MALLOC(A,(nbnd_coarse,nbnd_coarse))
1596  ABI_MALLOC(b,(nbnd_coarse))
1597  ABI_MALLOC(c,(nbnd_coarse))
1598  c = czero
1599 
1600  do ineighbour = 1,interpolator%nvert
1601    do ik_dense = 1,grid%nbz_dense
1602 
1603      do is1 = 1, Bsp%nsppol
1604        do iv1 = Bsp%lomo_spin(is1),Bsp%homo_spin(is1)
1605          do ic1 = BSp%lumo_spin(is1), Bsp%humo_spin(is1)
1606            ibnd_coarse = (iv1-BSp%lomo_spin(is1))*BSp%maxnbndc+(ic1-BSp%lumo_spin(is1)+1)
1607 
1608            ik_coarse = grid%dense_to_coarse(ik_dense)
1609            itt = BSp%vcks2t(iv1,ic1,ik_coarse,is1)
1610            b(ibnd_coarse) = test(interpolator%corresp(itt,ineighbour,is1))
1611          end do ! ic1
1612        end do ! iv1
1613      end do ! is1
1614 
1615      do is = 1, BSp%nsppol
1616        do iv = BSp%lomo_spin(is),Bsp%homo_spin(is)
1617          do ic = BSp%lumo_spin(is),BSp%humo_spin(is)
1618            ibnd_coarse = (iv-BSp%lomo_spin(is))*Bsp%maxnbndc+(ic-BSp%lumo_spin(is)+1)
1619            idense = BSp%vcks2t_interp(iv,ic,ik_dense,is)
1620 
1621            !A(ibnd_coarse,:) = CONJG(overlaps(idense,:,ineighbour))
1622 
1623            ! Should be optimized !!!
1624            do iv1 = BSp%lomo_spin(is),Bsp%homo_spin(is)
1625              do ic1 = BSp%lumo_spin(is),Bsp%humo_spin(is)
1626                ibnd_coarse1 = (iv1-BSp%lomo_spin(is))*BSp%maxnbndc+(ic1-BSp%lumo_spin(is)+1)
1627                A(ibnd_coarse,ibnd_coarse1) = (interpolator%overlaps(iv1,iv,ineighbour,ik_dense,is)) &
1628 &                                          *GWPC_CONJG(interpolator%overlaps(ic1,ic,ineighbour,ik_dense,is))
1629              end do !ic1
1630            end do !iv1
1631          end do ! ic
1632        end do !iv
1633      end do !is
1634 
1635      if(use_blas) then
1636        call xgemv('N',nbnd_coarse,nbnd_coarse,cone,A,nbnd_coarse,b,1,czero,c,1)
1637      else
1638        c = MATMUL(A,b)
1639      end if
1640 
1641      do is = 1, BSp%nsppol
1642        do iv = BSp%lomo_spin(is),Bsp%homo_spin(is)
1643          do ic = BSp%lumo_spin(is),BSp%humo_spin(is)
1644            ibnd_coarse = (iv-BSp%lomo_spin(is))*BSp%maxnbndc+(ic-BSp%lumo_spin(is)+1)
1645            idense = Bsp%vcks2t_interp(iv,ic,ik_dense,is)
1646            !ophi(ik_dense,ineighbour,ibnd_coarse) = c(idense)
1647            ophi(ik_dense,ineighbour,ibnd_coarse) = c(ibnd_coarse)
1648          end do
1649        end do
1650      end do
1651 
1652    end do ! ik_dense
1653  end do ! ineighbour
1654 
1655  !call xmpi_sum_(ophi,comm,ierr)
1656 
1657  ABI_FREE(A)
1658  ABI_FREE(b)
1659  ABI_FREE(c)
1660 
1661  ! Outer indices : it_dense
1662  ! Sum over neighbours
1663  ! Index of result : it_dense (ik,ic,iv dense)
1664 
1665  ABI_MALLOC(b,(interpolator%nvert))
1666  ABI_MALLOC(c,(interpolator%nvert))
1667 
1668  do is = 1, BSp%nsppol
1669    do itt = 1,BSp%nreh_interp(is)
1670     ! From itt -> ik_ibz,ic,iv
1671     ik_dense = BSp%Trans_interp(itt,is)%k
1672     ic = BSp%Trans_interp(itt,is)%c
1673     iv = BSp%Trans_interp(itt,is)%v
1674 
1675     ! From ik_ibz in the dense mesh -> indices_dense
1676     ik_coarse = grid%dense_to_coarse(ik_dense)
1677     it_coarse = BSp%vcks2t(iv,ic,ik_coarse,is)
1678 
1679     ibnd_coarse = (iv-BSp%lomo_spin(is))*BSp%maxnbndc+(ic-BSp%lumo_spin(is)+1)
1680 
1681     ikpt = kdense2div(ik_dense)
1682     !ikpt = -1
1683     !do ix = 1,grid%ndiv
1684     !  if (indices(it_coarse,ix) == ik_dense) then
1685     !    ikpt = ix
1686     !    exit
1687     !  end if
1688     !end do
1689     !ABI_CHECK(ikpt/=-1,"Cannot find ik_dense")
1690 
1691     !b = interp_factors(it_coarse,:,ikpt)
1692     b = interpolator%interp_factors(:,ikpt)
1693     c =  ophi(ik_dense,:,ibnd_coarse)
1694 
1695     hphi(itt) = hphi(itt) + xdotc(interpolator%nvert, b, 1, c, 1)
1696    end do
1697  end do
1698 
1699  ABI_FREE(b)
1700  ABI_FREE(c)
1701 
1702 end subroutine hexc_interp_matmul

m_hexc/hexc_interp_t [ Types ]

[ Top ] [ m_hexc ] [ Types ]

NAME

 hexc_interp_t

FUNCTION

  Store information about interpolation of excitonic hamiltonian

SOURCE

150  type,public :: hexc_interp_t
151 
152  !scalars
153     integer :: hsize_dense
154     ! Size of the dense hamiltonian
155 
156     real(dp) :: m3_width
157     ! Width of the region where M3 is applied instead of M1
158 
159     type(interpolator_t) :: interpolator
160     ! Interpolator containing overlaps and interpolation info
161 
162     ! Pointers to datatypes that are already in memory
163     type(kmesh_t),pointer :: kmesh_dense => null()
164     ! kmesh of the dense mesh
165 
166     type(vcoul_t),pointer :: vcp_dense => null()
167     ! coulomb interaction on the dense mesh
168 
169  !arrays
170     integer,allocatable :: kdense2div(:)
171     ! kdense2div(nbz_dense)
172     ! Index of kpoint -> Index of division
173 
174     integer,allocatable :: div2kdense(:,:)
175     ! div2kdense(nbz_coarse,ndiv)
176     ! Index of kpoint coarse + Index of division -> Index of kdense
177 
178     complex(dpc),allocatable :: diag_dense(:)
179     ! diag_dense(hsize_dense)
180     ! Diagonal part of the dense hamiltonian
181 
182     complex(dpc),allocatable :: hinterp(:,:)
183     ! hinterp(hsize_dense,hsize_dense)
184     ! Interpolated hamiltonian
185 
186     complex(dpc),allocatable :: all_hmat(:,:)
187     ! all_hmat,(hsize,hsize))
188     ! Coarse excitonic matrix in a format suitable for interpolation in k-space
189 
190     complex(dpc),allocatable :: all_acoeffs(:,:)
191     ! all_acoeffs(hsize,hsize))
192     ! a coefficients in a format suitable for interpolation in k-space
193 
194     complex(dpc),allocatable :: all_bcoeffs(:,:)
195     ! all_bcoeffs(hsize,hsize))
196     ! b coefficients in a format suitable for interpolation in k-space
197 
198     complex(dpc),allocatable :: all_ccoeffs(:,:)
199     ! all_ccoeffs(hsize,hsize))
200     ! c coefficients in a format suitable for interpolation in k-space
201 
202  end type hexc_interp_t

m_hexc/hexc_matmul_elphon [ Functions ]

[ Top ] [ m_hexc ] [ Functions ]

NAME

 hexc_matmul_elphon

FUNCTION

 Compute H | \psi > + E_{elphon} | \psi >

INPUTS

 hexc<hexc_t> = Excitonic hamiltonian
 phi = Input ket
 ep_renorm = vector with electron-phonon renorms
 op = 'N' for H | psi >, 'C' for H^\dagger | psi >

OUTPUT

 hphi = hreso * phi + ep_renorm * phi

PARENTS

      m_haydock

CHILDREN

      timab

SOURCE

1806 subroutine hexc_matmul_elphon(hexc, phi, hphi, op, ep_renorm)
1807 
1808 
1809 !This section has been created automatically by the script Abilint (TD).
1810 !Do not modify the following lines by hand.
1811 #undef ABI_FUNC
1812 #define ABI_FUNC 'hexc_matmul_elphon'
1813 !End of the abilint section
1814 
1815  implicit none
1816 
1817 !Arguments ---------------------------
1818  type(hexc_t),intent(in) :: hexc
1819  character,intent(in) :: op
1820  complex(dpc),intent(in) :: phi(hexc%my_nt)
1821  complex(dpc),intent(out) :: hphi(hexc%hsize)
1822  complex(dpc),intent(in) :: ep_renorm(hexc%hsize)
1823 
1824 !Local variables ---------------------
1825  integer :: ierr
1826  real(dp) :: tsec(2)
1827 
1828 !*****************************************************************************
1829 
1830  call timab(697,1,tsec)
1831 
1832  if(hexc%BSp%use_interp) then
1833    MSG_ERROR('Not yet implemented with interpolation !')
1834  else ! No interpolation
1835    ! As our matrix is hermitian (hreso), we should always use 'N' here (it is stored column-wise !)
1836    call xgemv('N',hexc%hsize,hexc%my_nt,cone,hexc%hreso,hexc%hsize,phi,1,czero,hphi,1)
1837 
1838    !!! ep_renorm is stored on each cpu
1839    if (op == 'N') then
1840      hphi(hexc%my_t1:hexc%my_t2) = hphi(hexc%my_t1:hexc%my_t2) + ep_renorm(hexc%my_t1:hexc%my_t2) * phi
1841    else if(op == 'C') then
1842      hphi(hexc%my_t1:hexc%my_t2) = hphi(hexc%my_t1:hexc%my_t2) + CONJG(ep_renorm(hexc%my_t1:hexc%my_t2)) * phi
1843    end if
1844    call xmpi_sum(hphi,hexc%comm,ierr)
1845  end if
1846 
1847  call timab(697,2,tsec)
1848 
1849 end subroutine hexc_matmul_elphon

m_hexc/hexc_matmul_full [ Functions ]

[ Top ] [ m_hexc ] [ Functions ]

NAME

 hexc_matmul_full

FUNCTION

 Compute H | \psi >

INPUTS

 hexc<hexc_t> = Excitonic hamiltonian
 hexc_i<hexc_interp_t> = Interpolated hamiltonian
 phi = Input ket
 parity = -1 or +1 parameter

OUTPUT

 hphi = hreso * phi + parity * hcoup * CONJ(phi)

PARENTS

      m_haydock

CHILDREN

      timab

SOURCE

1878 subroutine hexc_matmul_full(hexc, hexc_i, phi, hphi, parity)
1879 
1880 
1881 !This section has been created automatically by the script Abilint (TD).
1882 !Do not modify the following lines by hand.
1883 #undef ABI_FUNC
1884 #define ABI_FUNC 'hexc_matmul_full'
1885 !End of the abilint section
1886 
1887  implicit none
1888 
1889 !Arguments ---------------------------
1890  integer,intent(in) :: parity
1891  type(hexc_t),intent(in) :: hexc
1892  type(hexc_interp_t),intent(in) :: hexc_i
1893  complex(dpc),intent(in) :: phi(hexc%hsize)
1894  complex(dpc),intent(out) :: hphi(hexc%hsize)
1895 
1896 !Local variables ---------------------
1897  real(dp) :: tsec(2)
1898 
1899 !*****************************************************************************
1900 
1901  call timab(697,1,tsec)
1902 
1903  ABI_UNUSED(hexc_i%hsize_dense)
1904 
1905  if(hexc%BSp%use_interp) then
1906    MSG_ERROR("Coupling is not yet implemented with interpolation")
1907  else ! No interpolation
1908    hphi = MATMUL(hexc%hreso,phi) + parity * MATMUL(hexc%hcoup,CONJG(phi))
1909  end if
1910 
1911  call timab(697,2,tsec)
1912 
1913 end subroutine hexc_matmul_full

m_hexc/hexc_matmul_tda [ Functions ]

[ Top ] [ m_hexc ] [ Functions ]

NAME

 hexc_matmul_tda

FUNCTION

 Compute H | \psi >

INPUTS

 hexc<hexc_t> = Excitonic hamiltonian
 hexc_i<hexc_interp_t> = Interpolated excitonic hamiltonian
 phi = Input ket

OUTPUT

 hphi = hreso * phi

PARENTS

      m_haydock

CHILDREN

      timab

SOURCE

1730 subroutine hexc_matmul_tda(hexc, hexc_i, phi, hphi)
1731 
1732 
1733 !This section has been created automatically by the script Abilint (TD).
1734 !Do not modify the following lines by hand.
1735 #undef ABI_FUNC
1736 #define ABI_FUNC 'hexc_matmul_tda'
1737 !End of the abilint section
1738 
1739  implicit none
1740 
1741 !Arguments ---------------------------
1742  type(hexc_t),intent(in) :: hexc
1743  type(hexc_interp_t),intent(in) :: hexc_i
1744  complex(dpc),intent(in) :: phi(hexc%hsize)
1745  complex(dpc),intent(out) :: hphi(hexc%hsize)
1746 
1747 !Local variables ---------------------
1748  integer :: ierr
1749  real(dp) :: tsec(2)
1750 
1751 !*****************************************************************************
1752 
1753  call timab(697,1,tsec)
1754 
1755  if(hexc%BSp%use_interp) then
1756    hphi = hexc_i%diag_dense * phi
1757 
1758    if (any(hexc%BSp%interp_mode == [2,3,4])) then
1759      ! hphi = hphi + MATMUL(hinterp,phi)
1760      call xgemv('N',hexc_i%hsize_dense,hexc_i%hsize_dense,cone,hexc_i%hinterp,hexc_i%hsize_dense,phi,1,cone,hphi,1)
1761      if (any(hexc%BSp%interp_mode == [2,4])) then
1762        call timab(697,2,tsec)
1763        return ! We are done
1764      end if
1765    end if
1766 
1767    call hexc_interp_matmul(hexc%bsp, hexc%hsize_coarse, hexc_i%hsize_dense, hexc_i%all_hmat, phi, hphi, &
1768 &     hexc_i%interpolator%double_grid,hexc%nbnd_coarse, hexc_i%interpolator, hexc_i%div2kdense, hexc_i%kdense2div)
1769 
1770  else ! No interpolation
1771    call xgemv('N',hexc%hsize,hexc%my_nt,cone,hexc%hreso,hexc%hsize,phi,1,czero,hphi,1)
1772    call xmpi_sum(hphi,hexc%comm,ierr)
1773  end if
1774 
1775  call timab(697,2,tsec)
1776 
1777 end subroutine hexc_matmul_tda