MAPL_DistGridGet Subroutine

public subroutine MAPL_DistGridGet(distGrid, MinIndex, MaxIndex, rc)

Arguments

Type IntentOptional Attributes Name
type(ESMF_DistGrid), intent(inout) :: distGrid
integer, intent(inout) :: MinIndex(:,:)
integer, intent(inout) :: MaxIndex(:,:)
integer, intent(out), optional :: rc

Calls

proc~~mapl_distgridget~~CallsGraph proc~mapl_distgridget MAPL_DistGridGet ESMF_DistGridGet ESMF_DistGridGet proc~mapl_distgridget->ESMF_DistGridGet proc~mapl_verify MAPL_Verify proc~mapl_distgridget->proc~mapl_verify proc~mapl_throw_exception MAPL_throw_exception proc~mapl_verify->proc~mapl_throw_exception

Called by

proc~~mapl_distgridget~~CalledByGraph proc~mapl_distgridget MAPL_DistGridGet none~set~92 MaplGrid%set none~set~92->proc~mapl_distgridget proc~create_cubed_sphere_grid RegridSupport%create_cubed_sphere_grid proc~create_cubed_sphere_grid->proc~mapl_distgridget proc~get_esmf_grid_layout get_esmf_grid_layout proc~get_esmf_grid_layout->proc~mapl_distgridget proc~mapl_createrequest MAPL_CreateRequest proc~mapl_createrequest->proc~mapl_distgridget proc~mapl_genericinitialize MAPL_GenericInitialize proc~mapl_genericinitialize->proc~mapl_distgridget proc~mapl_gridget MAPL_GridGet proc~mapl_gridget->proc~mapl_distgridget proc~mapl_tilemaskget MAPL_TileMaskGet proc~mapl_tilemaskget->proc~mapl_distgridget

Source Code

  subroutine MAPL_DistGridGet(distGrid,minIndex,maxIndex,rc)

     type(ESMF_DistGrid), intent(inout) :: distGrid
     integer,             intent(inout) :: MinIndex(:,:)
     integer,             intent(inout) :: MaxIndex(:,:)
     integer, optional,   intent(out  ) :: rc

     integer :: status

     integer :: i,tileSize,tileCount,tile,deCount
     logical :: ESMFCubeSphere
     integer, allocatable  :: elementCountPTile(:)
     integer, allocatable :: deToTileMap(:)
     integer, allocatable :: oldMinIndex(:,:),oldMaxIndex(:,:)

     ESMFCubeSphere = .false.

     call ESMF_DistGridGet(distGrid,tileCount=tileCount,_RC)

     if (tileCount==6) ESMFCubeSphere = .true.

     if (ESMFCubeSphere) then
        allocate(elementCountPTile(tileCount),__STAT__)
        call ESMF_DistGridGet(distGrid,elementCountPTile=elementCountPTile,_RC)
        ! All tile should have same number of elements
        tileSize = elementCountPTile(1)
        tileSize = SQRT(real(tileSize))
        deallocate(elementCountPTile)

        deCount = size(minIndex,2)

        allocate(deToTileMap(deCount),__STAT__)
        allocate(oldMinIndex(2,deCount),oldMaxIndex(2,deCount),__STAT__)
        call ESMF_DistGridGet(distGrid,MaxIndexPDe=oldMaxIndex,MinIndexPDe=oldMinIndex, &
                              deToTileMap=deToTileMap,_RC)
        do i=1,deCount
           tile = deToTileMap(i)
           select case (tile)
              case (1)
                 minIndex(:,i)=oldMinIndex(:,i)
                 maxIndex(:,i)=oldMaxIndex(:,i)
              case (2)
                 minIndex(1,i)=oldMinIndex(1,i) -  tileSize
                 minIndex(2,i)=oldMinIndex(2,i) +  tileSize
                 maxIndex(1,i)=oldMaxIndex(1,i) -  tileSize
                 maxIndex(2,i)=oldMaxIndex(2,i) +  tileSize
              case (3)
                 minIndex(1,i)=oldMinIndex(1,i) -  tileSize
                 minIndex(2,i)=oldMinIndex(2,i) +  tileSize
                 maxIndex(1,i)=oldMaxIndex(1,i) -  tileSize
                 maxIndex(2,i)=oldMaxIndex(2,i) +  tileSize
              case (4)
                 minIndex(1,i)=oldMinIndex(1,i) -2*tileSize
                 minIndex(2,i)=oldMinIndex(2,i) +2*tileSize
                 maxIndex(1,i)=oldMaxIndex(1,i) -2*tileSize
                 maxIndex(2,i)=oldMaxIndex(2,i) +2*tileSize
              case (5)
                 minIndex(1,i)=oldMinIndex(1,i) -2*tileSize
                 minIndex(2,i)=oldMinIndex(2,i) +2*tileSize
                 maxIndex(1,i)=oldMaxIndex(1,i) -2*tileSize
                 maxIndex(2,i)=oldMaxIndex(2,i) +2*tileSize
              case (6)
                 minIndex(1,i)=oldMinIndex(1,i) -3*tileSize
                 minIndex(2,i)=oldMinIndex(2,i) +3*tileSize
                 maxIndex(1,i)=oldMaxIndex(1,i) -3*tileSize
                 maxIndex(2,i)=oldMaxIndex(2,i) +3*tileSize
           end select
        enddo
        deallocate(deToTileMap)
        deallocate(oldMaxIndex,oldMinIndex)

     else

        call ESMF_DistGridGet(distGrid,minIndexPDe=minIndex,maxIndexPDe=maxIndex,_RC)

     end if

 end subroutine MAPL_DistGridGet