1!------------------------------------------------------------------------------- 2 3! This file is part of Code_Saturne, a general-purpose CFD tool. 4! 5! Copyright (C) 1998-2021 EDF S.A. 6! 7! This program is free software; you can redistribute it and/or modify it under 8! the terms of the GNU General Public License as published by the Free Software 9! Foundation; either version 2 of the License, or (at your option) any later 10! version. 11! 12! This program is distributed in the hope that it will be useful, but WITHOUT 13! ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS 14! FOR A PARTICULAR PURPOSE. See the GNU General Public License for more 15! details. 16! 17! You should have received a copy of the GNU General Public License along with 18! this program; if not, write to the Free Software Foundation, Inc., 51 Franklin 19! Street, Fifth Floor, Boston, MA 02110-1301, USA. 20 21!------------------------------------------------------------------------------- 22 23!> \file parall.f90 24!> \brief Module for basic MPI and OpenMP parallelism-related values 25 26module parall 27 28 !============================================================================= 29 30 implicit none 31 32 !============================================================================= 33 34 !> \defgroup parall Module for basic MPI and OpenMP parallelism-related values 35 36 !> \addtogroup parall 37 !> \{ 38 39 !> thr_n_min : minimum number of elements for loops on threads 40 integer thr_n_min 41 parameter(thr_n_min = 128) 42 43 !> process rank 44 !> - -1 in sequential mode 45 !> - r (0 < r < n_processes) in distributed parallel run 46 integer, save :: irangp 47 48 !> number of processes (=1 if sequental) 49 integer, save :: nrangp 50 51 ! Global dimensions (i.e. independent of parallel partitioning) 52 53 !> global number of cells 54 integer(kind=8), save :: ncelgb 55 !> global number of interior faces 56 integer(kind=8), save :: nfacgb 57 !> global number of boundary faces 58 integer(kind=8), save :: nfbrgb 59 !> global number of vertices 60 integer(kind=8), save :: nsomgb 61 62 !> \} 63 64 !============================================================================= 65 66 interface 67 68 !--------------------------------------------------------------------------- 69 70 !> \brief Compute the global maximum of an integer in case of parellism. 71 72 !> \param[in, out] max local max in, global max out 73 74 subroutine parcmx(count) & 75 bind(C, name='cs_f_parall_max_i') 76 use, intrinsic :: iso_c_binding 77 implicit none 78 integer(c_int), intent(inout) :: count 79 end subroutine parcmx 80 81 !--------------------------------------------------------------------------- 82 83 !> \brief Compute the global maximum of a real number in case of parellism. 84 85 !> \param[in, out] max local max in, global max out 86 87 subroutine parmax(max) & 88 bind(C, name='cs_f_parall_max_r') 89 use, intrinsic :: iso_c_binding 90 implicit none 91 real(c_double), intent(inout) :: max 92 end subroutine parmax 93 94 !--------------------------------------------------------------------------- 95 96 !> \brief Compute the global minimum of an integer in case of parellism. 97 98 !> \param[in, out] min local min in, global min out 99 100 subroutine parcmn(count) & 101 bind(C, name='cs_f_parall_min_i') 102 use, intrinsic :: iso_c_binding 103 implicit none 104 integer(c_int), intent(inout) :: count 105 end subroutine parcmn 106 107 !--------------------------------------------------------------------------- 108 109 !> \brief Compute the global minimum of a real number in case of parellism. 110 111 !> \param[in, out] min local min in, global min out 112 113 subroutine parmin(min) & 114 bind(C, name='cs_f_parall_min_r') 115 use, intrinsic :: iso_c_binding 116 implicit none 117 real(c_double), intent(inout) :: min 118 end subroutine parmin 119 120 !--------------------------------------------------------------------------- 121 122 !> \brief Compute the global sum of an integer in case of parellism. 123 124 !> Note that for counters, on very large meshes, if the sum exceeds 125 !> 2**31, the result will be false on most machines. To avoid this, 126 !> using the C API (with counters as cs_gnum_t) is preferred. 127 128 !> \param[in, out] sum local sum in, global sum out 129 130 subroutine parcpt(count) & 131 bind(C, name='cs_f_parall_sum_i') 132 use, intrinsic :: iso_c_binding 133 implicit none 134 integer(c_int), intent(inout) :: count 135 end subroutine parcpt 136 137 !--------------------------------------------------------------------------- 138 139 !> \brief Compute the global sum of a real number in case of parellism. 140 141 !> \param[in, out] sum local sum in, global sum out 142 143 subroutine parsom(sum) & 144 bind(C, name='cs_f_parall_sum_r') 145 use, intrinsic :: iso_c_binding 146 implicit none 147 real(c_double), intent(inout) :: sum 148 end subroutine parsom 149 150 !--------------------------------------------------------------------------- 151 152 !> \brief Compute the global maxima of an array of integers 153 !> in case of parellism. 154 155 !> \param[in] n_elts size of array 156 !> \param[in, out] max local max in, global max out 157 158 subroutine parimx(n_elts, array) & 159 bind(C, name='cs_f_parall_max_n_i') 160 use, intrinsic :: iso_c_binding 161 implicit none 162 integer(c_int), value :: n_elts 163 integer(c_int), dimension(*), intent(inout) :: array 164 end subroutine parimx 165 166 !--------------------------------------------------------------------------- 167 168 !> \brief Compute the global maxima of an array of real numbers 169 !> in case of parellism. 170 171 !> \param[in] n_elts size of array 172 !> \param[in, out] max local max in, global max out 173 174 subroutine parrmx(n_elts, array) & 175 bind(C, name='cs_f_parall_max_n_r') 176 use, intrinsic :: iso_c_binding 177 implicit none 178 integer(c_int), value :: n_elts 179 real(c_double), dimension(*), intent(inout) :: array 180 end subroutine parrmx 181 182 !--------------------------------------------------------------------------- 183 184 !> \brief Compute the global minima of an array of integers 185 !> in case of parellism. 186 187 !> \param[in] n_elts size of array 188 !> \param[in, out] min local min in, global min out 189 190 subroutine parimn(n_elts, array) & 191 bind(C, name='cs_f_parall_min_n_i') 192 use, intrinsic :: iso_c_binding 193 implicit none 194 integer(c_int), value :: n_elts 195 integer(c_int), dimension(*), intent(inout) :: array 196 end subroutine parimn 197 198 !--------------------------------------------------------------------------- 199 200 !> \brief Compute the global minima of an array of real numbers 201 !> in case of parellism. 202 203 !> \param[in] n_elts size of array 204 !> \param[in, out] min local min in, global min out 205 206 subroutine parrmn(n_elts, array) & 207 bind(C, name='cs_f_parall_min_n_r') 208 use, intrinsic :: iso_c_binding 209 implicit none 210 integer(c_int), value :: n_elts 211 real(c_double), dimension(*), intent(inout) :: array 212 end subroutine parrmn 213 214 !--------------------------------------------------------------------------- 215 216 !> \brief Compute the global sums of an array of integers 217 !> in case of parellism. 218 219 !> Note that for counters, on very large meshes, if a sum exceeds 220 !> 2**31, the resuly will be false on most machines. To avoid this, 221 !> using the C API (with counters as cs_gnum_t) is preferred. 222 223 !> \param[in] n_elts size of array 224 !> \param[in, out] sum local sum in, global sum out 225 226 subroutine parism(n_elts, array) & 227 bind(C, name='cs_f_parall_sum_n_i') 228 use, intrinsic :: iso_c_binding 229 implicit none 230 integer(c_int), value :: n_elts 231 integer(c_int), dimension(*), intent(inout) :: array 232 end subroutine parism 233 234 !--------------------------------------------------------------------------- 235 236 !> \brief Compute the global sums of an array of real numbers 237 !> in case of parellism. 238 239 !> \param[in] n_elts size of array 240 !> \param[in, out] sum local sum in, global sum out 241 242 subroutine parrsm(n_elts, array) & 243 bind(C, name='cs_f_parall_sum_n_r') 244 use, intrinsic :: iso_c_binding 245 implicit none 246 integer(c_int), value :: n_elts 247 real(c_double), dimension(*), intent(inout) :: array 248 end subroutine parrsm 249 250 !--------------------------------------------------------------------------- 251 252 !> \brief Broadcast an integer in case of parellism. 253 254 !> \param[in] root_rank rank of the sending process 255 !> \param[in, out] val value to broadcast 256 !> (input on root_rank, output on others) 257 258 subroutine parall_bcast_i(root_rank, val) & 259 bind(C, name='cs_f_parall_bcast_i') 260 use, intrinsic :: iso_c_binding 261 implicit none 262 integer(c_int), value :: root_rank 263 integer(c_int), intent(inout) :: val 264 end subroutine parall_bcast_i 265 266 !--------------------------------------------------------------------------- 267 268 !> \brief Broadcast a real number in case of parellism. 269 270 !> \param[in] root_rank rank of the sending process 271 !> \param[in, out] val value to broadcast 272 !> (input on root_rank, output on others) 273 274 subroutine parall_bcast_r(root_rank, val) & 275 bind(C, name='cs_f_parall_bcast_r') 276 use, intrinsic :: iso_c_binding 277 implicit none 278 integer(c_int), value :: root_rank 279 real(c_double), intent(inout) :: val 280 end subroutine parall_bcast_r 281 282 !--------------------------------------------------------------------------- 283 284 !> \brief Broadcast an array of integers in case of parellism. 285 286 !> \param[in] root_rank rank of the sending process 287 !> \param[in] n_elts size of array 288 !> \param[in, out] array array to broadcast 289 !> (input on root_rank, output on others) 290 291 subroutine parbci(root_rank, n_elts, array) & 292 bind(C, name='cs_f_parall_bcast_n_i') 293 use, intrinsic :: iso_c_binding 294 implicit none 295 integer(c_int), value :: root_rank, n_elts 296 integer(c_int), dimension(*), intent(inout) :: array 297 end subroutine parbci 298 299 !--------------------------------------------------------------------------- 300 301 !> \brief Broadcast an array of real numbers in case of parellism. 302 303 !> \param[in] root_rank rank of the sending process 304 !> \param[in] n_elts size of array 305 !> \param[in, out] array array to broadcast 306 !> (input on root_rank, output on others) 307 308 subroutine parbcr(root_rank, n_elts, array) & 309 bind(C, name='cs_f_parall_bcast_n_r') 310 use, intrinsic :: iso_c_binding 311 implicit none 312 integer(c_int), value :: root_rank, n_elts 313 real(c_double), dimension(*), intent(inout) :: array 314 end subroutine parbcr 315 316 !--------------------------------------------------------------------------- 317 318 !> \brief Maximum value of a real and the value of related array on all 319 !> default communicator processes. 320 321 !> \param[in] n size of the related array 322 !> \param[in, out] max local max in, global max out 323 !> \param[in, out] max_loc_vals array values at location of local max in, 324 !> and at location of global max out 325 326 subroutine parmxl(n, max, max_loc_vals) & 327 bind(C, name='cs_parall_max_loc_vals') 328 use, intrinsic :: iso_c_binding 329 implicit none 330 integer(c_int), value :: n 331 real(c_double), intent(inout) :: max 332 real(c_double), dimension(*), intent(inout) :: max_loc_vals 333 end subroutine parmxl 334 335 !--------------------------------------------------------------------------- 336 337 !> \brief Minimum value of a real and the value of related array on all 338 !> default communicator processes. 339 340 !> \param[in] n size of the related array 341 !> \param[in, out] min local min in, global min out 342 !> \param[in, out] min_loc_vals array values at location of local min in, 343 !> and at location of global min out 344 345 subroutine parmnl(n, min, min_loc_vals) & 346 bind(C, name='cs_parall_min_loc_vals') 347 use, intrinsic :: iso_c_binding 348 implicit none 349 integer(c_int), value :: n 350 real(c_double), intent(inout) :: min 351 real(c_double), dimension(*), intent(inout) :: min_loc_vals 352 end subroutine parmnl 353 354 !--------------------------------------------------------------------------- 355 356 !> \brief Given an (id, rank, value) tuple, return the local id, rank, 357 !> and value corresponding to the global minimum value. 358 359 !> \param[in, out] elt_id element id for which the value is the smallest 360 !> (local in, global out) 361 !> \param[in, out] rank_id rank id for which the value is the smallest 362 !> (local in, global out) 363 !> \param[in] val associated local minimum value 364 365 subroutine parfpt(elt_id, rank_id, val) & 366 bind(C, name='cs_parall_min_id_rank_r') 367 use, intrinsic :: iso_c_binding 368 implicit none 369 integer(c_int), intent(inout) :: elt_id, rank_id 370 real(c_double), value :: val 371 end subroutine parfpt 372 373 !--------------------------------------------------------------------------- 374 375 !> \brief Build a global array from each local array in each domain. 376 377 !> Local arrays are appened in order of owning MPI rank. 378 !> The size of each local array may be different. 379 380 !> Use of this function may be quite practical, but should be limited 381 !> to user functions, as it may limit scalability (especially as regards 382 !> memory usage). 383 384 !> \param[in] n_elts size of the local array 385 !> \param[in] n_g_elts size of the global array 386 !> \param[in] array local array (size: n_elts) 387 !> \param[out] g_array global array (size: n_g_elts) 388 389 subroutine cs_parall_allgather_r(n_elts, n_g_elts, array, g_array) & 390 bind(C, name='cs_parall_allgather_r') 391 use, intrinsic :: iso_c_binding 392 implicit none 393 integer(c_int), value :: n_elts, n_g_elts 394 real(c_double), dimension(*), intent(in) :: array 395 real(c_double), dimension(*), intent(inout) :: g_array 396 end subroutine cs_parall_allgather_r 397 398 !--------------------------------------------------------------------------- 399 400 !> \brief Set a barrier on all default communicator processes. 401 402 !> The this function will exit only once all processes have called 403 !> it. This is not recommended for production code, but may be useful 404 !> when debugging. 405 406 subroutine parbar() & 407 bind(C, name='cs_f_parall_barrier') 408 use, intrinsic :: iso_c_binding 409 implicit none 410 end subroutine parbar 411 412 !--------------------------------------------------------------------------- 413 414 end interface 415 416 !============================================================================= 417 418contains 419 420 !============================================================================= 421 422 !> \brief Build a global array from each local array in each domain. 423 424 !> Local arrays are appened in order of owning MPI rank. 425 !> The size of each local array may be different. 426 427 !> Use of this function may be quite practical, but should be limited 428 !> to user functions, as it may limit scalability (especially as regards 429 !> memory usage). 430 431 !> \param[in] n_elts size of the local array 432 !> \param[in] n_g_elts size of the global array 433 !> \param[in] array local array (size: n_elts) 434 !> \param[out] g_array global array (size: n_g_elts) 435 436 subroutine paragv(n_elts, n_g_elts, array, g_array) 437 use, intrinsic :: iso_c_binding 438 implicit none 439 integer(c_int), value :: n_elts, n_g_elts 440 real(c_double), dimension(:), intent(in) :: array 441 real(c_double), dimension(:), intent(inout) :: g_array 442 call cs_parall_allgather_r(n_elts, n_g_elts, array(:), g_array(:)) 443 end subroutine paragv 444 445 !============================================================================= 446 447end module parall 448 449 450