Creating an "MPI-Scatter-like" routine - interface design

1

So I’m currently working on creating a routine that does something quite similar MPI_Scatter, in that a single task has some field that we want to send portions of to other processes within a communicator. Lets call the source fld global_fld, and the resulting, split field split_fld, and I want to call it like:

call split_global_fld(global_fld, split_fld)

noting that global_fld will really only have valid data on the root task, but the routine needs to be called by participating tasks.

This brings me to my question, what is a good way to setup this interface? My first idea was to set it up as follows:

subroutine split_global_fld(global_fld, split_fld)
   implicit none
   real, allocatable, dimension(:,:), intent(in   ) :: global_fld ! should only be allocated on root task
   real, allocatable, dimension(:,:), intent(  out) :: split_fld  ! will be allocated on all tasks

   ... 
   ! perform necessary MPI communication to split array

and then in the calling routine or program, it would look something like

real, allocatable, dimension(:,:) :: global_fld
real, allocatable, dimension(:,:) :: split_fld
...

if ( myrank == main ) then
    call populate_global_fld(global_fld) ! this returns an allocated array
endif
...
! split data across leaders
call split_global_fld(global_fld, split_fld) ! noting that global_fld isn't allocated on non root tasks

In some early tests, this seems to nominally work, but it feels odd and my colleague is worried that the behaviour here is undefined and the behaviour you actually get may be compiler dependent. Is this correct? Even if it isn’t, is this bad practice?

If so, I was wondering how others might try to set this up? Another solution would be to have something like this in the calling program:

if ( myrank == main ) then
   call populate_global_fld(global_fld)
else
   ! spoof field to satisfy interfaces
   allocate(global_fld(1,1))
   global_fld(:,:) = 0
end if

and getting rid of the allocatable key-word in the interface specification.

1 Like
2

First of all, welcome!

Sorry for a very stupid question, why not just using MPI_scatter or MPI_scatterv?

3

All good @CRquantum - its a valid question

Essentially we have some routines that act as an “interface” between two separate portions of our modelling system - I’ll refer to one as the coupler, and one as the atmosphere. One side of it handles receiving data from another executable running in MPMD mode, which then needs to get sent to the other side to be scattered in what ever way necessary. Effectively we are trying to separate the namespaces of the two sides. This lets us set up a neat transfer interface like this:

subroutine transfer_fields(step_count)
   use cpl_cap,   only : cpl_receive_fld, is_transfer_step
   use atm_utils, only : store_atmosphere_code
   implicit none

   ! INPUTS
   integer, intent(in) :: step_count
   
   ! LOCAL
   real, allocatable, dimension(:,:) :: global_field
   
   if ( is_transfer_step(step_count) ) then
       if ( myrank == atmosphere_root ) then
           call cpl_receive_field(global_field, "field_name") ! get global field from coupler exec (running in MPMD mode)
       end if
       ! store global field into atmosphere namespace to be picked up by later calculations
       call store_atmosphere_field(global_field, "field_name")
   end if
end subroutine transfer_fields

Which satisfies a few design ideas:

  • this code only brings in routines from eitherside, no variable imports
  • the coupler side can’t alter the agcm namespace
  • the agcm side can’t alter the coupler namespace

which allows developers on either side of the interface to change anything, provided they keep the arguments for the few imported routines the same.

We could definitely just import info on the communicators/ranks and other things from the other side, but they the coupler cap and the atmosphere code will become more tightly coupled and our developers will need to be more mindful of stepping on each other’s toes

4

@cseinen ,

Welcome to the forum!

A quick question for you first: have you tested COARRAYs in Fortran and tried out CO_BROADCAST`?

5

@FortranFan I personally have, yes! But its been a while since I’ve used them. I should note that I’m working on the coupler side of this interface, and interfacing with an atmosphere model code base that I don’t have heavy control over. On the atmosphere side of the split/scatter operations, they have might their own comm worlds with leaders that I don’t want to make assumptions about, which I why I just want to hand off the global array to their side and let their codes handle the splitting/scattering (note that in the current implementation, we first need to split the array onto two different root tasks within different comm worlds, before scattering to local fields within their respective comm worlds). I’ve wanted to look at applying co-arrays, in some way, in our system, but given the majority of our codes use a combination of MPI/OpenMP, I’ve been hesitant to muddy the waters with another system. That said, it is part of the standard now, so I’m open to ideas if I might be able to use it locally to some routines I have control over.

How are you thinking they could be applied here? Are you thinking I could just use CO_BROADCAST to broadcast the global array to all tasks? I am trying to avoid needing to populate the global array on every process, as it will be quite large