subroutine array_redist_test(PDS, test_failure, reportType, VM, rc)
!-----------------------------------------------------------------------------
!
!-----------------------------------------------------------------------------
! arguments
type(problem_descriptor_strings), intent(inout) :: PDS
character(THARN_MAXSTR), intent(in ) :: reportType
type(ESMF_VM), intent(in ) :: VM
integer, intent(inout) :: test_failure
integer, intent( out) :: rc
! local parameters
real(ESMF_KIND_R8), parameter :: initvalue = 0.0
! local ESMF Types
type(ESMF_Array) :: src_array, return_array, dst_array
type(ESMF_DistGrid) :: src_distgrid, dst_distgrid
type(ESMF_RouteHandle) :: redistHandle_forward, redistHandle_reverse
! local integers
integer :: localrc ! local error status
integer :: iDfile, iGfile, iD, iG
integer :: test_status
integer :: localPET
! local characters
character(THARN_MAXSTR) :: liG, liD
! debug
! real(ESMF_KIND_R8), pointer :: farrayPtr2(:,:)
! integer ::de, localDeCount, dimCount
! integer :: i1, i2
! integer, allocatable :: localDeToDeMap(:)
! type(ESMF_LocalArray), allocatable :: larrayList(:)
! integer, allocatable :: LBnd(:,:), UBnd(:,:)
! type(ESMF_Index_Flag) :: indexflag
! initialize return flag
localrc = ESMF_RC_NOT_IMPL
rc = ESMF_RC_NOT_IMPL
! initialize test counter
test_failure = 0
call ESMF_VMGet(VM, localPet=localPET, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"can not get local pet ", &
rcToReturn=rc)) return
!-----------------------------------------------------------------------------
! for a single problem descriptor string, loop through each specifier file
! combination
!-----------------------------------------------------------------------------
if (localPet == 0) then
print*,'-----------------======array redist test==========-----------------------'
end if
do iDfile=1,PDS%nDfiles ! distribution specifier files
do iGfile=1,PDS%nGfiles ! grid specifier files
do iD=1, PDS%Dfiles(iDfile)%nDspecs ! entries in distribution specifier
do iG=1, PDS%Gfiles(iGfile)%nGspecs ! entries in grid specifier file
write(liG,"(i5)") iG
write(liD,"(i5)") iD
if (localPet == 0) then
print '(A)',"array redist test - " &
// " with string " // trim(adjustL(PDS%pds)) // &
" with entry " // trim(adjustL(liD)) // " of file " // &
trim(adjustL(PDS%Dfiles(iDfile)%filename)) &
// " and entry " // trim(adjustL(liG)) // " of file " // &
trim(adjustL(PDS%Gfiles(iGfile)%filename))
end if
!---------------------------------------------------------------------
! create source and destination distributions
!---------------------------------------------------------------------
call create_distribution(PDS%SrcMem, PDS%Dfiles(iDfile)%src_dist(iD),&
PDS%Gfiles(iGfile)%src_grid(iG), src_distgrid, VM, localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error creating source distgrid " &
// " with string " // trim(adjustL(PDS%pds)) // &
" with entry " // trim(adjustL(liD)) // " of file " // &
trim(adjustL(PDS%Dfiles(iDfile)%filename)) &
// " and entry " // trim(adjustL(liG)) // " of file " // &
trim(adjustL(PDS%Gfiles(iGfile)%filename)), &
rcToReturn=rc)) return
!---------------------------------------------------------------------
! create source array
!---------------------------------------------------------------------
call create_array(src_array, src_distgrid, PDS%SrcMem, &
PDS%Gfiles(iGfile)%src_grid(iG), localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error creating source array " &
// " with string " // trim(adjustL(PDS%pds)) // &
" with entry " // trim(adjustL(liD)) // " of file " // &
trim(adjustL(PDS%Dfiles(iDfile)%filename)) &
// " and entry " // trim(adjustL(liG)) // " of file " // &
trim(adjustL(PDS%Gfiles(iGfile)%filename)), &
rcToReturn=rc)) return
!---------------------------------------------------------------------
! populate source array for redistribution test
!---------------------------------------------------------------------
call populate_redist_array(src_array, src_distgrid, PDS%SrcMem, &
PDS%Gfiles(iGfile)%src_grid(iG), localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error populating source array ", &
rcToReturn=rc)) return
!---------------------------------------------------------------------
! create return array
!---------------------------------------------------------------------
call create_array(return_array, src_distgrid, PDS%SrcMem, &
PDS%Gfiles(iGfile)%src_grid(iG), localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error creating return array " &
// " with string " // trim(adjustL(PDS%pds)) // &
" with entry " // trim(adjustL(liD)) // " of file " // &
trim(adjustL(PDS%Dfiles(iDfile)%filename)) &
// " and entry " // trim(adjustL(liG)) // " of file " // &
trim(adjustL(PDS%Gfiles(iGfile)%filename)), &
rcToReturn=rc)) return
!---------------------------------------------------------------------
! populate return array with zeros for redistribution test
!---------------------------------------------------------------------
call populate_array_value(return_array, initvalue, src_distgrid, &
PDS%SrcMem, PDS%Gfiles(iGfile)%src_grid(iG), localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error initializing value in " // &
"return array ", rcToReturn=rc)) return
!---------------------------------------------------------------------
! Create Destination distribution and array
!---------------------------------------------------------------------
call create_distribution(PDS%DstMem, PDS%Dfiles(iDfile)%dst_dist(iD),&
PDS%Gfiles(iGfile)%dst_grid(iG), dst_distgrid, VM, localrc)
write(liG,"(i5)") iG
write(liD,"(i5)") iD
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error creating source distgrid " &
// " with string " // trim(adjustL(PDS%pds)) // &
" with entry " // trim(adjustL(liD)) // " of file " // &
trim(adjustL(PDS%Dfiles(iDfile)%filename)) &
// " and entry " // trim(adjustL(liG)) // " of file " // &
trim(adjustL(PDS%Gfiles(iGfile)%filename)), &
rcToReturn=rc)) return
call create_array(dst_array, dst_distgrid, PDS%DstMem, &
PDS%Gfiles(iGfile)%dst_grid(iG), localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error creating destinationarray " &
// " with string " // trim(adjustL(PDS%pds)) // &
" with entry " // trim(adjustL(liD)) // " of file " // &
trim(adjustL(PDS%Dfiles(iDfile)%filename)) &
// " and entry " // trim(adjustL(liG)) // " of file " // &
trim(adjustL(PDS%Gfiles(iGfile)%filename)), &
rcToReturn=rc)) return
!---------------------------------------------------------------------
! populate destination array with zeros for redistribution test
!---------------------------------------------------------------------
call populate_array_value(dst_array, initvalue, dst_distgrid, &
PDS%DstMem, PDS%Gfiles(iGfile)%dst_grid(iG), localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error initializing value in " // &
"destination array ", rcToReturn=rc)) return
!-----------------------------------------------------------------------------
! Now conduct the redistribution test
!
! the test consists of a forward redistribution from the source
! distribution to the destination distribution, and a second backward
! redistribution from the destination back to the source distribution.
!-----------------------------------------------------------------------------
!---------------------------------------------------------------------
! redistribution store for forward direction
!---------------------------------------------------------------------
call ESMF_ArrayRedistStore(srcArray=src_array, dstArray=dst_array, &
routehandle=redistHandle_forward, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"Array redist store failed for" // &
" forward direction", rcToReturn=rc)) return
!---------------------------------------------------------------------
! forward redistribution run
!---------------------------------------------------------------------
call ESMF_ArrayRedist(srcArray=src_array, dstArray=dst_array, &
routehandle=redistHandle_forward, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"Array redist run failed for " // &
" forward failed ", rcToReturn=rc)) return
!---------------------------------------------------------------------
! release redistribution handles
!---------------------------------------------------------------------
call ESMF_ArrayRedistRelease(routehandle=redistHandle_forward, &
rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"redistribution release for" // &
" forward failed ", rcToReturn=rc)) return
!---------------------------------------------------------------------
! redistribution store for reverse direction
!---------------------------------------------------------------------
call ESMF_ArrayRedistStore(srcArray=dst_array, dstArray=return_array,&
routehandle=redistHandle_reverse, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"Array redist store failed for" // &
" reverse direction", rcToReturn=rc)) return
!---------------------------------------------------------------------
! forward redistribution run
!---------------------------------------------------------------------
call ESMF_ArrayRedist(srcArray=dst_array, dstArray=return_array, &
routehandle=redistHandle_reverse, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"Array redist run failed for " // &
" reverse failed ", rcToReturn=rc)) return
!---------------------------------------------------------------------
! release redistribution handles
!---------------------------------------------------------------------
call ESMF_ArrayRedistRelease(routehandle=redistHandle_reverse, &
rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"redistribution release for" // &
" reverse failed ", rcToReturn=rc)) return
!-----------------------------------------------------------------------------
! Check redistribution
!-----------------------------------------------------------------------------
!---------------------------------------------------------------------
! compare source array values with the return array values
!---------------------------------------------------------------------
call compare_redist_array(test_status, &
src_array, return_array, src_distgrid, src_distgrid, &
PDS%SrcMem, PDS%Gfiles(iGfile)%src_grid(iG), localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"redistribution array " // &
" comparison failed ", rcToReturn=rc)) return
! gblock - do we want to return on the compare failure? How about internal ESMF error?
! suspect already HarnessTest_FAILURE reported on compare error, ESMF_xxx on ESMF error
PDS%test_record(iDfile,iGfile)%test_status(iD,iG) = test_status
if( test_status == HarnessTest_FAILURE ) then
test_failure = test_failure + 1
endif
call report_descriptor_string(PDS, iG, iD, iGfile, iDfile, &
reportType, localPET, localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"redistribution array " // &
" test report failed ", rcToReturn=rc)) return
!-----------------------------------------------------------------------------
! Clean up!!!!!!
!-----------------------------------------------------------------------------
!---------------------------------------------------------------------
! Destroy Array objects before moving to next test
!---------------------------------------------------------------------
call ESMF_ArrayDestroy(src_array, rc=localrc) ! original source
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"unable to destroy src_array", &
rcToReturn=rc)) return
call ESMF_ArrayDestroy(return_array, rc=localrc) ! return to source
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"unable to destroy return_array", &
rcToReturn=rc)) return
call ESMF_ArrayDestroy(dst_array, rc=localrc) ! redistribution
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"unable to destroy dst_array ", &
rcToReturn=rc)) return
!---------------------------------------------------------------------
! Destroy DistGrid objects before running next test
!---------------------------------------------------------------------
call ESMF_DistGridDestroy(src_distgrid, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"unable to destroy src_distgrid", &
rcToReturn=rc)) return
call ESMF_DistGridDestroy(dst_distgrid, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"unable to destroy src_distgrid", &
rcToReturn=rc)) return
!---------------------------------------------------------------------
enddo ! iG
enddo ! iD
enddo ! iGfile
enddo ! iDfile
!-----------------------------------------------------------------------------
! if I've gotten this far without an error, then the routine has succeeded.
!-----------------------------------------------------------------------------
rc = ESMF_SUCCESS
if (localPet == 0) then
print*,'Array Redist Completed'
end if
!-----------------------------------------------------------------------------
end subroutine array_redist_test