1!{\src2tex{textfont=tt}}
2!!****f* ABINIT/xmpi_allgather
3!! NAME
4!!  xmpi_allgather
5!!
6!! FUNCTION
7!!  This module contains functions that calls MPI routine,
8!!  if we compile the code using the  MPI CPP flags.
9!!  xmpi_allgather is the generic function.
10!!
11!! COPYRIGHT
12!!  Copyright (C) 2001-2016 ABINIT group (AR,XG)
13!!  This file is distributed under the terms of the
14!!  GNU General Public License, see ~ABINIT/COPYING
15!!  or http://www.gnu.org/copyleft/gpl.txt .
16!!
17!! SOURCE
18
19!!***
20
21!!****f* ABINIT/xmpi_allgather_int
22!! NAME
23!!  xmpi_allgather_int
24!!
25!! FUNCTION
26!!  Gathers data from all tasks and distributes it to all.
27!!  Target: one-dimensional integer arrays.
28!!
29!! INPUTS
30!!  spaceComm= MPI communicator
31!!
32!! OUTPUT
33!!  ier= exit status, a non-zero value meaning there is an error
34!!
35!! SIDE EFFECTS
36!!  xval= buffer array
37!!  recvbuf= received elements
38!!
39!! PARENTS
40!!
41!! CHILDREN
42!!      mpi_allgather
43!!
44!! SOURCE
45
46subroutine xmpi_allgather_int(xval,recvbuf,spaceComm,ier)
47
48
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 'xmpi_allgather_int'
54!End of the abilint section
55
56 implicit none
57
58!Arguments-------------------------
59 integer,intent(inout) :: xval
60 integer, DEV_CONTARRD intent(inout) :: recvbuf(:)
61 integer, intent(in) :: spaceComm
62 integer,intent(out) :: ier
63
64!Local variables-------------------
65
66! *************************************************************************
67 ier=0
68#if defined HAVE_MPI
69 if (spaceComm /= MPI_COMM_SELF .and. spaceComm /= MPI_COMM_NULL) then
70!  allgather xval on all proc. in spaceComm
71   call MPI_ALLGATHER(xval,1,MPI_INTEGER,recvbuf,1,MPI_INTEGER,spaceComm,ier)
72 else if (spaceComm == MPI_COMM_SELF) then
73   recvbuf(1)=xval
74 end if
75#else
76 recvbuf(1)=xval
77#endif
78end subroutine xmpi_allgather_int
79!!***
80
81
82!!****f* ABINIT/xmpi_allgather_char
83!! NAME
84!!  xmpi_allgather_char
85!!
86!! FUNCTION
87!!  Gathers data from all tasks and distributes it to all.
88!!  Target: one-dimensional character(20) arrays.
89!!
90!! INPUTS
91!!  spaceComm= MPI communicator
92!!
93!! OUTPUT
94!!  ier= exit status, a non-zero value meaning there is an error
95!!
96!! SIDE EFFECTS
97!!  charval= buffer array
98!!  recvbuf= received elements
99!!
100!! PARENTS
101!!
102!! CHILDREN
103!!      mpi_allgather
104!!
105!! SOURCE
106
107subroutine xmpi_allgather_char(charval,recvbuf,spaceComm,ier)
108
109
110!This section has been created automatically by the script Abilint (TD).
111!Do not modify the following lines by hand.
112#undef ABI_FUNC
113#define ABI_FUNC 'xmpi_allgather_char'
114!End of the abilint section
115
116 implicit none
117
118!Arguments-------------------------
119 integer,intent(in)  :: spaceComm
120 integer,intent(out) :: ier
121 character(len=20),intent(inout) :: charval
122 character(len=20), DEV_CONTARRD intent(inout) :: recvbuf(:)
123
124! *************************************************************************
125 ier=0
126#if defined HAVE_MPI
127 if (spaceComm /= MPI_COMM_SELF .and. spaceComm /= MPI_COMM_NULL) then
128!  allgather xval on all proc. in spaceComm
129   call MPI_ALLGATHER(charval,20,MPI_CHARACTER,recvbuf,20,MPI_CHARACTER,spaceComm,ier)
130 else if (spaceComm == MPI_COMM_SELF) then
131   recvbuf=charval
132 end if
133#else
134 recvbuf=charval
135#endif
136
137end subroutine xmpi_allgather_char
138!!***
139
140
141!!****f* ABINIT/xmpi_allgather_int1d
142!! NAME
143!!  xmpi_allgather_int1d
144!!
145!! FUNCTION
146!!  Gathers data from all tasks and distributes it to all.
147!!  Target: one-dimensional integer arrays.
148!!
149!! INPUTS
150!!  xval= buffer array
151!!  nelem= number of elements
152!!  spaceComm= MPI communicator
153!!
154!! OUTPUT
155!!  ier= exit status, a non-zero value meaning there is an error
156!!
157!! SIDE EFFECTS
158!!  recvbuf= received elements
159!!
160!! PARENTS
161!!
162!! CHILDREN
163!!      mpi_allgather
164!!
165!! SOURCE
166
167subroutine xmpi_allgather_int1d(xval,nelem,recvbuf,spaceComm,ier)
168
169
170
171!This section has been created automatically by the script Abilint (TD).
172!Do not modify the following lines by hand.
173#undef ABI_FUNC
174#define ABI_FUNC 'xmpi_allgather_int1d'
175!End of the abilint section
176
177 implicit none
178
179!Arguments-------------------------
180 integer, DEV_CONTARRD intent(in) :: xval(:)
181 integer, DEV_CONTARRD intent(inout) :: recvbuf(:)
182 integer ,intent(in) :: nelem,spaceComm
183 integer ,intent(out) :: ier
184
185! *************************************************************************
186 ier=0
187#if defined HAVE_MPI
188 if (spaceComm /= MPI_COMM_SELF .and. spaceComm /= MPI_COMM_NULL) then
189!  allgather xval on all proc. in spaceComm
190   call MPI_ALLGATHER(xval,nelem,MPI_INTEGER,recvbuf,nelem,MPI_INTEGER,spaceComm,ier)
191 else if (spaceComm == MPI_COMM_SELF) then
192   recvbuf(1:nelem)=xval(1:nelem)
193 end if
194#else
195 recvbuf(1:nelem)=xval(1:nelem)
196#endif
197end subroutine xmpi_allgather_int1d
198!!***
199
200!!****f* ABINIT/xmpi_allgather_dp1d
201!! NAME
202!!  xmpi_allgather_dp1d
203!!
204!! FUNCTION
205!!  Gathers data from all tasks and distributes it to all.
206!!  Target: double precision one-dimensional arrays.
207!!
208!! INPUTS
209!!  xval= buffer array
210!!  nelem= number of elements
211!!  spaceComm= MPI communicator
212!!
213!! OUTPUT
214!!  ier= exit status, a non-zero value meaning there is an error
215!!
216!! SIDE EFFECTS
217!!  recvbuf= received elements
218!!
219!! PARENTS
220!!
221!! CHILDREN
222!!      mpi_allgather
223!!
224!! SOURCE
225
226subroutine xmpi_allgather_dp1d(xval,nelem,recvbuf,spaceComm,ier)
227
228
229
230!This section has been created automatically by the script Abilint (TD).
231!Do not modify the following lines by hand.
232#undef ABI_FUNC
233#define ABI_FUNC 'xmpi_allgather_dp1d'
234!End of the abilint section
235
236 implicit none
237
238!Arguments-------------------------
239 real(dp), DEV_CONTARRD intent(in) :: xval(:)
240 real(dp), DEV_CONTARRD intent(inout) :: recvbuf(:)
241 integer ,intent(in) :: nelem,spaceComm
242 integer ,intent(out) :: ier
243
244! *************************************************************************
245 ier=0
246#if defined HAVE_MPI
247 if (spaceComm /= MPI_COMM_SELF .and. spaceComm /= MPI_COMM_NULL) then
248!  allgather xval on all proc. in spaceComm
249   call MPI_ALLGATHER(xval,nelem,MPI_DOUBLE_PRECISION,recvbuf,nelem,MPI_DOUBLE_PRECISION,spaceComm,ier)
250 else if (spaceComm == MPI_COMM_SELF) then
251   recvbuf(1:nelem)=xval(1:nelem)
252 end if
253#else
254 recvbuf(1:nelem)=xval(1:nelem)
255#endif
256end subroutine xmpi_allgather_dp1d
257!!***
258
259!!****f* ABINIT/xmpi_allgather_dp2d
260!! NAME
261!!  xmpi_allgather_dp2d
262!!
263!! FUNCTION
264!!  Gathers data from all tasks and distributes it to all.
265!!  Target: double precision two-dimensional arrays.
266!!
267!! INPUTS
268!!  xval= buffer array
269!!  nelem= number of elements
270!!  spaceComm= MPI communicator
271!!
272!! OUTPUT
273!!  ier= exit status, a non-zero value meaning there is an error
274!!
275!! SIDE EFFECTS
276!!  recvbuf= received elements
277!!
278!! PARENTS
279!!
280!! CHILDREN
281!!      mpi_allgather
282!!
283!! SOURCE
284
285subroutine xmpi_allgather_dp2d(xval,nelem,recvbuf,spaceComm,ier)
286
287
288
289!This section has been created automatically by the script Abilint (TD).
290!Do not modify the following lines by hand.
291#undef ABI_FUNC
292#define ABI_FUNC 'xmpi_allgather_dp2d'
293!End of the abilint section
294
295 implicit none
296
297!Arguments-------------------------
298 real(dp), DEV_CONTARRD intent(in) :: xval(:,:)
299 real(dp), DEV_CONTARRD intent(inout) :: recvbuf(:,:)
300 integer ,intent(in) :: nelem,spaceComm
301 integer ,intent(out)   :: ier
302
303! *************************************************************************
304 ier=0
305#if defined HAVE_MPI
306 if (spaceComm /= MPI_COMM_SELF .and. spaceComm /= MPI_COMM_NULL) then
307!  allgather xval on all proc. in spaceComm
308   call MPI_ALLGATHER(xval,nelem,MPI_DOUBLE_PRECISION,recvbuf,nelem,MPI_DOUBLE_PRECISION,spaceComm,ier)
309 else if (spaceComm == MPI_COMM_SELF) then
310   recvbuf(:,:)=xval(:,:)
311 end if
312#else
313 recvbuf(:,:)=xval(:,:)
314#endif
315end subroutine xmpi_allgather_dp2d
316!!***
317
318!!****f* ABINIT/xmpi_allgather_dp3d
319!! NAME
320!!  xmpi_allgather_dp3d
321!!
322!! FUNCTION
323!!  Gathers data from all tasks and distributes it to all.
324!!  Target: double precision three-dimensional arrays.
325!!
326!! INPUTS
327!!  xval= buffer array
328!!  nelem= number of elements
329!!  spaceComm= MPI communicator
330!!
331!! OUTPUT
332!!  ier= exit status, a non-zero value meaning there is an error
333!!
334!! SIDE EFFECTS
335!!  recvbuf= received elements
336!!
337!! PARENTS
338!!
339!! CHILDREN
340!!      mpi_allgather
341!!
342!! SOURCE
343
344subroutine xmpi_allgather_dp3d(xval,nelem,recvbuf,spaceComm,ier)
345
346
347
348!This section has been created automatically by the script Abilint (TD).
349!Do not modify the following lines by hand.
350#undef ABI_FUNC
351#define ABI_FUNC 'xmpi_allgather_dp3d'
352!End of the abilint section
353
354 implicit none
355
356
357!Arguments-------------------------
358 real(dp), DEV_CONTARRD intent(in) :: xval(:,:,:)
359 real(dp), DEV_CONTARRD intent(inout) :: recvbuf(:,:,:)
360 integer ,intent(in) :: nelem,spaceComm
361 integer ,intent(out)   :: ier
362
363! *************************************************************************
364 ier=0
365#if defined HAVE_MPI
366 if (spaceComm /= MPI_COMM_SELF .and. spaceComm /= MPI_COMM_NULL) then
367!  allgather xval on all proc. in spaceComm
368   call MPI_ALLGATHER(xval,nelem,MPI_DOUBLE_PRECISION,recvbuf,nelem,MPI_DOUBLE_PRECISION,spaceComm,ier)
369 else if (spaceComm == MPI_COMM_SELF) then
370   recvbuf(:,:,:)=xval(:,:,:)
371 end if
372#else
373 recvbuf(:,:,:)=xval(:,:,:)
374#endif
375end subroutine xmpi_allgather_dp3d
376!!***
377
378!!****f* ABINIT/xmpi_allgather_dp4d
379!! NAME
380!!  xmpi_allgather_dp4d
381!!
382!! FUNCTION
383!!  Gathers data from all tasks and distributes it to all.
384!!  Target: double precision four-dimensional arrays.
385!!
386!! INPUTS
387!!  xval= buffer array
388!!  nelem= number of elements
389!!  spaceComm= MPI communicator
390!!
391!! OUTPUT
392!!  ier= exit status, a non-zero value meaning there is an error
393!!
394!! SIDE EFFECTS
395!!  recvbuf= received elements
396!!
397!! PARENTS
398!!
399!! CHILDREN
400!!      mpi_allgather
401!!
402!! SOURCE
403
404subroutine xmpi_allgather_dp4d(xval,nelem,recvbuf,spaceComm,ier)
405
406
407
408!This section has been created automatically by the script Abilint (TD).
409!Do not modify the following lines by hand.
410#undef ABI_FUNC
411#define ABI_FUNC 'xmpi_allgather_dp4d'
412!End of the abilint section
413
414 implicit none
415
416!Arguments-------------------------
417 real(dp), DEV_CONTARRD intent(in) :: xval(:,:,:,:)
418 real(dp), DEV_CONTARRD intent(inout) :: recvbuf(:,:,:,:)
419 integer ,intent(in) :: nelem,spaceComm
420 integer ,intent(out)   :: ier
421
422! *************************************************************************
423 ier=0
424#if defined HAVE_MPI
425 if (spaceComm /= MPI_COMM_SELF .and. spaceComm /= MPI_COMM_NULL) then
426!  allgather xval on all proc. in spaceComm
427   call MPI_ALLGATHER(xval,nelem,MPI_DOUBLE_PRECISION,recvbuf,nelem,MPI_DOUBLE_PRECISION,spaceComm,ier)
428 else if (spaceComm == MPI_COMM_SELF) then
429   recvbuf(:,:,:,:)=xval(:,:,:,:)
430 end if
431#else
432 recvbuf(:,:,:,:)=xval(:,:,:,:)
433#endif
434end subroutine xmpi_allgather_dp4d
435!!***
436