subroutine test_regrid3D_DTOS(rc)
integer, intent(out) :: rc
logical :: correct
integer :: localrc
type(ESMF_Grid) :: gridA
type(ESMF_Grid) :: gridB
type(ESMF_Field) :: srcFieldA
type(ESMF_Field) :: dstFieldA
type(ESMF_Field) :: fieldB
type(ESMF_Field) :: errorField
type(ESMF_Array) :: arrayB
type(ESMF_Array) :: errorArray
type(ESMF_Array) :: lonArrayA
type(ESMF_Array) :: srcArrayA, dstArrayA
type(ESMF_RouteHandle) :: routeHandle
type(ESMF_RouteHandle) :: routeHandle1
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 :: farrayPtrZC(:,:,:)
real(ESMF_KIND_R8), pointer :: farrayPtr(:,:,:),farrayPtr2(:,:,:),errorfarrayPtr(:,:,:)
integer :: clbnd(3),cubnd(3)
integer :: fclbnd(3),fcubnd(3)
integer :: i1,i2,i3, index(3)
integer :: lDE, localDECount
real(ESMF_KIND_R8) :: coord(3)
character(len=ESMF_MAXSTR) :: string
integer A_nx, A_ny, A_nz, B_nx, B_ny, B_nz
integer num_arrays
real(ESMF_KIND_R8) :: A_minx,A_miny,A_minz
real(ESMF_KIND_R8) :: A_maxx,A_maxy,A_maxz
real(ESMF_KIND_R8) :: B_minx,B_miny,B_minz
real(ESMF_KIND_R8) :: B_maxx,B_maxy,B_maxz
integer :: spherical_grid
integer, pointer :: larrayList(:)
integer :: localPet, petCount
! result code
integer :: finalrc
! 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
! Establish the resolution of the grids
! (Make the grids the same, because that's an easy way to test if nearest neighbor is correct)
! (Slightly different nearest neighbor is tested elsewhere)
B_nx = 10
B_ny = 10
B_nz = 10
A_nx = 10
A_ny = 10
A_nz = 10
! Establish the coordinates of the grids
B_minx = 0.0
B_miny = 0.0
B_minz = 0.0
B_maxx = 10.0
B_maxy = 10.0
B_maxz = 10.0
A_minx = 0.0
A_miny = 0.0
A_minz = 0.0
A_maxx = 10.0
A_maxy = 10.0
A_maxz = 10.0
! setup source grid
gridA=ESMF_GridCreateNoPeriDim(minIndex=(/1,1,1/),maxIndex=(/A_nx,A_ny,A_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
! setup dest. grid
gridB=ESMF_GridCreateNoPeriDim(minIndex=(/1,1,1/),maxIndex=(/B_nx,B_ny,B_nz/),regDecomp=(/1,1,petCount/), &
coordSys=ESMF_COORDSYS_CART, indexflag=ESMF_INDEX_GLOBAL, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Create source/destination fields
call ESMF_ArraySpecSet(arrayspec, 3, ESMF_TYPEKIND_R8, rc=rc)
srcFieldA = ESMF_FieldCreate(gridA, arrayspec, &
staggerloc=ESMF_STAGGERLOC_CENTER, name="source", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
dstFieldA = ESMF_FieldCreate(gridA, arrayspec, &
staggerloc=ESMF_STAGGERLOC_CENTER, name="source", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
errorField = ESMF_FieldCreate(gridA, arrayspec, &
staggerloc=ESMF_STAGGERLOC_CENTER, name="source", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
fieldB = ESMF_FieldCreate(gridB, arrayspec, &
staggerloc=ESMF_STAGGERLOC_CENTER, name="dest", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Allocate coordinates
call ESMF_GridAddCoord(gridA, staggerloc=ESMF_STAGGERLOC_CENTER, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_GridAddCoord(gridB, staggerloc=ESMF_STAGGERLOC_CENTER, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Get number of local DEs
call ESMF_GridGet(gridA, localDECount=localDECount, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Get arrays
! arrayB
call ESMF_FieldGet(fieldB, array=arrayB, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! srcArrayA
call ESMF_FieldGet(srcFieldA, array=srcArrayA, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! dstArrayA
call ESMF_FieldGet(dstFieldA, array=dstArrayA, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! errorArray
call ESMF_FieldGet(errorField, array=errorArray, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Write results to a mesh
num_arrays = 1
! Construct 3D Grid A
! (Get memory and set coords for src)
do lDE=0,localDECount-1
!! get coord 1
call ESMF_GridGetCoord(gridA, 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(gridA, 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(gridA, 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
! get src pointer
call ESMF_FieldGet(srcFieldA, lDE, farrayPtr, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! get src destination Field pointer
call ESMF_FieldGet(dstFieldA, lDE, farrayPtr2, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
!! set coords, interpolated function
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) = ((A_maxx-A_minx)*REAL(i1-1)/REAL(A_nx-1))+A_minx
farrayPtrYC(i1,i2,i3) = ((A_maxy-A_miny)*REAL(i2-1)/REAL(A_ny-1))+A_miny
farrayPtrZC(i1,i2,i3) = ((A_maxz-A_minz)*REAL(i3-1)/REAL(A_nz-1))+A_minz
! set src data
! (something smooth, that varies everywhere)
farrayPtr(i1,i2,i3) = farrayPtrXC(i1,i2,i3)+farrayPtrYC(i1,i2,i3)+farrayPtrZC(i1,i2,i3)+15.0
! initialize src destination field
farrayPtr2(i1,i2,i3)=0.0
enddo
enddo
enddo
enddo ! lDE
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Destination grid
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Get memory and set coords for dst
do lDE=0,localDECount-1
!! get coords
call ESMF_GridGetCoord(gridB, 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(gridB, 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(gridB, 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(fieldB, lDE, farrayPtr, computationalLBound=fclbnd, &
computationalUBound=fcubnd, 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) = ((B_maxx-B_minx)*REAL(i1-1)/REAL(B_nx-1))+B_minx
farrayPtrYC(i1,i2,i3) = ((B_maxy-B_miny)*REAL(i2-1)/REAL(B_ny-1))+B_miny
farrayPtrZC(i1,i2,i3) = ((B_maxz-B_minz)*REAL(i3-1)/REAL(B_nz-1))+B_minz
! initialize destination field
farrayPtr(i1,i2,i3)=0.0
enddo
enddo
enddo
enddo ! lDE
!!! Regrid forward from the A grid to the B grid
! Regrid store
call ESMF_FieldRegridStore(srcFieldA, dstField=fieldB, &
routeHandle=routeHandle, &
regridmethod=ESMF_REGRIDMETHOD_NEAREST_DTOS, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Do regrid
call ESMF_FieldRegrid(srcFieldA, fieldB, 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 back from the B grid to the A grid
! Regrid store
call ESMF_FieldRegridStore(fieldB, dstField=dstFieldA, &
routeHandle=routeHandle, &
regridmethod=ESMF_REGRIDMETHOD_NEAREST_DTOS, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Do regrid
call ESMF_FieldRegrid(fieldB, dstFieldA, 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 if the values are close
do lDE=0,localDECount-1
! get src Field
call ESMF_FieldGet(srcFieldA, lDE, farrayPtr, computationalLBound=clbnd, &
computationalUBound=cubnd, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! get src destination Field
call ESMF_FieldGet(dstFieldA, lDE, farrayPtr2, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! get src Field
call ESMF_FieldGet(errorField, lDE, errorfarrayPtr, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
!! check relative error
do i1=clbnd(1),cubnd(1)
do i2=clbnd(2),cubnd(2)
do i3=clbnd(3),cubnd(3)
if (farrayPtr(i1,i2,i3) .ne. 0.0) then
errorfarrayPtr(i1,i2,i3)=ABS((farrayPtr(i1,i2,i3) - farrayPtr2(i1,i2,i3))/farrayPtr(i1,i2,i3))
else
errorfarrayPtr(i1,i2,i3)=(farrayPtr(i1,i2,i3) - farrayPtr2(i1,i2,i3))
endif
if (ABS(errorfarrayPtr(i1,i2,i3)) .gt. 0.001) then
correct=.false.
endif
enddo
enddo
enddo
enddo ! lDE
! Uncomment these calls to see some actual regrid results
spherical_grid = 0
! call ESMF_MeshIO(vm, gridA, ESMF_STAGGERLOC_CENTER, &
! "srcmesh", srcArrayA, dstArrayA, errorArray, rc=localrc, &
! spherical=spherical_grid)
! call ESMF_MeshIO(vm, gridB, ESMF_STAGGERLOC_CENTER, &
! "dstmesh", arrayB, rc=localrc, &
! spherical=spherical_grid)
! Destroy the Fields
call ESMF_FieldDestroy(srcFieldA, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_FieldDestroy(dstFieldA, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_FieldDestroy(errorField, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_FieldDestroy(fieldB, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Free the grids
call ESMF_GridDestroy(gridA, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_GridDestroy(gridB, 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_regrid3D_DTOS