subroutine test_PatchMeshToLocStreamMask(rc)
integer, intent(out) :: rc
logical :: correct
integer :: localrc
type(ESMF_Mesh) :: srcMesh
type(ESMF_Field) :: srcField
type(ESMF_Field) :: dstField
type(ESMF_Array) :: dstArray
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 :: farrayPtr(:,:),farrayPtr2(:,:)
integer :: i1,i2,i3, index(2)
integer :: lDE, localDECount
integer :: cl,cu
real(ESMF_KIND_R8) :: x,y
integer :: localPet, petCount
integer, allocatable :: nodeIds(:),nodeOwners(:)
real(ESMF_KIND_R8), allocatable :: nodeCoords(:)
integer, allocatable :: elemIds(:),elemTypes(:),elemConn(:)
integer :: numNodes, numElems
integer :: numQuadElems,numTriElems, numTotElems
type(ESMF_LocStream) :: dstLocStream
integer :: numLocationsOnThisPet
real(ESMF_KIND_R8), pointer :: Xarray(:),Yarray(:)
integer(ESMF_KIND_I4) :: maskValues(2)
integer(ESMF_KIND_I4), pointer :: maskArray(:)
! 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
! Setup Src Mesh
if (petCount .eq. 1) then
! Set number of nodes
numNodes=9
! Allocate and fill the node id array.
allocate(nodeIds(numNodes))
nodeIds=(/100,20,30,40,50,60,70,80,90/)
! Allocate and fill node coordinate array.
! Since this is a 2D Mesh the size is 2x the
! number of nodes.
allocate(nodeCoords(2*numNodes))
nodeCoords=(/-0.1,-0.1, & ! node id 1
1.0,-0.1, & ! node id 2
2.1,-0.1, & ! node id 3
-0.1, 1.0, & ! node id 4
1.0, 1.0, & ! node id 5
2.1, 1.0, & ! node id 6
-0.1, 2.1, & ! node id 7
1.0, 2.1, & ! node id 8
2.1, 2.1 /) ! node id 9
! Allocate and fill the node owner array.
! Since this Mesh is all on PET 0, it's just set to all 0.
allocate(nodeOwners(numNodes))
nodeOwners=0 ! everything on PET 0
! Set the number of each type of element, plus the total number.
numQuadElems=3
numTriElems=2
numTotElems=numQuadElems+numTriElems
! Allocate and fill the element id array.
allocate(elemIds(numTotElems))
elemIds=(/1,2,3,4,5/)
! Allocate and fill the element topology type array.
allocate(elemTypes(numTotElems))
elemTypes=(/ESMF_MESHELEMTYPE_QUAD, & ! elem id 1
ESMF_MESHELEMTYPE_TRI, & ! elem id 2
ESMF_MESHELEMTYPE_TRI, & ! elem id 3
ESMF_MESHELEMTYPE_QUAD, & ! elem id 4
ESMF_MESHELEMTYPE_QUAD/) ! elem id 5
! Allocate and fill the element connection type array.
! Note that entries in this array refer to the
! positions in the nodeIds, etc. arrays and that
! the order and number of entries for each element
! reflects that given in the Mesh options
! section for the corresponding entry
! in the elemTypes array.
allocate(elemConn(4*numQuadElems+3*numTriElems))
elemConn=(/1,2,5,4, & ! elem id 1
2,3,5, & ! elem id 2
3,6,5, & ! elem id 3
4,5,8,7, & ! elem id 4
5,6,9,8/) ! elem id 5
else if (petCount .eq. 4) then
! Setup mesh data depending on PET
if (localPET .eq. 0) then !!! This part only for PET 0
! Set number of nodes
numNodes=4
! Allocate and fill the node id array.
allocate(nodeIds(numNodes))
nodeIds=(/100,20,40,50/)
! Allocate and fill node coordinate array.
! Since this is a 2D Mesh the size is 2x the
! number of nodes.
allocate(nodeCoords(2*numNodes))
nodeCoords=(/-0.1, -0.1, & ! node id 1
1.0, -0.1, & ! node id 2
-0.1, 1.0, & ! node id 4
1.0, 1.0 /) ! node id 5
! Allocate and fill the node owner array.
allocate(nodeOwners(numNodes))
nodeOwners=(/0, & ! node id 1
0, & ! node id 2
0, & ! node id 4
0/) ! node id 5
! Set the number of each type of element, plus the total number.
numQuadElems=1
numTriElems=0
numTotElems=numQuadElems+numTriElems
! Allocate and fill the element id array.
allocate(elemIds(numTotElems))
elemIds=(/1/)
! Allocate and fill the element topology type array.
allocate(elemTypes(numTotElems))
elemTypes=(/ESMF_MESHELEMTYPE_QUAD/) ! elem id 1
! Allocate and fill the element connection type array.
! Note that entry are local indices
allocate(elemConn(4*numQuadElems+3*numTriElems))
elemConn=(/1,2,4,3/) ! elem id 1
else if (localPET .eq. 1) then !!! This part only for PET 1
! Set number of nodes
numNodes=4
! Allocate and fill the node id array.
allocate(nodeIds(numNodes))
nodeIds=(/20,30,50,60/)
! Allocate and fill node coordinate array.
! Since this is a 2D Mesh the size is 2x the
! number of nodes.
allocate(nodeCoords(2*numNodes))
nodeCoords=(/1.0,-0.1, & ! node id 2
2.1,-0.1, & ! node id 3
1.0, 1.0, & ! node id 5
2.1, 1.0 /) ! node id 6
! Allocate and fill the node owner array.
allocate(nodeOwners(numNodes))
nodeOwners=(/0, & ! node id 2
1, & ! node id 3
0, & ! node id 5
1/) ! node id 6
! Set the number of each type of element, plus the total number.
numQuadElems=0
numTriElems=2
numTotElems=numQuadElems+numTriElems
! Allocate and fill the element id array.
allocate(elemIds(numTotElems))
elemIds=(/2,3/)
! Allocate and fill the element topology type array.
allocate(elemTypes(numTotElems))
elemTypes=(/ESMF_MESHELEMTYPE_TRI, & ! elem id 2
ESMF_MESHELEMTYPE_TRI/) ! elem id 3
! Allocate and fill the element connection type array.
allocate(elemConn(4*numQuadElems+3*numTriElems))
elemConn=(/1,2,3, & ! elem id 2
2,4,3/) ! elem id 3
else if (localPET .eq. 2) then !!! This part only for PET 2
! Set number of nodes
numNodes=4
! Allocate and fill the node id array.
allocate(nodeIds(numNodes))
nodeIds=(/40,50,70,80/)
! Allocate and fill node coordinate array.
! Since this is a 2D Mesh the size is 2x the
! number of nodes.
allocate(nodeCoords(2*numNodes))
nodeCoords=(/-0.1,1.0, & ! node id 4
1.0,1.0, & ! node id 5
-0.1,2.1, & ! node id 7
1.0,2.1 /) ! node id 8
! Allocate and fill the node owner array.
! Since this Mesh is all on PET 0, it's just set to all 0.
allocate(nodeOwners(numNodes))
nodeOwners=(/0, & ! node id 4
0, & ! node id 5
2, & ! node id 7
2/) ! node id 8
! Set the number of each type of element, plus the total number.
numQuadElems=1
numTriElems=0
numTotElems=numQuadElems+numTriElems
! Allocate and fill the element id array.
allocate(elemIds(numTotElems))
elemIds=(/4/)
! Allocate and fill the element topology type array.
allocate(elemTypes(numTotElems))
elemTypes=(/ESMF_MESHELEMTYPE_QUAD/) ! elem id 4
! Allocate and fill the element connection type array.
allocate(elemConn(4*numQuadElems+3*numTriElems))
elemConn=(/1,2,4,3/) ! elem id 4
else if (localPET .eq. 3) then !!! This part only for PET 3
! Set number of nodes
numNodes=4
! Allocate and fill the node id array.
allocate(nodeIds(numNodes))
nodeIds=(/50,60,80,90/)
! Allocate and fill node coordinate array.
! Since this is a 2D Mesh the size is 2x the
! number of nodes.
allocate(nodeCoords(2*numNodes))
nodeCoords=(/1.0,1.0, & ! node id 5
2.1,1.0, & ! node id 6
1.0,2.1, & ! node id 8
2.1,2.1 /) ! node id 9
! Allocate and fill the node owner array.
allocate(nodeOwners(numNodes))
nodeOwners=(/0, & ! node id 5
1, & ! node id 6
2, & ! node id 8
3/) ! node id 9
! Set the number of each type of element, plus the total number.
numQuadElems=1
numTriElems=0
numTotElems=numQuadElems+numTriElems
! Allocate and fill the element id array.
allocate(elemIds(numTotElems))
elemIds=(/5/)
! Allocate and fill the element topology type array.
allocate(elemTypes(numTotElems))
elemTypes=(/ESMF_MESHELEMTYPE_QUAD/) ! elem id 5
! Allocate and fill the element connection type array.
allocate(elemConn(4*numQuadElems+3*numTriElems))
elemConn=(/1,2,4,3/) ! elem id 5
endif
endif
! Create Mesh structure in 1 step
srcMesh=ESMF_MeshCreate(parametricDim=2,spatialDim=2, &
coordSys=ESMF_COORDSYS_CART, &
nodeIds=nodeIds, nodeCoords=nodeCoords, &
nodeOwners=nodeOwners, elementIds=elemIds,&
elementTypes=elemTypes, elementConn=elemConn, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Create source field
call ESMF_ArraySpecSet(arrayspec, 1, ESMF_TYPEKIND_R8, rc=rc)
srcField = ESMF_FieldCreate(srcMesh, arrayspec, &
name="source", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Load test data into the source Field
! Should only be 1 localDE
call ESMF_FieldGet(srcField, 0, farrayPtr1D, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! set interpolated function
i2=1
do i1=1,numNodes
if (nodeOwners(i1) .eq. localPet) then
! Get coordinates
x=nodeCoords(2*i1-1)
y=nodeCoords(2*i1)
! Set source function
farrayPtr1D(i2) = 20.0+x+y
! Advance to next owner
i2=i2+1
endif
enddo
! deallocate node data
deallocate(nodeIds)
deallocate(nodeCoords)
deallocate(nodeOwners)
! deallocate elem data
deallocate(elemIds)
deallocate(elemTypes)
deallocate(elemConn)
! Setup Dst LocStream
if (petCount .eq. 1) then
numLocationsOnThisPet=11
cl=1
cu=11
else
if (localpet .eq. 0) then
numLocationsOnThisPet=3
cl=1
cu=3
else if (localpet .eq. 1) then
numLocationsOnThisPet=3
cl=4
cu=6
else if (localpet .eq. 2) then
numLocationsOnThisPet=5
cl=7
cu=11
else if (localpet .eq. 3) then
numLocationsOnThisPet=0
cl=12
cu=11
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, &
indexflag=ESMF_INDEX_GLOBAL, &
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:Mask", &
KeyTypeKind=ESMF_TYPEKIND_I4, &
keyUnits="none", &
keyLongName="mask values", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble adding LocStream key for Mask'
rc=ESMF_FAILURE
return
endif
!-------------------------------------------------------------------
! Get key data.
!-------------------------------------------------------------------
call ESMF_LocStreamGetKey(dstLocStream, &
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, &
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, &
keyName="ESMF:Mask", &
farray=maskArray, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble getting LocStream key for Mask'
rc=ESMF_FAILURE
return
endif
!-------------------------------------------------------------------
! Set key data.
!-------------------------------------------------------------------
if (petCount .eq. 1) then
Xarray(1)=0.0
Xarray(2)=0.5
Xarray(3)=1.0
Xarray(4)=2.0
Xarray(5)=0.0
Xarray(6)=1.0
Xarray(7)=2.0
Xarray(8)=0.0
Xarray(9)=1.0
Xarray(10)=1.5
Xarray(11)=2.0
Yarray(1)=0.0
Yarray(2)=0.5
Yarray(3)=0.0
Yarray(4)=0.0
Yarray(5)=1.0
Yarray(6)=1.0
Yarray(7)=1.0
Yarray(8)=2.0
Yarray(9)=2.0
Yarray(10)=1.5
Yarray(11)=2.0
maskArray(1)=0
maskArray(2)=1
maskArray(3)=0
maskArray(4)=0
maskArray(5)=0
maskArray(6)=0
maskArray(7)=0
maskArray(8)=0
maskArray(9)=0
maskArray(10)=2
maskArray(11)=0
else
if (localpet .eq. 0) then
Xarray(1)=0.0
Xarray(2)=0.5
Xarray(3)=1.0
Yarray(1)=0.0
Yarray(2)=0.5
Yarray(3)=0.0
maskArray(1)=0
maskArray(2)=1
maskArray(3)=0
else if (localpet .eq. 1) then
Xarray(4)=2.0
Xarray(5)=0.0
Xarray(6)=1.0
Yarray(4)=0.0
Yarray(5)=1.0
Yarray(6)=1.0
maskArray(4)=0
maskArray(5)=0
maskArray(6)=0
else if (localpet .eq. 2) then
Xarray(7)=2.0
Xarray(8)=0.0
Xarray(9)=1.0
Xarray(10)=1.5
Xarray(11)=2.0
Yarray(7)=1.0
Yarray(8)=2.0
Yarray(9)=2.0
Yarray(10)=1.5
Yarray(11)=2.0
maskArray(7)=0
maskArray(8)=0
maskArray(9)=0
maskArray(10)=2
maskArray(11)=0
! else if (localpet .eq. 3) then
endif
endif
! Create dest field
call ESMF_ArraySpecSet(arrayspec, 1, ESMF_TYPEKIND_R8, rc=rc)
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 Field
! 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 mesh to the locstream
! Regrid store
call ESMF_FieldRegridStore( &
srcField, &
dstField=dstField, dstMaskValues=(/1,2/), &
routeHandle=routeHandle, &
regridmethod=ESMF_REGRIDMETHOD_PATCH, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=localrc
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=cl,cu
if (maskArray(i1) .gt. 0) then
if ( abs( farrayPtr1D(i1) ) > 0.0001) then
correct=.false.
endif
else
! Get coordinates
x=Xarray(i1)
y=Yarray(i1)
if ( abs( farrayPtr1D(i1)-(x+y+20.0) ) > 0.0001) then
print*,'ERROR: expecting ',x+y+20.0,' got ',abs(farrayPtr1d(i1))
correct=.false.
endif
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_MeshDestroy(srcMesh, 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_PatchMeshToLocStreamMask