function ESMF_GridCreateCubedSphereIReg(tileSize, &
countsPerDEDim1PTile, countsPerDEDim2PTile, &
keywordEnforcer, &
coordSys, coordTypeKind, &
deLabelList, staggerLocList, &
delayout, indexflag, name, transformArgs, rc)
!
! !RETURN VALUE:
type(ESMF_Grid) :: ESMF_GridCreateCubedSphereIReg
!
! !ARGUMENTS:
integer, intent(in) :: tilesize
integer, intent(in) :: countsPerDEDim1PTile(:,:)
integer, intent(in) :: countsPerDEDim2PTile(:,:)
type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below
type(ESMF_CoordSys_Flag), intent(in), optional :: coordSys
type(ESMF_TypeKind_Flag), intent(in), optional :: coordTypeKind
integer, intent(in), optional :: deLabelList(:)
type(ESMF_StaggerLoc), intent(in), optional :: staggerLocList(:)
type(ESMF_DELayout), intent(in), optional :: delayout
type(ESMF_Index_Flag), intent(in), optional :: indexflag
character(len=*), intent(in), optional :: name
type(ESMF_CubedSphereTransform_Args), intent(in), optional :: transformArgs
integer, intent(out), optional :: rc
!
! !DESCRIPTION:
! Create a six-tile {\tt ESMF\_Grid} for a Cubed Sphere grid using irregular decomposition. Each tile can
! have different decomposition. The grid coordinates are generated based on the algorithm used by GEOS-5,
! The tile resolution is defined by tileSize.
!
! The arguments are:
! \begin{description}
! \item[tilesize]
! The number of elements on each side of the tile of the Cubed Sphere grid.
! \item[countsPerDEDim1PTile]
! This array specifies the number of cells per DE for index dimension 1 for the
! center stagger location. The second index steps through the tiles. If each tile is
! decomposed into different number of DEs, the first dimension is the maximal DEs of
! all the tiles.
! \item[countsPerDEDim2PTile]
! This array specifies the number of cells per DE for index dimension 2 for the
! center stagger location. The second index steps through the tiles. If each tile is
! decomposed into different number of DEs, the first dimension is the maximal DEs of
! all the tiles.
! \item[{[coordSys]}]
! The coordinate system of the grid coordinate data.
! Only ESMF\_COORDSYS\_SPH\_DEG and ESMF\_COORDSYS\_SPH\_RAD are supported.
! If not specified then defaults to ESMF\_COORDSYS\_SPH\_DEG.
! \item[{[coordTypeKind]}]
! The type/kind of the grid coordinate data. Only ESMF\_TYPEKIND\_R4
! and ESMF\_TYPEKIND\_R8 are supported.
! If not specified then defaults to ESMF\_TYPEKIND\_R8.
! \item[{[deLabelList]}]
! List assigning DE labels to the default sequence of DEs. The default
! sequence is given by the column major order in the sequence as they appear
! in {\tt countsPerDEDim1PTile}, followed by {\tt countsPerDEDim2PTile}, then the
! tile index.
! \item[{[staggerLocList]}]
! The list of stagger locations to fill with coordinates. Only {\tt ESMF\_STAGGERLOC\_CENTER}
! and {\tt ESMF\_STAGGERLOC\_CORNER} are supported. If not present, no coordinates
! will be added or filled.
! \item[{[delayout]}]
! Optional ESMF\_DELayout object to be used. If a delayout object is specified,
! the number of DEs must match with the total DEs defined in {\tt countsPerDEDim1PTile}
! and {\tt countsPerDEDim2PTile}.
! \item[{[indexflag]}]
! Indicates the indexing scheme to be used in the new Grid. Please see
! Section~\ref{const:indexflag} for the list of options. If not present,
! defaults to ESMF\_INDEX\_DELOCAL.
! \item[{[name]}]
! {\tt ESMF\_Grid} name.
! \item[{[transformArgs]}]
! A data type containing the stretch factor, target longitude, and target latitude
! to perform a Schmidt transformation on the Cubed-Sphere grid. See section
! \ref{sec:usage:cubedspherewttransform} for details.
! \item[{[rc]}]
! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors.
! \end{description}
!
!EOP
type(ESMF_VM) :: vm
integer :: PetNo, PetCnt
type(ESMF_DELayout) :: defaultDELayout
type(ESMF_Grid) :: grid, newgrid
type(ESMF_DistGrid) :: distgrid, newdistgrid
integer :: localrc
integer :: i, j
integer :: nx, ny, nxy, bigFac, totalDE
integer :: localDeCount, localDe, DeNo, tile
integer, pointer :: minIndexPTile(:,:)
integer, pointer :: maxIndexPTile(:,:)
type(ESMF_DistGridConnection), pointer :: connectionList(:)
real(kind=ESMF_KIND_R8), pointer :: lonPtrR8(:,:), latPtrR8(:,:)
real(kind=ESMF_KIND_R8), pointer :: lonCornerPtrR8(:,:), latCornerPtrR8(:,:)
real(kind=ESMF_KIND_R4), pointer :: lonPtrR4(:,:), latPtrR4(:,:)
real(kind=ESMF_KIND_R4), pointer :: lonCornerPtrR4(:,:), latCornerPtrR4(:,:)
integer :: tileCount
integer :: start(2), count(2)
integer :: shapLon(2), shapLat(2)
integer :: decount
type(ESMF_Index_Flag) :: localIndexFlag
type(ESMF_CoordSys_Flag) :: coordSysLocal
type(ESMF_TypeKind_Flag) :: coordTypeKindLocal
integer :: s
logical :: docenter, docorner
integer, pointer :: deBlockList(:,:,:), deToTileMap(:)
integer, pointer :: DeDim1(:), DeDim2(:), demap(:)
integer :: k,t, minIndx, minIndy
integer :: myde, startde, endde
integer :: tiles, totalelmt
!real(ESMF_KIND_R8) :: starttime, endtime
if (present(rc)) rc=ESMF_SUCCESS
!------------------------------------------------------------------------
! get global vm information
!
call ESMF_VMGetCurrent(vm, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! set up local pet info
call ESMF_VMGet(vm, localPet=PetNo, petCount=PetCnt, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! calculate totalDE based on the decomposition
tiles=size(countsPerDEDim1PTile,2)
if (tiles /= 6) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- the second dimension of countsPerDEDim1PTile is not equal to 6", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
if (size(countsPerDEDim2PTile,2) /= 6) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- the second dimension of countsPerDEDim2PTile is not equal to 6", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
allocate(DeDim1(tiles), DeDim2(tiles))
do j=1,tiles
totalelmt = 0
DeDim1(j)=size(countsPerDEDim1Ptile,1)
DeDim2(j)=size(countsPerDEDim2Ptile,1)
do i=1,size(countsPerDEDim1PTile,1)
! check the total elements counts in dimension 1 is equal to tilesize
! count how many DEs for this tile
totalelmt = countsPerDEDim1PTile(i,j)+totalelmt
if (countsPerDEDim1PTile(i,j)==0) then
DeDim1(j)=i-1
exit
endif
enddo
if (totalelmt /= tilesize) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- the total number of elements in dimension 1 does not add up to tilesize", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
totalelmt = 0
do i=1,size(countsPerDEDim2PTile,1)
! check the total elements counts in dimension 1 is equal to tilesize
! count how many DEs for this tile
totalelmt = countsPerDEDim2PTile(i,j)+totalelmt
if (countsPerDEDim2PTile(i,j)==0) then
DeDim2(j)=i-1
exit
endif
enddo
if (totalelmt /= tilesize) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- the total number of elements in dimension 2 does not add up to tilesize", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
enddo
!if (PetNo == 0) then
! print *, 'DeDim: ', DeDim1(:), DeDim2(:)
!endif
! calculate totalDE
totalDE=0
do j=1,tiles
totalDE = totalDE+DeDim1(j)*DeDim2(j)
enddo
if (present(delayout)) then
!Check if delayout has the same number of DEs
call ESMF_DELayoutGet(delayout, deCount=deCount, rc=localrc)
if (deCount /= totalDE) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- the total number of DEs specified in delayout is inconsistent with the decomposition arguments", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
defaultDELayout = delayout
else
defaultDELayout = ESMF_DELayoutCreate(deCount = totalDE, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
if (present(indexflag)) then
localIndexFlag = indexflag
else
localIndexFlag = ESMF_INDEX_DELOCAL
endif
! Set Default coordSys
if (present(coordSys)) then
if (coordSys .eq. ESMF_COORDSYS_CART) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- only ESMF_TYPEKIND_CART is not supported", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
coordSysLocal=coordSys
else
coordSysLocal=ESMF_COORDSYS_SPH_DEG
endif
! Set Default coordTypeKind
if (present(coordTypeKind)) then
if (coordTypeKind .ne. ESMF_TYPEKIND_R4 .and. &
coordTypeKind .ne. ESMF_TYPEKIND_R8) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- only ESMF_TYPEKIND_R4 and ESMF_TYPEKIND_R8 are allowed", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
coordTypeKindLocal=coordTypeKind
else
coordTypeKindLocal=ESMF_TYPEKIND_R8
endif
! set defaults
docenter = .false.
docorner = .false.
tileCount = 6
allocate(minIndexPTile(2,tileCount))
allocate(maxIndexPTile(2,tileCount))
allocate(connectionList(12))
call CalculateConnection(tilesize, minIndexPTile, maxIndexPTile, connectionList, localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_DELayoutGet(defaultDelayout, DeCount = decount, localDeCount=localDeCount, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(demap(localdecount))
call ESMF_DELayoutGet(defaultDelayout, localDeToDeMap = demap, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! deBlockList and deToTileMap contains all the blocks for all the DEs (not
! just the localDEs
allocate(deBlockList(2,2,decount),deToTileMap(decount))
!print *, PetNo, 'total DE count ', decount
! minIndexPTile and maxIndexPTile are in ESMF_INDEX_GLOBAL, therefore, need
! to use global index in deBlockList as well
k=1
do t=1,tiles
do j=1,DeDim2(t)
do i=1,DeDim1(t)
minIndx = sum(countsPerDEDim1PTile(1:i-1,t))+minIndexPTile(1,t)
minIndy = sum(countsPerDEDim2PTile(1:j-1,t))+minIndexPTile(2,t)
deBlockList(1,1,k)=minIndx
deBlockList(2,1,k)=minIndy
deBlockList(1,2,k)=minIndx+countsPerDEDim1PTile(i,t)-1
deBlockList(2,2,k)=minIndy+countsPerDEDim2PTile(j,t)-1
deToTileMap(k)=t
k=k+1
enddo
enddo
enddo
! need to constrcut deBlockList and deToTileMap
distgrid = ESMF_DistGridCreate(&
minIndexPTile=minIndexPTile, maxIndexPTile=maxIndexPTile, &
deBlockList = deBlockList, deToTileMap = deToTileMap, &
indexflag=ESMF_INDEX_GLOBAL, connectionList=connectionList, &
deLabelList = deLabelList, &
delayout = defaultDelayout, &
rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! - create Grid
grid = ESMF_GridCreate(distgrid, coordSys=coordSysLocal, &
coordTypeKind=coordTypeKindLocal, &
indexflag=localIndexFlag, name=name, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (present(staggerLocList)) then
do s=1, size(staggerLocList)
if (staggerLocList(s) == ESMF_STAGGERLOC_EDGE1 .or. &
staggerLocList(s) == ESMF_STAGGERLOC_EDGE2) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="- only ESMF_STAGGERLOC_CENTER and ESMF_STAGGERLOC_CORNER are supported", &
ESMF_CONTEXT, rcToReturn=rc)
return
elseif (staggerLocList(s) == ESMF_STAGGERLOC_CENTER) then
docenter = .TRUE.
elseif (staggerLocList(s) == ESMF_STAGGERLOC_CORNER) then
docorner = .TRUE.
endif
call ESMF_GridAddCoord(grid, staggerloc=staggerLocList(s), rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
enddo
! calculate the actual cubed sphere coordiantes for each DE
do i = 1,localdecount
j = demap(i)+1
localDe = i-1
start(1)=deBlockList(1,1,j)-minIndexPTile(1,deToTileMap(j))+1
start(2)=deBlockList(2,1,j)-minIndexPTile(2,deToTileMap(j))+1
count(1)=deBlockList(1,2,j)-deBlockList(1,1,j)+1
count(2)=deBlockList(2,2,j)-deBlockList(2,1,j)+1
tile = deToTileMap(j)
if (coordTypeKindLocal .eq. ESMF_TYPEKIND_R8) then
if (docenter) then
call ESMF_GridGetCoord(grid, coordDim=1, localDe=localDe, &
farrayPtr=lonPtrR8, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridGetCoord(grid, coordDim=2, localDe=localDe, &
farrayPtr=latPtrR8, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
if (docorner) then
call ESMF_GridGetCoord(grid, coordDim=1, localDe=localDe, &
staggerloc=ESMF_STAGGERLOC_CORNER, farrayPtr=lonCornerPtrR8, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridGetCoord(grid, coordDim=2, localDe=localDe, &
staggerloc=ESMF_STAGGERLOC_CORNER, farrayPtr=latCornerPtrR8, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
else ! ESMF_TYPEKIND_R4
if (docenter) then
call ESMF_GridGetCoord(grid, coordDim=1, localDe=localDe, &
farrayPtr=lonPtrR4, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridGetCoord(grid, coordDim=2, localDe=localDe, &
farrayPtr=latPtrR4, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(lonPtrR8(count(1), count(2)), latPtrR8(count(1), count(2)))
endif
if (docorner) then
call ESMF_GridGetCoord(grid, coordDim=1, localDe=localDe, &
staggerloc=ESMF_STAGGERLOC_CORNER, farrayPtr=lonCornerPtrR4, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridGetCoord(grid, coordDim=2, localDe=localDe, &
staggerloc=ESMF_STAGGERLOC_CORNER, farrayPtr=latCornerPtrR4, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
shapLon=shape(lonCornerPtrR4) ! make sure lhs and rhs is same shape
shapLat=shape(latCornerPtrR4) ! make sure lhs and rhs is same shape
allocate(lonCornerPtrR8(shapLon(1), shapLon(2)), &
latCornerPtrR8(shapLat(1),shapLat(2)))
endif
endif
!call ESMF_VMWtime(starttime, rc=localrc)
! Generate glocal edge coordinates and local center coordinates
if (docenter .and. docorner) then
call ESMF_UtilCreateCSCoordsPar(tileSize, lonEdge=lonCornerPtrR8, &
latEdge=latCornerPtrR8, start=start, count=count, &
tile=tile, lonCenter=lonPtrR8, latCenter=latPtrR8, &
schmidtTransform=transformArgs)
elseif (docorner) then
call ESMF_UtilCreateCSCoordsPar(tileSize, lonEdge=lonCornerPtrR8, &
latEdge=latCornerPtrR8, start=start, count=count, tile=tile, &
schmidtTransform=transformArgs)
else
call ESMF_UtilCreateCSCoordsPar(tileSize, &
start=start, count=count, &
tile=tile, lonCenter=lonPtrR8, latCenter=latPtrR8, &
schmidtTransform=transformArgs)
endif
!call ESMF_VMWtime(endtime, rc=localrc)
if (coordTypeKindLocal .eq. ESMF_TYPEKIND_R8) then
if (docenter) then
if (coordSysLocal .eq. ESMF_COORDSYS_SPH_DEG) then
lonPtrR8 = lonPtrR8 * ESMF_COORDSYS_RAD2DEG
latPtrR8 = latPtrR8 * ESMF_COORDSYS_RAD2DEG
endif
endif
if (docorner) then
if (coordSysLocal .eq. ESMF_COORDSYS_SPH_DEG) then
lonCornerPtrR8 = lonCornerPtrR8 * ESMF_COORDSYS_RAD2DEG
latCornerPtrR8 = latCornerPtrR8 * ESMF_COORDSYS_RAD2DEG
endif
endif
else ! ESMF_TYPE_KIND_R4
if (docenter) then
if (coordSysLocal .eq. ESMF_COORDSYS_SPH_DEG) then
lonPtrR4 = lonPtrR8 * ESMF_COORDSYS_RAD2DEG
latPtrR4 = latPtrR8 * ESMF_COORDSYS_RAD2DEG
else
lonPtrR4 = lonPtrR8
latPtrR4 = latPtrR8
endif
deallocate(lonPtrR8, latPtrR8)
endif
if (docorner) then
if (coordSysLocal .eq. ESMF_COORDSYS_SPH_DEG) then
lonCornerPtrR4 = lonCornerPtrR8 * ESMF_COORDSYS_RAD2DEG
latCornerPtrR4 = latCornerPtrR8 * ESMF_COORDSYS_RAD2DEG
else
lonCornerPtrR4 = lonCornerPtrR8
latCornerPtrR4 = latCornerPtrR8
endif
deallocate(lonCornerPtrR8, latCornerPtrR8)
endif
endif
!print *, 'Create CS size ', tileSize, 'in', (endtime-starttime)*1000.0, ' msecs'
end do
endif
ESMF_GridCreateCubedSphereIReg = grid
! - deallocate connectionList
deallocate(minIndexPTile, maxIndexPTile)
deallocate(connectionList)
deallocate(deBlockList,deToTileMap)
deallocate(DeDim1, DeDim2, demap)
return
end function ESMF_GridCreateCubedSphereIReg