TABLE OF CONTENTS


ABINIT/initmpi_band [ Functions ]

[ Top ] [ Functions ]

NAME

  initmpi_band

FUNCTION

  Initializes the mpi informations for band parallelism (paralbd=1).

COPYRIGHT

  Copyright (C) 2008-2018 ABINIT group (MT)
  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

  mpi_enreg= informations about MPI parallelization
  nband(nkpt*nsppol)= number of bands per k point, for each spin
  nkpt= number of k-points
  nsppol= 1 for unpolarized, 2 for polarized

OUTPUT

  mpi_enreg=information about MPI parallelization
  mpi_enreg%comm_band=communicator of BAND set

PARENTS

      dfpt_looppert

CHILDREN

SOURCE

 34 #if defined HAVE_CONFIG_H
 35 #include "config.h"
 36 #endif
 37 
 38 #include "abi_common.h"
 39 
 40 subroutine initmpi_band(mpi_enreg,nband,nkpt,nsppol)
 41 
 42  use defs_basis
 43  use defs_abitypes
 44  use m_profiling_abi
 45  use m_errors
 46  use m_xmpi
 47 
 48 !This section has been created automatically by the script Abilint (TD).
 49 !Do not modify the following lines by hand.
 50 #undef ABI_FUNC
 51 #define ABI_FUNC 'initmpi_band'
 52 !End of the abilint section
 53 
 54  implicit none
 55 
 56 !Arguments ------------------------------------
 57 !scalars
 58  integer,intent(in) :: nkpt,nsppol
 59  integer,intent(in) :: nband(nkpt*nsppol)
 60  type(MPI_type),intent(inout) :: mpi_enreg
 61 
 62 !Local variables-------------------------------
 63 !scalars
 64  integer :: ii,ikpt,iproc_min,iproc_max,irank,isppol
 65  integer :: me,nband_k,nproc,nbsteps,nrank,nstates,spacecomm
 66  character(len=500) :: msg
 67 !arrays
 68  integer,allocatable :: ranks(:)
 69 
 70 ! ***********************************************************************
 71 
 72  mpi_enreg%comm_band=xmpi_comm_self
 73 
 74  if (mpi_enreg%paralbd==1.and.xmpi_paral==1) then
 75 
 76 !  Comm_kpt is supposed to treat spins, k-points and bands
 77    spacecomm=mpi_enreg%comm_kpt
 78    nproc=mpi_enreg%nproc_kpt
 79    me=mpi_enreg%me_kpt
 80 
 81    nstates=sum(nband(1:nkpt*nsppol))
 82    nbsteps=nstates/nproc
 83    if (mod(nstates,nproc)/=0) nbsteps=nbsteps+1
 84 
 85    if (nbsteps<maxval(nband(1:nkpt*nsppol))) then
 86 
 87      nrank=0
 88      do isppol=1,nsppol
 89        do ikpt=1,nkpt
 90          ii=ikpt+(isppol-1)*nkpt
 91          nband_k=nband(ii)
 92          if (nbsteps<nband_k) then
 93            iproc_min=minval(mpi_enreg%proc_distrb(ikpt,:,isppol))
 94            iproc_max=maxval(mpi_enreg%proc_distrb(ikpt,:,isppol))
 95            if ((me>=iproc_min).and.(me<=iproc_max)) then
 96              nrank=iproc_max-iproc_min+1
 97              if (.not.allocated(ranks)) then
 98                ABI_ALLOCATE(ranks,(nrank))
 99                if (nrank>0) ranks=(/((iproc_min+irank-1),irank=1,nrank)/)
100              else if (nrank/=size(ranks)) then
101                msg='Number of bands per proc should be the same for all k-points!'
102                MSG_BUG(msg)
103              end if
104            end if
105          end if
106        end do
107      end do
108      if (.not.allocated(ranks)) then
109        ABI_ALLOCATE(ranks,(0))
110      end if
111 
112      mpi_enreg%comm_band=xmpi_subcomm(spacecomm,nrank,ranks)
113 
114      ABI_DEALLOCATE(ranks)
115    end if
116  end if
117 
118 end subroutine initmpi_band