TABLE OF CONTENTS


ABINIT/m_xredistribute [ Modules ]

[ Top ] [ Modules ]

NAME

  m_xredistribute

FUNCTION

  This module contains a function to re-dristribute
  results on different procs.

COPYRIGHT

 Copyright (C) 2008-2018 ABINIT group (MMANCINI)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

OUTPUT

PARENTS

CHILDREN

SOURCE

26 #if defined HAVE_CONFIG_H
27 #include "config.h"
28 #endif
29 
30 #include "abi_common.h"
31 
32 
33 module m_xredistribute
34 
35  use defs_basis
36  use m_errors
37  use m_abicore
38 #if defined HAVE_MPI2
39  use mpi
40 #endif
41 
42  implicit none
43 
44 #if defined HAVE_MPI1
45  include 'mpif.h'
46 #endif
47 
48  private
49 
50  public :: xredistribute         ! Redistribute the work
51 
52 
53  interface xredistribute
54   module procedure xredistribute_mpi_dp
55   module procedure xredistribute_mpi_2d_dp
56  end interface xredistribute
57 
58 
59 CONTAINS  !===========================================================

ABINIT/xredistribute_mpi_2d_dp [ Functions ]

[ Top ] [ Functions ]

NAME

  xredistribute_mpi_2d_dp

FUNCTION

INPUTS

  xval= buffer array
  send_counts= number of sent elements (initial distribution)
  send_displs= postions of values sent by the processor (initial positions)
  rec_counts= number of received elements (final distribution)
  rec_displs= positions of values received by the processors (final position)
  nproc=number of processor
  me=proc me
  spaceComm= MPI communicator

OUTPUT

  ier= exit status, a non-zero value meaning there is an error

SIDE EFFECTS

  recvbuf= received buffer

PARENTS

CHILDREN

      mpi_allgatherv,mpi_scatterv

SOURCE

174 subroutine xredistribute_mpi_2d_dp(xval,send_counts,send_displs,recvbuf,&
175   &                              rec_counts,rec_displs,me,nproc,spaceComm,ier)
176 
177 
178 !This section has been created automatically by the script Abilint (TD).
179 !Do not modify the following lines by hand.
180 #undef ABI_FUNC
181 #define ABI_FUNC 'xredistribute_mpi_2d_dp'
182 !End of the abilint section
183 
184  implicit none
185 
186 
187 !Arguments-------------------------
188  integer ,intent(in) :: me,nproc
189  real(dp),intent(in) :: xval(:,:)
190  real(dp),intent(inout) :: recvbuf(:,:)
191  integer ,intent(in) :: send_displs(0:nproc-1),send_counts(0:nproc-1)
192  integer ,intent(in) :: rec_displs(0:nproc-1),rec_counts(0:nproc-1)
193  integer ,intent(in) :: spaceComm
194  integer ,intent(out):: ier
195 
196 !Local variables-------------------
197  integer :: size
198  real(dp),allocatable :: totbuff(:)
199  character(500) :: msg
200 ! *********************************************************************
201 
202  ier=0
203 #if defined HAVE_MPI
204  if (spaceComm /= MPI_COMM_SELF .and. spaceComm /= MPI_COMM_NULL) then
205  size = sum(send_counts)
206  if(size /=sum(rec_counts))then
207    msg = 'the total sizes of sent and receved msg are not equal'
208    MSG_ERROR(msg)
209  endif
210 
211  ABI_ALLOCATE(totbuff,(size))
212    !--join all the vector in to a single one
213    call MPI_ALLGATHERV(xval,send_counts(me),MPI_DOUBLE_PRECISION,totbuff,&
214 &                      send_counts,send_displs,MPI_DOUBLE_PRECISION,spaceComm,ier)
215 
216 
217    !--now distribute the total vector on the procs
218    call MPI_SCATTERV(totbuff,rec_counts,rec_displs,MPI_DOUBLE_PRECISION,&
219 &                    recvbuf,rec_counts(me),MPI_DOUBLE_PRECISION,&
220 &                    0,spaceComm,ier)
221    ABI_DEALLOCATE(totbuff)
222  end if
223 #endif
224 end subroutine xredistribute_mpi_2d_dp

ABINIT/xredistribute_mpi_dp [ Functions ]

[ Top ] [ Functions ]

NAME

  xredistribute_mpi_dp

FUNCTION

INPUTS

  xval= buffer array
  send_counts= number of sent elements (initial distribution)
  send_displs= postions of values sent by the processor (initial positions)
  rec_counts= number of received elements (final distribution)
  rec_displs= positions of values received by the processors (final position)
  nproc=number of processor
  me=proc me
  spaceComm= MPI communicator

OUTPUT

  ier= exit status, a non-zero value meaning there is an error

SIDE EFFECTS

  recvbuf= received buffer

PARENTS

CHILDREN

      mpi_allgatherv,mpi_scatterv

SOURCE

 92 subroutine xredistribute_mpi_dp(xval,send_counts,send_displs,recvbuf,&
 93   &                              rec_counts,rec_displs,me,nproc,spaceComm,ier)
 94 
 95 
 96 !This section has been created automatically by the script Abilint (TD).
 97 !Do not modify the following lines by hand.
 98 #undef ABI_FUNC
 99 #define ABI_FUNC 'xredistribute_mpi_dp'
100 !End of the abilint section
101 
102  implicit none
103 
104 !Arguments-------------------------
105  integer ,intent(in) :: me,nproc
106  real(dp),intent(in) :: xval(:)
107  real(dp),intent(inout) :: recvbuf(:)
108  integer ,intent(in) :: send_displs(0:nproc-1),send_counts(0:nproc-1)
109  integer ,intent(in) :: rec_displs(0:nproc-1),rec_counts(0:nproc-1)
110  integer ,intent(in) :: spaceComm
111  integer ,intent(out):: ier
112 
113 !Local variables-------------------
114  integer :: size
115  real(dp),allocatable :: totbuff(:)
116  character(500) :: msg
117 ! *********************************************************************
118 
119  ier=0
120 #if defined HAVE_MPI
121  if (spaceComm /= MPI_COMM_SELF .and. spaceComm /= MPI_COMM_NULL) then
122  size = sum(send_counts)
123  if(size /=sum(rec_counts))then
124    msg = 'the total sizes of sent and receved msg are not equal'
125    MSG_ERROR(msg)
126  endif
127 
128  ABI_ALLOCATE(totbuff,(size))
129    !--join all the vector in to a single one
130    call MPI_ALLGATHERV(xval,send_counts(me),MPI_DOUBLE_PRECISION,totbuff,&
131 &                      send_counts,send_displs,MPI_DOUBLE_PRECISION,spaceComm,ier)
132 
133 
134    !--now distribute the total vector on the procs
135    call MPI_SCATTERV(totbuff,rec_counts,rec_displs,MPI_DOUBLE_PRECISION,&
136 &                    recvbuf,rec_counts(me),MPI_DOUBLE_PRECISION,&
137 &                    0,spaceComm,ier)
138 
139    ABI_DEALLOCATE(totbuff)
140  end if
141 #endif
142 end subroutine xredistribute_mpi_dp