populate_redist_array Subroutine

private subroutine populate_redist_array(Array, DistGrid, Memory, Grid, rc)

Arguments

Type IntentOptional Attributes Name
type(ESMF_Array), intent(inout) :: Array
type(ESMF_DistGrid), intent(in) :: DistGrid
type(memory_config), intent(in) :: Memory
type(grid_specification_record), intent(in) :: Grid
integer, intent(inout) :: rc

Source Code

  subroutine populate_redist_array(Array, DistGrid, Memory, Grid, rc)
  !-----------------------------------------------------------------------------
  ! routine populates an esmf array to the values used for a redist test. These
  ! values are dependent on the global coordinates and the local DE number. For
  ! a rank 2 array the values are set to Value = localDE + 1000*i1 + 1/1000 *i2
  !
  !-----------------------------------------------------------------------------
  ! arguments
  type(ESMF_Array), intent(inout) :: Array
  type(ESMF_DistGrid), intent(in   ) :: DistGrid
  type(memory_config), intent(in   ) :: Memory
  type(grid_specification_record), intent(in   ) :: Grid
  integer, intent(inout) :: rc

  ! local ESMF types
  type(ESMF_LocalArray), allocatable :: larrayList(:)
  type(ESMF_Index_Flag) :: indexflag

  ! local parameters
  integer :: localrc ! local error status

  ! local integer variables
  integer :: de, localDeCount, dimCount
  integer, allocatable ::  localDeToDeMap(:)
  integer, allocatable :: LBnd(:,:), UBnd(:,:)
  integer :: i1, i2, i3, i4, i5, i6, i7
  ! integer :: irank, k, tensorsize, fsize(7)
  ! integer, allocatable :: haloL(:), haloR(:)
  ! integer, allocatable :: top(:), bottom(:)
  integer :: allocRcToTest

  ! local real variables
  real(ESMF_KIND_R8), pointer :: farrayPtr1(:), farrayPtr2(:,:)
  real(ESMF_KIND_R8), pointer :: farrayPtr3(:,:,:)
  real(ESMF_KIND_R8), pointer :: farrayPtr4(:,:,:,:)
  real(ESMF_KIND_R8), pointer :: farrayPtr5(:,:,:,:,:)
  real(ESMF_KIND_R8), pointer :: farrayPtr6(:,:,:,:,:,:)
  real(ESMF_KIND_R8), pointer :: farrayPtr7(:,:,:,:,:,:,:)

  ! initialize return flag
  localrc = ESMF_RC_NOT_IMPL
  rc = ESMF_RC_NOT_IMPL

  !-----------------------------------------------------------------------------
  ! get local array DE list
  !-----------------------------------------------------------------------------
  call ESMF_ArrayGet(array, localDeCount=localDeCount, rc=localrc)
  if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error getting local DE count from array", &
          rcToReturn=rc)) return

  allocate(localDeToDeMap(localDeCount), stat=allocRcToTest)
  if (ESMF_LogFoundAllocError(allocRcToTest, msg="integer variable "//          &
     " localDeToDeMap in populate_redist_array", rcToReturn=rc)) then
  endif
  call ESMF_ArrayGet(array, localDeToDeMap=localDeToDeMap, rc=localrc)
  if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error getting local DE list from array",  &
          rcToReturn=rc)) return

  allocate(larrayList(localDeCount), stat=allocRcToTest)
  if (ESMF_LogFoundAllocError(allocRcToTest, msg="type "//                      &
     " larrayList in populate_redist_array", rcToReturn=rc)) then
  endif
  call ESMF_ArrayGet(array, localarrayList=larrayList, rc=localrc)
  if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error getting local array list",          &
          rcToReturn=rc)) return

  !-----------------------------------------------------------------------------
  ! get dimcount to allocate bound arrays
  !-----------------------------------------------------------------------------
  call ESMF_DistGridGet(DistGrid, dimCount=dimCount, rc=localrc)
  if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error getting dimCount from distGrid",    &
          rcToReturn=rc)) return

  allocate(UBnd(dimCount, localDeCount), stat=allocRcToTest )
  if (ESMF_LogFoundAllocError(allocRcToTest, msg="integer variable "//          &
     " UBnd in populate_redist_array", rcToReturn=rc)) then
  endif
  allocate(LBnd(dimCount, localDeCount), stat=allocRcToTest )
  if (ESMF_LogFoundAllocError(allocRcToTest, msg="integer variable "//          &
     " LBnd in populate_redist_array", rcToReturn=rc)) then
  endif

  call ESMF_ArrayGet(array, indexflag=indexflag,                               &
           exclusiveLBound=LBnd, exclusiveUBound=UBnd, rc=localrc)
  if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error getting exclusive bound range",     &
          rcToReturn=rc)) return

  !-----------------------------------------------------------------------------
  ! associate the fortran pointer with the array object and populate the array
  !-----------------------------------------------------------------------------

  !-----------------------------------------------------------------------------
  ! Memory Rank = Grid Rank, then there are no tensor dimensions
  !-----------------------------------------------------------------------------
  if( Memory%memRank ==  Memory%GridRank ) then

     select case(dimCount)
     case(1)
     !--------------------------------------------------------------------------
     ! rank = 1
     !--------------------------------------------------------------------------
        do de=1, localDeCount
           call ESMF_LocalArrayGet(larrayList(de), farrayPtr=farrayPtr1, &
                                   datacopyflag=ESMF_DATACOPY_REFERENCE, rc=localrc)
           if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error connecting pointer to " // &
                   "array list", rcToReturn=rc)) return

           do i1=LBnd(1,de), UBnd(1,de)
              farrayPtr1(i1) = localDeToDeMap(de) + 1000.0d0*i1
           enddo    !   i1
        enddo    ! de
     case(2)
     !--------------------------------------------------------------------------
     ! rank = 2
     !--------------------------------------------------------------------------
        do de=1, localDeCount
           call ESMF_LocalArrayGet(larrayList(de), farrayPtr=farrayPtr2, &
                                   datacopyflag=ESMF_DATACOPY_REFERENCE, rc=localrc)
           if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error connecting pointer to " // &
                   "array list", rcToReturn=rc)) return

           do i1=LBnd(1,de), UBnd(1,de)
              do i2=LBnd(2,de), UBnd(2,de)
                 farrayPtr2(i1,i2) = localDeToDeMap(de) + 1000.0d0*i1 + 0.001d0*i2
              enddo   !   i2
           enddo    !   i1
        enddo    ! de
     case(3)
     !--------------------------------------------------------------------------
     ! rank = 3
     !--------------------------------------------------------------------------
        do de=1, localDeCount
           call ESMF_LocalArrayGet(larrayList(de), farrayPtr=farrayPtr3, &
                                   datacopyflag=ESMF_DATACOPY_REFERENCE, rc=localrc)
           if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error connecting pointer to " // &
                   "array list", rcToReturn=rc)) return

           do i1=LBnd(1,de), UBnd(1,de)
              do i2=LBnd(2,de), UBnd(2,de)
                 do i3=LBnd(3,de), UBnd(3,de)
                    farrayPtr3(i1,i2,i3) = localDeToDeMap(de) + 1.0d4*i1 + 10.0d2*i2 &
                          + 1.0d-2*i3
                 enddo   !   i3
              enddo   !   i2
           enddo    !   i1
        enddo    ! de
     case(4)
     !--------------------------------------------------------------------------
     ! rank = 4
     !--------------------------------------------------------------------------
        do de=1, localDeCount
           call ESMF_LocalArrayGet(larrayList(de), farrayPtr=farrayPtr4, &
                                   datacopyflag=ESMF_DATACOPY_REFERENCE, rc=localrc)
           if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error connecting pointer to " // &
                   "array list", rcToReturn=rc)) return

           do i1=LBnd(1,de), UBnd(1,de)
              do i2=LBnd(2,de), UBnd(2,de)
                 do i3=LBnd(3,de), UBnd(3,de)
                    do i4=LBnd(4,de), UBnd(4,de)
                       farrayPtr4(i1,i2,i3,i4) = localDeToDeMap(de) + 1.0d4*i1         &
                             + 1.0d2*i2 + 1.0d-2*i3 + 1.0d-4*i4
                    enddo   !   i4
                 enddo   !   i3
              enddo   !   i2
           enddo    !   i1
        enddo    ! de
#ifndef ESMF_NO_GREATER_THAN_4D
     case(5)
     !--------------------------------------------------------------------------
     ! rank = 5
     !--------------------------------------------------------------------------
        do de=1, localDeCount
           call ESMF_LocalArrayGet(larrayList(de), farrayPtr=farrayPtr5, &
                                   datacopyflag=ESMF_DATACOPY_REFERENCE, rc=localrc)
           if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error connecting pointer to " // &
                   "array list", rcToReturn=rc)) return

           do i1=LBnd(1,de), UBnd(1,de)
              do i2=LBnd(2,de), UBnd(2,de)
                 do i3=LBnd(3,de), UBnd(3,de)
                    do i4=LBnd(4,de), UBnd(4,de)
                       do i5=LBnd(5,de), UBnd(5,de)
                          farrayPtr5(i1,i2,i3,i4,i5) = localDeToDeMap(de) + 1.0d4*i1   &
                             + 1.0d2*i2 + 1.0d0*i3 + 1.0d-2*i4  + 1.0d-4*i5
                       enddo   !   i5
                    enddo   !   i4
                 enddo   !   i3
              enddo   !   i2
           enddo    !   i1
        enddo    ! de
     case(6)
     !--------------------------------------------------------------------------
     ! rank = 6
     !--------------------------------------------------------------------------
        do de=1, localDeCount
           call ESMF_LocalArrayGet(larrayList(de), farrayPtr=farrayPtr6, &
                                   datacopyflag=ESMF_DATACOPY_REFERENCE, rc=localrc)
           if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error connecting pointer to " // &
                   "array list", rcToReturn=rc)) return

           do i1=LBnd(1,de), UBnd(1,de)
              do i2=LBnd(2,de), UBnd(2,de)
                 do i3=LBnd(3,de), UBnd(3,de)
                    do i4=LBnd(4,de), UBnd(4,de)
                       do i5=LBnd(5,de), UBnd(5,de)
                       do i6=LBnd(6,de), UBnd(6,de)
                       farrayPtr6(i1,i2,i3,i4,i5,i6) = localDeToDeMap(de) +            &
                             1.0d5*i1 + 1.0d3*i2 + 1.0d1*i3 + 1.0d-1*i4        &
                             + 1.0d-3*i5 + 1.0d-5*i6
                       enddo   !   i6
                       enddo   !   i5
                    enddo   !   i4
                 enddo   !   i3
              enddo   !   i2
           enddo    !   i1
        enddo    ! de
     case(7)
     !--------------------------------------------------------------------------
     ! rank = 7
     !--------------------------------------------------------------------------
        do de=1, localDeCount
           call ESMF_LocalArrayGet(larrayList(de), farrayPtr=farrayPtr7, &
                                   datacopyflag=ESMF_DATACOPY_REFERENCE, rc=localrc)
           if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error connecting pointer to " // &
                   "array list", rcToReturn=rc)) return

           do i1=LBnd(1,de), UBnd(1,de)
              do i2=LBnd(2,de), UBnd(2,de)
                 do i3=LBnd(3,de), UBnd(3,de)
                    do i4=LBnd(4,de), UBnd(4,de)
                       do i5=LBnd(5,de), UBnd(5,de)
                       do i6=LBnd(6,de), UBnd(6,de)
                       do i7=LBnd(7,de), UBnd(7,de)
                          farrayPtr7(i1,i2,i3,i4,i5,i6,i7) = localDeToDeMap(de) +      &
                             1.0d5*i1 + 1.0d3*i2 + 1.0d1*i3 + 1.0d-1*i4        &
                             + 1.0d-3*i5 + 1.0d-5*i6 + 1.0d0*i7
                       enddo   !   i7
                       enddo   !   i6
                       enddo   !   i5
                    enddo   !   i4
                 enddo   !   i3
              enddo   !   i2
           enddo    !   i1
        enddo    ! de
#endif
     case default
     !--------------------------------------------------------------------------
     ! error
     !--------------------------------------------------------------------------
        localrc = ESMF_FAILURE
        call ESMF_LogSetError(ESMF_FAILURE, msg="DimCount out of range", &
                 rcToReturn=localrc)
        return
     end select

  !-----------------------------------------------------------------------------
  ! Memory Rank > Grid Rank, then there are MemRank-GridRank tensor dimensions
  !-----------------------------------------------------------------------------
  elseif( Memory%memRank >  Memory%GridRank ) then
! -----------
      localrc = ESMF_FAILURE
      call ESMF_LogSetError(ESMF_FAILURE, msg="tensor dimensions not supported",      &
                    rcToReturn=localrc)
      return

  endif

  !-----------------------------------------------------------------------------
  ! clean up allocated arrays
  !-----------------------------------------------------------------------------
  deallocate(localDeToDeMap)
  deallocate(LBnd, UBnd)
  deallocate(larrayList)

  !-----------------------------------------------------------------------------
  rc = ESMF_SUCCESS
  !-----------------------------------------------------------------------------

  !-----------------------------------------------------------------------------
  end subroutine populate_redist_array