function ESMF_XGridCreate(keywordEnforcer, &
sideAGrid, sideAMesh, &
sideBGrid, sideBMesh, &
sideAGridPriority, sideAMeshPriority, &
sideBGridPriority, sideBMeshPriority, &
sideAMaskValues, sideBMaskValues, &
storeOverlay, &
name, rc)
!
! !RETURN VALUE:
type(ESMF_XGrid) :: ESMF_XGridCreate
!
! !ARGUMENTS:
type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below
type(ESMF_Grid), intent(in), optional :: sideAGrid(:)
type(ESMF_Mesh), intent(in), optional :: sideAMesh(:)
type(ESMF_Grid), intent(in), optional :: sideBGrid(:)
type(ESMF_Mesh), intent(in), optional :: sideBMesh(:)
integer, intent(in), optional :: sideAGridPriority(:)
integer, intent(in), optional :: sideAMeshPriority(:)
integer, intent(in), optional :: sideBGridPriority(:)
integer, intent(in), optional :: sideBMeshPriority(:)
integer(ESMF_KIND_I4),intent(in), optional :: sideAMaskValues(:)
integer(ESMF_KIND_I4),intent(in), optional :: sideBMaskValues(:)
logical, intent(in), optional :: storeOverlay
character(len=*), intent(in), optional :: name
integer, intent(out),optional :: rc
!
! !DESCRIPTION:
! Create an XGrid from user supplied input: the list of Grids or Meshes on side A and side B,
! and other optional arguments. A user can supply both Grids and Meshes on one side to create
! the XGrid. By default, the Grids have a higher priority over Meshes but the order of priority
! can be adjusted by the optional GridPriority and MeshPriority arguments. The priority order
! of Grids and Meshes can also be interleaved by rearranging the optional
! GridPriority and MeshPriority arguments accordingly.
!
! Sparse matrix multiplication coefficients are internally computed and
! uniquely determined by the Grids or Meshes provided in {\tt sideA} and {\tt sideB}. User can supply
! a single {\tt ESMF\_Grid} or an array of {\tt ESMF\_Grid} on either side of the
! {\tt ESMF\_XGrid}. For an array of {\tt ESMF\_Grid} or {\tt ESMF\_Mesh} in {\tt sideA} or {\tt sideB},
! a merging process concatenates all the {\tt ESMF\_Grid}s and {\tt ESMF\_Mesh}es
! into a super mesh represented
! by {\tt ESMF\_Mesh}. The super mesh is then used to compute the XGrid.
! Grid or Mesh objects in {\tt sideA} and {\tt sideB} arguments must have coordinates defined for
! the corners of a Grid or Mesh cell. XGrid creation can be potentially memory expensive given the
! size of the input Grid and Mesh objects. By default, the super mesh is not stored
! to reduce memory usage.
! Once communication routehandles are computed using {\tt ESMF\_FieldRegridStore()} method through
! XGrid, all memory can be released by destroying the XGrid.
!
! If {\tt sideA} and {\tt sideB} have a single
! Grid or Mesh object, it's erroneous
! if the two Grids or Meshes are spatially disjoint.
! It is also erroneous to specify a Grid or Mesh object in {\tt sideA} or {\tt sideB}
! that is spatially disjoint from the {\tt ESMF\_XGrid}.
!
! This call is {\em collective} across the current VM. For more details please refer to the description
! \ref{sec:xgrid:desc} of the XGrid class. For an example and associated documentation using this method see section
! \ref{sec:xgrid:usage:xgrid_create}
!
! The arguments are:
! \begin{description}
! \item [{[sideAGrid]}]
! Parametric 2D Grids on side A, for example,
! these Grids can be either Cartesian 2D or Spherical.
! \item [{[sideAMesh]}]
! Parametric 2D Meshes on side A, for example,
! these Meshes can be either Cartesian 2D or Spherical.
! \item [{[sideBGrid]}]
! Parametric 2D Grids on side B, for example,
! these Grids can be either Cartesian 2D or Spherical.
! \item [{[sideBMesh]}]
! Parametric 2D Meshes on side B, for example,
! these Meshes can be either Cartesian 2D or Spherical.
! \item [{[sideAGridPriority]}]
! Priority array of Grids on sideA during overlay generation.
! The {\tt sideAGridPriority} array should be the same size as the {\tt sideAGrid} array. The values
! in the array should range from 1 to size(sideAGrid)+size(sideAMesh). A Grid whose corresponding
! value in this array is lower than another side A Grid or Mesh, will take precedence over that Grid or Mesh
! during side A merging. In other words, if both have parts in the same region, then the object with the lower value will win, and
! the other Grid or Mesh part will be clipped away.
! \item [{[sideAMeshPriority]}]
! Priority array of Meshes on sideA during overlay generation.
! The {\tt sideAMeshPriority} array should be the same size as the {\tt sideAMesh} array. The values
! in the array should range from 1 to size(sideAGrid)+size(sideAMesh). A Mesh whose corresponding
! value in this array is lower than another side A Grid or Mesh, will take precedence over that Grid or Mesh
! during side A merging. In other words, if both have parts in the same region, then the object with the lower value will win, and
! the other Grid or Mesh part will be clipped away.
! \item [{[sideBGridPriority]}]
! Priority array of Grids on sideB during overlay generation.
! The {\tt sideBGridPriority} array should be the same size as the {\tt sideBGrid} array. The values
! in the array should range from 1 to size(sideBGrid)+size(sideBMesh). A Grid whose corresponding
! value in this array is lower than another side B Grid or Mesh, will take precedence over that Grid or Mesh
! during side B merging. In other words, if both have parts in the same region, then the object with the lower value will win, and
! the other Grid or Mesh part will be clipped away.
! \item [{[sideBMeshPriority]}]
! Priority array of Meshes on sideB during overlay generation.
! The {\tt sideBMeshPriority} array should be the same size as the {\tt sideBMesh} array. The values
! in the array should range from 1 to size(sideBGrid)+size(sideBMesh). A Mesh whose corresponding
! value in this array is lower than another side B Grid or Mesh, will take precedence over that Grid or Mesh
! during side B merging. In other words, if both have parts in the same region, then the object with the lower value will win, and
! the other Grid or Mesh part will be clipped away.
! \item [{[sideAMaskValues]}]
! Mask information can be set in the Grid (see~\ref{sec:usage:items}) or Mesh (see~\ref{sec:mesh:mask})
! upon which the {\tt Field} is built. The {\tt sideAMaskValues} argument specifies the values in that
! mask information which indicate a point should be masked out. In other words, a location is masked if and only if the
! value for that location in the mask information matches one of the values listed in {\tt sideAMaskValues}.
! If {\tt sideAMaskValues} is not specified, no masking on side A will occur.
! \item [{[sideBMaskValues]}]
! Mask information can be set in the Grid (see~\ref{sec:usage:items}) or Mesh (see~\ref{sec:mesh:mask})
! upon which the {\tt Field} is built. The {\tt sideBMaskValues} argument specifies the values in that
! mask information which indicate a point should be masked out. In other words, a location is masked if and only if the
! value for that location in the mask information matches one of the values listed in {\tt sideBMaskValues}.
! If {\tt sideBMaskValues} is not specified, no masking on side B will occur.
! \item [{[storeOverlay]}]
! Setting the {\tt storeOverlay} optional argument to .false. (default)
! allows a user to bypass storage of the {\tt ESMF\_Mesh} used to represent the XGrid.
! Only a {\tt ESMF\_DistGrid} is stored to allow Field to be built on the XGrid.
! If the temporary mesh object is of interest, {\tt storeOverlay} can be set to .true.
! so a user can retrieve it for future use.
! \item [{[name]}]
! name of the xgrid object.
! \item [{[rc]}]
! Return code; equals {\tt ESMF\_SUCCESS} only if the {\tt ESMF\_XGrid}
! is created.
! \end{description}
!
!EOP
integer :: localrc, ngrid_a, ngrid_b, i, j
type(ESMF_XGridType), pointer :: xgtype
type(ESMF_Mesh) :: meshA, meshB, mesh, tmpmesh
type(ESMF_Mesh), allocatable :: meshAt(:), meshBt(:)
type(ESMF_Pointer) :: meshp
integer(ESMF_KIND_I4), pointer:: indicies(:,:)
real(ESMF_KIND_R8), pointer :: weights(:), sidemesharea(:)
integer :: nentries
type(ESMF_TempWeights) :: tweights
integer :: AisSphere, BisSphere
integer :: compute_midmesh
real(ESMF_KIND_R8) :: fraction = 1.0 ! all newly created Mesh has 1.0 frac2
type(ESMF_INDEX_FLAG) :: indexflag
type(ESMF_DistGrid) :: distgridTmp
!real(ESMF_KIND_R8), pointer :: fracFptr(:,:)
integer :: localElemCount, sdim, pdim, sdim2
type(ESMF_XGridGeomType_Flag), allocatable :: xggt_a(:), xggt_b(:)
integer :: tileCount
integer :: side
type(ESMF_CoordSys_Flag) :: coordSys
! Initialize
localrc = ESMF_RC_NOT_IMPL
! Initialize return code
if(present(rc)) rc = ESMF_RC_NOT_IMPL
! Presumably no longer useful
AisSphere = 0
BisSphere = 0
! check there are enough input to create the XGrid
if(.not. present(sideAGrid) .and. .not. present(sideAMesh)) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- Either Grid or Mesh must be provided on sideA", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
if(.not. present(sideBGrid) .and. .not. present(sideBMesh)) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- Either Grid or Mesh must be provided on sideB", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
if(present(sideAGrid) .and. present(sideAMesh)) then
if(size(sideAGrid, 1)+size(sideAMesh, 1) .le. 0) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- Either Grid or Mesh must be provided on sideA", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
endif
if(present(sideBGrid) .and. present(sideBMesh)) then
if(size(sideBGrid, 1)+size(sideBMesh, 1) .le. 0) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- Either Grid or Mesh must be provided on sideB", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
endif
! check init status of input Grids
if(present(sideAGrid)) then
ngrid_a = size(sideAGrid, 1)
do i = 1, ngrid_a
ESMF_INIT_CHECK_DEEP(ESMF_GridGetInit,sideAGrid(i),rc)
enddo
! Can only do conservative on 2D right now
! Make sure all Grids are 2 dimensional and
! has enough data points in each dimension for every de-element on a PET to clip
do i = 1, ngrid_a
call checkGrid(sideAGrid(i), ESMF_STAGGERLOC_CORNER, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
if(present(sideAGridPriority)) then
if(size(sideAGridPriority, 1) /= ngrid_a) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- Size of sideAGridPriority does not match size of sideAGrid", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
endif
endif
if(present(sideBGrid)) then
ngrid_b = size(sideBGrid, 1)
do i = 1, ngrid_b
ESMF_INIT_CHECK_DEEP(ESMF_GridGetInit,sideBGrid(i),rc)
enddo
do i = 1, ngrid_b
call checkGrid(sideBGrid(i), ESMF_STAGGERLOC_CORNER, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
if(present(sideBGridPriority)) then
if(size(sideBGridPriority, 1) /= ngrid_b) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- Size of sideBGridPriority does not match size of sideBGrid", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
endif
endif
!TODO: need to check Meshes are initialized properly too.
if(present(sideAMesh)) then
ngrid_a = size(sideAMesh, 1)
if(present(sideAMeshPriority)) then
if(size(sideAMeshPriority, 1) /= ngrid_a) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- Size of sideAMeshPriority does not match size of sideAMesh", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
endif
endif
if(present(sideBMesh)) then
ngrid_b = size(sideBMesh, 1)
if(present(sideBMeshPriority)) then
if(size(sideBMeshPriority, 1) /= ngrid_b) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- Size of sideBMeshPriority does not match size of sideBMesh", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
endif
endif
! Priority for Grid and Mesh must all be present simultaneously
if( (present(sideAGridPriority) .and. .not. present(sideAGrid)) .or. &
present(sideAMeshPriority) .and. .not. present(sideAMesh) .or. &
present(sideBGridPriority) .and. .not. present(sideBGrid) .or. &
present(sideBMeshPriority) .and. .not. present(sideBMesh)) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- Cannot specify Grid or Mesh Priority without actual list of Grids or Meshes", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
! Make sure that the input objects have a consisent coordSys
call error_check_coordsys(sideAGrid, sideAMesh, sideBGrid, sideBMesh, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Output Mask Values for debugging
#if 0
if (present(sideAMaskValues)) then
write(*,*) "sideAMaskValues=",sideAMaskValues
endif
if (present(sideBMaskValues)) then
write(*,*) "sideBMaskValues=",sideBMaskValues
endif
#endif
! At this point the inputs are consistently sized.
! Take care of ordering
ngrid_a = 0
if(present(sideAGrid)) ngrid_a = size(sideAGrid, 1)
if(present(sideAMesh)) ngrid_a = ngrid_a + size(sideAMesh, 1)
ngrid_b = 0
if(present(sideBGrid)) ngrid_b = size(sideBGrid, 1)
if(present(sideBMesh)) ngrid_b = ngrid_b + size(sideBMesh, 1)
! do some range checking on priority lists
if(present(sideAGridPriority) .and. present(sideAMeshPriority)) then
do i = 1, size(sideAGridPriority, 1)
if(sideAGridPriority(i) .le. 0 .or. sideAGridPriority(i) .gt. ngrid_a) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- sideAGridPriority value out of range", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
enddo
do i = 1, size(sideAMeshPriority, 1)
if(sideAMeshPriority(i) .le. 0 .or. sideAMeshPriority(i) .gt. ngrid_a) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- sideAMeshPriority value out of range", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
enddo
do i = 1, size(sideAGridPriority, 1)
do j = 1, size(sideAMeshPriority, 1)
if(sideAGridPriority(i) == sideAMeshPriority(j)) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- sideAGridPriority and sideAMeshPriority cannot have duplicate entry", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
enddo
enddo
endif
if(present(sideBGridPriority) .and. present(sideBMeshPriority)) then
do i = 1, size(sideBGridPriority, 1)
if(sideBGridPriority(i) .le. 0 .or. sideBGridPriority(i) .gt. ngrid_b) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- sideBGridPriority value out of range", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
enddo
do i = 1, size(sideBMeshPriority, 1)
if(sideBMeshPriority(i) .le. 0 .or. sideBMeshPriority(i) .gt. ngrid_b) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- sideBMeshPriority value out of range", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
enddo
do i = 1, size(sideBGridPriority, 1)
do j = 1, size(sideBMeshPriority, 1)
if(sideBGridPriority(i) == sideBMeshPriority(j)) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- sideBGridPriority and sideBMeshPriority cannot have duplicate entry", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
enddo
enddo
endif
! at this point, the input priority lists have valid entries within correct range.
! initialize XGridType object and its base object
nullify(xgtype)
nullify(ESMF_XGridCreate%xgtypep)
call ESMF_XGridConstructBaseObj(xgtype, name, localrc)
if (ESMF_LogFoundAllocError(localrc, &
msg="Constructing xgtype base object ", &
ESMF_CONTEXT, rcToReturn=rc)) return
! Choose a default coordSys for XGrid
! TODO: eventually also allow a user to pass one in
call choose_default_coordsys(sideAGrid, sideAMesh, sideBGrid, sideBMesh, xgtype%coordSys, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(xgtype%sideA(ngrid_a), xgtype%sideB(ngrid_b), stat=localrc)
if (ESMF_LogFoundAllocError(localrc, &
msg="- Allocating xgtype%sideA or xgtype%sideB ", &
ESMF_CONTEXT, rcToReturn=rc)) return
! Need to initialize xgtype%sideA based on sideAGrid and/or sideAMesh
if(present(sideAGrid)) then
if(present(sideAGridPriority)) then
do i = 1, size(sideAGrid, 1)
xgtype%sideA(sideAGridPriority(i)) = ESMF_XGridGeomBaseCreate(sideAGrid(i), &
ESMF_STAGGERLOC_CENTER, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
else
do i = 1, size(sideAGrid, 1)
xgtype%sideA(i) = ESMF_XGridGeomBaseCreate(sideAGrid(i), &
ESMF_STAGGERLOC_CENTER, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
endif
if(present(sideAMesh)) then
if(present(sideAMeshPriority)) then
do i = 1, size(sideAMesh, 1)
xgtype%sideA(sideAMeshPriority(i)) = ESMF_XGridGeomBaseCreate(sideAMesh(i), &
ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
else
do i = 1, size(sideAMesh, 1)
xgtype%sideA(i+size(sideAGrid, 1)) = ESMF_XGridGeomBaseCreate(sideAMesh(i), &
ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
endif
endif
else ! .not. present(sideAGrid)
if(present(sideAMesh)) then
if(present(sideAMeshPriority)) then
do i = 1, size(sideAMesh, 1)
xgtype%sideA(sideAMeshPriority(i)) = ESMF_XGridGeomBaseCreate(sideAMesh(i), &
ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
else
do i = 1, size(sideAMesh, 1)
xgtype%sideA(i) = ESMF_XGridGeomBaseCreate(sideAMesh(i), &
ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
endif
endif
endif
! Need to initialize xgtype%sideB based on sideBGrid and/or sideBMesh
if(present(sideBGrid)) then
if(present(sideBGridPriority)) then
do i = 1, size(sideBGrid, 1)
xgtype%sideB(sideBGridPriority(i)) = ESMF_XGridGeomBaseCreate(sideBGrid(i), &
ESMF_STAGGERLOC_CENTER, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
else
do i = 1, size(sideBGrid, 1)
xgtype%sideB(i) = ESMF_XGridGeomBaseCreate(sideBGrid(i), &
ESMF_STAGGERLOC_CENTER, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
endif
if(present(sideBMesh)) then
if(present(sideBMeshPriority)) then
do i = 1, size(sideBMesh, 1)
xgtype%sideB(sideBMeshPriority(i)) = ESMF_XGridGeomBaseCreate(sideBMesh(i), &
ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
else
do i = 1, size(sideBMesh, 1)
xgtype%sideB(i+size(sideBGrid, 1)) = ESMF_XGridGeomBaseCreate(sideBMesh(i), &
ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
endif
endif
else ! .not. present(sideBGrid)
if(present(sideBMesh)) then
if(present(sideBMeshPriority)) then
do i = 1, size(sideBMesh, 1)
xgtype%sideB(sideBMeshPriority(i)) = ESMF_XGridGeomBaseCreate(sideBMesh(i), &
ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
else
do i = 1, size(sideBMesh, 1)
xgtype%sideB(i) = ESMF_XGridGeomBaseCreate(sideBMesh(i), &
ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
endif
endif
endif
! allocate the temporary meshes
allocate(meshAt(ngrid_a), meshBt(ngrid_b), &
xggt_a(ngrid_a), xggt_b(ngrid_b), stat=localrc)
if (ESMF_LogFoundAllocError(localrc, &
msg="- Allocating temporary meshes for Xgrid creation", &
ESMF_CONTEXT, rcToReturn=rc)) return
do i = 1, ngrid_a
call ESMF_XGridGeomBaseGet(xgtype%sideA(i), geomtype=xggt_a(i), rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if(xggt_a(i) == ESMF_XGRIDGEOMTYPE_GRID) then
! Get number of tiles in Grid
call ESMF_GridGet(grid=xgtype%sideA(i)%gbcp%grid, &
tileCount=tileCount, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Use different methods for GridToMesh depending on tileCount.
! TODO: Get rid of this and just use new method for everything
if (tileCount .eq. 1) then
! Create Mesh from Grid
meshAt(i) = ESMF_GridToMesh(xgtype%sideA(i)%gbcp%grid, &
ESMF_STAGGERLOC_CORNER, AisSphere, &
maskValues=sideAMaskValues, &
regridConserve=ESMF_REGRID_CONSERVE_ON, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
! Create Mesh from Grid
meshAt(i) = ESMF_GridToMeshCell(xgtype%sideA(i)%gbcp%grid, &
rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Turn on masking
if (present(sideAMaskValues)) then
call ESMF_MeshTurnOnCellMask(meshAt(i), &
maskValues=sideAMaskValues, rc=localrc);
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
endif
else if(xggt_a(i) == ESMF_XGRIDGEOMTYPE_MESH) then
meshAt(i) = xgtype%sideA(i)%gbcp%mesh
if (present(sideAMaskValues)) then
call ESMF_MeshTurnOnCellMask(meshAt(i), maskValues=sideAMaskValues, rc=localrc);
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
call c_esmc_meshsetfraction(meshAt(i), fraction, localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- Invalid sideA xgridgeombase object", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
! Set temp xgrid info in mesh
side=1
call c_esmc_meshsetxgridinfo(meshAt(i), side, i, localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Merge meshes
if(i == 1) meshA = meshAt(i)
if(i .ge. 2) then
! call into mesh merge with priority taken into account
! meshAt is truncated(if necessary) and concatenated onto meshA
! and result stored in tmpmesh
call c_esmc_meshmerge(meshA, meshAt(i), tmpmesh, localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if(i .gt. 2) then
! the intermediate meshA is only a pointer type of mesh at this point,
! call the C api to destroy it
call C_ESMC_MeshDestroy(meshA, ESMF_TRUE, localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
meshA = tmpmesh
endif
enddo
do i = 1, ngrid_b
call ESMF_XGridGeomBaseGet(xgtype%sideB(i), geomtype=xggt_b(i), rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if(xggt_b(i) == ESMF_XGRIDGEOMTYPE_GRID) then
! Get number of tiles in Grid
call ESMF_GridGet(grid=xgtype%sideB(i)%gbcp%grid, &
tileCount=tileCount, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Use different methods for GridToMesh depending on tileCount.
! TODO: Get rid of this and just use new method for everything
if (tileCount .eq. 1) then
! Create Mesh from Grid
meshBt(i) = ESMF_GridToMesh(xgtype%sideB(i)%gbcp%grid, &
ESMF_STAGGERLOC_CORNER, BisSphere, &
maskValues=sideBMaskValues, &
regridConserve=ESMF_REGRID_CONSERVE_ON, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
! Create Mesh from Grid
meshBt(i) = ESMF_GridToMeshCell(xgtype%sideB(i)%gbcp%grid, &
rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Turn on masking
if (present(sideBMaskValues)) then
call ESMF_MeshTurnOnCellMask(meshBt(i), &
maskValues=sideBMaskValues, rc=localrc);
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
endif
else if(xggt_b(i) == ESMF_XGRIDGEOMTYPE_MESH) then
meshBt(i) = xgtype%sideB(i)%gbcp%mesh
if (present(sideBMaskValues)) then
call ESMF_MeshTurnOnCellMask(meshBt(i), maskValues=sideBMaskValues, rc=localrc);
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
call c_esmc_meshsetfraction(meshBt(i), fraction, localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- Invalid sideB xgridgeombase object", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
! Set temp xgrid info in mesh
side=2
call c_esmc_meshsetxgridinfo(meshBt(i), side, i, localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Merge meshes
if( i == 1) meshB = meshBt(i)
! call into mesh merge with priority taken into account
! meshBt is truncated(if necessary) and concatenated onto meshB
! and result stored in tmpmesh
if(i .ge. 2) then
call c_esmc_meshmerge(meshB, meshBt(i), tmpmesh, localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if(i .gt. 2) then
call C_ESMC_MeshDestroy(meshB, ESMF_TRUE, localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
meshB = tmpmesh
endif
enddo
allocate(xgtype%sparseMatA2X(ngrid_a), &
xgtype%sparseMatX2A(ngrid_a), &
xgtype%sparseMatB2X(ngrid_b), &
xgtype%sparseMatX2B(ngrid_b), stat=localrc)
if(localrc /= 0) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- Failed to allocate SMM parameters", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
! Call into a modifed conservative regrid call to create XGrid (super) Mesh
! input: MeshA: merged mesh on side A
! input: MeshB: merged mesh on side B
! output: Meshp: merged XGrid mesh in the middle (the super mesh)
compute_midmesh = 1
call c_esmc_xgridregrid_create(meshA, meshB, &
meshp, compute_midmesh, &
ESMF_REGRIDMETHOD_CONSERVE, &
ESMF_UNMAPPEDACTION_IGNORE, &
xgtype%coordSys, &
nentries, tweights, &
localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Wrap XGrid mesh pointer in proper F90 structure
mesh = ESMF_MeshCreate(meshp, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Get mesh info
call ESMF_MeshGet(mesh, numOwnedElements=localElemCount, &
spatialDim=sdim, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call C_ESMC_MeshGetDimensions(mesh, sdim2, pdim, coordSys, localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if(sdim /= sdim2) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_VAL_OUTOFRANGE, &
msg="spatial dim mismatch", ESMF_CONTEXT, rcToReturn=rc)
return
endif
allocate(xgtype%area(localElemCount), xgtype%centroid(localElemCount, sdim), stat=localrc)
if(localrc /= 0) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- Failed to allocate area or centroid", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
if(localElemCount .gt. 0) then
call C_ESMC_MeshGetArea(mesh, localElemCount, xgtype%area, localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call C_ESMC_MeshGetCentroid(mesh, localElemCount, xgtype%centroid, localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
! Create and Retrieve fraction Arrays for store use, only for online
allocate(xgtype%fracA2X(ngrid_a), xgtype%fracB2X(ngrid_b))
allocate(xgtype%fracX2A(ngrid_a), xgtype%fracX2B(ngrid_b))
allocate(xgtype%frac2A(ngrid_a), xgtype%frac2B(ngrid_b))
do i = 1, ngrid_A
if(xggt_a(i) == ESMF_XGRIDGEOMTYPE_GRID) then
call ESMF_GridGet(xgtype%sideA(i)%gbcp%grid, &
indexflag=indexflag, &
distgrid=distgridTmp, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_MeshGet(xgtype%sideA(i)%gbcp%mesh, elementDistgrid=distgridTmp, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! use global index for Mesh
indexflag=ESMF_INDEX_GLOBAL
endif
xgtype%fracA2X(i) = ESMF_ArrayCreate(distgridTmp, typekind=ESMF_TYPEKIND_R8, &
indexflag=indexflag, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
xgtype%fracX2A(i) = ESMF_ArrayCreate(distgridTmp, typekind=ESMF_TYPEKIND_R8, &
indexflag=indexflag, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
xgtype%frac2A(i) = ESMF_ArrayCreate(distgridTmp, typekind=ESMF_TYPEKIND_R8, &
indexflag=indexflag, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! retrieve frac2 Field
if(xggt_a(i) == ESMF_XGRIDGEOMTYPE_GRID) then
call ESMF_XGridGetFrac2Int(xgtype%sideA(i)%gbcp%grid, mesh=meshAt(i), array=xgtype%frac2A(i), &
staggerloc=ESMF_STAGGERLOC_CORNER, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_XGridGetFrac2Int(mesh=meshAt(i), array=xgtype%frac2A(i), &
meshloc=ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
enddo
do i = 1, ngrid_B
if(xggt_b(i) == ESMF_XGRIDGEOMTYPE_GRID) then
call ESMF_GridGet(xgtype%sideB(i)%gbcp%grid, &
indexflag=indexflag, &
distgrid=distgridTmp, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_MeshGet(xgtype%sideB(i)%gbcp%mesh, elementDistgrid=distgridTmp, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! use global index for Mesh
indexflag=ESMF_INDEX_GLOBAL
endif
xgtype%fracB2X(i) = ESMF_ArrayCreate(distgridTmp, typekind=ESMF_TYPEKIND_R8, &
indexflag=indexflag, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
xgtype%fracX2B(i) = ESMF_ArrayCreate(distgridTmp, typekind=ESMF_TYPEKIND_R8, &
indexflag=indexflag, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
xgtype%frac2B(i) = ESMF_ArrayCreate(distgridTmp, typekind=ESMF_TYPEKIND_R8, &
indexflag=indexflag, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! retrieve frac2 Field
if(xggt_b(i) == ESMF_XGRIDGEOMTYPE_GRID) then
call ESMF_XGridGetFrac2Int(xgtype%sideB(i)%gbcp%grid, mesh=meshBt(i), array=xgtype%frac2B(i), &
staggerloc=ESMF_STAGGERLOC_CORNER, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_XGridGetFrac2Int(mesh=meshBt(i), array=xgtype%frac2B(i), &
meshloc=ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
enddo
! TODO: investigate the possibility of optimization for the general case of multiple Grids.
! When there is only 1 grid per side, optimization can be done but it's not clear for multiple Grids.
compute_midmesh = 0
do i = 1, ngrid_a
#ifndef XGRID_WGT_CALC_OLDWAY
! Calculate weights from a side mesh to the xgrid
call c_esmc_xgrid_calc_wgts_from_mesh(meshAt(i), mesh, &
nentries, tweights, localrc)
#else
call c_esmc_xgridregrid_createP(meshAt(i), mesh, &
tmpmesh, compute_midmesh, &
ESMF_REGRIDMETHOD_CONSERVE, &
ESMF_UNMAPPEDACTION_IGNORE, &
xgtype%coordSys, &
nentries, tweights, &
localrc)
#endif
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(xgtype%sparseMatA2X(i)%factorIndexList(2,nentries))
allocate(xgtype%sparseMatA2X(i)%factorList(nentries))
if(nentries .ge. 1) then
call c_ESMC_Copy_TempWeights_xgrid(tweights, &
xgtype%sparseMatA2X(i)%factorIndexList(1,1), &
xgtype%sparseMatA2X(i)%factorList(1))
#if 0
! DEBUG OUTPUT
do j=1,size(xgtype%sparseMatA2X(i)%factorList,1)
write(*,*) "A2X: Grid=",i," src=",xgtype%sparseMatA2X(i)%factorIndexList(1,j)," wgt=",xgtype%sparseMatA2X(i)%factorList(j)," dst=",xgtype%sparseMatA2X(i)%factorIndexList(2,j)
enddo
#endif
endif
#ifdef XGRID_WGT_CALC_OLDWAY
if(xggt_a(i) == ESMF_XGRIDGEOMTYPE_GRID) then
call ESMF_XGridGetFracInt(xgtype%sideA(i)%gbcp%grid, mesh=meshAt(i), array=xgtype%fracA2X(i), &
staggerloc=ESMF_STAGGERLOC_CORNER, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_XGridGetFracInt(mesh=meshAt(i), array=xgtype%fracA2X(i), &
meshloc=ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
#endif
! Now the reverse direction
#ifndef BOB_XGRID_DEBUG
#ifndef XGRID_WGT_CALC_OLDWAY
! Calculate weights from the xgrid to a side mesh
! Also, set fraction information in meshAt(i) to be retrieved below
call c_esmc_xgrid_calc_wgts_to_mesh(mesh, meshAt(i), &
nentries, tweights, localrc)
#else
call c_esmc_xgridregrid_createP(mesh, meshAt(i), &
tmpmesh, compute_midmesh, &
ESMF_REGRIDMETHOD_CONSERVE, &
ESMF_UNMAPPEDACTION_IGNORE, &
xgtype%coordSys, &
nentries, tweights, &
localrc)
#endif
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
#else
nentries=0
#endif
allocate(xgtype%sparseMatX2A(i)%factorIndexList(2,nentries))
allocate(xgtype%sparseMatX2A(i)%factorList(nentries))
if(nentries .ge. 1) then
call c_ESMC_Copy_TempWeights_xgrid(tweights, &
xgtype%sparseMatX2A(i)%factorIndexList(1,1), &
xgtype%sparseMatX2A(i)%factorList(1))
#if 0
! DEBUG OUTPUT
do j=1,size(xgtype%sparseMatX2A(i)%factorList,1)
write(*,*) "X2A: Grid=",i," src=",xgtype%sparseMatX2A(i)%factorIndexList(1,j)," wgt=",xgtype%sparseMatX2A(i)%factorList(j)," dst=",xgtype%sparseMatX2A(i)%factorIndexList(2,j)
enddo
#endif
endif
#ifndef XGRID_WGT_CALC_OLDWAY
! If doing things the new way, the fractions are set after xgrid to side mesh weight calc, so get them here
if(xggt_a(i) == ESMF_XGRIDGEOMTYPE_GRID) then
call ESMF_XGridGetFracInt(xgtype%sideA(i)%gbcp%grid, mesh=meshAt(i), array=xgtype%fracA2X(i), &
staggerloc=ESMF_STAGGERLOC_CORNER, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_XGridGetFracInt(mesh=meshAt(i), array=xgtype%fracA2X(i), &
meshloc=ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
#endif
if(xggt_a(i) == ESMF_XGRIDGEOMTYPE_GRID) then
call ESMF_XGridGetFracInt(xgtype%sideA(i)%gbcp%grid, mesh=meshAt(i), array=xgtype%fracX2A(i), &
staggerloc=ESMF_STAGGERLOC_CORNER, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_XGridGetFracInt(mesh=meshAt(i), array=xgtype%fracX2A(i), &
meshloc=ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
enddo
! now do the B side
do i = 1, ngrid_b
#ifndef BOB_XGRID_DEBUG
#ifndef XGRID_WGT_CALC_OLDWAY
! Calculate weights from a side mesh to the xgrid
call c_esmc_xgrid_calc_wgts_from_mesh(meshBt(i), mesh, &
nentries, tweights, localrc)
#else
call c_esmc_xgridregrid_createP(meshBt(i), mesh, &
tmpmesh, compute_midmesh, &
ESMF_REGRIDMETHOD_CONSERVE, &
ESMF_UNMAPPEDACTION_IGNORE, &
xgtype%coordSys, &
nentries, tweights, &
localrc)
#endif
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
#else
nentries=0
#endif
allocate(xgtype%sparseMatB2X(i)%factorIndexList(2,nentries))
allocate(xgtype%sparseMatB2X(i)%factorList(nentries))
if(nentries .ge. 1) then
call c_ESMC_Copy_TempWeights_xgrid(tweights, &
xgtype%sparseMatB2X(i)%factorIndexList(1,1), &
xgtype%sparseMatB2X(i)%factorList(1))
endif
#ifdef XGRID_WGT_CALC_OLDWAY
if(xggt_b(i) == ESMF_XGRIDGEOMTYPE_GRID) then
call ESMF_XGridGetFracInt(xgtype%sideB(i)%gbcp%grid, mesh=meshBt(i), array=xgtype%fracB2X(i), &
staggerloc=ESMF_STAGGERLOC_CORNER, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_XGridGetFracInt(mesh=meshBt(i), array=xgtype%fracB2X(i), &
meshloc=ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
#endif
! Now the reverse direction
#ifndef BOB_XGRID_DEBUG
#ifndef XGRID_WGT_CALC_OLDWAY
! Calculate weights from the xgrid to a side mesh
! Also, set fraction information in meshBt(i) to be retrieved below
call c_esmc_xgrid_calc_wgts_to_mesh(mesh, meshBt(i), &
nentries, tweights, localrc)
#else
call c_esmc_xgridregrid_createP(mesh, meshBt(i), &
tmpmesh, compute_midmesh, &
ESMF_REGRIDMETHOD_CONSERVE, &
ESMF_UNMAPPEDACTION_IGNORE, &
xgtype%coordSys, &
nentries, tweights, &
localrc)
#endif
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
#else
nentries=0
#endif
allocate(xgtype%sparseMatX2B(i)%factorIndexList(2,nentries))
allocate(xgtype%sparseMatX2B(i)%factorList(nentries))
if(nentries .ge. 1) then
call c_ESMC_Copy_TempWeights_xgrid(tweights, &
xgtype%sparseMatX2B(i)%factorIndexList(1,1), &
xgtype%sparseMatX2B(i)%factorList(1))
endif
#ifndef XGRID_WGT_CALC_OLDWAY
! If doing things the new way, the fractions are set after xgrid to side mesh weight calc, so get them here
if(xggt_b(i) == ESMF_XGRIDGEOMTYPE_GRID) then
call ESMF_XGridGetFracInt(xgtype%sideB(i)%gbcp%grid, mesh=meshBt(i), array=xgtype%fracB2X(i), &
staggerloc=ESMF_STAGGERLOC_CORNER, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_XGridGetFracInt(mesh=meshBt(i), array=xgtype%fracB2X(i), &
meshloc=ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
#endif
if(xggt_b(i) == ESMF_XGRIDGEOMTYPE_GRID) then
call ESMF_XGridGetFracInt(xgtype%sideB(i)%gbcp%grid, mesh=meshBt(i), array=xgtype%fracX2B(i), &
staggerloc=ESMF_STAGGERLOC_CORNER, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_XGridGetFracInt(mesh=meshBt(i), array=xgtype%fracX2B(i), &
meshloc=ESMF_MESHLOC_ELEMENT, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
enddo
xgtype%storeOverlay = .false.
if(present(storeOverlay)) then
if(storeOverlay) xgtype%storeOverlay = .true.
endif
! call into offline xgrid create with the xgrid specs
call ESMF_XGridConstruct(xgtype, xgtype%sideA, xgtype%sideB, &
sparseMatA2X=xgtype%sparseMatA2X, sparseMatX2A=xgtype%sparseMatX2A, &
sparseMatB2X=xgtype%sparseMatB2X, sparseMatX2B=xgtype%sparseMatX2B, &
offline=.false., &
mesh=mesh, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! store the middle mesh if needed
! and clean up temporary memory used
if(xgtype%storeOverlay) then
xgtype%mesh = mesh
!! Debug output of xgrid mesh
#ifdef BOB_XGRID_DEBUG
call ESMF_MeshWrite(mesh, "xgrid_mid_mesh")
#endif
else
call ESMF_MeshDestroy(mesh, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
do i = 1, ngrid_a
if(present(sideAMaskValues)) then
call ESMF_MeshTurnOffCellMask(meshAt(i), rc=localrc);
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
if(xggt_a(i) == ESMF_XGRIDGEOMTYPE_GRID) then
call ESMF_MeshDestroy(meshAt(i), rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
enddo
do i = 1, ngrid_b
if(present(sideBMaskValues)) then
call ESMF_MeshTurnOffCellMask(meshBt(i), rc=localrc);
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
if(xggt_b(i) == ESMF_XGRIDGEOMTYPE_GRID) then
call ESMF_MeshDestroy(meshBt(i), rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
enddo
deallocate(meshAt, meshBt)
deallocate(xggt_a, xggt_b)
! Finalize XGrid Creation
xgtype%online = 1
xgtype%status = ESMF_STATUS_READY
ESMF_XGridCreate%xgtypep => xgtype
ESMF_INIT_SET_CREATED(ESMF_XGridCreate)
!call ESMF_XGridValidate(ESMF_XGridCreate, rc=localrc)
!if (ESMF_LogFoundError(localrc, &
! ESMF_ERR_PASSTHRU, &
! ESMF_CONTEXT, rcToReturn=rc)) return
if(present(rc)) rc = ESMF_SUCCESS
end function ESMF_XGridCreate