test_regrid2TileDG Subroutine

subroutine test_regrid2TileDG(itrp, csrv, rc)

Arguments

Type IntentOptional Attributes Name
logical, intent(out) :: itrp
logical, intent(out) :: csrv
integer, intent(out) :: rc

Source Code

 subroutine test_regrid2TileDG(itrp,csrv,rc)
  logical, intent(out)  :: itrp
  logical, intent(out)  :: csrv
  integer, intent(out)  :: rc
  logical :: correct
  integer :: localrc
  type(ESMF_Grid) :: srcGrid
  type(ESMF_Grid) :: dstGrid
  type(ESMF_Field) :: srcField, errorField
  type(ESMF_Field) :: dstField
  type(ESMF_Field) :: xdstField
  type(ESMF_Array) :: dstArray
  type(ESMF_Array) :: srcArray
  type(ESMF_Field) :: dstFracField
  type(ESMF_Field) :: srcFracField
  type(ESMF_Field) :: srcArea, dstArea
  type(ESMF_RouteHandle) :: routeHandle
  type(ESMF_ArraySpec) :: arrayspec
  type(ESMF_VM) :: vm
  type(ESMF_DistGrid) :: srcDistgrid
  type(ESMF_DELayout) :: delayout
  integer(ESMF_KIND_I4), pointer :: maskB(:,:), maskA(:,:)
  real(ESMF_KIND_R8), pointer :: srcAreaptr(:,:), dstAreaptr(:,:)
  real(ESMF_KIND_R8), pointer :: srcFracptr(:,:), dstFracptr(:,:)
  real(ESMF_KIND_R8), pointer :: farrayPtrXC(:,:)
  real(ESMF_KIND_R8), pointer :: farrayPtrYC(:,:)
  real(ESMF_KIND_R8), pointer :: farrayPtr(:,:), farrayPtr2(:,:)
  real(ESMF_KIND_R8), pointer :: xfarrayPtr(:,:)
  real(ESMF_KIND_R8), pointer :: fptr(:,:),xfptr(:,:),errorfptr(:,:)
  integer :: clbnd(2),cubnd(2)
  integer :: fclbnd(2),fcubnd(2)
  integer :: i1,i2,i3, index(2)
  integer :: lDE, srclocalDECount, dstlocalDECount
  real(ESMF_KIND_R8) :: coord(2)
   character(len=ESMF_MAXSTR) :: string

  integer :: src_cntr_nx(2), src_cntr_ny(2)
  real(ESMF_KIND_R8) :: src_cnr_minx(2), src_cnr_miny(2)
  real(ESMF_KIND_R8) :: src_cnr_maxx(2), src_cnr_maxy(2)

  integer :: dst_cntr_nx, dst_cntr_ny
  integer :: dst_cnr_nx, dst_cnr_ny
  real(ESMF_KIND_R8) :: dst_cnr_minx, dst_cnr_miny
  real(ESMF_KIND_R8) :: dst_cnr_maxx, dst_cnr_maxy

  integer :: tile
  integer :: tile_cnr_nx, tile_cnr_ny
  integer :: num_arrays

  real(ESMF_KIND_R8) :: tile_cnr_minx, tile_cnr_miny
  real(ESMF_KIND_R8) :: tile_cnr_maxx, tile_cnr_maxy
  real(ESMF_KIND_R8) :: cnr_x, cnr_y
  real(ESMF_KIND_R8) :: cnr_p1_x, cnr_p1_y

  real(ESMF_KIND_R8) :: dx,dy
  real(ESMF_KIND_R8) :: src_dx, src_dy
  real(ESMF_KIND_R8) :: dst_dx, dst_dy

  real(ESMF_KIND_R8) :: lon, lat, theta, phi, DEG2RAD
  real(ESMF_KIND_R8) :: half_dx, half_dy
   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

  integer :: spherical_grid
  integer :: localPet, petCount


  type(ESMF_DistGridConnection) :: connectionList(1)
!  type(ESMF_DistGridConnection) :: connectionList(2)
  integer :: minIndex(2,2), maxIndex(2,2)
  integer :: regDecomp(2,2)
  type(ESMF_Decomp_Flag) :: decomp(2,2)

  integer :: deCount
  integer, allocatable :: localDEtoDEMap(:), deToTileMap(:)


  ! result code
   integer :: finalrc
  
  ! 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

  ! Set src coordinates and resolution
  !! Src tile 1
  src_cntr_nx(1) = 20
  src_cntr_ny(1) = 10

  src_cnr_minx(1) = 0.0
  src_cnr_miny(1) = 0.0
  
  src_cnr_maxx(1) = 10.0
  src_cnr_maxy(1) = 10.0

  !! Src tile 2
  src_cntr_nx(2) = 20
  src_cntr_ny(2) = 10

  src_cnr_minx(2) = 10.0
  src_cnr_miny(2) = 0.0
  
  src_cnr_maxx(2) = 20.0
  src_cnr_maxy(2) = 10.0

  ! Set dst coordinates and resolution
  ! dst grid is set so that it fits entirely within src grid 
  dst_cntr_nx = 20
  dst_cntr_ny = 20

  dst_cnr_nx = dst_cntr_nx+1
  dst_cnr_ny = dst_cntr_ny+1

  dst_cnr_minx = 0.5
  dst_cnr_miny = 0.5
  
  dst_cnr_maxx = 19.5
  dst_cnr_maxy = 9.5


  ! Create connectionList
  ! periodicity
  call ESMF_DistgridConnectionSet(connection=connectionList(1), &
       tileIndexA=1,tileIndexB=2, &
       positionVector=(/-src_cntr_nx(1),0/), &
       rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return


  ! Setup index space
  minIndex(:,1)=(/1,1/)
  maxIndex(:,1)=(/src_cntr_nx(1),src_cntr_ny(1)/)
  regDecomp(:,1)=(/petCount,1/)

  minIndex(:,2)=(/1,1/)
  maxIndex(:,2)=(/src_cntr_nx(2),src_cntr_ny(2)/)
  regDecomp(:,2)=(/petCount,1/)
  
  decomp(:,:)=ESMF_DECOMP_BALANCED

  ! Create source distgrid
  srcDistgrid=ESMF_DistgridCreate(minIndexPTile=minIndex, maxIndexPTile=maxIndex, regDecompPTile=regDecomp, &
       decompflagPTile=decomp,indexflag=ESMF_INDEX_GLOBAL, &
       connectionList=connectionList,                 &
       rc=localrc)
  if (ESMF_LogFoundError(localrc, &
       ESMF_ERR_PASSTHRU, &
       ESMF_CONTEXT, rcToReturn=rc)) return
  
  
  ! setup source grid
  srcGrid=ESMF_GridCreate(distgrid=srcDistgrid, indexflag=ESMF_INDEX_GLOBAL, &
       coordSys=ESMF_COORDSYS_CART, &
       rc=localrc)
  if (ESMF_LogFoundError(localrc, &
       ESMF_ERR_PASSTHRU, &
       ESMF_CONTEXT, rcToReturn=rc)) return
  
  
  ! setup dest. grid
  dstGrid=ESMF_GridCreateNoPeriDim(minIndex=(/1,1/),maxIndex=(/dst_cntr_nx,dst_cntr_ny/),regDecomp=(/1,petCount/), &
       coordSys=ESMF_COORDSYS_CART, indexflag=ESMF_INDEX_GLOBAL, &
       rc=localrc)
  if (ESMF_LogFoundError(localrc, &
       ESMF_ERR_PASSTHRU, &
       ESMF_CONTEXT, rcToReturn=rc)) return
  
  ! Create source/destination fields
  call ESMF_ArraySpecSet(arrayspec, 2, ESMF_TYPEKIND_R8, rc=localrc)
  if (ESMF_LogFoundError(localrc, &
       ESMF_ERR_PASSTHRU, &
       ESMF_CONTEXT, rcToReturn=rc)) return
  
  
  srcField = ESMF_FieldCreate(srcGrid, arrayspec, &
       staggerloc=ESMF_STAGGERLOC_CENTER, name="source", rc=localrc)
  if (ESMF_LogFoundError(localrc, &
       ESMF_ERR_PASSTHRU, &
       ESMF_CONTEXT, rcToReturn=rc)) return
  
  
  dstField = ESMF_FieldCreate(dstGrid, arrayspec, &
       staggerloc=ESMF_STAGGERLOC_CENTER, name="dest", rc=localrc)
  if (ESMF_LogFoundError(localrc, &
       ESMF_ERR_PASSTHRU, &
       ESMF_CONTEXT, rcToReturn=rc)) return
  
  xdstField = ESMF_FieldCreate(dstGrid, arrayspec, &
                         staggerloc=ESMF_STAGGERLOC_CENTER, name="xdest", rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return

   errorField = ESMF_FieldCreate(dstGrid, arrayspec, &
                         staggerloc=ESMF_STAGGERLOC_CENTER, name="dest", rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return

   srcFracField = ESMF_FieldCreate(srcGrid, arrayspec, &
                         staggerloc=ESMF_STAGGERLOC_CENTER, name="dest", rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return


   srcArea = ESMF_FieldCreate(srcGrid, arrayspec, &
                         staggerloc=ESMF_STAGGERLOC_CENTER, name="source", rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return


   dstFracField = ESMF_FieldCreate(dstGrid, arrayspec, &
                         staggerloc=ESMF_STAGGERLOC_CENTER, name="dest", rc=localrc)
   if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return

   dstArea = ESMF_FieldCreate(dstGrid, arrayspec, &
                         staggerloc=ESMF_STAGGERLOC_CENTER, name="dest", rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return


  ! Allocate coordinates
  call ESMF_GridAddCoord(srcGrid, staggerloc=ESMF_STAGGERLOC_CENTER, rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return

 call ESMF_GridAddCoord(srcGrid, staggerloc=ESMF_STAGGERLOC_CORNER, rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return


  call ESMF_GridAddCoord(dstGrid, staggerloc=ESMF_STAGGERLOC_CENTER, rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return

  call ESMF_GridAddCoord(dstGrid, staggerloc=ESMF_STAGGERLOC_CORNER, rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return


  ! Get arrays
  ! dstArray
  call ESMF_FieldGet(dstField, array=dstArray, rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return


  ! srcArray
  call ESMF_FieldGet(srcField, array=srcArray, rc=localrc)
  if (ESMF_LogFoundError(localrc, &
       ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return
 

  ! Get information for going from localDE to Tile
  call ESMF_GridGet(srcGrid, localDECount=srclocalDECount, rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return
 
  call ESMF_DistgridGet(srcDistgrid, delayout=delayout, rc=localrc)
   if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return

  allocate(localDEtoDEMap(srclocalDECount))  

  call ESMF_DELayoutGet(delayout, deCount=deCount, localDeToDeMap=localDeToDEMap, &
                        rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return

  allocate(deToTileMap(deCount))  


  call ESMF_DistgridGet(srcDistgrid, deToTileMap=detoTileMap, rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return

  ! Construct Src Grid
  !!!!  CORNERS !!!!
  do lDE=0,srclocalDECount-1
 
     !! get coord 1
     call ESMF_GridGetCoord(srcGrid, localDE=lDE, staggerLoc=ESMF_STAGGERLOC_CORNER, coordDim=1, &
                            computationalLBound=clbnd, computationalUBound=cubnd, farrayPtr=farrayPtrXC, rc=localrc)
     if (ESMF_LogFoundError(localrc, &
         ESMF_ERR_PASSTHRU, &
         ESMF_CONTEXT, rcToReturn=rc)) return

     call ESMF_GridGetCoord(srcGrid, localDE=lDE, staggerLoc=ESMF_STAGGERLOC_CORNER, coordDim=2, &
                            computationalLBound=clbnd, computationalUBound=cubnd, farrayPtr=farrayPtrYC, rc=localrc)
     if (ESMF_LogFoundError(localrc, &
         ESMF_ERR_PASSTHRU, &
         ESMF_CONTEXT, rcToReturn=rc)) return

     !! Get Tile from localDE
     tile=deToTileMap(localDEtoDEMap(lDE+1)+1)  

     !! Set values based on tile
     tile_cnr_nx=src_cntr_nx(tile)+1 ! corners are 1 bigger than centers
     tile_cnr_ny=src_cntr_ny(tile)+1 ! corners are 1 bigger than centers

     tile_cnr_minx=src_cnr_minx(tile)
     tile_cnr_maxx=src_cnr_maxx(tile)

     tile_cnr_miny=src_cnr_miny(tile)
     tile_cnr_maxy=src_cnr_maxy(tile)

    ! write(*,*) "t=",tile," cnr bnds l=",clbnd(1),clbnd(2)," u=",cubnd(1),cubnd(2)
    ! write(*,*) "t=",tile," nx=",tile_cnr_nx," ny=",tile_cnr_ny


     !! set coords, interpolated function
     do i1=clbnd(1),cubnd(1)
     do i2=clbnd(2),cubnd(2)

        ! Set source coordinates
        farrayPtrXC(i1,i2) = ((tile_cnr_maxx-tile_cnr_minx)*REAL(i1-1)/REAL(tile_cnr_nx-1))+tile_cnr_minx
        farrayPtrYC(i1,i2) = ((tile_cnr_maxy-tile_cnr_miny)*REAL(i2-1)/REAL(tile_cnr_ny-1))+tile_cnr_miny

     enddo
     enddo

  enddo    ! lDE


  ! Construct Src Grid
  !!!!  CENTERS !!!!
  do lDE=0,srclocalDECount-1
 
     !! get coord 1
     call ESMF_GridGetCoord(srcGrid, localDE=lDE, staggerLoc=ESMF_STAGGERLOC_CENTER, coordDim=1, &
                            computationalLBound=clbnd, computationalUBound=cubnd, farrayPtr=farrayPtrXC, rc=localrc)
     if (ESMF_LogFoundError(localrc, &
         ESMF_ERR_PASSTHRU, &
         ESMF_CONTEXT, rcToReturn=rc)) return

     call ESMF_GridGetCoord(srcGrid, localDE=lDE, staggerLoc=ESMF_STAGGERLOC_CENTER, coordDim=2, &
                            computationalLBound=clbnd, computationalUBound=cubnd, farrayPtr=farrayPtrYC, rc=localrc)
     if (ESMF_LogFoundError(localrc, &
         ESMF_ERR_PASSTHRU, &
         ESMF_CONTEXT, rcToReturn=rc)) return


     ! get src pointer
     call ESMF_FieldGet(srcField, lDE, farrayPtr,  rc=localrc)
     if (ESMF_LogFoundError(localrc, &
         ESMF_ERR_PASSTHRU, &
         ESMF_CONTEXT, rcToReturn=rc)) return

 
      !! Get Tile from localDE
      tile=deToTileMap(localDEtoDEMap(lDE+1)+1)  


     !! Set values based on tile
     tile_cnr_nx=src_cntr_nx(tile)+1 ! corners 1 bigger than center
     tile_cnr_ny=src_cntr_ny(tile)+1 ! corners 1 bigger than center

     tile_cnr_minx=src_cnr_minx(tile)
     tile_cnr_maxx=src_cnr_maxx(tile)
 
     tile_cnr_miny=src_cnr_miny(tile)
     tile_cnr_maxy=src_cnr_maxy(tile)

 
     !! set coords, interpolated function
     do i1=clbnd(1),cubnd(1)
     do i2=clbnd(2),cubnd(2)

        ! corner coordinate
        cnr_x = ((tile_cnr_maxx-tile_cnr_minx)*REAL(i1-1)/REAL(tile_cnr_nx-1))+tile_cnr_minx
        cnr_y = ((tile_cnr_maxy-tile_cnr_miny)*REAL(i2-1)/REAL(tile_cnr_ny-1))+tile_cnr_miny

        ! corner +1 coordinate
        cnr_p1_x = ((tile_cnr_maxx-tile_cnr_minx)*REAL(i1+1-1)/REAL(tile_cnr_nx-1))+tile_cnr_minx
        cnr_p1_y = ((tile_cnr_maxy-tile_cnr_miny)*REAL(i2+1-1)/REAL(tile_cnr_ny-1))+tile_cnr_miny


        ! Set source coordinates
        farrayPtrXC(i1,i2) = 0.5*(cnr_x+cnr_p1_x)
        farrayPtrYC(i1,i2) = 0.5*(cnr_y+cnr_p1_y)

        ! set src data
        farrayPtr(i1,i2) = farrayPtrXC(i1,i2) + farrayPtrYC(i1,i2) + 20.0
     enddo
     enddo

  enddo    ! lDE




  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
  ! Destination grid
  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

  ! Get number of local DEs
  call ESMF_GridGet(dstGrid, localDECount=dstlocalDECount, rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return

  ! Get memory and set coords for dst
  !! CORNERS !!
  do lDE=0,dstlocalDECount-1
 
     !! get coords
     call ESMF_GridGetCoord(dstGrid, localDE=lDE, staggerLoc=ESMF_STAGGERLOC_CORNER, coordDim=1, &
                            computationalLBound=clbnd, computationalUBound=cubnd, farrayPtr=farrayPtrXC, rc=localrc)
     if (ESMF_LogFoundError(localrc, &
         ESMF_ERR_PASSTHRU, &
         ESMF_CONTEXT, rcToReturn=rc)) return

     call ESMF_GridGetCoord(dstGrid, localDE=lDE, staggerLoc=ESMF_STAGGERLOC_CORNER, coordDim=2, &
                            computationalLBound=clbnd, computationalUBound=cubnd, farrayPtr=farrayPtrYC, rc=localrc)
     if (ESMF_LogFoundError(localrc, &
         ESMF_ERR_PASSTHRU, &
         ESMF_CONTEXT, rcToReturn=rc)) return

    ! write(*,*) "D: cnr bnds l=",clbnd(1),clbnd(2)," u=",cubnd(1),cubnd(2), " dcn=",dst_cnr_nx,dst_cnr_ny

     !! set coords
     do i1=clbnd(1),cubnd(1)
     do i2=clbnd(2),cubnd(2)
 
        ! set source coordinates
        farrayPtrXC(i1,i2) = ((dst_cnr_maxx-dst_cnr_minx)*REAL(i1-1)/REAL(dst_cnr_nx-1))+dst_cnr_minx
        farrayPtrYC(i1,i2) = ((dst_cnr_maxy-dst_cnr_miny)*REAL(i2-1)/REAL(dst_cnr_ny-1))+dst_cnr_miny

     enddo
     enddo

  enddo    ! lDE

  ! Get memory and set coords for dst
  !! CENTERS
  do lDE=0,dstlocalDECount-1
 
     !! get coords
     call ESMF_GridGetCoord(dstGrid, localDE=lDE, staggerLoc=ESMF_STAGGERLOC_CENTER, coordDim=1, &
                            computationalLBound=clbnd, computationalUBound=cubnd, farrayPtr=farrayPtrXC, rc=localrc)
     if (ESMF_LogFoundError(localrc, &
         ESMF_ERR_PASSTHRU, &
         ESMF_CONTEXT, rcToReturn=rc)) return

     call ESMF_GridGetCoord(dstGrid, localDE=lDE, staggerLoc=ESMF_STAGGERLOC_CENTER, coordDim=2, &
                            computationalLBound=clbnd, computationalUBound=cubnd, farrayPtr=farrayPtrYC, rc=localrc)
     if (ESMF_LogFoundError(localrc, &
         ESMF_ERR_PASSTHRU, &
         ESMF_CONTEXT, rcToReturn=rc)) return


     call ESMF_FieldGet(dstField, lDE, farrayPtr, computationalLBound=fclbnd, &
                             computationalUBound=fcubnd,  rc=localrc)
     if (ESMF_LogFoundError(localrc, &
         ESMF_ERR_PASSTHRU, &
         ESMF_CONTEXT, rcToReturn=rc)) return


     call ESMF_FieldGet(xdstField, lDE, xfarrayPtr,  rc=localrc)
     if (ESMF_LogFoundError(localrc, &
          ESMF_ERR_PASSTHRU, &
         ESMF_CONTEXT, rcToReturn=rc)) return


     !! set coords
     do i1=clbnd(1),cubnd(1)
     do i2=clbnd(2),cubnd(2)

        ! corner coordinate
        cnr_x = ((dst_cnr_maxx-dst_cnr_minx)*REAL(i1-1)/REAL(dst_cnr_nx-1))+dst_cnr_minx
        cnr_y = ((dst_cnr_maxy-dst_cnr_miny)*REAL(i2-1)/REAL(dst_cnr_ny-1))+dst_cnr_miny

        ! corner coordinate
        cnr_p1_x = ((dst_cnr_maxx-dst_cnr_minx)*REAL(i1+1-1)/REAL(dst_cnr_nx-1))+dst_cnr_minx
        cnr_p1_y = ((dst_cnr_maxy-dst_cnr_miny)*REAL(i2+1-1)/REAL(dst_cnr_ny-1))+dst_cnr_miny


        ! Set destination center coordinates
        farrayPtrXC(i1,i2) = 0.5*(cnr_x+cnr_p1_x)
        farrayPtrYC(i1,i2) = 0.5*(cnr_y+cnr_p1_y)

        ! set expected result field
        xfarrayPtr(i1,i2) = farrayPtrXC(i1,i2) + farrayPtrYC(i1,i2) + 20.0

        ! initialize destination field
        farrayPtr(i1,i2)=0.0

     enddo
     enddo
  enddo    ! lDE




#if 0
  call ESMF_GridWriteVTK(srcGrid,staggerloc=ESMF_STAGGERLOC_CENTER, &
        filename="srcGrid", &
       rc=localrc)
  if (ESMF_LogFoundError(localrc, &
       ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return

  call ESMF_GridWriteVTK(srcGrid,staggerloc=ESMF_STAGGERLOC_CORNER, &
        filename="srcGridCnr", &
       rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return

  call ESMF_GridWriteVTK(dstGrid,staggerloc=ESMF_STAGGERLOC_CENTER, &
        filename="dstGrid", &
       rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return

  call ESMF_GridWriteVTK(dstGrid,staggerloc=ESMF_STAGGERLOC_CORNER, &
        filename="dstGridCnr", &
       rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return
#endif


  !!! Regrid forward from the A grid to the B grid
  ! Regrid store
  call ESMF_FieldRegridStore( &
          srcField, &
          dstField=dstField, &
          dstFracField=dstFracField, &
          srcFracField=srcFracField, &
          unmappedaction=ESMF_UNMAPPEDACTION_IGNORE, &
          regridmethod=ESMF_REGRIDMETHOD_CONSERVE, &
          normType=ESMF_NORMTYPE_FRACAREA, &
          routeHandle=routeHandle, &
          rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
       ESMF_CONTEXT, rcToReturn=rc)) return

   ! Do regrid
  call ESMF_FieldRegrid(srcField, dstField, routeHandle, rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return

  call ESMF_FieldRegridRelease(routeHandle, rc=localrc)
  if (ESMF_LogFoundError(localrc, &
      ESMF_ERR_PASSTHRU, &
      ESMF_CONTEXT, rcToReturn=rc)) return


  ! Get the integration weights
   call ESMF_FieldRegridGetArea(srcArea, &
           rc=localrc)
  if (localrc /=ESMF_SUCCESS) then
      rc=ESMF_FAILURE
      return
  endif

  ! Get the integration weights
  call ESMF_FieldRegridGetArea(dstArea, &
          rc=localrc)
  if (localrc /=ESMF_SUCCESS) then
      rc=ESMF_FAILURE
      return
  endif


  ! Init
  minerror(1) = 100000.
  maxerror(1) = 0.
  error = 0.
  dstmass = 0.
 
  ! Check if the values are close
  do lDE=0,dstLocalDECount-1

     ! get dst Field
     call ESMF_FieldGet(dstField, lDE, fptr, computationalLBound=clbnd, &
                              computationalUBound=cubnd,  rc=localrc)
     if (localrc /=ESMF_SUCCESS) then
        rc=ESMF_FAILURE
        return
     endif

     ! get exact destination Field
     call ESMF_FieldGet(xdstField, lDE, xfptr,  rc=localrc)
     if (localrc /=ESMF_SUCCESS) then
        rc=ESMF_FAILURE
        return
     endif

     ! get error Field
      call ESMF_FieldGet(errorField, lDE, errorfptr,  rc=localrc)
     if (localrc /=ESMF_SUCCESS) then
        rc=ESMF_FAILURE
        return
     endif

     ! get dst area Field
     call ESMF_FieldGet(dstArea, lDE, dstAreaptr,  rc=localrc)
     if (localrc /=ESMF_SUCCESS) then
        rc=ESMF_FAILURE
        return
      endif

     ! get frac Field
     call ESMF_FieldGet(dstFracField, lDE, dstFracptr,  rc=localrc)
      if (localrc /=ESMF_SUCCESS) then
         rc=ESMF_FAILURE
        return
     endif

     ! destination grid
     !! check relative error
     do i1=clbnd(1),cubnd(1)
     do i2=clbnd(2),cubnd(2)

        ! This is WRONG, shouldn't include Frac
        dstmass = dstmass + dstFracptr(i1,i2)*dstAreaptr(i1,i2)*fptr(i1,i2)

        ! 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,i2) .lt. 0.999) cycle

        if (xfptr(i1,i2) .ne. 0.0) then
           errorfptr(i1,i2)=ABS(fptr(i1,i2) - xfptr(i1,i2))/ABS(xfptr(i1,i2))
           error = error + errorfptr(i1,i2)
           if (errorfptr(i1,i2) > maxerror(1)) then
             maxerror(1) = errorfptr(i1,i2)
           endif
           if (errorfptr(i1,i2) < minerror(1)) then
             minerror(1) = errorfptr(i1,i2)
           endif
        else
           errorfptr(i1,i2)=ABS(fptr(i1,i2) - xfptr(i1,i2))
           error = error + errorfptr(i1,i2)
           if (errorfptr(i1,i2) > maxerror(1)) then
             maxerror(1) = errorfptr(i1,i2)
           endif
           if (errorfptr(i1,i2) < minerror(1)) then
             minerror(1) = errorfptr(i1,i2)
           endif
        endif

     enddo
      enddo

  enddo    ! lDE


  srcmass(1) = 0.
  do lDE=0,srcLocalDECount-1

     ! get src pointer
     call ESMF_FieldGet(srcField, lDE, fptr, computationalLBound=clbnd, &
                             computationalUBound=cubnd,  rc=localrc)
     if (localrc /=ESMF_SUCCESS) then
        rc=ESMF_FAILURE
        return
     endif

     ! get src Field
     call ESMF_FieldGet(srcArea, lDE, srcAreaptr,  rc=localrc)
     if (localrc /=ESMF_SUCCESS) then
        rc=ESMF_FAILURE
        return
     endif

     ! get frac Field
     call ESMF_FieldGet(srcFracField, lDE, srcFracptr,  rc=localrc)
     if (localrc /=ESMF_SUCCESS) then
        rc=ESMF_FAILURE
        return
     endif

     do i1=clbnd(1),cubnd(1)
     do i2=clbnd(2),cubnd(2)
        srcmass(1) = srcmass(1) + srcFracptr(i1,i2)*srcAreaptr(i1,i2)*fptr(i1,i2)
     enddo
     enddo

  enddo    ! lDE


  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.

  ! Output results
  if (localPet == 0) then

    write(*,*) "=== 2 Tile ==="
    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

#if 0
  call ESMF_GridWriteVTK(srcGrid,staggerloc=ESMF_STAGGERLOC_CENTER, &
       filename="srcGridT", array1=srcArray, &
       rc=localrc)
  if (localrc /=ESMF_SUCCESS) then
     rc=ESMF_FAILURE
     return
  endif

  call ESMF_GridWriteVTK(dstGrid,staggerloc=ESMF_STAGGERLOC_CENTER, &
       filename="dstGridT", array1=dstArray, &
       rc=localrc)
  if (localrc /=ESMF_SUCCESS) then
     rc=ESMF_FAILURE
     return
  endif
#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(srcArea, 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(errorField, 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

   call ESMF_FieldDestroy(dstArea, 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


  ! Free the grids
  call ESMF_GridDestroy(srcGrid, rc=localrc)
  if (localrc /=ESMF_SUCCESS) then
      rc=ESMF_FAILURE
      return
   endif

  ! Free the srcDistgrid
   call ESMF_DistgridDestroy(srcDistgrid, rc=localrc)
  if (localrc /=ESMF_SUCCESS) then
      rc=ESMF_FAILURE
      return
   endif


  call ESMF_GridDestroy(dstGrid, 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_regrid2TileDG