subroutine test_regrid2xgSph(rc)
integer, intent(out) :: rc
type(ESMF_Grid) :: grid_atm, grid_ocn
type(ESMF_Field) :: f_atm, f_ocn
type(ESMF_XGrid) :: xgrid
real(ESMF_KIND_R8) :: atm_dx, atm_dy, ocn_dx, ocn_dy, startx, starty
real(ESMF_KIND_R8) :: atm_sx, atm_sy, ocn_sx, ocn_sy
integer :: atm_nx, atm_ny, ocn_nx, ocn_ny
integer :: localrc, npet, i, j, lpet
real(ESMF_KIND_R8), pointer :: weights(:)
integer(ESMF_KIND_I4), pointer :: indices(:,:)
real(ESMF_KIND_R8), pointer :: coordX(:,:), coordY(:,:)
type(ESMF_Grid) :: sideA(1), sideB(1)
type(ESMF_XGridSpec) :: A2X(1)
type(ESMF_VM) :: vm
localrc = ESMF_SUCCESS
call ESMF_VMGetCurrent(vm, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_VMGet(vm, petCount=npet, localPet=lpet, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Grid constants
! Atmosphere covers the area (-165,30) - (-15, 50)
! North Pacific Ocean covers (-165,30) - (-120, 50)
! Running the Grids from NW corner to SE corner,
! change the starting Y coord and dy
!
! 1 degree atmosphere, 1 degree ocean
atm_nx = 5
atm_ny = 20
atm_dx = 2.
atm_dy = 1.
ocn_nx = 5
ocn_ny = 20
ocn_dx = 1.
ocn_dy = 1.
atm_sx = -165.
atm_sy = 30.
ocn_sx = -165.
ocn_sy = 30.
!------------- ATM --------------
! atm grid, horizontally decomposed
grid_atm = ESMF_GridCreateNoPeriDim(maxIndex=(/atm_nx, atm_ny/), &
indexflag=ESMF_INDEX_GLOBAL, &
!gridEdgeLWidth=(/0,0/), gridEdgeUWidth=(/0,1/), &
regDecomp=(/npet, 1/), &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridAddCoord(grid_atm, staggerloc=ESMF_STAGGERLOC_CENTER, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridAddCoord(grid_atm, 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(grid_atm, 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(grid_atm, 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), 'coordX: ', coordX(:,1), coordX(1,:)
! X corner
call ESMF_GridGetCoord(grid_atm, 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(grid_atm, 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
!------------- OCN --------------
! ocn grid, horizontally decomposed
grid_ocn = ESMF_GridCreateNoPeriDim(maxIndex=(/ocn_nx, ocn_ny/), &
indexflag=ESMF_INDEX_GLOBAL, &
!gridEdgeLWidth=(/0,0/), gridEdgeUWidth=(/0,1/), &
regDecomp=(/npet, 1/), &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridAddCoord(grid_ocn, staggerloc=ESMF_STAGGERLOC_CENTER, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridAddCoord(grid_ocn, staggerloc=ESMF_STAGGERLOC_CORNER, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! ocn grid is not decomposed in the y direction
!startx = lpet*ocn_nx/npet*ocn_dx
startx = ocn_sx
starty = ocn_sy
! compute coord
! X center
call ESMF_GridGetCoord(grid_ocn, 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(grid_ocn, 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 + ocn_dx/2. + (i-1)*ocn_dx
coordY(i,j) = starty + ocn_dy/2. + (j-1)*ocn_dy
enddo
enddo
print *, 'startx: ', startx, lbound(coordX, 1), 'coordX: ', coordX(:,1), coordX(1,:)
! X corner
call ESMF_GridGetCoord(grid_ocn, 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(grid_ocn, 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)*ocn_dx
coordY(i,j) = starty + (j-1)*ocn_dy
enddo
enddo
! build Fields on the Grids
f_atm = ESMF_FieldCreate(grid_atm, typekind=ESMF_TYPEKIND_R8, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
f_ocn = ESMF_FieldCreate(grid_ocn, typekind=ESMF_TYPEKIND_R8, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Call RegridStore here to compute SMM weights and indices
call ESMF_FieldRegridStore(srcField=f_atm, dstField=f_ocn, &
regridmethod=ESMF_REGRIDMETHOD_CONSERVE, &
unmappedaction = ESMF_UNMAPPEDACTION_IGNORE, &
factorIndexList=indices, factorList=weights, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! clean up
call ESMF_FieldDestroy(f_atm, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_FieldDestroy(f_ocn, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridDestroy(grid_atm, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridDestroy(grid_ocn, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
print *, lpet, '-', associated(weights), '-', size(weights), '-', weights
print *, lpet, '-', associated(indices), '-', size(indices),'-', indices
sideA(1) = grid_ocn
sideB(1) = grid_atm
A2X(1)%factorIndexList = indices
A2X(1)%factorList = weights
xgrid = ESMF_XGridCreateFromSparseMat(sideAGrid=sideA, sideBGrid=sideB, &
sparseMatA2X=A2X, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_XGridDestroy(xgrid, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
rc = ESMF_SUCCESS
end subroutine test_regrid2xgSph