subroutine test_regridGridToLocStreamRegDist(rc)
integer, intent(out) :: rc
logical :: correct
integer :: localrc
type(ESMF_Grid) :: srcGrid
type(ESMF_LocStream) :: dstLocStream
type(ESMF_Field) :: srcField,dstField,dstFieldPatch
type(ESMF_RouteHandle) :: routeHandle,routeHandlePatch
type(ESMF_ArraySpec) :: arrayspec
type(ESMF_VM) :: vm
real(ESMF_KIND_R8), pointer :: farrayPtrXC(:,:)
real(ESMF_KIND_R8), pointer :: farrayPtrYC(:,:)
real(ESMF_KIND_R8), pointer :: farrayPtr(:,:), farrayPtr1D(:)
real(ESMF_KIND_R8), pointer :: latArray(:),lonArray(:)
real(ESMF_KIND_R8) :: src_dx, src_dy, dst_dx, dst_dy
real(ESMF_KIND_R8) :: RAD2DEG,DEG2RAD,theta,phi
real(ESMF_KIND_R8) :: lat,lon
real(ESMF_KIND_R8) :: x,y,z,expected
integer :: clbnd(2),cubnd(2)
integer :: fclbnd(2),fcubnd(2)
integer :: dclbnd(1),dcubnd(1)
integer :: i1,i2
integer :: lDE, localDECount
integer :: src_nx, src_ny, dst_nx, dst_ny
integer :: cl,cu,idx
integer :: localPet, petCount
real(ESMF_KIND_R8) :: beg_time, end_time
#if defined (ESMF_LAPACK)
logical, external :: LSAME
logical :: tf
#endif
! 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 with failure
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 = 100
src_ny = 50
src_dx = 360./src_nx
src_dy = 180./src_ny
DEG2RAD = 3.14159265/180.0
RAD2DEG = 1./DEG2RAD
! setup source grid
srcGrid=ESMF_GridCreate1PeriDim(minIndex=(/1,1/),maxIndex=(/src_nx,src_ny/),regDecomp=(/petCount,1/), &
coordSys=ESMF_COORDSYS_SPH_DEG, indexflag=ESMF_INDEX_GLOBAL, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Create source field
call ESMF_ArraySpecSet(arrayspec, 2, 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
! Allocate coordinates
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
! Get memory and set coords for src
do lDE=0,localDECount-1
!! get coord 1
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
! get src pointer
call ESMF_FieldGet(srcField, lDE, farrayPtr, computationalLBound=fclbnd, &
computationalUBound=fcubnd, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
if (clbnd(1) .ne. fclbnd(1)) print *, 'Error clbnd != fclbnd'
if (clbnd(2) .ne. fclbnd(2)) print *, 'Error clbnd != fclbnd'
if (cubnd(1) .ne. fcubnd(1)) print *, 'Error cubnd != fcubnd'
if (cubnd(2) .ne. fcubnd(2)) print *, 'Error cubnd != fcubnd'
!! set coords, interpolated function
do i1=clbnd(1),cubnd(1)
do i2=clbnd(2),cubnd(2)
! Set source coordinates as 0 to 360
farrayPtrXC(i1,i2) = REAL(i1-1)*src_dx
farrayPtrYC(i1,i2) = -90. + (REAL(i2-1)*src_dy + 0.5*src_dy)
lon = farrayPtrXC(i1,i2)
lat = farrayPtrYC(i1,i2)
! Set the source to be a function of the x,y,z coordinate
theta = DEG2RAD*(lon)
phi = DEG2RAD*(90.-lat)
x = cos(theta)*sin(phi)
y = sin(theta)*sin(phi)
z = cos(phi)
! set src data
! (something relatively smooth, that varies everywhere)
farrayPtr(i1,i2) = x+y+z+15.0
enddo
enddo
enddo ! lDE
! Setup Dst LocStream
dst_nx = 90
dst_ny = 40
dst_dx = 360./dst_nx
dst_dy = 180./dst_ny
dstLocStream=ESMF_LocStreamCreate(minIndex=1, maxIndex=dst_nx*dst_ny, regDecomp=petCount, &
indexflag=ESMF_INDEX_GLOBAL, &
coordSys=ESMF_COORDSYS_SPH_DEG, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble creating locStream'
rc=ESMF_FAILURE
return
endif
call ESMF_LocStreamGet(dstLocStream, localDECount=localDECount, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble with locStreamGet'
rc=ESMF_FAILURE
return
endif
!-------------------------------------------------------------------
! Add key data (internally allocating memory).
!-------------------------------------------------------------------
call ESMF_LocStreamAddKey(dstLocStream, &
keyName="ESMF:Lat", &
KeyTypeKind=ESMF_TYPEKIND_R8, &
keyUnits="degrees", &
keyLongName="Latitude", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble adding LocStream key for latitude'
rc=ESMF_FAILURE
return
endif
call ESMF_LocStreamAddKey(dstLocStream, &
keyName="ESMF:Lon", &
KeyTypeKind=ESMF_TYPEKIND_R8, &
keyUnits="degrees", &
keyLongName="Longitude", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble adding LocStream key for longitude'
rc=ESMF_FAILURE
return
endif
!-------------------------------------------------------------------
! Get key data.
!-------------------------------------------------------------------
do lDE=0,localDECount-1
call ESMF_LocStreamGetBounds(dstLocStream, localDE=lDE, &
computationalLBound=cl, computationalUBound=cu, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble with LocStreamGet'
rc=ESMF_FAILURE
return
endif
!-------------------------------------------------------------------
! Get key data.
!-------------------------------------------------------------------
call ESMF_LocStreamGetKey(dstLocStream, &
localDE=lDE, &
keyName="ESMF:Lat", &
farray=latArray, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble getting LocStream key for latitude'
rc=ESMF_FAILURE
return
endif
call ESMF_LocStreamGetKey(dstLocStream, &
localDE=lDE, &
keyName="ESMF:Lon", &
farray=lonArray, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble getting LocStream key for longitude'
rc=ESMF_FAILURE
return
endif
do idx=cl,cu
i1=(idx-1)/dst_nx + 1
i2=mod((idx-1),dst_nx) + 1
! Set source coordinates as 0 to 360
lonArray(idx) = REAL(i2-1)*dst_dx
latArray(idx) = -90. + (REAL(i1-1)*dst_dy + 0.5*dst_dy)
enddo
enddo
! 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
dstFieldPatch = ESMF_FieldCreate(dstLocStream, arrayspec, &
name="destPatch", 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
call ESMF_FieldGet(dstFieldPatch, 0, farrayPtr1D, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
farrayPtr1D=0.0
!!! Regrid forward from the 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
!!! Regrid forward from the Src grid to the LocStream - this time with PATCH
! Make sure LAPACK lib is mapped in before the timing block
#if defined (ESMF_LAPACK)
tf = LSAME ('a', 'A')
#endif
! Get start time
call ESMF_VMWtime(beg_time)
! Regrid store
call ESMF_FieldRegridStore( &
srcField, &
dstField=dstFieldPatch, &
routeHandle=routeHandlePatch, &
regridmethod=ESMF_REGRIDMETHOD_PATCH, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=localrc
return
endif
! Get end time
call ESMF_VMWtime(end_time)
write(*,*) localPet," Time to do Patch FieldRegridStore()=",end_time-beg_time
! Do regrid
call ESMF_FieldRegrid(srcField, dstFieldPatch, routeHandlePatch, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_FieldRegridRelease(routeHandlePatch, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Check results
do lDE=0,localDECount-1
call ESMF_FieldGet(dstField, lDE, farrayPtr1D, computationalLBound=dclbnd, &
computationalUBound=dcubnd, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Make sure everthing looks ok
do i1=dclbnd(1),dcubnd(1)
! Get coordinates
lon=lonArray(i1)
lat=latArray(i1)
! get the x,y,z coordinates
theta = DEG2RAD*(lon)
phi = DEG2RAD*(90.-lat)
x = cos(theta)*sin(phi)
y = sin(theta)*sin(phi)
z = cos(phi)
! determine validation data
expected = x+y+z+15.0
!! if error is too big report an error
if ( abs( farrayPtr1D(i1)-(expected) )/expected > 0.001) then
print*,'ERROR: larger than expected difference, expected ',expected, &
' got ',farrayPtr1D(i1),' diff= ',abs(farrayPtr1D(i1)-expected)
correct=.false.
endif
enddo
! now for patch
call ESMF_FieldGet(dstFieldPatch, lDE, farrayPtr1D, computationalLBound=dclbnd, &
computationalUBound=dcubnd, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Make sure everthing looks ok
do i1=dclbnd(1),dcubnd(1)
! Get coordinates
lon=lonArray(i1)
lat=latArray(i1)
! get the x,y,z coordinates
theta = DEG2RAD*(lon)
phi = DEG2RAD*(90.-lat)
x = cos(theta)*sin(phi)
y = sin(theta)*sin(phi)
z = cos(phi)
! determine validation data
expected = x+y+z+15.0
!! if error is too big report an error
if ( abs( farrayPtr1D(i1)-(expected) )/expected > 0.001) then
print*,'ERROR: larger than expected difference, expected ',expected, &
' got ',farrayPtr1D(i1),' diff= ',abs(farrayPtr1D(i1)-expected)
correct=.false.
endif
enddo
enddo ! lDE
! 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
call ESMF_FieldDestroy(dstFieldPatch, 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_regridGridToLocStreamRegDist