MAPL_CollectiveScatter3D Subroutine

public subroutine MAPL_CollectiveScatter3D(Grid, GlobArray, LocArray, hw, rc)

Arguments

Type IntentOptional Attributes Name
type(ESMF_Grid), intent(in) :: Grid
real, intent(in) :: GlobArray(:,:,:)
real, intent(inout), target :: LocArray(:,:,:)
integer, intent(in), optional :: hw
integer, intent(out), optional :: rc

Calls

proc~~mapl_collectivescatter3d~~CallsGraph proc~mapl_collectivescatter3d MAPL_CollectiveScatter3D ESMF_VMGet ESMF_VMGet proc~mapl_collectivescatter3d->ESMF_VMGet ESMF_VMGetCurrent ESMF_VMGetCurrent proc~mapl_collectivescatter3d->ESMF_VMGetCurrent interface~mapl_arrayiscatter MAPL_ArrayIScatter proc~mapl_collectivescatter3d->interface~mapl_arrayiscatter interface~mapl_assert MAPL_Assert proc~mapl_collectivescatter3d->interface~mapl_assert proc~mapl_collectivewait MAPL_CollectiveWait proc~mapl_collectivescatter3d->proc~mapl_collectivewait proc~mapl_createrequest MAPL_CreateRequest proc~mapl_collectivescatter3d->proc~mapl_createrequest proc~mapl_return MAPL_Return proc~mapl_collectivescatter3d->proc~mapl_return proc~mapl_roundrobinpelist MAPL_RoundRobinPEList proc~mapl_collectivescatter3d->proc~mapl_roundrobinpelist proc~mapl_verify MAPL_Verify proc~mapl_collectivescatter3d->proc~mapl_verify proc~mapl_collectivewait->proc~mapl_return proc~mapl_collectivewait->proc~mapl_verify mpi_recv mpi_recv proc~mapl_collectivewait->mpi_recv mpi_wait mpi_wait proc~mapl_collectivewait->mpi_wait proc~mapl_createrequest->ESMF_VMGet proc~mapl_createrequest->ESMF_VMGetCurrent proc~mapl_createrequest->interface~mapl_assert proc~mapl_createrequest->proc~mapl_return proc~mapl_createrequest->proc~mapl_verify ESMF_GridGet ESMF_GridGet proc~mapl_createrequest->ESMF_GridGet mpi_irecv mpi_irecv proc~mapl_createrequest->mpi_irecv proc~mapl_distgridget MAPL_DistGridGet proc~mapl_createrequest->proc~mapl_distgridget at at proc~mapl_return->at insert insert proc~mapl_return->insert proc~mapl_throw_exception MAPL_throw_exception proc~mapl_return->proc~mapl_throw_exception proc~mapl_roundrobinpelist->proc~mapl_return proc~mapl_roundrobinpelist->proc~mapl_verify interface~mapl_getnewrank MAPL_GetNewRank proc~mapl_roundrobinpelist->interface~mapl_getnewrank proc~mapl_verify->proc~mapl_throw_exception proc~mapl_distgridget->proc~mapl_verify ESMF_DistGridGet ESMF_DistGridGet proc~mapl_distgridget->ESMF_DistGridGet

Source Code

  subroutine MAPL_CollectiveScatter3D(Grid, GlobArray, LocArray, hw, rc)

    type (ESMF_Grid),        intent(IN   ) :: Grid
    real, target,            intent(INOUT) :: LocArray(:,:,:)
    real,                    intent(IN   ) :: GlobArray(:,:,:)
    integer, optional,       intent(IN   ) :: hw
    integer, optional,       intent(  OUT) :: rc

! Locals
!-------

    integer                       :: status


    type (MAPL_CommRequest)       :: reqs(size(LocArray,3))
    integer                       :: root(size(LocArray,3))
    integer                       :: nNodes
    integer                       :: LM, L, nc, npes, mype
    integer                       :: nn
    type(ESMF_VM)                 :: VM
    logical                       :: HaveGlobal
    integer                       :: comm
    integer                       :: hw_

! Begin
!------

    call ESMF_VMGetCurrent(VM,         RC=STATUS)
    _VERIFY(STATUS)
    call ESMF_VMGet(VM, petcount=npes, localpet=MYPE, mpiCommunicator=comm, RC=STATUS)
    _VERIFY(STATUS)

    if(present(hw)) then
       hw_ = hw
    else
       hw_ = 0
    endif

    nNodes = size(MAPL_NodeRankList)
    call MAPL_RoundRobinPEList(Root, nNodes, RC=STATUS)
    _VERIFY(STATUS)

    LM = size(LocArray,3)
    NC = count(Root==mype)

    HaveGlobal = NC>0

    do L=1,LM
       call MAPL_CreateRequest(GRID, Root(L), reqs(L), tag=L,       &
                               RequestType=MAPL_IsScatter,          &
                               DstArray=LocArray(:,:,L),            &
                               PrePost=.true., hw=hw_, RC=STATUS)
       _VERIFY(STATUS)
    enddo

    if(HaveGlobal) then
       _ASSERT(size(GlobArray,3)==NC, 'inconsisntent rank')

       nn = 0
       do L=1,LM
          if(Root(L)==mype) then

             nn = nn + 1
             call MAPL_ArrayIScatter (GlobArray(:,:,nn), reqs(L), hw=hw_, RC=STATUS)
             _VERIFY(STATUS)
             if(nn==NC) exit
          endif
       enddo
    end if

    do L=1,LM
       call MAPL_CollectiveWait(reqs(L), rc=status)
       _VERIFY(STATUS)
    end do

    _RETURN(ESMF_SUCCESS)
  end subroutine MAPL_CollectiveScatter3D