subroutine test_RegridSph4ConcaveMesh(itrp, csrv, rc)
logical, intent(out) :: itrp
logical, intent(out) :: csrv
integer, intent(out) :: rc
integer :: localrc
type(ESMF_Mesh) :: srcMesh
type(ESMF_Mesh) :: dstMesh
type(ESMF_Field) :: srcField
type(ESMF_Field) :: dstField
type(ESMF_Field) :: xdstField
type(ESMF_Field) :: srcAreaField, dstAreaField
type(ESMF_Field) :: srcFracField, dstFracField
type(ESMF_RouteHandle) :: routeHandle
type(ESMF_ArraySpec) :: arrayspec
type(ESMF_VM) :: vm
real(ESMF_KIND_R8), pointer :: srcFarrayPtr(:), dstFarrayPtr(:), xdstFarrayPtr(:)
real(ESMF_KIND_R8), pointer :: srcAreaPtr(:), dstAreaPtr(:)
real(ESMF_KIND_R8), pointer :: srcFracPtr(:), dstFracPtr(:)
integer :: clbnd(1),cubnd(1)
integer :: i1,i2,i3
real(ESMF_KIND_R8) :: x,y,z
integer :: localPet, petCount
real(ESMF_KIND_R8) :: srcmass(1), dstmass(1), srcmassg(1), dstmassg(1)
real(ESMF_KIND_R8) :: maxerror(1), minerror(1), error
real(ESMF_KIND_R8) :: maxerrorg(1), minerrorg(1), errorg
real(ESMF_KIND_R8) :: errorTot, errorTotG
integer, pointer :: nodeIds(:),nodeOwners(:)
real(ESMF_KIND_R8), pointer :: nodeCoords(:)
integer, pointer :: elemIds(:),elemTypes(:),elemConn(:),elemMask(:)
integer :: numNodes
integer :: iconn,inode
integer :: numQuadElems,numTriElems
integer :: numPentElems,numHexElems,numTotElems
integer :: numElemConn
! result code
integer :: finalrc
! Init to success
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
! XMRKX
!!!! Setup source mesh !!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
! Creates the following mesh on
! 1 or 4 PETs. Returns an error
! if run on other than 1 or 4 PETs
!
! Mesh Ids
!
! 20 7 ------- 8
! | |\
! 13 | 3 | 9
! | |4 \
! 10 4 ------- 5 -------- 6
! | | |
! | 1 | 2 |
! | | |
! 0.0 1 ------- 2 -------- 3
!
! 0.0 10 13 20
!
! Node Ids at corners
! Element Ids in centers
!
!!!!!
!
! The owners for 1 PET are all Pet 0.
! The owners for 4 PETs are as follows:
!
! Mesh Owners
!
! 20 2 ------- 2
! | |\
! 12 | 2 | 3
! | |3 \
! 10 0 ------- 0 -------- 1
! | | |
! | 0 | 1 |
! | | |
! 0.0 0 ------- 0 -------- 1
!
! 0.0 10 12 20
!
! Node Owners at corners
! Element Owners in centers
!
! Setup mesh info depending on the
! number of PETs
if (petCount .eq. 1) then
! Fill in node data
numNodes=9
!! node ids
allocate(nodeIds(numNodes))
nodeIds=(/1,2,3,4,5,6,7,8,9/)
!! node Coords
allocate(nodeCoords(numNodes*2))
nodeCoords=(/0.0,0.0, &
10.0,0.0, &
20.0,0.0, &
0.0,10.0, &
10.0,10.0, &
20.0,10.0, &
0.0,20.0, &
10.0,20.0, &
12.0,12.0 /)
!! node owners
allocate(nodeOwners(numNodes))
nodeOwners=0 ! everything on proc 0
! Fill in elem data
numTotElems=4
!! elem ids
allocate(elemIds(numTotElems))
elemIds=(/1,2,3,4/)
!! elem types
allocate(elemTypes(numTotElems))
elemTypes=ESMF_MESHELEMTYPE_QUAD
!! elem conn
allocate(elemConn(numTotElems*4))
elemConn=(/1,2,5,4, &
2,3,6,5, &
4,5,8,7, &
5,6,9,8/)
else if (petCount .eq. 4) then
! Setup mesh data depending on PET
if (localPet .eq. 0) then
! Fill in node data
numNodes=4
!! node ids
allocate(nodeIds(numNodes))
nodeIds=(/1,2,4,5/)
!! node Coords
allocate(nodeCoords(numNodes*2))
nodeCoords=(/0.0,0.0, &
10.0,0.0, &
0.0,10.0, &
10.0,10.0/)
!! node owners
allocate(nodeOwners(numNodes))
nodeOwners=(/0,0,0,0/) ! everything on proc 0
! Fill in elem data
numTotElems=1
!! elem ids
allocate(elemIds(numTotElems))
elemIds=(/1/)
!! elem type
allocate(elemTypes(numTotElems))
elemTypes=ESMF_MESHELEMTYPE_QUAD
!! elem conn
allocate(elemConn(numTotElems*4))
elemConn=(/1,2,4,3/)
else if (localPet .eq. 1) then
! Fill in node data
numNodes=4
!! node ids
allocate(nodeIds(numNodes))
nodeIds=(/2,3,5,6/)
!! node Coords
allocate(nodeCoords(numNodes*2))
nodeCoords=(/10.0,0.0, &
20.0,0.0, &
10.0,10.0, &
20.0,10.0/)
!! node owners
allocate(nodeOwners(numNodes))
nodeOwners=(/0,1,0,1/)
! Fill in elem data
numTotElems=1
!! elem ids
allocate(elemIds(numTotElems))
elemIds=(/2/)
!! elem type
allocate(elemTypes(numTotElems))
elemTypes=ESMF_MESHELEMTYPE_QUAD
!! elem conn
allocate(elemConn(numTotElems*4))
elemConn=(/1,2,4,3/)
else if (localPet .eq. 2) then
! Fill in node data
numNodes=4
!! node ids
allocate(nodeIds(numNodes))
nodeIds=(/4,5,7,8/)
!! node Coords
allocate(nodeCoords(numNodes*2))
nodeCoords=(/0.0,10.0, &
10.0,10.0, &
0.0,20.0, &
10.0,20.0/)
!! node owners
allocate(nodeOwners(numNodes))
nodeOwners=(/0,0,2,2/)
! Fill in elem data
numTotElems=1
!! elem ids
allocate(elemIds(numTotElems))
elemIds=(/3/)
!! elem type
allocate(elemTypes(numTotElems))
elemTypes=ESMF_MESHELEMTYPE_QUAD
!! elem conn
allocate(elemConn(numTotElems*4))
elemConn=(/1,2,4,3/)
else
! Fill in node data
numNodes=4
!! node ids
allocate(nodeIds(numNodes))
nodeIds=(/5,6,8,9/)
!! node Coords
allocate(nodeCoords(numNodes*2))
nodeCoords=(/10.0,10.0, &
20.0,10.0, &
10.0,20.0, &
12.0,12.0/)
!! node owners
allocate(nodeOwners(numNodes))
nodeOwners=(/0,1,2,3/)
! Fill in elem data
numTotElems=1
!! elem ids
allocate(elemIds(numTotElems))
elemIds=(/4/)
!! elem type
allocate(elemTypes(numTotElems))
elemTypes=ESMF_MESHELEMTYPE_QUAD
!! elem conn
allocate(elemConn(numTotElems*4))
elemConn=(/1,2,4,3/)
endif
endif
! Create Mesh structure in 1 step
srcMesh=ESMF_MeshCreate(parametricDim=2,spatialDim=2, &
coordSys=ESMF_COORDSYS_SPH_DEG, &
nodeIds=nodeIds, nodeCoords=nodeCoords, &
nodeOwners=nodeOwners, elementIds=elemIds,&
elementTypes=elemTypes, elementConn=elemConn, &
rc=rc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Array spec for fields
call ESMF_ArraySpecSet(arrayspec, 1, ESMF_TYPEKIND_R8, rc=rc)
! Create source field
srcField = ESMF_FieldCreate(srcMesh, arrayspec, meshloc=ESMF_MESHLOC_ELEMENT, &
name="source", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Create source area field
srcAreaField = ESMF_FieldCreate(srcMesh, arrayspec, meshloc=ESMF_MESHLOC_ELEMENT, &
name="source_area", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Create source frac field
srcFracField = ESMF_FieldCreate(srcMesh, arrayspec, meshloc=ESMF_MESHLOC_ELEMENT, &
name="source_frac", 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, srcFarrayPtr, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! set interpolated function
iconn=1
do i1=1,numTotElems
! Loop through nodes in elem
! to compute point in center
x=0.0
y=0.0
do i2=1,elemTypes(i1)
inode=2*(elemConn(iconn)-1)
x=x+nodeCoords(inode+1)
y=y+nodeCoords(inode+2)
iconn=iconn+1
enddo
x=x*(1.0/REAL(elemTypes(i1),ESMF_KIND_R8))
y=y*(1.0/REAL(elemTypes(i1),ESMF_KIND_R8))
! Set source function
srcFarrayPtr(i1) = 20.0+x+y
enddo
! deallocate node data
deallocate(nodeIds)
deallocate(nodeCoords)
deallocate(nodeOwners)
! deallocate elem data
deallocate(elemIds)
deallocate(elemTypes)
deallocate(elemConn)
!!!!!!!!!!!!!!! Setup Destination Mesh !!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
! Creates the following mesh on
! 1 or 4 PETs. Returns an error
! if run on other than 1 or 4 PETs
!
! Mesh Ids
!
! 20.0 7 ------- 8
! | |\
! 14.0 | 3 | 9
! | |4 \
! 10.0 4 ------- 5 -------- 6
! | | |
! | 1 | 2 |
! | | |
! 0.0 1 ------- 2 -------- 3
!
! 0.0 10.0 14.0 20.0
!
! Node Ids at corners
! Element Ids in centers
!
!!!!!
!
! The owners for 1 PET are all Pet 0.
! The owners for 4 PETs are as follows:
!
! Mesh Owners
!
! 20.0 2 ------- 2
! | |\
! 14.0 | 2 | 3
! | |3 \
! 10.0 0 ------- 0 -------- 1
! | | |
! | 0 | 1 |
! | | |
! 0.0 0 ------- 0 -------- 1
!
! 0.0 10.0 14.0 20.0
!
! Node Owners at corners
! Element Owners in centers
!
! Setup mesh info depending on the
! number of PETs
if (petCount .eq. 1) then
! Fill in node data
numNodes=9
!! node ids
allocate(nodeIds(numNodes))
nodeIds=(/1,2,3,4,5,6,7,8,9/)
!! node Coords
allocate(nodeCoords(numNodes*2))
nodeCoords=(/0.0,0.0, &
10.0,0.0, &
20.0,0.0, &
0.0,10.0, &
10.0,10.0, &
20.0,10.0, &
0.0,20.0, &
10.0,20.0, &
14.0,14.0 /)
!! node owners
allocate(nodeOwners(numNodes))
nodeOwners=0 ! everything on proc 0
! Fill in elem data
numTotElems=4
!! elem ids
allocate(elemIds(numTotElems))
elemIds=(/1,2,3,4/)
!! elem types
allocate(elemTypes(numTotElems))
elemTypes=ESMF_MESHELEMTYPE_QUAD
!! elem conn
allocate(elemConn(numTotElems*4))
elemConn=(/1,2,5,4, &
2,3,6,5, &
4,5,8,7, &
5,6,9,8/)
else if (petCount .eq. 4) then
! Setup mesh data depending on PET
if (localPet .eq. 0) then
! Fill in node data
numNodes=4
!! node ids
allocate(nodeIds(numNodes))
nodeIds=(/1,2,4,5/)
!! node Coords
allocate(nodeCoords(numNodes*2))
nodeCoords=(/0.0,0.0, &
10.0,0.0, &
0.0,10.0, &
10.0,10.0/)
!! node owners
allocate(nodeOwners(numNodes))
nodeOwners=(/0,0,0,0/) ! everything on proc 0
! Fill in elem data
numTotElems=1
!! elem ids
allocate(elemIds(numTotElems))
elemIds=(/1/)
!! elem type
allocate(elemTypes(numTotElems))
elemTypes=ESMF_MESHELEMTYPE_QUAD
!! elem conn
allocate(elemConn(numTotElems*4))
elemConn=(/1,2,4,3/)
else if (localPet .eq. 1) then
! Fill in node data
numNodes=4
!! node ids
allocate(nodeIds(numNodes))
nodeIds=(/2,3,5,6/)
!! node Coords
allocate(nodeCoords(numNodes*2))
nodeCoords=(/10.0,0.0, &
20.0,0.0, &
10.0,10.0, &
20.0,10.0/)
!! node owners
allocate(nodeOwners(numNodes))
nodeOwners=(/0,1,0,1/)
! Fill in elem data
numTotElems=1
!! elem ids
allocate(elemIds(numTotElems))
elemIds=(/2/)
!! elem type
allocate(elemTypes(numTotElems))
elemTypes=ESMF_MESHELEMTYPE_QUAD
!! elem conn
allocate(elemConn(numTotElems*4))
elemConn=(/1,2,4,3/)
else if (localPet .eq. 2) then
! Fill in node data
numNodes=4
!! node ids
allocate(nodeIds(numNodes))
nodeIds=(/4,5,7,8/)
!! node Coords
allocate(nodeCoords(numNodes*2))
nodeCoords=(/0.0,10.0, &
10.0,10.0, &
0.0,20.0, &
10.0,20.0/)
!! node owners
allocate(nodeOwners(numNodes))
nodeOwners=(/0,0,2,2/)
! Fill in elem data
numTotElems=1
!! elem ids
allocate(elemIds(numTotElems))
elemIds=(/3/)
!! elem type
allocate(elemTypes(numTotElems))
elemTypes=ESMF_MESHELEMTYPE_QUAD
!! elem conn
allocate(elemConn(numTotElems*4))
elemConn=(/1,2,4,3/)
else
! Fill in node data
numNodes=4
!! node ids
allocate(nodeIds(numNodes))
nodeIds=(/5,6,8,9/)
!! node Coords
allocate(nodeCoords(numNodes*2))
nodeCoords=(/10.0,10.0, &
20.0,10.0, &
10.0,20.0, &
14.0,14.0/)
!! node owners
allocate(nodeOwners(numNodes))
nodeOwners=(/0,1,2,3/)
! Fill in elem data
numTotElems=1
!! elem ids
allocate(elemIds(numTotElems))
elemIds=(/4/)
!! elem type
allocate(elemTypes(numTotElems))
elemTypes=ESMF_MESHELEMTYPE_QUAD
!! elem conn
allocate(elemConn(numTotElems*4))
elemConn=(/1,2,4,3/)
endif
endif
! Create Mesh structure in 1 step
dstMesh=ESMF_MeshCreate(parametricDim=2,spatialDim=2, &
coordSys=ESMF_COORDSYS_SPH_DEG, &
nodeIds=nodeIds, nodeCoords=nodeCoords, &
nodeOwners=nodeOwners, elementIds=elemIds,&
elementTypes=elemTypes, elementConn=elemConn, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Array spec
call ESMF_ArraySpecSet(arrayspec, 1, ESMF_TYPEKIND_R8, rc=rc)
! Create dest. field
dstField = ESMF_FieldCreate(dstMesh, arrayspec, meshloc=ESMF_MESHLOC_ELEMENT, &
name="dest", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Create dest. area field
dstAreaField = ESMF_FieldCreate(dstMesh, arrayspec, meshloc=ESMF_MESHLOC_ELEMENT, &
name="dest_area", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Create dest. frac field
dstFracField = ESMF_FieldCreate(dstMesh, arrayspec, meshloc=ESMF_MESHLOC_ELEMENT, &
name="dest_frac", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Create exact dest. field
xdstField = ESMF_FieldCreate(dstMesh, arrayspec, meshloc=ESMF_MESHLOC_ELEMENT, &
name="xdest", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Init destination field to 0.0
! Should only be 1 localDE
call ESMF_FieldGet(dstField, 0, dstFarrayPtr, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Init destination field to 0.0
dstFarrayPtr=0.0
! Init exact destination field
! Should only be 1 localDE
call ESMF_FieldGet(xdstField, 0, xdstFarrayPtr, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! set interpolated function
iconn=1
do i1=1,numTotElems
! Loop through nodes in elem
! to compute point in center
x=0.0
y=0.0
do i2=1,elemTypes(i1)
inode=2*(elemConn(iconn)-1)
x=x+nodeCoords(inode+1)
y=y+nodeCoords(inode+2)
iconn=iconn+1
enddo
x=x*(1.0/REAL(elemTypes(i1),ESMF_KIND_R8))
y=y*(1.0/REAL(elemTypes(i1),ESMF_KIND_R8))
! Set source function
xdstFarrayPtr(i1) = 20.0+x+y
enddo
! For now, Easy set interpolated function
!xdstFarrayPtr=1.0
! deallocate node data
deallocate(nodeIds)
deallocate(nodeCoords)
deallocate(nodeOwners)
! deallocate elem data
deallocate(elemIds)
deallocate(elemTypes)
deallocate(elemConn)
#if 0
call ESMF_MeshWrite(srcMesh,"srcMesh")
call ESMF_MeshWrite(dstMesh,"dstMesh")
#endif
!!! Regrid forward from the A grid to the B grid
! Regrid store
call ESMF_FieldRegridStore( &
srcField, &
srcMaskValues=(/1/), &
dstField=dstField, &
routeHandle=routeHandle, &
regridmethod=ESMF_REGRIDMETHOD_CONSERVE, &
dstFracField=dstFracField, &
srcFracField=srcFracField, &
unmappedaction=ESMF_UNMAPPEDACTION_IGNORE, &
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
! Get the integration weights
call ESMF_FieldRegridGetArea(srcAreaField, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Get the integration weights
call ESMF_FieldRegridGetArea(dstAreaField, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Check if the values are close
minerror(1) = 100000.
maxerror(1) = 0.
error = 0.
errorTot=0.0
dstmass = 0.
! get dst Field
call ESMF_FieldGet(dstField, 0, dstFarrayPtr, computationalLBound=clbnd, &
computationalUBound=cubnd, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! get exact destination Field
call ESMF_FieldGet(xdstField, 0, xdstFarrayPtr, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! get dst area Field
call ESMF_FieldGet(dstAreaField, 0, dstAreaPtr, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! get frac Field
call ESMF_FieldGet(dstFracField, 0, dstFracptr, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! destination grid
!! check relative error
do i1=clbnd(1),cubnd(1)
! This is WRONG, shouldn't include Frac
! dstmass = dstmass + dstFracptr(i1,i2)*dstAreaptr(i1)*fptr(i1)
! Instead do this
dstmass = dstmass + dstAreaptr(i1)*dstFarrayPtr(i1)
! If this destination cell isn't covered by a sig. amount of source, then compute error on it.
! (Note that this is what SCRIP does)
if (dstFracptr(i1) .lt. 0.999) cycle
! write(*,*) i1,"::",dstFarrayPtr(i1),xdstFarrayPtr(i1)
if (xdstFarrayPtr(i1) .ne. 0.0) then
error=ABS(dstFarrayPtr(i1) - xdstFarrayPtr(i1))/ABS(xdstFarrayPtr(i1))
errorTot=errorTot+error
if (error > maxerror(1)) then
maxerror(1) = error
endif
if (error < minerror(1)) then
minerror(1) = error
endif
else
error=ABS(dstFarrayPtr(i1) - xdstFarrayPtr(i1))/ABS(xdstFarrayPtr(i1))
errorTot=errorTot+error
if (error > maxerror(1)) then
maxerror(1) = error
endif
if (error < minerror(1)) then
minerror(1) = error
endif
endif
enddo
srcmass(1) = 0.
! get src pointer
call ESMF_FieldGet(srcField, 0, srcFarrayPtr, computationalLBound=clbnd, &
computationalUBound=cubnd, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! get src Field
call ESMF_FieldGet(srcAreaField, 0, srcAreaptr, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! get frac Field
call ESMF_FieldGet(srcFracField, 0, srcFracptr, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
do i1=clbnd(1),cubnd(1)
srcmass(1) = srcmass(1) + srcFracptr(i1)*srcAreaptr(i1)*srcFarrayPtr(i1)
enddo
srcmassg(1) = 0.
dstmassg(1) = 0.
call ESMF_VMAllReduce(vm, srcmass, srcmassg, 1, ESMF_REDUCE_SUM, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_VMAllReduce(vm, dstmass, dstmassg, 1, ESMF_REDUCE_SUM, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_VMAllReduce(vm, maxerror, maxerrorg, 1, ESMF_REDUCE_MAX, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_VMAllReduce(vm, minerror, minerrorg, 1, ESMF_REDUCE_MIN, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! return answer based on correct flags
csrv = .false.
if (ABS(dstmassg(1)-srcmassg(1))/srcmassg(1) < 10E-10) csrv = .true.
itrp = .false.
if (maxerrorg(1) < 10E-2) itrp = .true.
! Uncomment these calls to see some actual regrid results
if (localPet == 0) then
write(*,*) "=== Mesh with concave quads ==="
write(*,*) "Conservation:"
write(*,*) "Rel Error = ", ABS(dstmassg(1)-srcmassg(1))/srcmassg(1)
write(*,*) "SRC mass = ", srcmassg(1)
write(*,*) "DST mass = ", dstmassg(1)
write(*,*) " "
write(*,*) "Interpolation:"
write(*,*) "Max Error = ", maxerrorg(1)
write(*,*) "Min Error = ", minerrorg(1)
write(*,*) "Avg Error = ", (maxerrorg(1) + minerrorg(1))/2
write(*,*) " "
endif
! 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(srcAreaField, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_FieldDestroy(dstAreaField, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_FieldDestroy(srcFracField, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_FieldDestroy(dstFracField, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_FieldDestroy(xdstField, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Free the meshes
call ESMF_MeshDestroy(srcMesh, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_MeshDestroy(dstMesh, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! return success if we've gotten this far
rc=ESMF_SUCCESS
end subroutine test_RegridSph4ConcaveMesh