TABLE OF CONTENTS
ABINIT/cgq_builder [ Functions ]
NAME
cgq_builder
FUNCTION
This routine locates cgq for efield calculations, especially for // case
COPYRIGHT
Copyright (C) 2003-2018 ABINIT group This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt.
INPUTS
berryflag = logical flag determining use of electric field variables cg(2,mcg)=planewave coefficients of wavefunctions. dtset <type(dataset_type)>=all input variables for this dataset ikpt=index of current k kpt ikpt_loc=index of k point on current processor (see vtorho.F90) isspol=value of spin polarization currently treated me_distrb=current value from spaceComm_distrb (see vtorho.F90) mcg=size of wave-functions array (cg) =mpw*nspinor*mband*mkmem*nsppol mcgq=size of cgq array (see vtorho.F90) mkgq=size of pwnsfacq array (see vtorho.F90) my_nspinor=nspinor value determined by current // set up nband_k=number of bands at each k point nproc_distrb=nproc from spaceComm_distrb (see vtorho.F90) npwarr(nkpt)=number of planewaves in basis at this k point pwnsfac(2,pwind_alloc) = phase factors for non-symmorphic translations (see initberry.f) pwind_alloc = first dimension of pwind spaceComm_distrb=comm_cell from mpi_enreg
OUTPUT
cgq(2,mcgq)=planewave coefficients of wavenfunctions adjacent to cg at ikpt pwnsfacq(2,mkgq)=phase factors for non-symmorphic translations for cg's adjacent to cg(ikpt)
SIDE EFFECTS
Input/Output dtefield <type(efield_type)> = efield variables mpi_enreg=information about MPI parallelization
TODO
NOTES
PARENTS
vtorho
CHILDREN
timab,xmpi_recv,xmpi_send
SOURCE
58 #if defined HAVE_CONFIG_H 59 #include "config.h" 60 #endif 61 62 #include "abi_common.h" 63 64 subroutine cgq_builder(berryflag,cg,cgq,dtefield,dtset,ikpt,ikpt_loc,isppol,mcg,mcgq,& 65 & me_distrb,mkgq,mpi_enreg,my_nspinor,nband_k,nproc_distrb,& 66 & npwarr,pwnsfac,pwnsfacq,pwind_alloc,spaceComm_distrb) 67 68 use defs_basis 69 use defs_abitypes 70 use m_xmpi 71 use m_errors 72 use m_efield 73 use m_profiling_abi 74 75 !This section has been created automatically by the script Abilint (TD). 76 !Do not modify the following lines by hand. 77 #undef ABI_FUNC 78 #define ABI_FUNC 'cgq_builder' 79 use interfaces_18_timing 80 !End of the abilint section 81 82 implicit none 83 84 !Arguments ------------------------------------ 85 integer,intent(in) :: ikpt,ikpt_loc,isppol,me_distrb,mcg,mcgq,mkgq,my_nspinor,nband_k 86 integer,intent(in) :: nproc_distrb,pwind_alloc,spaceComm_distrb 87 logical,intent(in) :: berryflag 88 type(dataset_type), intent(in) :: dtset 89 type(efield_type), intent(inout) :: dtefield 90 type(MPI_type), intent(in) :: mpi_enreg 91 !arrays 92 integer,intent(in) :: npwarr(dtset%nkpt) 93 real(dp),intent(in) :: cg(2,mcg),pwnsfac(2,pwind_alloc) 94 real(dp),intent(out) :: cgq(2,mcgq),pwnsfacq(2,mkgq) 95 96 !Local variables ------------------------- 97 !scalars 98 integer :: count,count1,icg1,icg2,dest,his_source 99 integer :: idir,ierr,ifor,ikg1,ikg2,ikptf,ikpt1f,ikpt1i 100 integer :: jkpt,jkpt1i,jkptf,jkpt1f,jsppol,my_source,npw_k1,tag 101 !arrays 102 integer,allocatable :: flag_send(:,:), flag_receive(:) 103 real(dp) :: tsec(2) 104 real(dp),allocatable :: buffer(:,:) 105 106 ! ************************************************************************* 107 108 !DEBUG 109 !write(std_out,'(a)')'cgq_builder enter' 110 !DEBUG 111 112 if (mcgq==0.or.mkgq==0) return 113 114 call timab(983,1,tsec) 115 116 !Test compatbility of berryflag 117 if (berryflag) then 118 ABI_ALLOCATE(flag_send,(0:nproc_distrb-1,dtefield%fnkpt)) 119 end if 120 ABI_ALLOCATE(flag_receive,(dtset%nkpt)) 121 flag_send(:,:) = 0 122 flag_receive(:) = 0 123 124 if (berryflag) ikptf = dtefield%i2fbz(ikpt) 125 126 do idir = 1, 3 127 128 ! skip idir values for which efield_dot(idir) = 0 129 if (berryflag .and. abs(dtefield%efield_dot(idir)) < tol12 ) cycle 130 131 do ifor = 1, 2 132 133 if(berryflag) then 134 dtefield%sflag(:,ikpt + dtset%nkpt*(isppol - 1),ifor,idir) = 0 135 ikpt1f = dtefield%ikpt_dk(ikptf,ifor,idir) 136 ikpt1i = dtefield%indkk_f2ibz(ikpt1f,1) 137 end if 138 139 npw_k1 = npwarr(ikpt1i) 140 count = npw_k1*my_nspinor*nband_k 141 my_source = mpi_enreg%proc_distrb(ikpt1i,1,isppol) 142 143 do dest = 0, nproc_distrb-1 144 145 if ((dest==me_distrb).and.(ikpt_loc <= dtset%mkmem)) then 146 ! I am dest and have something to do 147 148 if ( my_source == me_distrb ) then 149 ! I am destination and source 150 151 if(berryflag) then 152 ikg1 = dtefield%fkgindex(ikpt1f) 153 ikg2 = dtefield%cgqindex(3,ifor+2*(idir-1),ikpt+(isppol-1)*dtset%nkpt) 154 icg1 = dtefield%cgindex(ikpt1i,isppol) 155 icg2 = dtefield%cgqindex(2,ifor+2*(idir-1),ikpt+(isppol-1)*dtset%nkpt) 156 end if 157 158 pwnsfacq(:,ikg2 + 1:ikg2 + npw_k1) = pwnsfac(:,ikg1 + 1:ikg1 + npw_k1) 159 cgq(:,icg2 + 1:icg2 + count) = cg(:,icg1 + 1:icg1 + count) 160 161 else ! I am the destination but not the source -> receive 162 ! receive pwnsfacq 163 if(berryflag) then 164 tag = ikpt1f + (isppol - 1)*dtefield%fnkpt 165 ikg1 = dtefield%cgqindex(3,ifor+2*(idir-1),ikpt+(isppol-1)*dtset%nkpt) 166 end if 167 ABI_ALLOCATE(buffer,(2,npw_k1)) 168 call xmpi_recv(buffer,my_source,tag,spaceComm_distrb,ierr) 169 pwnsfacq(:,ikg1+1:ikg1+npw_k1) = buffer(:,1:npw_k1) 170 ABI_DEALLOCATE(buffer) 171 172 ! receive cgq if necessary 173 if(flag_receive(ikpt1i) == 0) then 174 ABI_ALLOCATE(buffer,(2,count)) 175 tag = ikpt1i + (isppol - 1)*dtset%nkpt 176 call xmpi_recv(buffer,my_source,tag,spaceComm_distrb,ierr) 177 if(berryflag) icg1 = dtefield%cgqindex(2,ifor+2*(idir-1),ikpt+(isppol-1)*dtset%nkpt) 178 cgq(:,icg1+1:icg1+count) = buffer(:,1:count) 179 ABI_DEALLOCATE(buffer) 180 flag_receive(ikpt1i) = 1 181 end if ! end if flag_receive == 0 182 end if ! end tasks if I am the destination 183 184 else if (ikpt_loc <= mpi_enreg%mkmem(dest)) then ! dest != me and the dest has a k-point to treat 185 186 ! jkpt is the kpt which is being treated by dest (in ibz) 187 ! jsppol is his isppol 188 jkpt = mpi_enreg%kpt_loc2ibz_sp(dest, ikpt_loc,1) 189 jsppol = mpi_enreg%kpt_loc2ibz_sp(dest, ikpt_loc,2) 190 191 if(jkpt > 0 .and. jsppol > 0) then 192 193 if(berryflag) then 194 jkptf = dtefield%i2fbz(jkpt) 195 jkpt1f = dtefield%ikpt_dk(jkptf,ifor,idir) 196 jkpt1i = dtefield%indkk_f2ibz(jkpt1f,1) 197 end if 198 his_source = mpi_enreg%proc_distrb(jkpt1i,1,jsppol) 199 200 if (his_source == me_distrb) then 201 202 ! send 203 ! pwnsfacq 204 if(berryflag) then 205 ikg1 = dtefield%fkgindex(jkpt1f) 206 tag = jkpt1f + (jsppol - 1)*dtefield%fnkpt 207 end if 208 count1 = npwarr(jkpt1i) 209 ABI_ALLOCATE(buffer,(2,count1)) 210 buffer(:,1:count1) = pwnsfac(:,ikg1+1:ikg1+count1) 211 call xmpi_send(buffer,dest,tag,spaceComm_distrb,ierr) 212 ABI_DEALLOCATE(buffer) 213 214 ! send cgq if necessary 215 if(flag_send(dest, jkpt1i)==0) then 216 if(berryflag) icg1 = dtefield%cgindex(jkpt1i,jsppol) 217 tag = jkpt1i + (jsppol - 1)*dtset%nkpt 218 count1 = npwarr(jkpt1i)*nband_k*my_nspinor 219 ABI_ALLOCATE(buffer,(2,count1)) 220 buffer(:,1:count1) = cg(:,icg1+1:icg1+count1) 221 call xmpi_send(buffer,dest,tag,spaceComm_distrb,ierr) 222 ABI_DEALLOCATE(buffer) 223 flag_send(dest, jkpt1i)=1 224 end if ! if send cgq 225 226 end if ! end check that his_source == me 227 228 end if ! end check on jkpt > 0 and jsppol > 0 229 230 end if ! end check on me = dest else if me != dest 231 232 end do ! end loop over dest = 0, nproc-1 233 234 end do !end loop over ifor 235 236 end do !end loop over idir 237 238 call timab(983,2,tsec) 239 240 ABI_DEALLOCATE(flag_send) 241 ABI_DEALLOCATE(flag_receive) 242 243 !DEBUG 244 !write(std_out,'(a)')'cgq_builder exit' 245 !END_DEBUG 246 247 end subroutine cgq_builder