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-2022 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

22 #include "libpaw.h"
23 
24 MODULE m_paral_atom
25 
26  USE_DEFS
27  USE_MSG_HANDLING
28  USE_MPI_WRAPPERS
29  USE_MEMORY_PROFILING
30 
31  implicit none
32 
33  private
34 
35 !public procedures.
36  public :: get_my_natom
37  public :: get_my_atmtab
38  public :: free_my_atmtab
39  public :: get_proc_atmtab
40  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

SOURCE

204 subroutine free_my_atmtab(my_atmtab,my_atmtab_allocated)
205 
206 !Arguments ---------------------------------------------
207 !scalars
208  logical,intent(inout) :: my_atmtab_allocated
209 !arrays
210  integer,pointer :: my_atmtab(:)
211 !Local variables ---------------------------------------
212 !scalars
213 !arrays
214 
215 ! *************************************************************************
216 
217  if (my_atmtab_allocated) then
218    LIBPAW_POINTER_DEALLOCATE(my_atmtab)
219    nullify(my_atmtab)
220    my_atmtab_allocated=.false.
221  end if
222 
223 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-2022 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.

SOURCE

328  subroutine get_atm_proc(atom_list,natom,nproc,proc_list)
329 
330 !Arguments ---------------------------------------------
331 !scalars
332  integer, intent(in) :: natom,nproc
333 !arrays
334  integer, intent(in) :: atom_list(:)
335  integer, intent(out) :: proc_list(:)
336 
337 !Local variables ---------------------------------------
338 !scalars
339  integer :: nb_atom,dn,dn1,iatom,natomlim,iatm,jproclim,nmod
340 !arrays
341 
342 ! *************************************************************************
343 
344  nmod=mod(natom,nproc);nb_atom=size(atom_list)
345 
346  if (nmod==0) then
347    dn=natom/nproc
348    do iatm =1, nb_atom
349      iatom=atom_list(iatm)
350      proc_list(iatm)=(iatom-1)/dn
351    end do
352  else
353    dn=natom/nproc
354    dn1=natom/nproc + 1
355 !  Under (jproclim+1), 1 atome by proc is added
356 !  The rest nmod is distributed among jproclim+1 first procs
357    jproclim=nmod -1
358    natomlim=dn1*(jproclim+1)
359    do iatm=1,nb_atom
360      iatom=atom_list(iatm)
361      if (iatom<=natomlim) then
362        proc_list(iatm)=(iatom -1 )/dn1
363      else
364        proc_list(iatm)=jproclim + 1 + (iatom - 1 -(natomlim))/dn
365      end if
366    enddo
367  end if
368 
369 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

SOURCE

118 subroutine get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,&
119 &                        my_natom_ref) ! optional argument
120 
121 !Arguments ---------------------------------------------
122 !scalars
123  integer,intent(in) :: comm_atom,natom
124  integer,intent(in),optional :: my_natom_ref
125  logical,intent(inout) :: my_atmtab_allocated,paral_atom
126 !arrays
127  integer,pointer :: my_atmtab(:)
128 !Local variables ---------------------------------------
129 !scalars
130  integer :: iatom,me,my_natom,natom_bef,nmod,nproc
131  character(len=500) :: msg
132 !arrays
133 
134 ! *************************************************************************
135 
136  my_atmtab_allocated=.false.
137  if (.not.paral_atom) return
138 
139  if (comm_atom==xmpi_comm_self.or.comm_atom==xmpi_comm_null) paral_atom=.false.
140  if (paral_atom)  then
141    nproc=xmpi_comm_size(comm_atom)
142    paral_atom=(nproc>1)
143    if (paral_atom) then
144      if (.not.associated(my_atmtab)) then
145 !      Get local number of atoms
146        me=xmpi_comm_rank(comm_atom)
147        my_natom=natom/nproc
148        if (me<=(mod(natom,nproc)-1)) my_natom=natom/nproc + 1
149 !      Get table of atoms
150        if (my_natom>0) then
151          LIBPAW_POINTER_ALLOCATE(my_atmtab,(my_natom))
152          my_atmtab_allocated=.true.
153          if (my_natom==natom) then
154            my_atmtab(1:my_natom)=(/(iatom,iatom=1,natom)/)
155          else
156 !          The atoms are distributed contigously by egal part
157 !          (the rest is distributed on all the procs)
158            nmod=mod(natom,nproc)
159            if (me<=(nmod-1)) then
160              natom_bef=me*(natom/nproc)+me
161            else
162              natom_bef=me*(natom/nproc)+nmod
163            endif
164            do iatom=1,my_natom
165              my_atmtab(iatom)=iatom+natom_bef
166            enddo
167          end if
168        end if
169      else
170        my_natom=size(my_atmtab)
171      end if
172      if (present(my_natom_ref).and.(my_natom>0)) then
173        if (my_natom_ref/=size(my_atmtab)) then
174          msg='my_atmtab should have a size equal to my_natom !'
175          LIBPAW_BUG(msg)
176        end if
177      end if
178    end if
179  endif
180 
181 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

SOURCE

65 subroutine get_my_natom(comm_atom,my_natom,natom)
66 
67 !Arguments ---------------------------------------------
68 !scalars
69  integer,intent(in) :: comm_atom,natom
70  integer,intent(out) :: my_natom
71 !arrays
72 
73 !Local variables ---------------------------------------
74 !scalars
75  integer ::  me,nproc
76 !arrays
77 
78 ! *************************************************************************
79 
80  my_natom=natom
81  if (comm_atom/=xmpi_comm_self.and.comm_atom/=xmpi_comm_null)  then
82    nproc=xmpi_comm_size(comm_atom)
83    me=xmpi_comm_rank(comm_atom)
84    my_natom=natom/nproc
85    if (me<=(mod(natom,nproc)-1)) my_natom=natom/nproc + 1
86  endif
87 
88 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

SOURCE

251  subroutine get_proc_atmtab(iproc,atmtab,natom_out,natom,comm_atom_size)
252 
253 !Arguments ---------------------------------------------
254 !scalars
255  integer, intent(in) :: comm_atom_size,iproc,natom
256 !arrays
257  integer,intent(out) :: natom_out
258  integer,allocatable, intent(out):: atmtab(:)
259 
260 !Local variables ---------------------------------------
261 !scalars
262  integer :: iatom,natom_bef,nmod,nproc
263 !arrays
264 
265 ! *************************************************************************
266 
267  nproc=comm_atom_size
268 
269  natom_out=natom/nproc ; if (iproc<=(mod(natom,nproc)-1)) natom_out=natom/nproc+1
270 
271 ! Get table of atoms
272  if (natom_out>0) then
273    LIBPAW_ALLOCATE(atmtab,(natom_out))
274 !  The atoms are distributed contigously by egal part
275 !  The rest is distributed on all the procs
276 !  (see get_my_atmtab)
277    nmod=mod(natom,nproc)
278    if (iproc<=(nmod-1)) then
279      natom_bef=iproc*(natom/nproc)+iproc
280    else
281      natom_bef=iproc*(natom/nproc)+nmod
282    end if
283    do iatom=1,natom_out
284      atmtab(iatom)=iatom+natom_bef
285    end do
286 
287  else
288    natom_out=0
289    LIBPAW_ALLOCATE(atmtab,(0))
290  end if
291 
292 end subroutine get_proc_atmtab