Understanding MPI_Gather

Hi there,

I’m relatively new to Fortran and to this community. I’m trying to write some code that performs some calculations on a multidimensional array using MPI. And for the life of me I can’t get MPI_Gather to work, and I’m sure it’s something very obvious that I’m missing here.

My code is:

program test

    USE mpi
    USE netcdf

    IMPLICIT NONE

    INTEGER,ALLOCATABLE,DIMENSION(:,:,:) :: array, global_array
    INTEGER :: j_start, j_end, i_start, i_end, J, I, t

    ! MPI
    INTEGER :: mpirank, mpisize, mpierr

    ! Initialize MPI
    call MPI_INIT(mpierr)
    call MPI_COMM_SIZE(MPI_COMM_WORLD, mpisize, mpierr)
    call MPI_COMM_RANK(MPI_COMM_WORLD, mpirank, mpierr)

    ALLOCATE(array(4,4,2))
    ALLOCATE(global_array(4,4,2))
    
    j_start = 1
    j_end   = 4
    i_start = 1 + (mpirank*4/mpisize)
    i_end   = (mpirank + 1)*4/mpisize
    
    DO t=1,2

        DO J=j_start,j_end

            DO I=i_start,i_end

                array(J,I,t) = J+I+t
            
            END DO

        END DO ! J

    END DO ! t

    IF (mpirank==0) WRITE(*,*) 'before MPI_Gather...'
    WRITE(*,*) mpirank, array(1,:,1)
    
    CALL MPI_Barrier(MPI_COMM_WORLD, mpierr)

    CALL MPI_Gather(array(:,i_start:i_end,:), 4*2*2, MPI_INTEGER, & 
                    global_array(:,i_start:i_end,:), 4*2*2, MPI_INTEGER, &
                    0, MPI_COMM_WORLD, mpierr)

    IF (mpirank == 0) THEN
        WRITE(*,*) 'AFTER MPI_Gather...'
        DO J=1,4
            WRITE(*,*) (global_array(J,I,t), I=1,4,1)
        END DO
    END IF
    
    DEALLOCATE(array)
    DEALLOCATE(global_array)

    ! Finalise MPI
    CALL MPI_FINALIZE(mpierr)

end program test

And I get these results, which are obviously wrong:

 before MPI_Gather...
           1           0           0           5           6
           0           3           4           0           0
 AFTER MPI_Gather...
           0     5289440           0           0
           0           0           0           0
         225           0           0           0
           0           0           0           0

I’m using Intel Fortran 2022, and I compile the program with mpiifort

Any help/comment/advice will be truly appreciated!

1 Like

Welcome!

You can find mpi_gather here,
https://mpitutorial.com/tutorials/mpi-scatter-gather-and-allgather/

The below illustration is helpful,

Note that the ‘‘gather’’ array on the rank 0 process (your global_array), needs to be the ‘‘total’’ size of all the small arrays (your array) on each of the processes.
In your case, your global_array seems to have the same size of your array (both have 4*4*2=32 elements). If not careful enough, it probably will not work as you expected.

Since you have i_start and i_end on each process, perhaps you are right. But the 4*4*2 looks wrong to me, because it does not look like you want to gather all the 4*4*2 elements of array on each process to the global_array on rank 0 process, right?
If you want to gather all the 4*4*2 elements of array on each process to the global_array on rank 0 process, your `global_array’ needs to be the size of 4*4*2*(total number of process).

I personally suggest that, you may begin with mpi_gather for 1D array first. You know, gather all the 1D arrays on all the process to the 1D array on the rank 0 process. Once you figure that out how the 1D array mpi_gather works, you could do 2D array gathering.

Below is my wrapper code for 1D array gathering.

  subroutine gatherr1d(r,rgather)
  real(kind=r8) :: r(:),rgather(:)
  integer :: ierror
  call mpi_gather(r,size(r),mpir8,rgather,size(r),mpir8,0, &
    mpi_comm_world,ierror)
  return
  end subroutine gatherr1d 

It gathers r array on each process to the rgather array on the rank 0 process. Note that, the size of rgather array needs to be size(r) * number of processes.

You can define a wrapper of mpi_gather by overloading like the below module shows,

  module mympi
  use mpi
  implicit none
  !include 'mpif.h'
  integer, private, parameter :: i4=selected_int_kind(9)
  integer, private, parameter :: i8=selected_int_kind(15)
  integer, private, parameter :: r8=selected_real_kind(15,9)
  integer, private, save :: mpii4,mpii8,mpir8,mpic8 ! c8 not defined now, use mpi_double_complex
  integer(kind=i4), private, save :: irank,iproc,irank_tot

interface gather ! gather to process 0
   module procedure gatheri1,gatheri1d,gatheri1_8
   module procedure gatherr1,gatherr1d
   module procedure gatherr2d_1d,gatherr3d_1d,gatherr4d_1d,gatherr5d_1d
end interface gather

contains
  subroutine init0 ! call this before anything else
  integer :: ierror,isize,ir,ip
  integer(kind=i4) :: itest4
  integer(kind=i8) :: itest8
  real(kind=r8) :: rtest8
  complex(kind=r8) :: ctest8
  call mpi_init(ierror)
  call mpi_comm_rank(mpi_comm_world,ir,ierror)
  irank=ir
  call mpi_comm_size(mpi_comm_world,ip,ierror)
  iproc=ip
  irank_tot=iproc-1
  call mpi_sizeof(itest4,isize,ierror)
  call mpi_type_match_size(mpi_typeclass_integer,isize,mpii4,ierror)
  call mpi_sizeof(itest8,isize,ierror)
  call mpi_type_match_size(mpi_typeclass_integer,isize,mpii8,ierror)
  call mpi_sizeof(rtest8,isize,ierror)
  call mpi_type_match_size(mpi_typeclass_real,isize,mpir8,ierror)
  !call mpi_sizeof(ctest8,isize,ierror)
  !call mpi_type_match_size(mpi_typeclass_complex,isize,mpic8,ierror)
  return
  end subroutine init0

  subroutine done ! wrapper for finalize routine
  integer :: ierror
  call mpi_finalize(ierror)
  return
  end subroutine done

  subroutine gatheri1(i,igather)
  integer(kind=i4) :: i,igather(:)
  integer :: ierror
  call mpi_gather(i,1,mpii4,igather,1,mpii4,0,mpi_comm_world,ierror)
  return
  end subroutine gatheri1

  subroutine gatheri1d(i,igather)
  integer(kind=i4) :: i(:),igather(:)
  integer :: ierror
  call mpi_gather(i,size(i),mpii4,igather,size(i),mpii4,0, &
    mpi_comm_world,ierror)
  return
  end subroutine gatheri1d

  subroutine gatheri1_8(i,igather)
  integer(kind=i8) :: i,igather(:)
  integer :: ierror
  call mpi_gather(i,1,mpii8,igather,1,mpii8,0, &
    mpi_comm_world,ierror)
  return
  end subroutine gatheri1_8   

  subroutine gatherr1(r,rgather)
  real(kind=r8) :: r,rgather(:)
  integer :: ierror
  call mpi_gather(r,1,mpir8,rgather,1,mpir8,0,mpi_comm_world,ierror)
  return
  end subroutine gatherr1

  subroutine gatherr1d(r,rgather)
  real(kind=r8) :: r(:),rgather(:)
  integer :: ierror
  call mpi_gather(r,size(r),mpir8,rgather,size(r),mpir8,0, &
    mpi_comm_world,ierror)
  return
  end subroutine gatherr1d  
   
  subroutine gatherr2d_1d(r,rgather)
  real(kind=r8) :: r(:,:),rgather(:)
  integer :: ierror
  call mpi_gather(r,size(r),mpir8,rgather,size(r),mpir8,0, &
    mpi_comm_world,ierror)
  return
  end subroutine gatherr2d_1d 
   
  subroutine gatherr3d_1d(r,rgather)
  real(kind=r8) :: r(:,:,:),rgather(:)
  integer :: ierror
  call mpi_gather(r,size(r),mpir8,rgather,size(r),mpir8,0, &
    mpi_comm_world,ierror)
  return
  end subroutine gatherr3d_1d
   
  subroutine gatherr4d_1d(r,rgather)
  real(kind=r8) :: r(:,:,:,:),rgather(:)
  integer :: ierror
  call mpi_gather(r,size(r),mpir8,rgather,size(r),mpir8,0, &
    mpi_comm_world,ierror)
  return
  end subroutine gatherr4d_1d
   
  subroutine gatherr5d_1d(r,rgather)
  real(kind=r8) :: r(:,:,:,:,:),rgather(:)
  integer :: ierror
  call mpi_gather(r,size(r),mpir8,rgather,size(r),mpir8,0, &
    mpi_comm_world,ierror)
  return
  end subroutine gatherr5d_1d

Then in the code, if you wanted to gather something, just type things like

call gather(aaa, bbb )

Hi @CRquantum , first of all happy birthday! :slight_smile: And thank you so much for your detailed examples! I grasp the concept of MPI_Gather and MPI_Reduce a bit better now…

I will play around with this to see if I can get it to work, I think I’ll have to do a significant re-structure of my code (the one I posted was just an example to understand the concept).

Thanks again!

1 Like