8.0
general documentation
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-2023 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  * General types and macros used throughout code_saturne
42  *============================================================================*/
43 
44 /*----------------------------------------------------------------------------
45  * Variable value type.
46  *----------------------------------------------------------------------------*/
47 
51 typedef enum {
52 
68 
69 /*============================================================================
70  * Global variables
71  *============================================================================*/
72 
73 /* Preferred indexed sum option, adapted to shared-memory parallelism */
74 
76 
77 /*=============================================================================
78  * Public function prototypes
79  *============================================================================*/
80 
81 /*----------------------------------------------------------------------------*/
88 /*----------------------------------------------------------------------------*/
89 
90 #if defined(HAVE_MPI_IN_PLACE)
91 
92 inline static void
94  const int n)
95 {
96  if (cs_glob_n_ranks > 1) {
97  MPI_Allreduce(MPI_IN_PLACE, cpt, n, CS_MPI_GNUM, MPI_SUM,
99  }
100 }
101 
102 #elif defined(HAVE_MPI)
103 
104 void
106  const int n);
107 
108 #else
109 
110 #define cs_parall_counter(_cpt, _n)
111 
112 #endif
113 
114 /*----------------------------------------------------------------------------*/
121 /*----------------------------------------------------------------------------*/
122 
123 #if defined(HAVE_MPI_IN_PLACE)
124 
125 inline static void
127  const int n)
128 {
129  if (cs_glob_n_ranks > 1) {
130  MPI_Allreduce(MPI_IN_PLACE, cpt, n, CS_MPI_LNUM, MPI_MAX,
132  }
133 }
134 
135 #elif defined(HAVE_MPI)
136 
137 void
139  const int n);
140 
141 #else
142 
143 #define cs_parall_counter_max(_cpt, _n)
144 
145 #endif
146 
147 /*----------------------------------------------------------------------------*/
155 /*----------------------------------------------------------------------------*/
156 
157 #if defined(HAVE_MPI_IN_PLACE)
158 
159 inline static void
161  cs_datatype_t datatype,
162  void *val)
163 {
164  if (cs_glob_n_ranks > 1) {
165  MPI_Allreduce(MPI_IN_PLACE, val, n, cs_datatype_to_mpi[datatype], MPI_SUM,
167  }
168 }
169 
170 #elif defined(HAVE_MPI)
171 
172 void
173 cs_parall_sum(int n,
174  cs_datatype_t datatype,
175  void *val);
176 
177 #else
178 
179 #define cs_parall_sum(_n, _datatype, _val) { };
180 
181 #endif
182 
183 /*----------------------------------------------------------------------------*/
192 /*----------------------------------------------------------------------------*/
193 
194 #if defined(HAVE_MPI_IN_PLACE)
195 
196 inline static void
198  cs_datatype_t datatype,
199  void *val)
200 {
201  if (cs_glob_n_ranks > 1) {
202  MPI_Allreduce(MPI_IN_PLACE, val, n, cs_datatype_to_mpi[datatype], MPI_MAX,
204  }
205 }
206 
207 #elif defined(HAVE_MPI)
208 
209 void
210 cs_parall_max(int n,
211  cs_datatype_t datatype,
212  void *val);
213 
214 #else
215 
216 #define cs_parall_max(_n, _datatype, _val);
217 
218 #endif
219 
220 /*----------------------------------------------------------------------------*/
229 /*----------------------------------------------------------------------------*/
230 
231 #if defined(HAVE_MPI_IN_PLACE)
232 
233 inline static void
235  cs_datatype_t datatype,
236  void *val)
237 {
238  if (cs_glob_n_ranks > 1) {
239  MPI_Allreduce(MPI_IN_PLACE, val, n, cs_datatype_to_mpi[datatype], MPI_MIN,
241  }
242 }
243 
244 #elif defined(HAVE_MPI)
245 
246 void
247 cs_parall_min(int n,
248  cs_datatype_t datatype,
249  void *val);
250 
251 #else
252 
253 #define cs_parall_min(_n, _datatype, _val);
254 
255 #endif
256 
257 /*----------------------------------------------------------------------------*/
268 /*----------------------------------------------------------------------------*/
269 
270 #if defined(HAVE_MPI)
271 
272 inline static void
273 cs_parall_bcast(int root_rank,
274  int n,
275  cs_datatype_t datatype,
276  void *val)
277 {
278  if (cs_glob_n_ranks > 1)
279  MPI_Bcast(val, n, cs_datatype_to_mpi[datatype], root_rank,
281 }
282 
283 #else
284 
285 #define cs_parall_bcast(_root_rank, _n, _datatype, _val);
286 
287 #endif
288 
289 /*----------------------------------------------------------------------------*/
305 /*----------------------------------------------------------------------------*/
306 
307 void
308 cs_parall_allgather_r(int n_elts,
309  int n_g_elts,
310  cs_real_t array[],
311  cs_real_t g_array[]);
312 
313 /*----------------------------------------------------------------------------*/
333 /*----------------------------------------------------------------------------*/
334 
335 void
337  int n_g_elts,
338  int stride,
339  cs_real_t o_key[],
340  cs_real_t array[],
341  cs_real_t g_array[]);
342 
343 /*----------------------------------------------------------------------------*/
361 /*----------------------------------------------------------------------------*/
362 
363 void
364 cs_parall_gather_r(int root_rank,
365  int n_elts,
366  int n_g_elts,
367  const cs_real_t array[],
368  cs_real_t g_array[]);
369 
370 /*----------------------------------------------------------------------------*/
392 /*----------------------------------------------------------------------------*/
393 
394 void
395 cs_parall_gather_ordered_r(int root_rank,
396  int n_elts,
397  int n_g_elts,
398  int stride,
399  cs_real_t o_key[],
400  cs_real_t array[],
401  cs_real_t g_array[]);
402 
403 /*----------------------------------------------------------------------------*/
421 /*----------------------------------------------------------------------------*/
422 
423 void
424 cs_parall_scatter_r(int root_rank,
425  int n_elts,
426  int n_g_elts,
427  const cs_real_t g_array[],
428  cs_real_t array[]);
429 
430 /*----------------------------------------------------------------------------*/
449 /*----------------------------------------------------------------------------*/
450 
451 void
452 cs_parall_gather_f(int root_rank,
453  int n_elts,
454  int n_g_elts,
455  const float array[],
456  float g_array[]);
457 
458 /*----------------------------------------------------------------------------*/
477 /*----------------------------------------------------------------------------*/
478 
479 void
480 cs_parall_scatter_f(int root_rank,
481  int n_elts,
482  int n_g_elts,
483  const float g_array[],
484  float array[]);
485 
486 /*----------------------------------------------------------------------------*/
496 /*----------------------------------------------------------------------------*/
497 
498 void
500  cs_real_t *max,
501  cs_real_t max_loc_vals[]);
502 
503 /*----------------------------------------------------------------------------*/
513 /*----------------------------------------------------------------------------*/
514 
515 void
517  cs_real_t *min,
518  cs_real_t min_loc_vals[]);
519 
520 /*----------------------------------------------------------------------------*/
531 /*----------------------------------------------------------------------------*/
532 
533 void
535  int *rank_id,
536  cs_real_t val);
537 
538 /*----------------------------------------------------------------------------*/
547 /*----------------------------------------------------------------------------*/
548 
549 size_t
551 
552 /*----------------------------------------------------------------------------*/
562 /*----------------------------------------------------------------------------*/
563 
564 void
565 cs_parall_set_min_coll_buf_size(size_t buffer_size);
566 
567 /*----------------------------------------------------------------------------*/
580 /*----------------------------------------------------------------------------*/
581 
582 inline static void
584  size_t type_size,
585  cs_lnum_t *s_id,
586  cs_lnum_t *e_id)
587 {
588 #if defined(HAVE_OPENMP)
589  const int t_id = omp_get_thread_num();
590  const int n_t = omp_get_num_threads();
591  const cs_lnum_t t_n = (n + n_t - 1) / n_t;
592  const cs_lnum_t cl_m = CS_CL_SIZE / type_size; /* Cache line multiple */
593 
594  *s_id = t_id * t_n;
595  *e_id = (t_id+1) * t_n;
596  *s_id = cs_align(*s_id, cl_m);
597  *e_id = cs_align(*e_id, cl_m);
598  if (*e_id > n) *e_id = n;
599 #else
600  CS_UNUSED(type_size); /* avoid compiler warning */
601  *s_id = 0;
602  *e_id = n;
603 #endif
604 }
605 
606 /*----------------------------------------------------------------------------*/
607 
609 
610 #endif /* __CS_PARALL_H__ */
int cs_glob_n_ranks
Definition: cs_defs.c:175
MPI_Datatype cs_datatype_to_mpi[]
Definition: cs_defs.c:157
MPI_Comm cs_glob_mpi_comm
Definition: cs_defs.c:183
cs_datatype_t
Definition: cs_defs.h:272
#define BEGIN_C_DECLS
Definition: cs_defs.h:509
double cs_real_t
Floating-point value.
Definition: cs_defs.h:319
#define CS_MPI_LNUM
Definition: cs_defs.h:405
#define CS_MPI_GNUM
Definition: cs_defs.h:385
unsigned long cs_gnum_t
global mesh entity number
Definition: cs_defs.h:298
static cs_lnum_t cs_align(cs_lnum_t i, cs_lnum_t m)
Given a base index i, return the next index aligned with a size m.
Definition: cs_defs.h:591
#define CS_UNUSED(x)
Definition: cs_defs.h:495
#define END_C_DECLS
Definition: cs_defs.h:510
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:313
#define CS_CL_SIZE
Definition: cs_defs.h:465
void cs_parall_gather_r(int root_rank, int n_elts, int n_g_elts, const cs_real_t array[], cs_real_t g_array[])
Build a global array on the given root rank from all local arrays.
Definition: cs_parall.c:1016
static void cs_parall_bcast(int root_rank, int n, cs_datatype_t datatype, void *val)
Broadcast values of a given datatype to all default communicator processes.
Definition: cs_parall.h:273
void cs_parall_set_min_coll_buf_size(size_t buffer_size)
Define minimum recommended scatter or gather buffer size.
Definition: cs_parall.c:1339
void cs_parall_gather_ordered_r(int root_rank, int n_elts, int n_g_elts, int stride, cs_real_t o_key[], cs_real_t array[], cs_real_t g_array[])
Build an ordered global array on the given root rank from all local arrays.
Definition: cs_parall.c:1083
void cs_parall_min_id_rank_r(cs_lnum_t *elt_id, int *rank_id, cs_real_t val)
Given an (id, rank, value) tuple, return the local id and rank corresponding to the global minimum va...
Definition: cs_parall.c:841
static void cs_parall_max(int n, cs_datatype_t datatype, void *val)
Maximum values of a given datatype on all default communicator processes.
Definition: cs_parall.h:197
static void cs_parall_counter_max(cs_lnum_t cpt[], const int n)
Maximum values of a counter on all default communicator processes.
Definition: cs_parall.h:126
void cs_parall_allgather_r(int n_elts, int n_g_elts, cs_real_t array[], cs_real_t g_array[])
Build a global array from each local array in each domain.
Definition: cs_parall.c:895
static void cs_parall_counter(cs_gnum_t cpt[], const int n)
Sum values of a counter on all default communicator processes.
Definition: cs_parall.h:93
void cs_parall_scatter_r(int root_rank, int n_elts, int n_g_elts, const cs_real_t g_array[], cs_real_t array[])
Distribute a global array from a given root rank over all ranks. Each rank receive the part related t...
Definition: cs_parall.c:1132
static void cs_parall_sum(int n, cs_datatype_t datatype, void *val)
Sum values of a given datatype on all default communicator processes.
Definition: cs_parall.h:160
void cs_parall_allgather_ordered_r(int n_elts, int n_g_elts, int stride, cs_real_t o_key[], cs_real_t array[], cs_real_t g_array[])
Build an ordered global array from each local array in each domain.
Definition: cs_parall.c:972
void cs_parall_scatter_f(int root_rank, int n_elts, int n_g_elts, const float g_array[], float array[])
Distribute a global array from a given root rank over all ranks. Each rank receive the part related t...
Definition: cs_parall.c:1262
static void cs_parall_thread_range(cs_lnum_t n, size_t type_size, cs_lnum_t *s_id, cs_lnum_t *e_id)
Compute array index bounds for a local thread. When called inside an OpenMP parallel section,...
Definition: cs_parall.h:583
void cs_parall_min_loc_vals(int n, cs_real_t *min, cs_real_t min_loc_vals[])
Minimum value of a real and the value of related array on all default communicator processes.
Definition: cs_parall.c:802
static void cs_parall_min(int n, cs_datatype_t datatype, void *val)
Minimum values of a given datatype on all default communicator processes.
Definition: cs_parall.h:234
void cs_parall_gather_f(int root_rank, int n_elts, int n_g_elts, const float array[], float g_array[])
Build a global array on the given root rank from all local arrays. Function dealing with single-preci...
Definition: cs_parall.c:1197
void cs_parall_max_loc_vals(int n, cs_real_t *max, cs_real_t max_loc_vals[])
Maximum value of a real and the value of related array on all default communicator processes.
Definition: cs_parall.c:764
size_t cs_parall_get_min_coll_buf_size(void)
Return minimum recommended scatter or gather buffer size.
Definition: cs_parall.c:1317
cs_e2n_sum_t
Definition: cs_parall.h:51
@ CS_E2N_SUM_SCATTER_ATOMIC
Definition: cs_parall.h:56
@ CS_E2N_SUM_SCATTER
Definition: cs_parall.h:53
@ CS_E2N_SUM_STORE_THEN_GATHER
Definition: cs_parall.h:60
@ CS_E2N_SUM_GATHER
Definition: cs_parall.h:58
cs_e2n_sum_t cs_glob_e2n_sum_type