TABLE OF CONTENTS


ABINIT/m_paral_atom [ Modules ]

[ Top ] [ Modules ]

NAME

  m_paral_atom

FUNCTION

  This module provides routines and methods used to manage the parallelisation/distribution
  of PAW data over atomic sites

COPYRIGHT

 Copyright (C) 2012-2018 ABINIT group (MT, MD)
 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.

NOTES

  FOR DEVELOPPERS: in order to preserve the portability of libPAW library,
  please consult ~abinit/src/??_libpaw/libpaw-coding-rules.txt

SOURCE

23 #include "libpaw.h"
24 
25 MODULE m_paral_atom
26 
27  USE_DEFS
28  USE_MSG_HANDLING
29  USE_MPI_WRAPPERS
30  USE_MEMORY_PROFILING
31 
32  implicit none
33 
34  private
35 
36 !public procedures.
37  public :: get_my_natom
38  public :: get_my_atmtab
39  public :: free_my_atmtab
40  public :: get_proc_atmtab
41  public :: get_atm_proc

m_paral_atom/free_my_atmtab [ Functions ]

[ Top ] [ m_paral_atom ] [ Functions ]

NAME

 free_my_atmtab

FUNCTION

 Cleanly deallocate a table of atom indexes (my_atmtab)

INPUTS

OUTPUT

SIDE EFFECTS

  my_atmtab_allocated=true if my_atmtab is allocated
  my_atmtab(:)=indexes of atoms treated by current process
               nothing is done if my_atmtab(:) is already associated to a target

PARENTS

      calc_efg,calc_fc,denfgr,dfpt_ewald,elt_ewald,eltxccore,initrhoij
      m_hamiltonian,m_paw_an,m_paw_ij,m_paw_pwaves_lmn,m_pawdij,m_pawfgrtab
      m_pawrhoij,make_efg_onsite,make_fc_paw,newfermie1,nhatgrid,outscfcv
      paw_mknewh0,pawaccrhoij,pawdenpot,pawdfptenergy,pawgrnl,pawmkaewf
      pawmknhat,pawmknhat_psipsi,pawnhatfr,pawprt,pawsushat,pawuj_red
      setnoccmmp,setrhoijpbe0

CHILDREN

SOURCE

248 subroutine free_my_atmtab(my_atmtab,my_atmtab_allocated)
249 
250 
251 !This section has been created automatically by the script Abilint (TD).
252 !Do not modify the following lines by hand.
253 #undef ABI_FUNC
254 #define ABI_FUNC 'free_my_atmtab'
255 !End of the abilint section
256 
257  implicit none
258 
259 !Arguments ---------------------------------------------
260 !scalars
261  logical,intent(inout) :: my_atmtab_allocated
262 !arrays
263  integer,pointer :: my_atmtab(:)
264 !Local variables ---------------------------------------
265 !scalars
266 !arrays
267 
268 ! *************************************************************************
269 
270  if (my_atmtab_allocated) then
271    LIBPAW_POINTER_DEALLOCATE(my_atmtab)
272    nullify(my_atmtab)
273    my_atmtab_allocated=.false.
274  end if
275 
276 end subroutine free_my_atmtab

m_paral_atom/get_atm_proc [ Functions ]

[ Top ] [ m_paral_atom ] [ Functions ]

NAME

 get_atm_proc

FUNCTION

  Given a list of atoms and a MPI communicator size, return a table
  containing the corresponding processor indexes.

COPYRIGHT

 Copyright (C) 2012-2018 ABINIT group (MD)
 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

  atom_list(:)= index of atoms
  nproc=size of communicator over atoms
  natom=total number of atoms

OUTPUT

 proc_list(:) = index of procs

NOTES

  The atoms are distributed contigously by egal part; the rest is distributed
  on all the procs (see get_my_atmtab).
  In case of modification of the distribution of atom over proc, this routine
  must be modified accordingly.

PARENTS

      m_paral_pert

CHILDREN

SOURCE

399  subroutine get_atm_proc(atom_list,natom,nproc,proc_list)
400 
401 !Arguments ---------------------------------------------
402 !scalars
403 
404 !This section has been created automatically by the script Abilint (TD).
405 !Do not modify the following lines by hand.
406 #undef ABI_FUNC
407 #define ABI_FUNC 'get_atm_proc'
408 !End of the abilint section
409 
410  integer, intent(in) :: natom,nproc
411 !arrays
412  integer, intent(in) :: atom_list(:)
413  integer, intent(out) :: proc_list(:)
414 
415 !Local variables ---------------------------------------
416 !scalars
417  integer :: nb_atom,dn,dn1,iatom,natomlim,iatm,jproclim,nmod
418 !arrays
419 
420 ! *************************************************************************
421 
422  nmod=mod(natom,nproc);nb_atom=size(atom_list)
423 
424  if (nmod==0) then
425    dn=natom/nproc
426    do iatm =1, nb_atom
427      iatom=atom_list(iatm)
428      proc_list(iatm)=(iatom-1)/dn
429    end do
430  else
431    dn=natom/nproc
432    dn1=natom/nproc + 1
433 !  Under (jproclim+1), 1 atome by proc is added
434 !  The rest nmod is distributed among jproclim+1 first procs
435    jproclim=nmod -1
436    natomlim=dn1*(jproclim+1)
437    do iatm=1,nb_atom
438      iatom=atom_list(iatm)
439      if (iatom<=natomlim) then
440        proc_list(iatm)=(iatom -1 )/dn1
441      else
442        proc_list(iatm)=jproclim + 1 + (iatom - 1 -(natomlim))/dn
443      end if
444    enddo
445  end if
446 
447 end subroutine get_atm_proc

m_paral_atom/get_my_atmtab [ Functions ]

[ Top ] [ m_paral_atom ] [ Functions ]

NAME

 get_my_atmtab

FUNCTION

 Given the total number of atoms and a MPI communicator return a table
 containing the indexes of atoms treated by current processor.

INPUTS

  comm_atom=communicator over atoms
  my_natom_ref= --optional-- a reference value for the local number of atoms
                (just for checking purposes)
  natom=total number of atoms

OUTPUT

SIDE EFFECTS

  my_atmtab(:)=indexes of atoms treated by current process
               nothing is done if my_atmtab(:) is already associated to a target
  my_atmtab_allocated=true if my_atmtab is allocated
  paral_atom=flag controlling parallelism over atoms

PARENTS

      calc_efg,calc_fc,denfgr,dfpt_accrho,dfpt_ewald,elt_ewald,eltxccore
      initmpi_atom,initrhoij,m_hamiltonian,m_paw_an,m_paw_ij,m_paw_pwaves_lmn
      m_pawdij,m_pawfgrtab,m_pawrhoij,make_efg_onsite,make_fc_paw,newfermie1
      nhatgrid,outscfcv,paw_mknewh0,pawaccrhoij,pawdenpot,pawdfptenergy
      pawgrnl,pawmkaewf,pawmknhat,pawmknhat_psipsi,pawnhatfr,pawprt,pawsushat
      pawuj_red,setnoccmmp,setrhoijpbe0

CHILDREN

SOURCE

143 subroutine get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,&
144 &                        my_natom_ref) ! optional argument
145 
146 
147 !This section has been created automatically by the script Abilint (TD).
148 !Do not modify the following lines by hand.
149 #undef ABI_FUNC
150 #define ABI_FUNC 'get_my_atmtab'
151 !End of the abilint section
152 
153  implicit none
154 
155 !Arguments ---------------------------------------------
156 !scalars
157  integer,intent(in) :: comm_atom,natom
158  integer,intent(in),optional :: my_natom_ref
159  logical,intent(inout) :: my_atmtab_allocated,paral_atom
160 !arrays
161  integer,pointer :: my_atmtab(:)
162 !Local variables ---------------------------------------
163 !scalars
164  integer :: iatom,me,my_natom,natom_bef,nmod,nproc
165  character(len=500) :: msg
166 !arrays
167 
168 ! *************************************************************************
169 
170  my_atmtab_allocated=.false.
171  if (.not.paral_atom) return
172 
173  if (comm_atom==xmpi_comm_self.or.comm_atom==xmpi_comm_null) paral_atom=.false.
174  if (paral_atom)  then
175    nproc=xmpi_comm_size(comm_atom)
176    paral_atom=(nproc>1)
177    if (paral_atom) then
178      if (.not.associated(my_atmtab)) then
179 !      Get local number of atoms
180        me=xmpi_comm_rank(comm_atom)
181        my_natom=natom/nproc
182        if (me<=(mod(natom,nproc)-1)) my_natom=natom/nproc + 1
183 !      Get table of atoms
184        if (my_natom>0) then
185          LIBPAW_POINTER_ALLOCATE(my_atmtab,(my_natom))
186          my_atmtab_allocated=.true.
187          if (my_natom==natom) then
188            my_atmtab(1:my_natom)=(/(iatom,iatom=1,natom)/)
189          else
190 !          The atoms are distributed contigously by egal part
191 !          (the rest is distributed on all the procs)
192            nmod=mod(natom,nproc)
193            if (me<=(nmod-1)) then
194              natom_bef=me*(natom/nproc)+me
195            else
196              natom_bef=me*(natom/nproc)+nmod
197            endif
198            do iatom=1,my_natom
199              my_atmtab(iatom)=iatom+natom_bef
200            enddo
201          end if
202        end if
203      else
204        my_natom=size(my_atmtab)
205      end if
206      if (present(my_natom_ref).and.(my_natom>0)) then
207        if (my_natom_ref/=size(my_atmtab)) then
208          msg='my_atmtab should have a size equal to my_natom !'
209          MSG_BUG(msg)
210        end if
211      end if
212    end if
213  endif
214 
215 end subroutine get_my_atmtab

m_paral_atom/get_my_natom [ Functions ]

[ Top ] [ m_paral_atom ] [ Functions ]

NAME

 get_my_natom

FUNCTION

 Given the total number of atoms, return the number of atoms treated by current process

INPUTS

  comm_atom=communicator over atoms
  natom=total number of atoms

OUTPUT

  my_natom=number of atoms treated by current process

PARENTS

      initmpi_atom,m_paw_an,m_paw_ij,m_pawfgrtab,m_pawrhoij

CHILDREN

SOURCE

 71 subroutine get_my_natom(comm_atom,my_natom,natom)
 72 
 73 
 74 !This section has been created automatically by the script Abilint (TD).
 75 !Do not modify the following lines by hand.
 76 #undef ABI_FUNC
 77 #define ABI_FUNC 'get_my_natom'
 78 !End of the abilint section
 79 
 80  implicit none
 81 
 82 !Arguments ---------------------------------------------
 83 !scalars
 84  integer,intent(in) :: comm_atom,natom
 85  integer,intent(out) :: my_natom
 86 !arrays
 87 
 88 !Local variables ---------------------------------------
 89 !scalars
 90  integer ::  me,nproc
 91 !arrays
 92 
 93 ! *************************************************************************
 94 
 95  my_natom=natom
 96  if (comm_atom/=xmpi_comm_self.and.comm_atom/=xmpi_comm_null)  then
 97    nproc=xmpi_comm_size(comm_atom)
 98    me=xmpi_comm_rank(comm_atom)
 99    my_natom=natom/nproc
100    if (me<=(mod(natom,nproc)-1)) my_natom=natom/nproc + 1
101  endif
102 
103 end subroutine get_my_natom

m_paral_atom/get_proc_atmtab [ Functions ]

[ Top ] [ m_paral_atom ] [ Functions ]

NAME

 get_proc_atmtab

FUNCTION

  Given the total number of atoms and the size of a communicator,
  return a table containing the indexes of atoms treated by a processor (with index iproc)

INPUTS

  comm_atom_size= size of communicator (over atoms)
  iproc= rank of the processor
  natom= total number of atoms

OUTPUT

  atmtab(natom_out)= indexes of atoms treated by process iproc
  natom_out= number of atoms treated by process iproc

NOTES

 In case of modification of the distribution of atom over proc,
   get_atmtab must be modified accordingly

PARENTS

CHILDREN

SOURCE

308  subroutine get_proc_atmtab(iproc,atmtab,natom_out,natom,comm_atom_size)
309 
310 
311 !This section has been created automatically by the script Abilint (TD).
312 !Do not modify the following lines by hand.
313 #undef ABI_FUNC
314 #define ABI_FUNC 'get_proc_atmtab'
315 !End of the abilint section
316 
317  implicit none
318 
319 !Arguments ---------------------------------------------
320 !scalars
321  integer, intent(in) :: comm_atom_size,iproc,natom
322 !arrays
323  integer,intent(out) :: natom_out
324  integer,allocatable, intent(out):: atmtab(:)
325 
326 !Local variables ---------------------------------------
327 !scalars
328  integer :: iatom,natom_bef,nmod,nproc
329 !arrays
330 
331 ! *************************************************************************
332 
333  nproc=comm_atom_size
334 
335  natom_out=natom/nproc ; if (iproc<=(mod(natom,nproc)-1)) natom_out=natom/nproc+1
336 
337 ! Get table of atoms
338  if (natom_out>0) then
339    LIBPAW_ALLOCATE(atmtab,(natom_out))
340 !  The atoms are distributed contigously by egal part
341 !  The rest is distributed on all the procs
342 !  (see get_my_atmtab)
343    nmod=mod(natom,nproc)
344    if (iproc<=(nmod-1)) then
345      natom_bef=iproc*(natom/nproc)+iproc
346    else
347      natom_bef=iproc*(natom/nproc)+nmod
348    end if
349    do iatom=1,natom_out
350      atmtab(iatom)=iatom+natom_bef
351    end do
352 
353  else
354    natom_out=0
355    LIBPAW_ALLOCATE(atmtab,(0))
356  end if
357 
358 end subroutine get_proc_atmtab