TABLE OF CONTENTS
ABINIT/initmpi_atom [ Functions ]
NAME
initmpi_atom
FUNCTION
Initializes the mpi informations for parallelism over atoms (PAW).
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
dtset <type(dataset_type)>=all input variables for this dataset mpi_enreg= informations about MPI parallelization
OUTPUT
mpi_enreg= informations about MPI parallelization comm_atom =communicator over atoms nproc_atom =size of the communicator over atoms my_natom =number of atoms treated by current proc my_atmtab(mpi_enreg%natom)=indexes of the atoms treated by current processor
PARENTS
m_paral_pert,mpi_setup
CHILDREN
get_my_atmtab,get_my_natom
SOURCE
36 #if defined HAVE_CONFIG_H 37 #include "config.h" 38 #endif 39 40 #include "abi_common.h" 41 42 subroutine initmpi_atom(dtset,mpi_enreg) 43 44 use defs_basis 45 use defs_abitypes 46 use m_errors 47 use m_profiling_abi 48 use m_xmpi 49 50 use m_paral_atom, only : get_my_natom, get_my_atmtab 51 52 !This section has been created automatically by the script Abilint (TD). 53 !Do not modify the following lines by hand. 54 #undef ABI_FUNC 55 #define ABI_FUNC 'initmpi_atom' 56 !End of the abilint section 57 58 implicit none 59 60 !Arguments ------------------------------------ 61 !scalars 62 type(dataset_type),intent(in) :: dtset 63 type(MPI_type),intent(inout) :: mpi_enreg 64 !arrays 65 66 !Local variables------------------------------- 67 !scalars 68 logical :: my_atmtab_allocated,paral_atom 69 character(len=500) :: msg 70 integer :: iatom 71 72 ! *********************************************************************** 73 74 DBG_ENTER("COLL") 75 76 mpi_enreg%nproc_atom=1 77 mpi_enreg%comm_atom=xmpi_comm_self 78 mpi_enreg%my_natom=dtset%natom 79 if (associated(mpi_enreg%my_atmtab))then 80 ABI_DEALLOCATE(mpi_enreg%my_atmtab) 81 end if 82 nullify(mpi_enreg%my_atmtab) 83 84 if (xmpi_paral==0) then 85 mpi_enreg%nproc_atom=0 86 ABI_ALLOCATE(mpi_enreg%my_atmtab,(0)) 87 return 88 end if 89 90 91 !Check compatibility 92 if (dtset%paral_atom>0) then 93 msg='' 94 if (dtset%usepaw==0) msg= 'Parallelisation over atoms not compatible with usepaw=0 !' 95 if (dtset%usedmft==1) msg=' Parallelisation over atoms not compatible with usedmft=1 !' 96 if (dtset%usewvl==1) msg= 'Parallelisation over atoms not compatible with usewvl=1 !' 97 if (dtset%prtden>1.and.dtset%paral_kgb<=0) & 98 & msg= 'Parallelisation over atoms not compatible with prtden>1 (PAW AE densities) !' 99 if (dtset%optdriver/=RUNL_GSTATE.and.dtset%optdriver/=RUNL_RESPFN) & 100 & msg=' Parallelisation over atoms only compatible with GS or RF !' 101 if (dtset%macro_uj/=0)msg=' Parallelisation over atoms not compatible with macro_uj!=0 !' 102 if (msg/='') then 103 MSG_ERROR(msg) 104 end if 105 end if 106 107 if (mpi_enreg%comm_atom==xmpi_comm_null) then 108 mpi_enreg%nproc_atom=0;mpi_enreg%my_natom=0 109 ABI_ALLOCATE(mpi_enreg%my_atmtab,(0)) 110 return 111 end if 112 113 if (dtset%paral_atom>0) then 114 115 ! Build correct atom communicator 116 if (dtset%optdriver==RUNL_GSTATE.and.dtset%paral_kgb==1) then 117 mpi_enreg%comm_atom=mpi_enreg%comm_kptband 118 else 119 mpi_enreg%comm_atom=mpi_enreg%comm_cell 120 end if 121 122 ! Get number of processors sharing the atomic data distribution 123 mpi_enreg%nproc_atom=xmpi_comm_size(mpi_enreg%comm_atom) 124 125 ! Get local number of atoms 126 call get_my_natom(mpi_enreg%comm_atom,mpi_enreg%my_natom,dtset%natom) 127 paral_atom=(mpi_enreg%my_natom/=dtset%natom) 128 129 ! Build atom table 130 if (mpi_enreg%my_natom>0.and.paral_atom) then 131 my_atmtab_allocated=.false. 132 call get_my_atmtab(mpi_enreg%comm_atom,mpi_enreg%my_atmtab,my_atmtab_allocated, & 133 & paral_atom,dtset%natom) 134 else if (.not.paral_atom) then 135 ABI_ALLOCATE(mpi_enreg%my_atmtab,(dtset%natom)) 136 mpi_enreg%my_atmtab(1:dtset%natom)=(/(iatom, iatom=1,dtset%natom)/) 137 else if (mpi_enreg%my_natom==0) then 138 ABI_ALLOCATE(mpi_enreg%my_atmtab,(0)) 139 end if 140 141 end if 142 143 DBG_EXIT("COLL") 144 145 end subroutine initmpi_atom