TABLE OF CONTENTS
ABINIT/m_wfutils [ Modules ]
NAME
m_wfutils
FUNCTION
parameters and function for wave functions copy
COPYRIGHT
Copyright (C) 2001-2024 ABINIT group (CS,GZ,FB) This file is distributed under the terms of the GNU General Public License, see ~ABINIT/Infos/copyright or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_wfutils 23 24 use defs_basis 25 use m_abicore 26 use m_errors 27 28 use m_time, only : timab 29 30 implicit none 31 32 private 33 34 public :: setWFParameter ! Definition of wave functions parameters for copy functions. 35 public :: wfcopy ! Copy of wave functions arrays. 36 37 ! Global variables 38 integer,save,public ABI_PROTECTED :: x_cplx ! fortran data type for wave functions arrays real or complex. 39 integer,save,public ABI_PROTECTED :: x_me_g0 ! processor number 40 integer,save,public ABI_PROTECTED :: x_npw_k ! number of plane waves at this k point 41 integer,save,public ABI_PROTECTED :: x_nspinor! number of spinorial components of the wavefunctions on current proc 42 integer,save,public ABI_PROTECTED :: x_icg ! shift to be applied on the location of data in the array cg 43 integer,save,public ABI_PROTECTED :: x_igsc ! shift to be applied on the location of data in the array gsc 44 integer,save,public ABI_PROTECTED :: x_blocksize 45 46 contains
m_wfutils/setWFParameter [ Functions ]
[ Top ] [ m_wfutils ] [ Functions ]
NAME
setWFParameter
FUNCTION
Initialize wave functions parameters for copy functions
INPUTS
cplx=fortran data type for wave functions arrays real or complex me_g0=1 if this node treats G=0. npw_k=number of plane waves at this k point nspinor=number of spinorial components of the wavefunctions on current proc icg=shift to be applied on the location of data in the array cg igsc=shift to be applied on the location of data in the array gsc blocksize=size of blocks
SOURCE
68 subroutine setWFParameter(cplx,me_g0,npw_k,nspinor,icg,igsc,blocksize) 69 70 implicit none 71 72 !Arguments ------------------------------------ 73 integer, intent(in) :: cplx,me_g0,npw_k,nspinor 74 integer, intent(in) :: icg,igsc,blocksize 75 76 ! ********************************************************************* 77 78 ! Copy values in global variables 79 x_cplx=cplx 80 x_me_g0=me_g0 81 x_npw_k=npw_k 82 x_nspinor=nspinor 83 x_icg=icg 84 x_igsc=igsc 85 x_blocksize=blocksize 86 87 end subroutine setWFParameter
m_wfutils/wfcopy [ Functions ]
[ Top ] [ m_wfutils ] [ Functions ]
NAME
wfcopy
FUNCTION
copy of wave functions arrays. called from lobpcg in REAL or COMPLEX wave functions.
INPUTS
direction=copy direction 'D' direct (global array to local) 'I' indirect (local to global) size=number of elements tsrc=source array incsrc=size increment for tsrc array tdest=destination array incdest=increment for tdest array blockiter=number of block iteration in case REAL iblock=block index indtype=indexation type in array withbbloc=apply block on band for each block
TODO
Split the two cases so that we can avoid the array descriptors.
SOURCE
172 subroutine wfcopy(direction,size,tsrc,incsrc,tdest,incdest,blockiter,iblock,indtype,& 173 & withbbloc,timopt,tim_wfcopy) ! optional arguments 174 175 implicit none 176 177 !Arguments ------------------------------------ 178 character(len=1), intent(in) :: direction 179 integer, intent(in) :: size,incsrc,incdest 180 integer, intent(in) :: blockiter,iblock 181 character(len=1), intent(in) :: indtype 182 logical, optional, intent(in) :: withbbloc 183 integer, intent(in), optional :: timopt,tim_wfcopy 184 !arrays 185 real(dp), DEV_CONTARRD intent(in) :: tsrc(:,:) 186 real(dp), DEV_CONTARRD intent(inout) :: tdest(:,:) 187 188 !Local variables ------------------------------------ 189 integer,parameter :: ndat1=1 190 logical :: bblock=.false. 191 integer :: lig,g1,g2,vectsize,rvectsize 192 integer :: blocksize,bblocksize,iblocksize,iband 193 real(dp) :: factor 194 real(dp) :: tsec(2) 195 196 ! ********************************************************************* 197 198 if (present(tim_wfcopy).and.present(timopt)) then 199 if(abs(timopt)==3) then 200 call timab(tim_wfcopy,1,tsec) 201 end if 202 end if 203 204 if (present(withbbloc)) bblock=withbbloc 205 206 if (indtype == 'C') then 207 lig=x_icg 208 else if (indtype == 'S') then 209 lig=x_igsc 210 else 211 lig=0 212 endif 213 214 rvectsize=x_npw_k*x_nspinor 215 if (x_me_g0 == 1) then 216 vectsize=2*rvectsize-1 217 else 218 vectsize=2*rvectsize 219 endif 220 if (x_cplx == 2) vectsize=x_npw_k*x_nspinor 221 222 blocksize = x_blocksize 223 bblocksize=(iblock-1)*blocksize 224 225 if (direction == 'D') then 226 227 if (x_cplx == 1) then 228 ! Pack real and imag part. 229 factor=sqrt(two) 230 do iblocksize=1,blockiter 231 iband=iblocksize 232 if (bblock) then 233 iband=iblocksize+bblocksize 234 endif 235 if (x_me_g0 == 1) then 236 tdest(1 ,iblocksize)=tsrc(1,wfindex(iband,indtype)) 237 g1 = wfindex(iband, indtype)+1 238 g2 = wfindex(iband+1,indtype)-1 239 call dcopy(rvectsize-1,tsrc(1,g1:g2),1,tdest(2:rvectsize,iblocksize),1) 240 tdest(2:rvectsize,iblocksize) = factor * tdest(2:rvectsize,iblocksize) 241 242 call dcopy(rvectsize-1,tsrc(2,g1:g2),1,tdest(rvectsize+1:vectsize,iblocksize),1) 243 tdest(rvectsize+1:vectsize,iblocksize) = factor * tdest(rvectsize+1:vectsize,iblocksize) 244 ! MG FIXME: Here gfortran4.9 allocates temporary arrays due to factor. 245 else 246 g1 = wfindex(iband, indtype) 247 g2 = wfindex(iband+1,indtype)-1 248 call dcopy(rvectsize,tsrc(1,g1:g2),1,tdest(1:rvectsize,iblocksize),1) 249 tdest(1:rvectsize,iblocksize) = factor * tdest(1:rvectsize,iblocksize) 250 251 call dcopy(rvectsize,tsrc(2,g1:g2),1,tdest(rvectsize+1:vectsize,iblocksize),1) 252 tdest(rvectsize+1:vectsize,iblocksize) = factor * tdest(rvectsize+1:vectsize,iblocksize) 253 end if 254 end do 255 else 256 if (indtype == 'C') then 257 if (bblock) then 258 g1 = vectsize*((iblock-1)*blocksize)+lig+1 259 g2 = vectsize*(iblock*blocksize)+lig 260 call zcopy(size,tsrc(:,g1:g2),incsrc,tdest(:,1:blocksize),incdest) 261 else 262 g1 = lig+1 263 g2 = vectsize*((iblock-1)*blocksize)+lig 264 call zcopy(size,tsrc(:,g1:g2),incsrc,tdest(:,1:bblocksize),incdest) 265 end if 266 else if (indtype == 'S') then 267 g1 = lig+1 268 g2 = vectsize*((iblock-1)*blocksize)+lig 269 call zcopy(size,tsrc(:,g1:g2),incsrc,tdest(:,1:bblocksize),incdest) 270 else 271 call zcopy(size,tsrc,incsrc,tdest,incdest) 272 endif 273 end if 274 275 else if (direction == 'I') then 276 if (x_cplx == 1) then 277 ! Unpack real and imag blocks. 278 factor=one/sqrt(two) 279 do iblocksize=1,blockiter 280 iband=iblocksize 281 if ( bblock ) then 282 iband=iblocksize+(iblock-1)*blocksize 283 end if 284 if (x_me_g0 == 1) then 285 tdest(1,wfindex(iband,indtype))=tsrc(1,iblocksize) 286 tdest(2,wfindex(iband,indtype))=zero 287 288 g1 = wfindex(iband,indtype)+1 289 g2 = wfindex(iband+1,indtype)-1 290 call dcopy(rvectsize-1,tsrc(2:rvectsize,iblocksize),1,tdest(1,g1:g2),1) 291 tdest(1,g1:g2) = factor * tdest(1,g1:g2) 292 293 call dcopy(rvectsize-1,tsrc(rvectsize+1:vectsize,iblocksize),1,tdest(2,g1:g2),1) 294 tdest(2,g1:g2) = factor * tdest(2,g1:g2) 295 else 296 g1 = wfindex(iband,indtype) 297 g2 = wfindex(iband+1,indtype)-1 298 call dcopy(rvectsize,tsrc(1:rvectsize,iblocksize),1,tdest(1,g1:g2),1) 299 tdest(1,g1:g2) = factor * tdest(1,g1:g2) 300 call dcopy(rvectsize,tsrc(rvectsize+1:vectsize,iblocksize),1,tdest(2,g1:g2),1) 301 tdest(2,g1:g2) = factor * tdest(2,g1:g2) 302 end if 303 end do 304 else 305 if (indtype == 'C') then 306 g1 = vectsize*((iblock-1)*blocksize)+lig+1 307 g2 = vectsize*(iblock*blocksize)+lig 308 call zcopy(size,tsrc,1,tdest(:,g1:g2),1) 309 else if ( indtype == 'S' ) then 310 g1 = vectsize*((iblock-1)*blocksize)+lig+1 311 g2 = vectsize*(iblock*blocksize)+lig 312 call zcopy(size,tsrc,1,tdest(:,g1:g2),1) 313 else 314 call zcopy(size,tsrc,1,tdest,1) 315 end if 316 end if 317 else 318 ABI_ERROR("Wrong direction: "//trim(direction)) 319 endif 320 321 if (present(tim_wfcopy).and.present(timopt)) then 322 if(abs(timopt)==3) then 323 call timab(tim_wfcopy,2,tsec) 324 end if 325 end if 326 327 end subroutine wfcopy