From 66f5fc5c58f30fc4d7e871abfdfc23bdfbbb012f Mon Sep 17 00:00:00 2001 From: Grzegorz Kowal Date: Sun, 17 Dec 2023 21:47:38 -0300 Subject: [PATCH] PARAMETERS: Read parameters on master and distribute them. This commit restores the previous way of processing the parameter file. It is read and processed on the MPI master process, and then the list of parameters is distributed to other MPI processes. This way only one process accesses the parameter file, reducing the number of I/O operations, which can be significant in the case of multiprocess MPI jobs. Signed-off-by: Grzegorz Kowal --- sources/mpitools.F90 | 77 +++++++++++++++++++++++++++++++++++++ sources/parameters.F90 | 87 ++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 164 insertions(+) diff --git a/sources/mpitools.F90 b/sources/mpitools.F90 index e04915e..bf40e3c 100644 --- a/sources/mpitools.F90 +++ b/sources/mpitools.F90 @@ -41,6 +41,10 @@ module mpitools ! subroutine interfaces ! #ifdef MPI + interface broadcast + module procedure broadcast_integer + module procedure broadcast_string + end interface interface reduce_minimum module procedure reduce_minimum_double_array end interface @@ -83,6 +87,7 @@ module mpitools public :: initialize_mpitools, finalize_mpitools public :: check_status #ifdef MPI + public :: broadcast public :: reduce_minimum, reduce_maximum, reduce_sum public :: send_array, receive_array public :: exchange_arrays @@ -353,6 +358,78 @@ module mpitools ! !=============================================================================== ! +! subroutine BROADCAST_INTEGER: +! ---------------------------- +! +! Subroutine broadcasts an integer buffer. +! +!=============================================================================== +! + subroutine broadcast_integer(buf) + + use helpers, only : print_message + + implicit none + + integer, dimension(..), intent(inout) :: buf + + integer :: ierror + + character(len=*), parameter :: loc = 'MPITOOLS::broadcast_integer()' + +!------------------------------------------------------------------------------- +! + call start_timer(imc) + + call mpi_bcast(buf, size(buf), MPI_INTEGER, 0, MPI_COMM_WORLD, ierror) + + if (ierror /= MPI_SUCCESS) & + call print_message(loc, 'Could not broadcast an integer buffer.') + + call stop_timer(imc) + +!------------------------------------------------------------------------------- +! + end subroutine broadcast_integer +! +!=============================================================================== +! +! subroutine BROADCAST_STRING: +! --------------------------- +! +! Subroutine broadcasts a string buffer. +! +!=============================================================================== +! + subroutine broadcast_string(buf) + + use helpers, only : print_message + + implicit none + + character(len=*), intent(inout) :: buf + + integer :: ierror + + character(len=*), parameter :: loc = 'MPITOOLS::broadcast_string()' + +!------------------------------------------------------------------------------- +! + call start_timer(imc) + + call mpi_bcast(buf, len(buf), MPI_CHARACTER, 0, MPI_COMM_WORLD, ierror) + + if (ierror /= MPI_SUCCESS) & + call print_message(loc, 'Could not broadcast a string buffer.') + + call stop_timer(imc) + +!------------------------------------------------------------------------------- +! + end subroutine broadcast_string +! +!=============================================================================== +! ! subroutine SEND_ARRAY: ! --------------------- ! diff --git a/sources/parameters.F90 b/sources/parameters.F90 index 86fcb2e..a64bea8 100644 --- a/sources/parameters.F90 +++ b/sources/parameters.F90 @@ -139,6 +139,10 @@ module parameters status = 111 end if +#ifdef MPI + call distribute_parameters() +#endif /* MPI */ + !------------------------------------------------------------------------------- ! end subroutine read_parameters @@ -196,6 +200,7 @@ module parameters subroutine parse_parameters(status) use iso_fortran_env, only : error_unit + use mpitools , only : master implicit none @@ -215,6 +220,7 @@ module parameters ! status = 0 + if (.not. master) return n = 0 j = 1024 @@ -457,6 +463,87 @@ module parameters !------------------------------------------------------------------------------- ! end subroutine get_parameter_string +#ifdef MPI +!=============================================================================== +! +! subroutine DISTRIBUTE_PARAMETERS: +! -------------------------------- +! +! Description: +! +! Subroutine distributes parameters among the MPI processes. +! +!=============================================================================== +! + subroutine distribute_parameters() + + use mpitools, only : master, broadcast + + implicit none + + type(item), pointer :: item_ptr + + character(len=:), allocatable :: str + + integer, dimension(2) :: counters ! 1: nitems, 2: maxlen + integer :: n, i + +!------------------------------------------------------------------------------- +! + counters = 0 + + if (master) then + item_ptr => parameter_list + do while(associated(item_ptr)) + counters(1) = counters(1) + 1 + counters(2) = max(counters(2), len(item_ptr%key // '|' // item_ptr%val)) + item_ptr => item_ptr%next + end do + end if + +! broadcast the number of items and the maximum item length +! + call broadcast(counters) + +! allocate string buffer +! + allocate(character(len=counters(2)) :: str) + +! iterate over all items in the list and broadcast them +! + if (master) then + item_ptr => parameter_list + do while(associated(item_ptr)) + write(str,"(a)") item_ptr%key // '|' // item_ptr%val + + call broadcast(str) + + item_ptr => item_ptr%next + end do + + else + + do n = 1, counters(1) + + call broadcast(str) + + i = index(str, '|') + + allocate(item_ptr) + item_ptr%key = trim(adjustl(str(:i-1))) + item_ptr%val = trim(adjustl(str(i+1:))) + item_ptr%next => parameter_list + parameter_list => item_ptr + end do + + end if + + deallocate(str) + +!------------------------------------------------------------------------------- +! + end subroutine distribute_parameters +#endif /* MPI */ !=============================================================================== !