subroutine check_field(test_status, Field, Grid, Grid_info, TestFunction, rc)
!-----------------------------------------------------------------------------
! routine checks destination field against analytical solution evaluated at
! the destination grid points
!-----------------------------------------------------------------------------
! arguments
integer, intent(inout) :: test_status
type(ESMF_Field), intent(inout) :: Field
type(ESMF_Grid), intent(in ) :: Grid
type(grid_specification_record), intent(in ) :: Grid_info
type(test_function_record), intent(in ) :: TestFunction
integer, intent(inout) :: rc
! local ESMF types
! local integer variables
integer :: i,j,k
! integer :: l,m
integer, allocatable :: lbnd(:), ubnd(:)
integer :: localDECount, lDE, GridRank
integer :: localrc ! local error status
! local real variables
real(ESMF_KIND_R8) :: a, b, kx, ly, lenk, lenl, exact
real(ESMF_KIND_R8), pointer :: coordX2D(:,:), coordY2D(:,:)
real(ESMF_KIND_R8), pointer :: interp2D(:,:)
real(ESMF_KIND_R8), pointer :: coordX3D(:,:,:), coordY3D(:,:,:)
real(ESMF_KIND_R8), pointer :: coordZ3D(:,:,:)
real(ESMF_KIND_R8), pointer :: interp3D(:,:,:)
! initialize return flag
test_status = HarnessTest_SUCCESS
localrc = ESMF_RC_NOT_IMPL
rc = ESMF_RC_NOT_IMPL
!-----------------------------------------------------------------------------
! Get the number of local DEs from the Grid
!-----------------------------------------------------------------------------
call ESMF_GridGet(grid=Grid, localDECount=localDECount, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error getting local DE count from grid", &
rcToReturn=rc)) return
!-----------------------------------------------------------------------------
! Get the grid rank
!-----------------------------------------------------------------------------
call ESMF_GridGet(grid=Grid, dimCount=GridRank, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error getting Grid rank from grid", &
rcToReturn=rc)) return
!-----------------------------------------------------------------------------
! allocate arrays to hold the local array bounds
!-----------------------------------------------------------------------------
allocate( lbnd(GridRank) )
allocate( ubnd(GridRank) )
!-----------------------------------------------------------------------------
! check Field
!-----------------------------------------------------------------------------
do lDE=0,localDECount-1
select case(GridRank)
case(1)
!-------------------------------------------------------------------------
! grid rank = 1
!-------------------------------------------------------------------------
localrc = ESMF_FAILURE
call ESMF_LogSetError(ESMF_FAILURE, msg="Grid rank=1 not supported ", &
rcToReturn=localrc)
return
case(2)
!-------------------------------------------------------------------------
! grid rank = 2
!-------------------------------------------------------------------------
call ESMF_GridGetCoord(Grid, localDE=lDE, &
staggerLoc=ESMF_STAGGERLOC_CENTER, coordDim=1, &
computationalLBound=lbnd, computationalUBound=ubnd, &
farrayPtr=coordX2D, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error getting grid=1 coordinates", &
rcToReturn=rc)) return
call ESMF_GridGetCoord(Grid, localDE=lDE, &
staggerLoc=ESMF_STAGGERLOC_CENTER, coordDim=2, &
computationalLBound=lbnd, computationalUBound=ubnd, &
farrayPtr=coordY2D, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error getting grid=2 coordinates", &
rcToReturn=rc)) return
!-------------------------------------------------------------------------
! Get pointers to field array
!-------------------------------------------------------------------------
call ESMF_FieldGet(Field, lDE, interp2D, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error field get", rcToReturn=rc)) return
!-------------------------------------------------------------------------
! compute relative error
!-------------------------------------------------------------------------
do j=lbnd(2),ubnd(2)
do i=lbnd(1),ubnd(1)
!-------------------------------------------------------------------
select case( trim(TestFunction%string) )
case("CONSTANT")
if (.not. check_value (TestFunction%param(1), interp2D(i,j), TestFunction%param(2))) then
test_status = HarnessTest_FAILURE
print *,' fields disagree ',i,j,interp2D(i,j), &
TestFunction%param(1)
endif
case("COORDINATEX")
if (.not. check_value (TestFunction%param(1)*coordX2D(i,j), interp2D(i,j), &
TestFunction%param(2))) then
test_status = HarnessTest_FAILURE
print*,' fields disagree ',i,j,interp2D(i,j), &
TestFunction%param(1)*coordX2D(i,j)
endif
case("COORDINATEY")
if (.not. check_value (TestFunction%param(1)*coordY2D(i,j), interp2D(i,j), &
TestFunction%param(2))) then
test_status = HarnessTest_FAILURE
print*,' fields disagree ',i,j,interp2D(i,j), &
TestFunction%param(1)*coordY2D(i,j)
endif
case("SPHERICAL_HARMONIC")
! (a+b) + acos(k*2pi*x/Lx) + b*sin(l*2pi*y/Ly)
a = TestFunction%param(1)
kx = TestFunction%param(2)
b = TestFunction%param(3)
ly = TestFunction%param(4)
lenk = Grid_info%grange(1,2) - Grid_info%grange(1,1)
lenl = Grid_info%grange(2,2) - Grid_info%grange(2,1)
exact = abs(a+b) +a*cos(pi2*kx*coordX2D(i,j)/lenk) + &
b*sin(pi2*ly*coordY2D(i,j)/lenl)
if (.not. check_value (exact, interp2D(i,j), TestFunction%param(5))) then
test_status = HarnessTest_FAILURE
print*,' fields disagree ',i,j,interp2D(i,j), exact
endif
case("PEAK_VALLEY")
! (1.-x*y)*sin(k*pi*x/Lx)*cos(l*pi*y)+2.
test_status = HarnessTest_FAILURE
print *, 'test function PEAK_VALLEY not implemented'
case default
test_status = HarnessTest_FAILURE
print *, 'undefined test function ',trim(TestFunction%string)
end select
!-------------------------------------------------------------------
enddo ! i loop
enddo ! j loop
case(3)
!-------------------------------------------------------------------------
! grid rank = 3
!-------------------------------------------------------------------------
call ESMF_GridGetCoord(Grid, localDE=lDE, &
staggerLoc=ESMF_STAGGERLOC_CENTER, coordDim=1, &
computationalLBound=lbnd, computationalUBound=ubnd, &
farrayPtr=coordX3D, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error getting grid=1 coordinates", &
rcToReturn=rc)) return
call ESMF_GridGetCoord(Grid, localDE=lDE, &
staggerLoc=ESMF_STAGGERLOC_CENTER, coordDim=2, &
computationalLBound=lbnd, computationalUBound=ubnd, &
farrayPtr=coordY3D, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error getting grid=2 coordinates", &
rcToReturn=rc)) return
call ESMF_GridGetCoord(Grid, localDE=lDE, &
staggerLoc=ESMF_STAGGERLOC_CENTER, coordDim=3, &
computationalLBound=lbnd, computationalUBound=ubnd, &
farrayPtr=coordZ3D, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error getting grid=3 coordinates", &
rcToReturn=rc)) return
!-------------------------------------------------------------------------
! Get pointers to field array
!-------------------------------------------------------------------------
call ESMF_FieldGet(Field, lDE, interp3D, rc=localrc)
if (CheckError(checkpoint, __LINE__, __FILE__, localrc,"error field get", rcToReturn=rc)) return
!-------------------------------------------------------------------------
! check all values
!-------------------------------------------------------------------------
do k=lbnd(3),ubnd(3)
do j=lbnd(2),ubnd(2)
do i=lbnd(1),ubnd(1)
!-------------------------------------------------------------------
select case( trim(TestFunction%string) )
case("CONSTANT")
if (.not. check_value (TestFunction%param(1), interp3D(i,j,k), TestFunction%param(2))) then
test_status = HarnessTest_FAILURE
print*,' fields disagree ',i,j,k,interp3D(i,j,k), &
TestFunction%param(1)
endif
case("COORDINATEX")
if (.not. check_value (TestFunction%param(1)*coordX3D(i,j,k), interp3D(i,j,k), &
TestFunction%param(2))) then
test_status = HarnessTest_FAILURE
print*,' fields disagree ',i,j,k,interp3D(i,j,k), &
TestFunction%param(1)*coordX3D(i,j,k)
endif
case("COORDINATEY")
if (.not. check_value (TestFunction%param(1)*coordY3D(i,j,k), interp3D(i,j,k), &
TestFunction%param(2))) then
test_status = HarnessTest_FAILURE
print*,' fields disagree ',i,j,k,interp3D(i,j,k), &
TestFunction%param(1)*coordY3D(i,j,k)
endif
case("COORDINATEZ")
if (.not. check_value (TestFunction%param(1)*coordZ3D(i,j,k), interp3D(i,j,k), &
TestFunction%param(2))) then
test_status = HarnessTest_FAILURE
print*,' fields disagree ',i,j,k,interp3D(i,j,k), &
TestFunction%param(1)*coordZ3D(i,j,k)
endif
case default
test_status = HarnessTest_FAILURE
print *, 'undefined test function ',trim(TestFunction%string)
end select
!-------------------------------------------------------------------
enddo ! i loop
enddo ! j loop
enddo ! k loop
case(4)
!-------------------------------------------------------------------------
! grid rank = 4
!-------------------------------------------------------------------------
localrc = ESMF_FAILURE
call ESMF_LogSetError(ESMF_FAILURE, msg="Grid rank=4 not supported ", &
rcToReturn=localrc)
return
case(5)
!-------------------------------------------------------------------------
! grid rank = 5
!-------------------------------------------------------------------------
localrc = ESMF_FAILURE
call ESMF_LogSetError(ESMF_FAILURE, msg="Grid rank=5 not supported ", &
rcToReturn=localrc)
return
case(6)
!-------------------------------------------------------------------------
! grid rank = 6
!-------------------------------------------------------------------------
localrc = ESMF_FAILURE
call ESMF_LogSetError(ESMF_FAILURE, msg="Grid rank=6 not supported ", &
rcToReturn=localrc)
return
case(7)
!-------------------------------------------------------------------------
! grid rank = 7
!-------------------------------------------------------------------------
localrc = ESMF_FAILURE
call ESMF_LogSetError(ESMF_FAILURE, msg="Grid rank=7 not supported ", &
rcToReturn=localrc)
return
case default
!-------------------------------------------------------------------------
! error
!-------------------------------------------------------------------------
localrc = ESMF_FAILURE
call ESMF_LogSetError(ESMF_FAILURE, msg="Grid rank not between 1 & 7",&
rcToReturn=localrc)
return
end select
enddo ! lDE
!-----------------------------------------------------------------------------
!-----------------------------------------------------------------------------
! clean up
!-----------------------------------------------------------------------------
deallocate( lbnd, ubnd)
!-----------------------------------------------------------------------------
rc = ESMF_SUCCESS
!-----------------------------------------------------------------------------
!-----------------------------------------------------------------------------
end subroutine check_field