TABLE OF CONTENTS
ABINIT/pareigocc [ Functions ]
NAME
pareigocc
FUNCTION
This subroutine transmit to all processors, using MPI : - the eigenvalues and, - if ground-state, the occupation numbers (In fact, in the present status of the routine, occupation numbers are NOT transmitted) transmit_occ = 2 is used in case the occ should be transmitted. Yet the code is not already written.
COPYRIGHT
Copyright (C) 2000-2017 ABINIT group (XG, AR) 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
formeig=format of eigenvalues (0 for GS, 1 for RF) localrdwf=(for parallel case) if 1, the eig and occ initial values are local to each machine, if 0, they are on proc me=0. mband=maximum number of bands of the output wavefunctions mpi_enreg=information about MPI parallelization nband(nkpt*nsppol)=desired number of bands at each k point nkpt=number of k points nsppol=1 for unpolarized, 2 for spin-polarized, output wf file processors, Warning : defined only when paralbd=1 transmit_occ/=2 transmit only eigenvalues, =2 for transmission of occ also (yet transmit_occ=2 is not safe or finished at all)
OUTPUT
(see side effects)
SIDE EFFECTS
eigen(mband*nkpt*nsppol)=eigenvalues (input or init to large number), (Ha) occ(mband*nkpt*nsppol)=occupation (input or init to 0.0) NOT USED NOW
NOTES
* The case paralbd=1 with formeig=0 is implemented, but not yet used. * The transmission of occ is not activated yet ! * The routine takes the eigenvalues in the eigen array on one of the processors that possess the wavefunctions, and transmit it to all procs. If localrdwf==0, me=0 has the full array at start, If localrdwf==1, the transfer might be more complex. * This routine should not be used for RF wavefunctions, since it does not treat the eigenvalues as a matrix.
PARENTS
newkpt,wfsinp
CHILDREN
timab,xmpi_bcast,xmpi_sum
SOURCE
63 #if defined HAVE_CONFIG_H 64 #include "config.h" 65 #endif 66 67 #include "abi_common.h" 68 69 subroutine pareigocc(eigen,formeig,localrdwf,mpi_enreg,mband,nband,nkpt,nsppol,occ,transmit_occ) 70 71 use defs_basis 72 use defs_abitypes 73 use m_profiling_abi 74 use m_xmpi 75 76 !This section has been created automatically by the script Abilint (TD). 77 !Do not modify the following lines by hand. 78 #undef ABI_FUNC 79 #define ABI_FUNC 'pareigocc' 80 use interfaces_18_timing 81 use interfaces_32_util 82 !End of the abilint section 83 84 implicit none 85 86 !Arguments ------------------------------------ 87 !scalars 88 integer,intent(in) :: formeig,localrdwf,mband,nkpt,nsppol,transmit_occ 89 type(MPI_type),intent(in) :: mpi_enreg 90 !arrays 91 integer,intent(in) :: nband(nkpt*nsppol) 92 real(dp),intent(inout) :: eigen(mband*(2*mband)**formeig*nkpt*nsppol) 93 real(dp),intent(inout) :: occ(mband*nkpt*nsppol) 94 95 !Local variables------------------------------- 96 !scalars 97 integer :: band_index,iband,ierr,ikpt,isppol,me,nbks,spaceComm 98 !character(len=500) :: message 99 !arrays 100 real(dp) :: tsec(2) 101 real(dp),allocatable :: buffer1(:),buffer2(:) 102 103 ! ************************************************************************* 104 105 if(xmpi_paral==1)then 106 107 ! Init mpi_comm 108 spaceComm=mpi_enreg%comm_cell 109 if(mpi_enreg%paral_kgb==1) spaceComm=mpi_enreg%comm_kpt 110 if(mpi_enreg%paral_hf==1) spaceComm=mpi_enreg%comm_kpt 111 ! Init me 112 me=mpi_enreg%me_kpt 113 114 if(localrdwf==0)then 115 call xmpi_bcast(eigen,0,spaceComm,ierr) 116 117 else if(localrdwf==1)then 118 119 ! Prepare transmission of eigen (and occ) 120 ABI_ALLOCATE(buffer1,(2*mband**(formeig+1)*nkpt*nsppol)) 121 ABI_ALLOCATE(buffer2,(2*mband**(formeig+1)*nkpt*nsppol)) 122 buffer1(:)=zero 123 buffer2(:)=zero 124 125 band_index=0 126 do isppol=1,nsppol 127 do ikpt=1,nkpt 128 nbks=nband(ikpt+(isppol-1)*nkpt) 129 130 if(mpi_enreg%paralbd==0)then 131 132 if(formeig==0)then 133 buffer1(2*band_index+1:2*band_index+nbks)=& 134 & eigen(band_index+1:band_index+nbks) 135 if(transmit_occ==2) then 136 buffer1(2*band_index+nbks+1:2*band_index+2*nbks)=& 137 & occ(band_index+1:band_index+nbks) 138 end if 139 band_index=band_index+nbks 140 else if(formeig==1)then 141 buffer1(band_index+1:band_index+2*nbks**2)=& 142 & eigen(band_index+1:band_index+2*nbks**2) 143 band_index=band_index+2*nbks**2 144 end if 145 146 else if(mpi_enreg%paralbd==1)then 147 148 ! Skip this k-point if not the proper processor 149 if(proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nbks,isppol,me)) then 150 if(formeig==0) then 151 band_index=band_index+nbks 152 else 153 band_index=band_index+2*nbks**2 154 end if 155 cycle 156 end if 157 ! Loop on bands 158 do iband=1,nbks 159 if(mpi_enreg%proc_distrb(ikpt, iband,isppol) /= me)cycle 160 if(formeig==0)then 161 buffer1(2*band_index+iband)=eigen(band_index+iband) 162 ! if(transmit_occ==2) then 163 ! buffer1(2*band_index+iband+nbdks)=occ(band_index+iband) 164 ! end if 165 else if (formeig==1)then 166 buffer1(band_index+(iband-1)*2*nbks+1: & 167 & band_index+(iband-1)*2*nbks+2*nbks)=& 168 & eigen(band_index+(iband-1)*2*nbks+1: & 169 & band_index+(iband-1)*2*nbks+2*nbks) 170 end if 171 end do 172 if(formeig==0)then 173 band_index=band_index+nbks 174 else 175 band_index=band_index+2*nbks**2 176 end if 177 end if 178 179 end do 180 end do 181 182 ! Build sum of everything 183 call timab(48,1,tsec) 184 ! call wrtout(std_out,' pareigocc : MPI_ALLREDUCE','COLL') 185 if(formeig==0)band_index=band_index*2 186 call xmpi_sum(buffer1,buffer2,band_index,spaceComm,ierr) 187 call timab(48,2,tsec) 188 189 band_index=0 190 do isppol=1,nsppol 191 do ikpt=1,nkpt 192 nbks=nband(ikpt+(isppol-1)*nkpt) 193 if(formeig==0)then 194 eigen(band_index+1:band_index+nbks)=& 195 & buffer2(2*band_index+1:2*band_index+nbks) 196 if(transmit_occ==2) then 197 occ(band_index+1:band_index+nbks)=& 198 & buffer2(2*band_index+nbks+1:2*band_index+2*nbks) 199 end if 200 band_index=band_index+nbks 201 else if(formeig==1)then 202 eigen(band_index+1:band_index+2*nbks**2)=& 203 & buffer1(band_index+1:band_index+2*nbks**2) 204 band_index=band_index+2*nbks**2 205 end if 206 end do 207 end do 208 209 ABI_DEALLOCATE(buffer1) 210 ABI_DEALLOCATE(buffer2) 211 212 end if 213 214 end if 215 216 end subroutine pareigocc