TABLE OF CONTENTS
ABINIT/m_wfutils [ Modules ]
NAME
m_wfutils
FUNCTION
parameters and function for wave functions copy
COPYRIGHT
Copyright (C) 2001-2018 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
17 #if defined HAVE_CONFIG_H 18 #include "config.h" 19 #endif 20 21 #include "abi_common.h" 22 23 module m_wfutils 24 25 use defs_basis 26 use m_profiling_abi 27 use m_errors 28 29 use m_cgtools, only : cg_to_reim, cg_from_reim 30 31 implicit none 32 33 private 34 35 public :: setWFParameter ! Definition of wave functions parameters for copy functions. 36 public :: wfcopy ! Copy of wave functions arrays. 37 38 ! Global variables 39 integer,save,public ABI_PROTECTED :: x_cplx ! fortran data type for wave functions arrays real or complex. 40 integer,save,public ABI_PROTECTED :: x_me_g0 ! processor number 41 integer,save,public ABI_PROTECTED :: x_npw_k ! number of plane waves at this k point 42 integer,save,public ABI_PROTECTED :: x_nspinor! number of spinorial components of the wavefunctions on current proc 43 integer,save,public ABI_PROTECTED :: x_icg ! shift to be applied on the location of data in the array cg 44 integer,save,public ABI_PROTECTED :: x_igsc ! shift to be applied on the location of data in the array gsc 45 integer,save,public ABI_PROTECTED :: x_blocksize 46 47 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
PARENTS
lobpcgwf
CHILDREN
SOURCE
74 subroutine setWFParameter(cplx,me_g0,npw_k,nspinor,icg,igsc,blocksize) 75 76 77 !This section has been created automatically by the script Abilint (TD). 78 !Do not modify the following lines by hand. 79 #undef ABI_FUNC 80 #define ABI_FUNC 'setWFParameter' 81 !End of the abilint section 82 83 implicit none 84 85 !Arguments ------------------------------------ 86 integer, intent(in) :: cplx,me_g0,npw_k,nspinor 87 integer, intent(in) :: icg,igsc,blocksize 88 89 ! ********************************************************************* 90 91 ! Copy values in global variables 92 x_cplx=cplx 93 x_me_g0=me_g0 94 x_npw_k=npw_k 95 x_nspinor=nspinor 96 x_icg=icg 97 x_igsc=igsc 98 x_blocksize=blocksize 99 100 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.
PARENTS
lobpcgwf
CHILDREN
dcopy, zcopy
SOURCE
219 subroutine wfcopy(direction,size,tsrc,incsrc,tdest,incdest,blockiter,iblock,indtype,& 220 & withbbloc,timopt,tim_wfcopy) ! optional arguments 221 222 223 !This section has been created automatically by the script Abilint (TD). 224 !Do not modify the following lines by hand. 225 #undef ABI_FUNC 226 #define ABI_FUNC 'wfcopy' 227 use interfaces_18_timing 228 !End of the abilint section 229 230 implicit none 231 232 !Arguments ------------------------------------ 233 character(len=1), intent(in) :: direction 234 integer, intent(in) :: size,incsrc,incdest 235 integer, intent(in) :: blockiter,iblock 236 character(len=1), intent(in) :: indtype 237 logical, optional, intent(in) :: withbbloc 238 integer, intent(in), optional :: timopt,tim_wfcopy 239 !arrays 240 real(dp), DEV_CONTARRD intent(in) :: tsrc(:,:) 241 real(dp), DEV_CONTARRD intent(inout) :: tdest(:,:) 242 243 !Local variables ------------------------------------ 244 integer,parameter :: ndat1=1 245 logical :: bblock=.false. 246 integer :: lig,g1,g2,vectsize,rvectsize 247 integer :: blocksize,bblocksize,iblocksize,iband 248 real(dp) :: factor 249 real(dp) :: tsec(2) 250 251 ! ********************************************************************* 252 253 if (present(tim_wfcopy).and.present(timopt)) then 254 if(abs(timopt)==3) then 255 call timab(tim_wfcopy,1,tsec) 256 end if 257 end if 258 259 if (present(withbbloc)) bblock=withbbloc 260 261 if (indtype == 'C') then 262 lig=x_icg 263 else if (indtype == 'S') then 264 lig=x_igsc 265 else 266 lig=0 267 endif 268 269 rvectsize=x_npw_k*x_nspinor 270 if (x_me_g0 == 1) then 271 vectsize=2*rvectsize-1 272 else 273 vectsize=2*rvectsize 274 endif 275 if (x_cplx == 2) vectsize=x_npw_k*x_nspinor 276 277 blocksize = x_blocksize 278 bblocksize=(iblock-1)*blocksize 279 280 if (direction == 'D') then 281 282 if (x_cplx == 1) then 283 ! Pack real and imag part. 284 factor=sqrt(two) 285 do iblocksize=1,blockiter 286 iband=iblocksize 287 if (bblock) then 288 iband=iblocksize+bblocksize 289 endif 290 if (x_me_g0 == 1) then 291 tdest(1 ,iblocksize)=tsrc(1,wfindex(iband,indtype)) 292 g1 = wfindex(iband, indtype)+1 293 g2 = wfindex(iband+1,indtype)-1 294 call dcopy(rvectsize-1,tsrc(1,g1:g2),1,tdest(2:rvectsize,iblocksize),1) 295 tdest(2:rvectsize,iblocksize) = factor * tdest(2:rvectsize,iblocksize) 296 297 call dcopy(rvectsize-1,tsrc(2,g1:g2),1,tdest(rvectsize+1:vectsize,iblocksize),1) 298 tdest(rvectsize+1:vectsize,iblocksize) = factor * tdest(rvectsize+1:vectsize,iblocksize) 299 ! MG FIXME: Here gfortran4.9 allocates temporary arrays due to factor. 300 !call cg_to_reim(rvectsize-1,ndat1,tsrc(1:,g1:),factor,tdest(2:,iblocksize)) 301 else 302 g1 = wfindex(iband, indtype) 303 g2 = wfindex(iband+1,indtype)-1 304 call dcopy(rvectsize,tsrc(1,g1:g2),1,tdest(1:rvectsize,iblocksize),1) 305 tdest(1:rvectsize,iblocksize) = factor * tdest(1:rvectsize,iblocksize) 306 307 call dcopy(rvectsize,tsrc(2,g1:g2),1,tdest(rvectsize+1:vectsize,iblocksize),1) 308 tdest(rvectsize+1:vectsize,iblocksize) = factor * tdest(rvectsize+1:vectsize,iblocksize) 309 !call cg_to_reim(rvectsize,ndat1,tsrc(1:,g1:),factor,tdest(1:,iblocksize)) 310 end if 311 end do 312 else 313 if (indtype == 'C') then 314 if (bblock) then 315 g1 = vectsize*((iblock-1)*blocksize)+lig+1 316 g2 = vectsize*(iblock*blocksize)+lig 317 call zcopy(size,tsrc(:,g1:g2),incsrc,tdest(:,1:blocksize),incdest) 318 else 319 g1 = lig+1 320 g2 = vectsize*((iblock-1)*blocksize)+lig 321 call zcopy(size,tsrc(:,g1:g2),incsrc,tdest(:,1:bblocksize),incdest) 322 end if 323 else if (indtype == 'S') then 324 g1 = lig+1 325 g2 = vectsize*((iblock-1)*blocksize)+lig 326 call zcopy(size,tsrc(:,g1:g2),incsrc,tdest(:,1:bblocksize),incdest) 327 else 328 call zcopy(size,tsrc,incsrc,tdest,incdest) 329 endif 330 end if 331 332 else if (direction == 'I') then 333 if (x_cplx == 1) then 334 ! Unpack real and imag blocks. 335 factor=one/sqrt(two) 336 do iblocksize=1,blockiter 337 iband=iblocksize 338 if ( bblock ) then 339 iband=iblocksize+(iblock-1)*blocksize 340 end if 341 if (x_me_g0 == 1) then 342 tdest(1,wfindex(iband,indtype))=tsrc(1,iblocksize) 343 tdest(2,wfindex(iband,indtype))=zero 344 345 g1 = wfindex(iband,indtype)+1 346 g2 = wfindex(iband+1,indtype)-1 347 call dcopy(rvectsize-1,tsrc(2:rvectsize,iblocksize),1,tdest(1,g1:g2),1) 348 tdest(1,g1:g2) = factor * tdest(1,g1:g2) 349 350 call dcopy(rvectsize-1,tsrc(rvectsize+1:vectsize,iblocksize),1,tdest(2,g1:g2),1) 351 tdest(2,g1:g2) = factor * tdest(2,g1:g2) 352 !call cg_from_reim(rvectsize-1,ndat1,tsrc(2:,iblocksize),factor,tdest(1,g1:)) 353 else 354 g1 = wfindex(iband,indtype) 355 g2 = wfindex(iband+1,indtype)-1 356 call dcopy(rvectsize,tsrc(1:rvectsize,iblocksize),1,tdest(1,g1:g2),1) 357 tdest(1,g1:g2) = factor * tdest(1,g1:g2) 358 call dcopy(rvectsize,tsrc(rvectsize+1:vectsize,iblocksize),1,tdest(2,g1:g2),1) 359 tdest(2,g1:g2) = factor * tdest(2,g1:g2) 360 !call cg_from_reim(rvectsize,ndat1,tsrc(1:,iblocksize),factor,tdest(1:,g1)) 361 end if 362 end do 363 else 364 if (indtype == 'C') then 365 g1 = vectsize*((iblock-1)*blocksize)+lig+1 366 g2 = vectsize*(iblock*blocksize)+lig 367 call zcopy(size,tsrc,1,tdest(:,g1:g2),1) 368 else if ( indtype == 'S' ) then 369 g1 = vectsize*((iblock-1)*blocksize)+lig+1 370 g2 = vectsize*(iblock*blocksize)+lig 371 call zcopy(size,tsrc,1,tdest(:,g1:g2),1) 372 else 373 call zcopy(size,tsrc,1,tdest,1) 374 end if 375 end if 376 else 377 MSG_ERROR("Wrong direction: "//trim(direction)) 378 endif 379 380 if (present(tim_wfcopy).and.present(timopt)) then 381 if(abs(timopt)==3) then 382 call timab(tim_wfcopy,2,tsec) 383 end if 384 end if 385 386 end subroutine wfcopy