Code_Saturne
CFD tool
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
cs_parall.h
Go to the documentation of this file.
1 #ifndef __CS_PARALL_H__
2 #define __CS_PARALL_H__
3 
4 /*============================================================================
5  * Functions dealing with parallelism
6  *============================================================================*/
7 
8 /*
9  This file is part of Code_Saturne, a general-purpose CFD tool.
10 
11  Copyright (C) 1998-2012 EDF S.A.
12 
13  This program is free software; you can redistribute it and/or modify it under
14  the terms of the GNU General Public License as published by the Free Software
15  Foundation; either version 2 of the License, or (at your option) any later
16  version.
17 
18  This program is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
20  FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
21  details.
22 
23  You should have received a copy of the GNU General Public License along with
24  this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
25  Street, Fifth Floor, Boston, MA 02110-1301, USA.
26 */
27 
28 /*----------------------------------------------------------------------------*/
29 
30 /*----------------------------------------------------------------------------
31  * Local headers
32  *----------------------------------------------------------------------------*/
33 
34 #include "cs_defs.h"
35 
36 /*----------------------------------------------------------------------------*/
37 
39 
40 /*============================================================================
41  * Public function prototypes for Fortran API
42  *============================================================================*/
43 
44 /*----------------------------------------------------------------------------
45  * Compute the maximum value of a counter (int) for the entire domain in
46  * case of parallelism.
47  *
48  * Fortran Interface
49  *
50  * subroutine parcmx (counter)
51  * *****************
52  *
53  * integer counter <-> : input = local counter
54  * output = global max counter
55  *----------------------------------------------------------------------------*/
56 
57 void
58 CS_PROCF (parcmx, PARCMX)(cs_int_t *counter);
59 
60 /*----------------------------------------------------------------------------
61  * Compute the minimum value of a counter (int) for the entire domain in
62  * case of parallelism.
63  *
64  * Fortran Interface
65  *
66  * subroutine parcmn (counter)
67  * *****************
68  *
69  * integer counter <-> : input = local counter
70  * output = global min counter
71  *----------------------------------------------------------------------------*/
72 
73 void
74 CS_PROCF (parcmn, PARCMN)(cs_int_t *counter);
75 
76 /*----------------------------------------------------------------------------
77  * Compute the global sum of a counter (int) for the entire domain in case
78  * of parallelism.
79  *
80  * Fortran Interface :
81  *
82  * subroutine parcpt (counter)
83  * *****************
84  *
85  * integer counter : <-> : input = counter to sum
86  * output = global sum
87  *----------------------------------------------------------------------------*/
88 
89 void
90 CS_PROCF (parcpt, PARCPT)(cs_int_t *counter);
91 
92 /*----------------------------------------------------------------------------
93  * Compute the global sum of a real for the entire domain in case of parellism
94  *
95  * Fortran Interface :
96  *
97  * subroutine parsom (var)
98  * *****************
99  *
100  * double precision var : <-> : input = value to sum
101  * output = global sum
102  *----------------------------------------------------------------------------*/
103 
104 void
105 CS_PROCF (parsom, PARSOM)(cs_real_t *var);
106 
107 /*----------------------------------------------------------------------------
108  * Compute the maximum value of a real variable for the entire domain in case
109  * of parallelism.
110  *
111  * Fortran Interface :
112  *
113  * subroutine parmax (var)
114  * *****************
115  *
116  * double precision var : <-> : input = local maximum value
117  * output = global maximum value
118  *----------------------------------------------------------------------------*/
119 
120 void
121 CS_PROCF (parmax, PARMAX)(cs_real_t *var);
122 
123 /*----------------------------------------------------------------------------
124  * Compute the minimum value of a real variable for the entire domain in case
125  * of parallelism.
126  *
127  * Fortran Interface :
128  *
129  * subroutine parmin (var)
130  * *****************
131  *
132  * double precision var : <-> : input = local minimum value
133  * output = global minimum value
134  *----------------------------------------------------------------------------*/
135 
136 void
137 CS_PROCF (parmin, PARMIN)(cs_real_t *var);
138 
139 /*----------------------------------------------------------------------------
140  * Maximum value of a real and its related rank for the entire domain in
141  * case of parallelism.
142  *
143  * Fortran Interface
144  *
145  * subroutine parmxl (nbr, var, xyzvar)
146  * *****************
147  *
148  * integer nbr : <-- : size of the related variable
149  * double precision var : <-> : input: local max. value
150  * output: global max. value
151  * double precision xyzvar(nbr) : <-> : input: value related to local max.
152  * output: value related to global max.
153  *----------------------------------------------------------------------------*/
154 
155 void
156 CS_PROCF (parmxl, PARMXL)(cs_int_t *nbr,
157  cs_real_t *var,
158  cs_real_t xyzvar[]);
159 
160 /*----------------------------------------------------------------------------
161  * Minimum value of a real and its related rank for the entire domain in
162  * case of parallelism.
163  *
164  * Fortran Interface
165  *
166  * Interface Fortran :
167  *
168  * subroutine parmnl (nbr, var, xyzvar)
169  * *****************
170  *
171  * integer nbr : <-- : size of the related variable
172  * double precision var : <-> : input: local max. value
173  * output: global max. value
174  * double precision xyzvar(nbr) : <-> : input: value related to local max.
175  * output: value related to global max.
176  *----------------------------------------------------------------------------*/
177 
178 void
179 CS_PROCF (parmnl, PARMNL)(cs_int_t *nbr,
180  cs_real_t *var,
181  cs_real_t xyzvar[]);
182 
183 /*----------------------------------------------------------------------------
184  * Compute the global sum for each element of an array of int in case of
185  * parallelism.
186  *
187  * Fortran Interface
188  *
189  * subroutine parism (n_elts, array)
190  * *****************
191  *
192  * integer n_elts : <-- : size of the array.
193  * integer array(*) : <-> : input = local array
194  * output = array of global sum values.
195  *----------------------------------------------------------------------------*/
196 
197 void
198 CS_PROCF (parism, PARISM)(cs_int_t *n_elts,
199  cs_int_t array[]);
200 
201 /*----------------------------------------------------------------------------
202  * Compute the global maximum value for each element of an array of int in
203  * case of parallelism.
204  *
205  * Fortran Interface :
206  *
207  * subroutine parimx (n_elts, array)
208  * *****************
209  *
210  * integer n_elts : <-- : size of the array.
211  * integer array(*) : <-> : input = local array
212  * output = array of global max. values.
213  *----------------------------------------------------------------------------*/
214 
215 void
216 CS_PROCF (parimx, PARIMX)(cs_int_t *n_elts,
217  cs_int_t array[]);
218 
219 /*----------------------------------------------------------------------------
220  * Compute the global minimum value for each element of an array of int in
221  * case of parallelism.
222  *
223  * Fortran Interface :
224  *
225  * subroutine parimn (n_elts, array)
226  * *****************
227  *
228  * integer n_elts : <-- : size of the array.
229  * integer array(*) : <-> : input = local array
230  * output = array of global min. values.
231  *----------------------------------------------------------------------------*/
232 
233 void
234 CS_PROCF (parimn, PARIMN)(cs_int_t *n_elts,
235  cs_int_t array[]);
236 
237 /*----------------------------------------------------------------------------
238  * Compute the global sum for each element of an array of real in case of
239  * parallelism.
240  *
241  * Fortran Interface
242  *
243  * subroutine parrsm (n_elts, array)
244  * *****************
245  *
246  * integer n_elts : <-- : size of the array.
247  * double precision array(*) : <-> : input = local array
248  * output = array of global sum values.
249  *----------------------------------------------------------------------------*/
250 
251 void
252 CS_PROCF (parrsm, PARRSM)(cs_int_t *n_elts,
253  cs_real_t array[]);
254 
255 /*----------------------------------------------------------------------------
256  * Compute the global maximum value for each element of an array of real in
257  * case of parallelism.
258  *
259  * Fortran Interface :
260  *
261  * subroutine parrmx (n_elts, array)
262  * *****************
263  *
264  * integer n_elts : <-- : size of the array
265  * double precision array(*) : <-> : input = local array
266  * output = array of global max. values.
267  *----------------------------------------------------------------------------*/
268 
269 void
270 CS_PROCF (parrmx, PARRMX)(cs_int_t *n_elts,
271  cs_real_t array[]);
272 
273 /*----------------------------------------------------------------------------
274  * Compute the global minimum value for each element of an array of real in
275  * case of parallelism.
276  *
277  * Fortran Interface :
278  *
279  * subroutine parrmn (n_elts, array)
280  * *****************
281  *
282  * integer n_elts : <-- : size of the array
283  * double precision array(*) : <-> : input = local array
284  * output = array of global min. values.
285  *----------------------------------------------------------------------------*/
286 
287 void
288 CS_PROCF (parrmn, PARRMN)(cs_int_t *n_elts,
289  cs_real_t array[]);
290 
291 /*----------------------------------------------------------------------------
292  * Broadcast to all the ranks the value of each element of an array of int.
293  * (encapsulation of MPI_Bcast())
294  *
295  * Fortran Interface :
296  *
297  * subroutine parbci (irank, n_elts, array)
298  * *****************
299  *
300  * integer irank : <-- : rank related to the sending process
301  * integer n_elts : <-- : size of the array
302  * integer array(*) : <-> : array to broadcast
303  *----------------------------------------------------------------------------*/
304 
305 void
306 CS_PROCF (parbci, PARBCI)(cs_int_t *irank,
307  cs_int_t *n_elts,
308  cs_int_t array[]);
309 
310 /*----------------------------------------------------------------------------
311  * Broadcast to all the ranks the value of each element of an array of real.
312  * (encapsulation of MPI_Bcast())
313  *
314  * Fortran Interface :
315  *
316  * subroutine parbcr (irank, n_elts, array)
317  * *****************
318  *
319  * integer irank : <-- : rank related to the sending process
320  * integer n_elts : <-- : size of the array
321  * double precision array(*) : <-> : array to broadcast
322  *----------------------------------------------------------------------------*/
323 
324 void
325 CS_PROCF (parbcr, PARBCR)(cs_int_t *irank,
326  cs_int_t *n_elts,
327  cs_real_t array[]);
328 
329 /*----------------------------------------------------------------------------
330  * Build a global array from each local array in each domain. The size of
331  * each local array can be different.
332  *
333  * Fortran Interface :
334  *
335  * subroutine paragv (nvar, nvargb, var, vargb)
336  * *****************
337  *
338  * integer n_elts : <-- : size of the local array
339  * integer n_g_elts : <-- : size of the global array
340  * double precision array(*) : <-- : local array
341  * double precision g_array(*) : --> : global array
342  *----------------------------------------------------------------------------*/
343 
344 void
345 CS_PROCF (paragv, PARAGV)(cs_int_t *n_elts,
346  cs_int_t *n_g_elts,
347  cs_real_t array[],
348  cs_real_t *g_array);
349 
350 /*----------------------------------------------------------------------------
351  * Find a node which minimizes a given distance and its related rank.
352  * May be used to locate a node among several domains.
353  *
354  * Fortran Interface :
355  *
356  * subroutine parfpt (node, ndrang, dis2mn)
357  * *****************
358  *
359  * integer node : <-> : local number of the closest node
360  * integer ndrang : --> : rank id for which the distance is the
361  * smallest
362  * double precision dis2mn : <-- : square distance between the closest node
363  * and the wanted node.
364  *----------------------------------------------------------------------------*/
365 
366 void
367 CS_PROCF (parfpt, PARFPT)(cs_int_t *node,
368  cs_int_t *ndrang,
369  cs_real_t *dis2mn);
370 
371 /*----------------------------------------------------------------------------
372  * Return the value associated to a probe.
373  *
374  * Fortran Interface :
375  *
376  * subroutine parhis (node, ndrang, var, varcap)
377  * *****************
378  *
379  * integer node : <-- : local number of the element related to
380  * a measure node
381  * integer ndrang : <-- : rank of the process owning the closest
382  * node from the measure node
383  * double precision var(*) : <-- : values of the variable on local elements
384  * double precision varcap : --> : value of the variable for the element
385  * related to the measure node
386  *----------------------------------------------------------------------------*/
387 
388 void
389 CS_PROCF (parhis, PARHIS)(cs_int_t *node,
390  cs_int_t *ndrang,
391  cs_real_t var[],
392  cs_real_t *varcap);
393 
394 /*----------------------------------------------------------------------------
395  * Call a barrier in case of parallelism
396  *
397  * This function should not be necessary in production code,
398  * but it may be useful for debugging purposes.
399  *
400  * Fortran interface :
401  *
402  * subroutine parbar
403  * *****************
404  *----------------------------------------------------------------------------*/
405 
406 void
407 CS_PROCF (parbar, PARBAR)(void);
408 
409 /*=============================================================================
410  * Public function prototypes
411  *============================================================================*/
412 
413 /*----------------------------------------------------------------------------
414  * Sum counters on all default communicator processes.
415  *
416  * parameters:
417  * cpt <-> local counter value input, global counter value output (array)
418  * n <-- number of counter array values
419  *----------------------------------------------------------------------------*/
420 
421 #if defined(HAVE_MPI_IN_PLACE)
422 
423 inline static void
425  const int n)
426 {
427  if (cs_glob_n_ranks > 1) {
428  MPI_Allreduce(MPI_IN_PLACE, cpt, n, CS_MPI_GNUM, MPI_SUM,
429  cs_glob_mpi_comm);
430  }
431 }
432 
433 #elif defined(HAVE_MPI)
434 
435 void
437  const int n);
438 
439 #else
440 
441 #define cs_parall_counter(_cpt, _n)
442 
443 #endif
444 
445 /*----------------------------------------------------------------------------
446  * Maximum values of a counter on all default communicator processes.
447  *
448  * parameters:
449  * cpt <-> local counter value input, global counter value output (array)
450  * n <-- number of counter array values
451  *----------------------------------------------------------------------------*/
452 
453 #if defined(HAVE_MPI_IN_PLACE)
454 
455 inline static void
457  const int n)
458 {
459  if (cs_glob_n_ranks > 1) {
460  MPI_Allreduce(MPI_IN_PLACE, cpt, n, CS_MPI_LNUM, MPI_MAX,
461  cs_glob_mpi_comm);
462  }
463 }
464 
465 #elif defined(HAVE_MPI)
466 
467 void
469  const int n);
470 
471 #else
472 
473 #define cs_parall_counter_max(_cpt, _n)
474 
475 #endif
476 
477 /*----------------------------------------------------------------------------
478  * Sum values of a given datatype on all default communicator processes.
479  *
480  * parameters:
481  * n <-- number of values
482  * datatype <-- matching Code_Saturne datatype
483  * val <-> local value input, global value output (array)
484  *----------------------------------------------------------------------------*/
485 
486 #if defined(HAVE_MPI_IN_PLACE)
487 
488 inline static void
489 cs_parall_sum(int n,
490  cs_datatype_t datatype,
491  void *val)
492 {
493  if (cs_glob_n_ranks > 1) {
494  MPI_Allreduce(MPI_IN_PLACE, val, n, cs_datatype_to_mpi[datatype], MPI_SUM,
495  cs_glob_mpi_comm);
496  }
497 }
498 
499 #elif defined(HAVE_MPI)
500 
501 void
502 cs_parall_sum(int n,
503  cs_datatype_t datatype,
504  void *val);
505 
506 #else
507 
508 #define cs_parall_sum(_n, _datatype, _val);
509 
510 #endif
511 
512 /*----------------------------------------------------------------------------
513  * Maximum values of a given datatype on all default communicator processes.
514  *
515  * parameters:
516  * n <-- number of values
517  * datatype <-- matching Code_Saturne datatype
518  * val <-> local value input, global value output (array)
519  *----------------------------------------------------------------------------*/
520 
521 #if defined(HAVE_MPI_IN_PLACE)
522 
523 inline static void
524 cs_parall_max(int n,
525  cs_datatype_t datatype,
526  void *val)
527 {
528  if (cs_glob_n_ranks > 1) {
529  MPI_Allreduce(MPI_IN_PLACE, val, n, cs_datatype_to_mpi[datatype], MPI_MAX,
530  cs_glob_mpi_comm);
531  }
532 }
533 
534 #elif defined(HAVE_MPI)
535 
536 void
537 cs_parall_max(int n,
538  cs_datatype_t datatype,
539  void *val);
540 
541 #else
542 
543 #define cs_parall_max(_n, _datatype, _val);
544 
545 #endif
546 
547 /*----------------------------------------------------------------------------
548  * Minimum values of a given datatype on all default communicator processes.
549  *
550  * parameters:
551  * n <-- number of values
552  * datatype <-- matching Code_Saturne datatype
553  * val <-> local value input, global value output (array)
554  *----------------------------------------------------------------------------*/
555 
556 #if defined(HAVE_MPI_IN_PLACE)
557 
558 inline static void
559 cs_parall_min(int n,
560  cs_datatype_t datatype,
561  void *val)
562 {
563  if (cs_glob_n_ranks > 1) {
564  MPI_Allreduce(MPI_IN_PLACE, val, n, cs_datatype_to_mpi[datatype], MPI_MIN,
565  cs_glob_mpi_comm);
566  }
567 }
568 
569 #elif defined(HAVE_MPI)
570 
571 void
572 cs_parall_min(int n,
573  cs_datatype_t datatype,
574  void *val);
575 
576 #else
577 
578 #define cs_parall_min(_n, _datatype, _val);
579 
580 #endif
581 
582 /*----------------------------------------------------------------------------
583  * Return minimum recommended scatter or gather buffer size.
584  *
585  * This is used by some internal part to block or scatter/gather algorithms,
586  * so as to allow I/O buffer size tuning.
587  *
588  * returns:
589  * minimum recommended part to block or gather buffer size (in bytes)
590  *----------------------------------------------------------------------------*/
591 
592 size_t
594 
595 /*----------------------------------------------------------------------------
596  * Define minimum recommended scatter or gather buffer size.
597  *
598  * This is used by some internal part to block or scatter/gather algorithms,
599  * so as to allow I/O buffer size tuning.
600  *
601  * parameters:
602  * minimum recommended part to block or gather buffer size (in bytes)
603  *----------------------------------------------------------------------------*/
604 
605 void
606 cs_parall_set_min_coll_buf_size(size_t buffer_size);
607 
608 /*----------------------------------------------------------------------------*/
609 
611 
612 #endif /* __CS_PARALL_H__ */
void parbcr(cs_int_t *irank, cs_int_t *n_elts, cs_real_t array[])
Definition: cs_parall.c:742
cs_datatype_t
Definition: cs_defs.h:223
void parmin(cs_real_t *var)
Definition: cs_parall.c:287
void parmnl(cs_int_t *nbr, cs_real_t *var, cs_real_t xyzvar[])
Definition: cs_parall.c:360
#define BEGIN_C_DECLS
Definition: cs_defs.h:365
#define cs_parall_counter(_cpt, _n)
Definition: cs_parall.h:441
BEGIN_C_DECLS void parcmx(cs_int_t *counter)
Definition: cs_parall.c:138
#define cs_parall_max(_n, _datatype, _val)
Definition: cs_parall.h:543
void parrmx(cs_int_t *n_elts, cs_real_t array[])
Definition: cs_parall.c:610
void paragv(cs_int_t *n_elts, cs_int_t *n_g_elts, cs_real_t array[], cs_real_t *g_array)
Definition: cs_parall.c:769
#define cs_parall_counter_max(_cpt, _n)
Definition: cs_parall.h:473
int cs_glob_n_ranks
Definition: cs_defs.c:82
int cs_int_t
Definition: cs_defs.h:263
void parbar(void)
Definition: cs_parall.c:903
void parsom(cs_real_t *var)
Definition: cs_parall.c:227
void parrmn(cs_int_t *n_elts, cs_real_t array[])
Definition: cs_parall.c:663
void parimx(cs_int_t *n_elts, cs_int_t array[])
Definition: cs_parall.c:451
void parcpt(cs_int_t *counter)
Definition: cs_parall.c:198
void parmxl(cs_int_t *nbr, cs_real_t *var, cs_real_t xyzvar[])
Definition: cs_parall.c:320
void parfpt(cs_int_t *node, cs_int_t *ndrang, cs_real_t *dis2mn)
Definition: cs_parall.c:822
int cs_lnum_t
Definition: cs_defs.h:260
size_t cs_parall_get_min_coll_buf_size(void)
Definition: cs_parall.c:1024
unsigned cs_gnum_t
Definition: cs_defs.h:255
void cs_parall_set_min_coll_buf_size(size_t buffer_size)
Definition: cs_parall.c:1044
void parism(cs_int_t *n_elts, cs_int_t array[])
Definition: cs_parall.c:398
void parhis(cs_int_t *node, cs_int_t *ndrang, cs_real_t var[], cs_real_t *varcap)
Definition: cs_parall.c:871
void parbci(cs_int_t *irank, cs_int_t *n_elts, cs_int_t array[])
Definition: cs_parall.c:716
#define END_C_DECLS
Definition: cs_defs.h:366
double cs_real_t
Definition: cs_defs.h:264
#define CS_PROCF(x, y)
Definition: cs_defs.h:379
void parcmn(cs_int_t *counter)
Definition: cs_parall.c:168
void parmax(cs_real_t *var)
Definition: cs_parall.c:257
#define cs_parall_min(_n, _datatype, _val)
Definition: cs_parall.h:578
void parrsm(cs_int_t *n_elts, cs_real_t array[])
Definition: cs_parall.c:557
void parimn(cs_int_t *n_elts, cs_int_t array[])
Definition: cs_parall.c:504
#define cs_parall_sum(_n, _datatype, _val)
Definition: cs_parall.h:508