subroutine field_regrid_test(PDS, test_failure, reportType, VM, rc)
!-----------------------------------------------------------------------------
! routine conducts the field regrid test
!-----------------------------------------------------------------------------
! 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_Grid) :: gridSrc
type(ESMF_Grid) :: gridDst
type(ESMF_Field) :: srcField
type(ESMF_Field) :: dstField
type(ESMF_ArraySpec) :: SrcArraySpec
type(ESMF_ArraySpec) :: DstArraySpec
type(ESMF_DistGrid) :: src_distgrid, dst_distgrid
type(ESMF_RouteHandle) :: routeHandle
! local integers
integer :: localrc ! local error status
integer :: iDfile, iGfile, iD, iG, ix
integer :: test_status
integer :: localPET
integer :: libflag
! local characters
character(THARN_MAXSTR) :: liG, liD
! debug
! real(ESMF_KIND_R8), pointer :: farrayPtr2(:,:)
! integer :: i1, i2, de, localDeCount, dimCount
! 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
! Create source and destination distributions, Fields and conduct regrid
!-----------------------------------------------------------------------------
if (localPet == 0) then
print*,'-----------------======field regrid 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)', "field regrid" &
// " 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 distribution
!---------------------------------------------------------------------
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 grid from from distribution
!---------------------------------------------------------------------
call create_grid_from_distgrid(gridSrc, src_distgrid, PDS%SrcMem, &
PDS%Gfiles(iGfile)%src_grid(iG), &
PDS%Dfiles(iDfile)%src_dist(iD), 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
!---------------------------------------------------------------------
! create array spec
!---------------------------------------------------------------------
!---------------------------------------------------------------------
! set the dimensionality of actual data storage to the memory size
! specified by the problem descriptor string
!---------------------------------------------------------------------
call ESMF_ArraySpecSet(SrcArraySpec, typekind=ESMF_TYPEKIND_R8, &
rank=PDS%SrcMem%memRank, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error creating ArraySpecSet", &
rcToReturn=rc)) return
!---------------------------------------------------------------------
! create source field from grid and arrayspec
!---------------------------------------------------------------------
srcField = ESMF_FieldCreate(grid=gridSrc, arrayspec=SrcArraySpec, &
staggerloc=ESMF_STAGGERLOC_CENTER, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error creating source field", &
rcToReturn=rc)) return
!---------------------------------------------------------------------
! populate src field with test function for regridding test
!---------------------------------------------------------------------
call populate_field(srcField, gridSrc, PDS%SrcMem, &
PDS%Gfiles(iGfile)%src_grid(iG), &
PDS%Gfiles(iGfile)%testfunction(iG), localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error initializing value in " // &
"source array ", rcToReturn=rc)) return
!-------------------------------------------------------------------------------
! Destination
!-------------------------------------------------------------------------------
!print*,'field regrid create destination'
!---------------------------------------------------------------------
! Create Destination distribution
!---------------------------------------------------------------------
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
!---------------------------------------------------------------------
! create destination grid from from distribution
!---------------------------------------------------------------------
call create_grid_from_distgrid(gridDst, dst_distgrid, PDS%DstMem, &
PDS%Gfiles(iGfile)%dst_grid(iG), &
PDS%Dfiles(iDfile)%dst_dist(iD), 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
!---------------------------------------------------------------------
! create array spec
!---------------------------------------------------------------------
!---------------------------------------------------------------------
! set the dimensionality of actual data storage to the memory size
! specified by the problem descriptor string
!---------------------------------------------------------------------
call ESMF_ArraySpecSet(DstArraySpec, typekind=ESMF_TYPEKIND_R8, &
rank=PDS%DstMem%memRank, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error creating dst ArraySpecSet", &
rcToReturn=rc)) return
!---------------------------------------------------------------------
! create source field from grid and arrayspec
!---------------------------------------------------------------------
dstField = ESMF_FieldCreate(gridDst, DstArraySpec, &
staggerloc=ESMF_STAGGERLOC_CENTER, name="dest", rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error creating dst field", &
rcToReturn=rc)) return
!-----------------------------------------------------------------------------
!-----------------------------------------------------------------------------
!!!!!! NEED TO FIX THIS WHEN GRID AND MESH RECOGNIZE UNITS !!!!!!!
!-----------------------------------------------------------------------------
!-----------------------------------------------------------------------------
! check that the grid units agree. If they don't log an error. Eventually the
! mesh class should be able to reconcile different units.
!-----------------------------------------------------------------------------
do ix=1, PDS%Gfiles(iGfile)%src_grid(iG)%grank
if( trim(PDS%Gfiles(iGfile)%src_grid(iG)%gunits(ix)%string) .ne. &
trim(PDS%Gfiles(iGfile)%dst_grid(iG)%gunits(ix)%string) ) then
print*,'ERROR: source and destination grid units do not agree'
print*,'Source units are: ', &
trim(PDS%Gfiles(iGfile)%src_grid(iG)%gunits(ix)%string), &
' while destination units are: ', &
trim(PDS%Gfiles(iGfile)%dst_grid(iG)%gunits(ix)%string)
call ESMF_LogSetError( ESMF_FAILURE, msg="Incompatible grid" // &
" units.",rcToReturn=localrc)
return
endif
enddo ! ix
!-----------------------------------------------------------------------------
! Now conduct the interpolation
!
! the test consists of interpolating an analytical test function evaluated
! on a source grid onto a destination grid, and checking that the result
! agrees with the analytical solution to within a set tolerance.
!-----------------------------------------------------------------------------
!print*,'field regrid store'
!---------------------------------------------------------------------
! select the correct regridding method and do a Field Regrid store
!---------------------------------------------------------------------
select case( PDS%process%tag )
case( Harness_BilinearRegrid )
call ESMF_FieldRegridStore(srcField, dstField=dstField, routeHandle=routeHandle, &
regridmethod=ESMF_REGRIDMETHOD_BILINEAR, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"Field Bilinear Regrid " // &
"store failed", rcToReturn=rc)) return
case( Harness_PatchRegrid )
libflag = 0 ! flag to catch if framework built w/ LAPACK and BLAS libs
#ifdef ESMF_LAPACK
libflag = 1 !
call ESMF_FieldRegridStore(srcField, dstField=dstField, routeHandle=routeHandle, &
regridmethod=ESMF_REGRIDMETHOD_PATCH, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"Field Patch Regrid " // &
"store failed", rcToReturn=rc)) return
#endif
if(libflag==0) call ESMF_LogSetError( ESMF_FAILURE, msg="Patch " // &
"regridding requires LAPACK & BLAS libraries. These" // &
" libraries appear not to have been set. The "// &
" environment variable ESMF_LAPACK must be set ON," // &
" and ESMF_LAPACK_LIBS set to the LAPACK and BLAS "// &
" library paths",rcToReturn=localrc)
case( Harness_ConservRegrid )
call ESMF_LogSetError( ESMF_FAILURE, msg="Conservative " // &
"regridding not currently supported",rcToReturn=localrc)
return
case( Harness_2ndConservRegrid )
call ESMF_LogSetError( ESMF_FAILURE, msg="Conservative " // &
"regridding not currently supported",rcToReturn=localrc)
return
case( Harness_NearNeighRegrid )
call ESMF_LogSetError( ESMF_FAILURE, msg="Nearest Neighbor " // &
"regridding not currently supported",rcToReturn=localrc)
return
case( Harness_Error )
! error - invalid regrid method
call ESMF_LogSetError( ESMF_FAILURE, msg="Invalid regrid " // &
"method.",rcToReturn=localrc)
case default
! error
call ESMF_LogSetError( ESMF_FAILURE, msg="Invalid regrid " // &
"method.",rcToReturn=localrc)
end select
!---------------------------------------------------------------------
! regrid run
!---------------------------------------------------------------------
!print*,'field regrid run'
call ESMF_FieldRegrid(srcField, dstField, routeHandle=routeHandle, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"Field Regrid run failed for " // &
" forward failed ", rcToReturn=rc)) return
!-----------------------------------------------------------------------------
! Check regrid
!-----------------------------------------------------------------------------
!print*,'check field regrid'
!---------------------------------------------------------------------
! compare interpolated array values with the exact solution
!---------------------------------------------------------------------
call check_field(test_status, dstField, gridDst, &
PDS%Gfiles(iGfile)%dst_grid(iG), &
PDS%Gfiles(iGfile)%testfunction(iG), localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"redistribution array " // &
" comparison failed ", rcToReturn=rc)) return
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!!!
!-----------------------------------------------------------------------------
!---------------------------------------------------------------------
! release handles
!---------------------------------------------------------------------
call ESMF_FieldRegridRelease(routeHandle, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"Regrid routehandle Release failed",&
rcToReturn=rc)) return
!---------------------------------------------------------------------
! release Fields
!---------------------------------------------------------------------
call ESMF_FieldDestroy(srcField, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"SRC Field Regrid Release failed", &
rcToReturn=rc)) return
call ESMF_FieldDestroy(dstField, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"DST Field Regrid Release failed", &
rcToReturn=rc)) return
!---------------------------------------------------------------------
! release Grids
!---------------------------------------------------------------------
call ESMF_GridDestroy(gridSrc, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"SRC grid Release failed", &
rcToReturn=rc)) return
call ESMF_GridDestroy(gridDst, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"DST Grid Release failed", &
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*,'Regrid Completed'
end if
!-----------------------------------------------------------------------------
end subroutine field_regrid_test