TABLE OF CONTENTS
ABINIT/initmpi_pert [ Functions ]
NAME
initmpi_pert
FUNCTION
Creates group for Parallelization over Perturbations.
COPYRIGHT
Copyright (C) 2005-2018 ABINIT group (FJ,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
dtset <type(dataset_type)>=all input variables in this dataset
OUTPUT
SIDE EFFECTS
mpi_enreg=information about MPI parallelization
PARENTS
mpi_setup
CHILDREN
get_npert_rbz,xmpi_comm_free
SOURCE
33 #if defined HAVE_CONFIG_H 34 #include "config.h" 35 #endif 36 37 #include "abi_common.h" 38 39 40 subroutine initmpi_pert(dtset,mpi_enreg) 41 42 use defs_basis 43 use defs_abitypes 44 use m_errors 45 use m_profiling_abi 46 use m_xmpi 47 48 use m_dtset, only : get_npert_rbz 49 50 !This section has been created automatically by the script Abilint (TD). 51 !Do not modify the following lines by hand. 52 #undef ABI_FUNC 53 #define ABI_FUNC 'initmpi_pert' 54 !End of the abilint section 55 56 implicit none 57 58 !Arguments ------------------------------------ 59 !scalars 60 type(MPI_type),intent(inout) :: mpi_enreg 61 type(dataset_type),intent(in) :: dtset 62 63 !Local variables------------------------------- 64 !scalars 65 integer:: iprocmin,irank,npert,nproc_per_cell,nrank,numproc 66 integer,allocatable :: ranks(:) 67 character(len=500) :: msg 68 !arrays 69 integer,pointer :: nkpt_rbz(:) 70 real(dp),pointer :: nband_rbz(:,:) 71 72 ! *********************************************************************** 73 74 if (mpi_enreg%me_pert<0) then 75 msg='Error in MPI distribution! Change your proc(s) distribution or use autoparal>0.' 76 MSG_ERROR(msg) 77 end if 78 79 call get_npert_rbz(dtset,nband_rbz,nkpt_rbz,npert) 80 81 if (dtset%nppert>=1) then 82 if (mpi_enreg%comm_cell/=mpi_enreg%comm_world) then 83 call xmpi_comm_free(mpi_enreg%comm_cell) 84 end if 85 mpi_enreg%comm_cell=mpi_enreg%comm_world 86 87 ! These values will be properly set in set_pert_comm 88 mpi_enreg%me_cell=mpi_enreg%me 89 mpi_enreg%nproc_cell=mpi_enreg%nproc 90 91 if (mpi_enreg%me>=0) then 92 nproc_per_cell=mpi_enreg%nproc/dtset%nppert 93 ABI_ALLOCATE(ranks,(dtset%nppert)) 94 iprocmin=mod(mpi_enreg%me,nproc_per_cell) 95 ranks=(/((iprocmin+(irank-1)*nproc_per_cell),irank=1,dtset%nppert)/) 96 mpi_enreg%comm_pert=xmpi_subcomm(mpi_enreg%comm_world,dtset%nppert,ranks) 97 ABI_DEALLOCATE(ranks) 98 mpi_enreg%me_pert=xmpi_comm_rank(mpi_enreg%comm_pert) 99 mpi_enreg%nproc_pert=dtset%nppert 100 if (iprocmin==0.and.mpi_enreg%me_pert==0.and.mpi_enreg%me/=0) then 101 MSG_BUG('Error on me_pert!') 102 end if 103 ! Define mpi_enreg%distrb_pert 104 ABI_ALLOCATE(mpi_enreg%distrb_pert,(npert)) 105 nrank=0 106 do irank=1,npert 107 nrank=nrank+1 108 mpi_enreg%distrb_pert(irank)=mod(nrank,dtset%nppert)-1 109 if (mpi_enreg%distrb_pert(irank)==-1) mpi_enreg%distrb_pert(irank)=dtset%nppert-1 110 end do 111 ! Make sure that subrank 0 is working on the last perturbation 112 ! Swap the ranks if necessary 113 numproc=mpi_enreg%distrb_pert(npert) 114 if(numproc/=0) then 115 do irank=1,npert 116 if (mpi_enreg%distrb_pert(irank)==numproc) mpi_enreg%distrb_pert(irank)=-2 117 if (mpi_enreg%distrb_pert(irank)==0) mpi_enreg%distrb_pert(irank)=-3 118 end do 119 do irank=1,npert 120 if (mpi_enreg%distrb_pert(irank)==-2) mpi_enreg%distrb_pert(irank)=0 121 if (mpi_enreg%distrb_pert(irank)==-3) mpi_enreg%distrb_pert(irank)=numproc 122 end do 123 end if 124 ! Communicator over one cell 125 ABI_ALLOCATE(ranks,(nproc_per_cell)) 126 iprocmin=(mpi_enreg%me/nproc_per_cell)*nproc_per_cell 127 ranks=(/((iprocmin+irank-1),irank=1,nproc_per_cell)/) 128 mpi_enreg%comm_cell_pert=xmpi_subcomm(mpi_enreg%comm_world,nproc_per_cell,ranks) 129 ABI_DEALLOCATE(ranks) 130 end if 131 132 else !nppert<=1 133 mpi_enreg%nproc_pert=1 134 mpi_enreg%comm_pert=xmpi_comm_self 135 mpi_enreg%me_pert=0 136 ABI_ALLOCATE(mpi_enreg%distrb_pert,(npert)) 137 mpi_enreg%distrb_pert(:)=0 138 end if 139 140 ABI_DEALLOCATE(nband_rbz) 141 ABI_DEALLOCATE(nkpt_rbz) 142 143 end subroutine initmpi_pert