program ESMF_MeshOpUTest
!------------------------------------------------------------------------------
#include "ESMF.h"
!==============================================================================
!BOPI
! !PROGRAM: ESMF_MeshOpUTest - Unit tests for Field Create and Get methods
!
! !DESCRIPTION:
!
! The code in this file drives F90 Field Create and Get unit tests.
! The companion folder Fieldsrc contains the definitions for the
! Field methods.
!EOPI
!-----------------------------------------------------------------------------
! !USES:
use ESMF_TestMod ! test methods
use ESMF
implicit none
!------------------------------------------------------------------------------
! The following line turns the CVS identifier string into a printable variable.
character(*), parameter :: version = &
'$Id$'
! cumulative result: count failures; no failures equals "all pass"
integer :: result = 0
! individual test result code
integer :: rc = 1
! individual test failure message
character(ESMF_MAXSTR) :: failMsg
character(512) :: name
call ESMF_TestStart(ESMF_SRCLINE, rc=rc)
if (rc /= ESMF_SUCCESS) call ESMF_Finalize(endflag=ESMF_END_ABORT)
!------------------------------------------------------------------------
!NEX_UTest
! exercise difference operation
print *, 'Starting test3'
call test3(rc)
write(failMsg, *) ""
write(name, *) "exercise difference operation"
call ESMF_Test((rc.eq.ESMF_SUCCESS), name, failMsg, result, ESMF_SRCLINE)
call ESMF_TestEnd(ESMF_SRCLINE)
contains
#define ESMF_METHOD "ESMF_TESTS"
!------------------------------------------------------------------------
subroutine test3(rc)
integer, intent(out) :: rc
integer :: localrc
type(ESMF_Mesh) :: meshA, meshB, mesh
rc = ESMF_SUCCESS
localrc = ESMF_SUCCESS
meshA = make_mesh_sph(40,20,1.,1.,30.,-5.,rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
meshB = make_mesh_sph(20,10,0.5,1.,50.,0.,rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
mesh = ESMF_MeshCreate(meshA, meshB, ESMF_MESHOP_DIFFERENCE, &
areaThreshold=1.e-7_ESMF_KIND_R8, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
!call ESMF_MeshWrite(mesh, "meshdiff", rc=localrc)
!if (ESMF_LogFoundError(localrc, &
! ESMF_ERR_PASSTHRU, &
! ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_MeshDestroy(mesh, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
end subroutine test3
!------------------------------------------------------------------------
! Utility methods
!------------------------------------------------------------------------
function make_grid(atm_nx, atm_ny, atm_dx, atm_dy, atm_sx, atm_sy, field, rc)
! return value
type(ESMF_Grid) :: make_grid
! arguments
integer, intent(in) :: atm_nx, atm_ny
real(ESMF_KIND_R4), intent(in) :: atm_dx, atm_dy
real(ESMF_KIND_R4), intent(in) :: atm_sx, atm_sy
type(ESMF_Field), intent(out), optional :: field
integer, intent(out), optional :: rc
! local variables
integer :: localrc, i, j
real(ESMF_KIND_R8), pointer :: coordX(:), coordY(:)
real(ESMF_KIND_R8) :: startx, starty
make_grid = ESMF_GridCreateNoPeriDim(maxIndex=(/atm_nx, atm_ny/), &
coordSys=ESMF_COORDSYS_CART, &
indexflag=ESMF_INDEX_GLOBAL, &
coordDep1=(/1/), &
coordDep2=(/2/), &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridAddCoord(make_grid, staggerloc=ESMF_STAGGERLOC_CENTER, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridAddCoord(make_grid, staggerloc=ESMF_STAGGERLOC_CORNER, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! global indexing
! atm grid is not decomposed in the y direction
!startx = lpet*atm_nx/npet*atm_dx
startx = atm_sx
starty = atm_sy
! compute coord
! X center
!call ESMF_GridGetCoord(make_grid, localDE=0, staggerLoc=ESMF_STAGGERLOC_CENTER, &
! coordDim=1, farrayPtr=coordX, rc=localrc)
!if (ESMF_LogFoundError(localrc, &
! ESMF_ERR_PASSTHRU, &
! ESMF_CONTEXT, rcToReturn=rc)) return
!do i = lbound(coordX,1), ubound(coordX,1)
! coordX(i) = startx + atm_dx/2. + (i-1)*atm_dx
!enddo
!print *, 'coordX: ', coordX
! X corner
call ESMF_GridGetCoord(make_grid, localDE=0, staggerLoc=ESMF_STAGGERLOC_CORNER, &
coordDim=1, farrayPtr=coordX, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
do i = lbound(coordX,1), ubound(coordX,1)
coordX(i) = startx + (i-1)*atm_dx
enddo
!print *, 'startx: ', startx, lbound(coordX, 1), 'coordX: ', coordX
! Y center
!call ESMF_GridGetCoord(make_grid, localDE=0, staggerLoc=ESMF_STAGGERLOC_CENTER, &
! coordDim=2, farrayPtr=coordY, rc=localrc)
!if (ESMF_LogFoundError(localrc, &
! ESMF_ERR_PASSTHRU, &
! ESMF_CONTEXT, rcToReturn=rc)) return
!do i = lbound(coordY,1), ubound(coordY,1)
! coordY(i) = starty + atm_dy/2. + (i-1)*atm_dy
!enddo
! Y corner
call ESMF_GridGetCoord(make_grid, localDE=0, staggerLoc=ESMF_STAGGERLOC_CORNER, &
coordDim=2, farrayPtr=coordY, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
do i = lbound(coordY,1), ubound(coordY,1)
coordY(i) = starty + (i-1)*atm_dy
enddo
if(present(field)) then
field = ESMF_FieldCreate(make_grid, typekind=ESMF_TYPEKIND_R8, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
if(present(rc)) rc = ESMF_SUCCESS
end function make_grid
function make_grid_sph(atm_nx, atm_ny, atm_dx, atm_dy, atm_sx, atm_sy, area_adj, tag, scheme, rc)
! return value
type(ESMF_Grid) :: make_grid_sph
! arguments
integer, intent(in) :: atm_nx, atm_ny
real(ESMF_KIND_R4), intent(in) :: atm_dx, atm_dy
real(ESMF_KIND_R4), intent(in) :: atm_sx, atm_sy
real(ESMF_KIND_R4), intent(in), optional :: area_adj
character(len=*), intent(in), optional :: tag
integer, intent(in) , optional :: scheme
integer, intent(out), optional :: rc
! local variables
integer :: localrc, i, j
real(ESMF_KIND_R8), pointer :: coordX(:,:), coordY(:,:)
real(ESMF_KIND_R8), pointer :: f_area(:,:), f_area_m(:), o_area(:,:)
real(ESMF_KIND_R8) :: startx, starty
integer :: l_scheme
type(ESMF_Mesh) :: mesh
type(ESMF_Field) :: field
l_scheme = ESMF_REGRID_SCHEME_REGION3D
if(present(scheme)) l_scheme = scheme
if(l_scheme == ESMF_REGRID_SCHEME_FULL3D) then
make_grid_sph = ESMF_GridCreate1PeriDim(maxIndex=(/atm_nx, atm_ny/), &
indexflag=ESMF_INDEX_GLOBAL, &
gridEdgeLWidth=(/0,0/), gridEdgeUWidth=(/0,1/), &
!regDecomp=(/npet, 1/), &
rc=localrc)
else
make_grid_sph = ESMF_GridCreateNoPeriDim(maxIndex=(/atm_nx, atm_ny/), &
indexflag=ESMF_INDEX_GLOBAL, &
gridEdgeLWidth=(/0,0/), gridEdgeUWidth=(/1,1/), &
!regDecomp=(/npet, 1/), &
rc=localrc)
endif
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridAddCoord(make_grid_sph, staggerloc=ESMF_STAGGERLOC_CENTER, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridAddCoord(make_grid_sph, staggerloc=ESMF_STAGGERLOC_CORNER, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! global indexing
! atm grid is not decomposed in the y direction
!startx = lpet*atm_nx/npet*atm_dx
startx = atm_sx
starty = atm_sy
! compute coord
! X center
call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CENTER, &
coordDim=1, farrayPtr=coordX, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Y center
call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CENTER, &
coordDim=2, farrayPtr=coordY, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
do i = lbound(coordX,1), ubound(coordX,1)
do j = lbound(coordX, 2), ubound(coordX, 2)
coordX(i,j) = startx + atm_dx/2. + (i-1)*atm_dx
coordY(i,j) = starty + atm_dy/2. + (j-1)*atm_dy
enddo
enddo
!print *, 'startx: ', startx, lbound(coordX, 1), ubound(coordX, 1), 'coordX: ', coordX(:,1)
! X corner
call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CORNER, &
coordDim=1, farrayPtr=coordX, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Y corner
call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CORNER, &
coordDim=2, farrayPtr=coordY, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
do i = lbound(coordX,1), ubound(coordX,1)
do j = lbound(coordX, 2), ubound(coordX, 2)
coordX(i,j) = startx + (i-1)*atm_dx
coordY(i,j) = starty + (j-1)*atm_dy
enddo
enddo
if(present(area_adj)) then
! retrieve area
!mesh = ESMF_GridToMesh(make_grid_sph, &
! ESMF_STAGGERLOC_CORNER, 0, &
! regridConserve=ESMF_REGRID_CONSERVE_ON, rc=localrc)
!if (ESMF_LogFoundError(localrc, &
! ESMF_ERR_PASSTHRU, &
! ESMF_CONTEXT, rcToReturn=rc)) return
!allocate(f_area_m(mesh%NumOwnedElements))
!call ESMF_MeshGetElemArea(mesh, arealist=f_area_m, rc=localrc)
!if (ESMF_LogFoundError(localrc, &
! ESMF_ERR_PASSTHRU, &
! ESMF_CONTEXT, rcToReturn=rc)) return
!deallocate(f_area_m)
! find out original Grid cell area
field = ESMF_FieldCreate(make_grid_sph, typekind=ESMF_TYPEKIND_R8, &
staggerloc=ESMF_STAGGERLOC_CENTER, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_FieldRegridGetArea(field, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_FieldGet(field, farrayPtr=o_area, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! add area to Grid
call ESMF_GridAddItem(make_grid_sph, ESMF_GRIDITEM_AREA, &
staggerloc=ESMF_STAGGERLOC_CENTER, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridGetItem(make_grid_sph, ESMF_GRIDITEM_AREA, &
staggerloc=ESMF_STAGGERLOC_CENTER, farrayptr=f_area, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! adjust Grid area
f_area = area_adj*o_area
endif
if(present(rc)) rc = ESMF_SUCCESS
end function make_grid_sph
function make_mesh_sph(atm_nx, atm_ny, atm_dx, atm_dy, atm_sx, atm_sy, tag, scheme, rc)
! return value
type(ESMF_Mesh) :: make_mesh_sph
! arguments
integer, intent(in) :: atm_nx, atm_ny
real(ESMF_KIND_R4), intent(in) :: atm_dx, atm_dy
real(ESMF_KIND_R4), intent(in) :: atm_sx, atm_sy
character(len=*), intent(in), optional :: tag
integer, intent(in) , optional :: scheme
integer, intent(out), optional :: rc
! local variables
integer :: localrc, i, j
real(ESMF_KIND_R8), pointer :: coordX(:,:), coordY(:,:)
real(ESMF_KIND_R8) :: startx, starty
integer :: l_scheme
type(ESMF_Grid) :: make_grid_sph
l_scheme = ESMF_REGRID_SCHEME_REGION3D
if(present(scheme)) l_scheme = scheme
if(l_scheme == ESMF_REGRID_SCHEME_FULL3D) then
make_grid_sph = ESMF_GridCreate1PeriDim(maxIndex=(/atm_nx, atm_ny/), &
indexflag=ESMF_INDEX_GLOBAL, &
gridEdgeLWidth=(/0,0/), gridEdgeUWidth=(/0,1/), &
!regDecomp=(/npet, 1/), &
rc=localrc)
else
make_grid_sph = ESMF_GridCreateNoPeriDim(maxIndex=(/atm_nx, atm_ny/), &
indexflag=ESMF_INDEX_GLOBAL, &
gridEdgeLWidth=(/0,0/), gridEdgeUWidth=(/1,1/), &
!regDecomp=(/npet, 1/), &
rc=localrc)
endif
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridAddCoord(make_grid_sph, staggerloc=ESMF_STAGGERLOC_CENTER, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridAddCoord(make_grid_sph, staggerloc=ESMF_STAGGERLOC_CORNER, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! global indexing
! atm grid is not decomposed in the y direction
!startx = lpet*atm_nx/npet*atm_dx
startx = atm_sx
starty = atm_sy
! compute coord
! X center
call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CENTER, &
coordDim=1, farrayPtr=coordX, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Y center
call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CENTER, &
coordDim=2, farrayPtr=coordY, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
do i = lbound(coordX,1), ubound(coordX,1)
do j = lbound(coordX, 2), ubound(coordX, 2)
coordX(i,j) = startx + atm_dx/2. + (i-1)*atm_dx
coordY(i,j) = starty + atm_dy/2. + (j-1)*atm_dy
enddo
enddo
!print *, 'startx: ', startx, lbound(coordX, 1), ubound(coordX, 1), 'coordX: ', coordX(:,1)
! X corner
call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CORNER, &
coordDim=1, farrayPtr=coordX, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Y corner
call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CORNER, &
coordDim=2, farrayPtr=coordY, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
do i = lbound(coordX,1), ubound(coordX,1)
do j = lbound(coordX, 2), ubound(coordX, 2)
coordX(i,j) = startx + (i-1)*atm_dx
coordY(i,j) = starty + (j-1)*atm_dy
enddo
enddo
make_mesh_sph = ESMF_GridToMesh(make_grid_sph, &
ESMF_STAGGERLOC_CORNER, 0, &
regridConserve=ESMF_REGRID_CONSERVE_ON, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if(present(rc)) rc = ESMF_SUCCESS
end function make_mesh_sph
end program ESMF_MeshOpUTest