subroutine test_regridGridToLocStream3d(rc)
integer, intent(out) :: rc
logical :: correct
integer :: localrc
type(ESMF_Grid) :: srcGrid
type(ESMF_Field) :: srcField
type(ESMF_Field) :: dstField
type(ESMF_Array) :: dstArray
type(ESMF_Array) :: srcArrayA
type(ESMF_RouteHandle) :: routeHandle
type(ESMF_ArraySpec) :: arrayspec
type(ESMF_VM) :: vm
real(ESMF_KIND_R8), pointer :: farrayPtrXC(:,:,:), farrayPtr1D(:)
real(ESMF_KIND_R8), pointer :: farrayPtrYC(:,:,:)
real(ESMF_KIND_R8), pointer :: farrayPtrZC(:,:,:)
real(ESMF_KIND_R8), pointer :: farrayPtr(:,:,:)
integer :: clbnd(3),cubnd(3)
integer :: fclbnd(2),fcubnd(2)
integer :: i1,i2,i3, index(2)
integer :: lDE, localDECount
real(ESMF_KIND_R8) :: coord(2)
integer src_nx, src_ny, src_nz
integer num_arrays
real(ESMF_KIND_R8) :: dx,dy
real(ESMF_KIND_R8) :: src_minx,src_miny,src_minz
real(ESMF_KIND_R8) :: src_maxx,src_maxy,src_maxz
real(ESMF_KIND_R8) :: x,y,z,expected
integer :: localPet, petCount
! result code
type(ESMF_LocStream) :: dstLocStream
integer :: numLocationsOnThisPet,i
real(ESMF_KIND_R8), pointer :: Xarray(:),Yarray(:),Zarray(:)
! init success flag
correct=.true.
rc=ESMF_SUCCESS
! get pet info
call ESMF_VMGetGlobal(vm, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_VMGet(vm, petCount=petCount, localPet=localpet, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! If we don't have 1 or 4 PETS then exit successfully
if ((petCount .ne. 1) .and. (petCount .ne. 4)) then
print*,'ERROR: test must be run using exactly 1 or 4 PETS - detected ',petCount
rc=ESMF_FAILURE
return
endif
! Establish the resolution of the grids
src_nx = 10
src_ny = 10
src_nz = 10
! Establish the coordinates of the grids
src_minx = -0.1
src_miny = -0.1
src_minz = -0.1
src_maxx = 2.1
src_maxy = 2.1
src_maxz = 2.1
! setup src grid
srcGrid=ESMF_GridCreateNoPeriDim(minIndex=(/1,1,1/),maxIndex=(/src_nx,src_ny,src_nz/), &
regDecomp=(/petCount,1,1/), &
coordSys=ESMF_COORDSYS_CART, indexflag=ESMF_INDEX_GLOBAL, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Create source fields
call ESMF_ArraySpecSet(arrayspec, 3, ESMF_TYPEKIND_R8, rc=rc)
srcField = ESMF_FieldCreate(srcGrid, arrayspec, &
staggerloc=ESMF_STAGGERLOC_CENTER, name="source", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_GridAddCoord(srcGrid, staggerloc=ESMF_STAGGERLOC_CENTER, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Get number of local DEs
call ESMF_GridGet(srcGrid, localDECount=localDECount, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Source grid
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Get memory and set coords for dst
do lDE=0,localDECount-1
!! get coords
call ESMF_GridGetCoord(srcGrid, localDE=lDE, staggerLoc=ESMF_STAGGERLOC_CENTER, &
coordDim=1, &
computationalLBound=clbnd, computationalUBound=cubnd, &
farrayPtr=farrayPtrXC, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_GridGetCoord(srcGrid, localDE=lDE, staggerLoc=ESMF_STAGGERLOC_CENTER, &
coordDim=2, &
computationalLBound=clbnd, computationalUBound=cubnd, &
farrayPtr=farrayPtrYC, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_GridGetCoord(srcGrid, localDE=lDE, staggerLoc=ESMF_STAGGERLOC_CENTER, &
coordDim=3, &
computationalLBound=clbnd, computationalUBound=cubnd, &
farrayPtr=farrayPtrZC, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_FieldGet(srcField, lDE, farrayPtr, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
!! set coords
do i1=clbnd(1),cubnd(1)
do i2=clbnd(2),cubnd(2)
do i3=clbnd(3),cubnd(3)
! Set source coordinates
farrayPtrXC(i1,i2,i3) = ((src_maxx-src_minx)*REAL(i1-1)/REAL(src_nx-1))+src_minx
farrayPtrYC(i1,i2,i3) = ((src_maxy-src_miny)*REAL(i2-1)/REAL(src_ny-1))+src_miny
farrayPtrZC(i1,i2,i3) = ((src_maxz-src_minz)*REAL(i3-1)/REAL(src_nz-1))+src_minz
! initialize source field
farrayPtr(i1,i2,i3) = farrayPtrXC(i1,i2,i3)+farrayPtrYC(i1,i2,i3)+farrayPtrZC(i1,i2,i3)+20.0
enddo
enddo
enddo
enddo ! lDE
! Setup Dst LocStream
if (petCount .eq. 1) then
numLocationsOnThisPet=10
else
if (localpet .eq. 0) then
numLocationsOnThisPet=3
else if (localpet .eq. 1) then
numLocationsOnThisPet=3
else if (localpet .eq. 2) then
numLocationsOnThisPet=2
else if (localpet .eq. 3) then
numLocationsOnThisPet=2
endif
endif
!-------------------------------------------------------------------
! Create the LocStream: Allocate space for the LocStream object,
! define the number and distribution of the locations.
!-------------------------------------------------------------------
dstLocStream=ESMF_LocStreamCreate(name="Equatorial Measurements", &
localCount=numLocationsOnThisPet, &
coordSys=ESMF_COORDSYS_CART, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble creating locStream'
rc=ESMF_FAILURE
return
endif
!-------------------------------------------------------------------
! Add key data (internally allocating memory).
!-------------------------------------------------------------------
call ESMF_LocStreamAddKey(dstLocStream, &
keyName="ESMF:Y", &
KeyTypeKind=ESMF_TYPEKIND_R8, &
keyUnits="Units", &
keyLongName="Ydimension", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble adding LocStream key for Y'
rc=ESMF_FAILURE
return
endif
call ESMF_LocStreamAddKey(dstLocStream, &
keyName="ESMF:X", &
KeyTypeKind=ESMF_TYPEKIND_R8, &
keyUnits="Units", &
keyLongName="Xdimension", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble adding LocStream key for X'
rc=ESMF_FAILURE
return
endif
call ESMF_LocStreamAddKey(dstLocStream, &
keyName="ESMF:Z", &
KeyTypeKind=ESMF_TYPEKIND_R8, &
keyUnits="Units", &
keyLongName="Zdimension", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble adding LocStream key for Z'
rc=ESMF_FAILURE
return
endif
!-------------------------------------------------------------------
! Get key data.
!-------------------------------------------------------------------
call ESMF_LocStreamGetKey(dstLocStream, &
localDE=0, &
keyName="ESMF:Y", &
farray=Yarray, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble getting LocStream key for Y coordinate'
rc=ESMF_FAILURE
return
endif
call ESMF_LocStreamGetKey(dstLocStream, &
localDE=0, &
keyName="ESMF:X", &
farray=Xarray, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble getting LocStream key for X coordinate'
rc=ESMF_FAILURE
return
endif
call ESMF_LocStreamGetKey(dstLocStream, &
localDE=0, &
keyName="ESMF:Z", &
farray=Zarray, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble getting LocStream key for Z coordinate'
rc=ESMF_FAILURE
return
endif
!-------------------------------------------------------------------
! Set key data.
!-------------------------------------------------------------------
if (petCount .eq. 1) then
!test multiple points at same location
Xarray(1)=0.0
Xarray(2)=1.0
Xarray(3)=2.0
Xarray(4)=0.0
Xarray(5)=1.0
Xarray(6)=2.0
Xarray(7)=0.0
Xarray(8)=1.0
Xarray(9)=2.0
Xarray(10)=1.0
Yarray(1)=0.0
Yarray(2)=0.0
Yarray(3)=0.0
Yarray(4)=1.0
Yarray(5)=1.0
Yarray(6)=1.0
Yarray(7)=2.0
Yarray(8)=2.0
Yarray(9)=2.0
Yarray(10)=1.0
Zarray(1)=0.0
Zarray(2)=0.0
Zarray(3)=0.0
Zarray(4)=1.0
Zarray(5)=1.0
Zarray(6)=1.0
Zarray(7)=2.0
Zarray(8)=2.0
Zarray(9)=2.0
Zarray(10)=1.0
else
if (localpet .eq. 0) then
Xarray(1)=0.0
Xarray(2)=1.0
Xarray(3)=2.0
Yarray(1)=0.0
Yarray(2)=0.0
Yarray(3)=0.0
Zarray(1)=0.0
Zarray(2)=0.0
Zarray(3)=0.0
else if (localpet .eq. 1) then
Xarray(1)=0.0
Xarray(2)=1.0
Xarray(3)=2.0
Yarray(1)=1.0
Yarray(2)=1.0
Yarray(3)=1.0
Zarray(1)=1.0
Zarray(2)=1.0
Zarray(3)=1.0
else if (localpet .eq. 2) then
Xarray(1)=0.0
Xarray(2)=1.0
Yarray(1)=2.0
Yarray(2)=2.0
Zarray(1)=2.0
Zarray(2)=2.0
else if (localpet .eq. 3) then
!test multiple points at same location
Xarray(1)=2.0
Xarray(2)=1.0
Yarray(1)=2.0
Yarray(2)=1.0
Zarray(1)=2.0
Zarray(2)=1.0
endif
endif
! Create dest fields
call ESMF_ArraySpecSet(arrayspec, 1, ESMF_TYPEKIND_R8, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble calling ArraySpecSet'
rc=ESMF_FAILURE
return
endif
dstField = ESMF_FieldCreate(dstLocStream, arrayspec, &
name="dest", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble creating field on locStream'
rc=ESMF_FAILURE
return
endif
! clear destination Fields
! Should only be 1 localDE
call ESMF_FieldGet(dstField, 0, farrayPtr1D, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
farrayPtr1D=0.0
!!! Regrid forward from the Src grid to the LocStream
! Regrid store
call ESMF_FieldRegridStore( &
srcField, &
dstField=dstField, &
routeHandle=routeHandle, &
regridmethod=ESMF_REGRIDMETHOD_BILINEAR, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Do regrid
call ESMF_FieldRegrid(srcField, dstField, routeHandle, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_FieldRegridRelease(routeHandle, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Check destination field
! Should only be 1 localDE
call ESMF_FieldGet(dstField, 0, farrayPtr1D, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! loop through nodes and make sure interpolated values are reasonable
do i1=1,numLocationsOnThisPet
! Get coordinates
x=Xarray(i1)
y=Yarray(i1)
z=Zarray(i1)
expected = x+y+z+20.0
!! if error is too big report an error
if ( abs( farrayPtr1D(i1)-expected )/expected > 0.001) then
print*,'ERROR: larger than expected error, expected ',expected, &
' got ',farrayPtr1D(i1)
correct=.false.
endif
enddo
! Destroy the Fields
call ESMF_FieldDestroy(srcField, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_FieldDestroy(dstField, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Free the grids
call ESMF_LocStreamDestroy(dstLocStream, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_GridDestroy(srcGrid, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! return answer based on correct flag
if (correct) then
rc=ESMF_SUCCESS
else
rc=ESMF_FAILURE
endif
end subroutine test_regridGridToLocStream3d