LCOV - code coverage report
Current view: | top level - dynamics/se/dycore - edge_mod.F90 (source / functions) |  |  | Hit | Total | Coverage
Test: | coverage.info |  | Lines: | 600 | 1145 | 52.4 %
Date: | 2024-12-17 17:57:11 |  | Functions: | 14 | 28 | 50.0 %

              Line data    Source code
           1             : module edge_mod
           2             :
           3             :   use shr_kind_mod,           only: r8=>shr_kind_r8, i8=>shr_kind_i8
           4             :   use dimensions_mod,         only: max_neigh_edges, nelemd
           5             :   use perf_mod,               only: t_startf, t_stopf, t_adj_detailf ! _EXTERNAL
           6             :   use thread_mod,             only: max_num_threads, omp_get_num_threads, omp_get_thread_num
           7             :   use coordinate_systems_mod, only: cartesian3D_t
           8             :   use schedtype_mod,          only: cycle_t, schedule_t, pgindex_t, schedule, HME_Ordinal,HME_Cardinal
           9             :   use cam_abortutils,         only: endrun
          10             :   use cam_logfile,            only: iulog
          11             :   use parallel_mod,           only: parallel_t, &
          12             :       MAX_ACTIVE_MSG, HME_status_size, BNDRY_TAG_BASE, HME_BNDRY_A2A, HME_BNDRY_P2P, &
          13             :       HME_BNDRY_A2AO
          14             :   use edgetype_mod,           only: edgedescriptor_t, edgebuffer_t,           &
          15             :       Longedgebuffer_t, initedgebuffer_callid, Ghostbuffer3D_t
          16             :   use element_mod,            only: element_t
          17             :   use gbarrier_mod,           only: gbarrier_init, gbarrier_delete
          18             :   use spmd_utils,             only: mpi_real8, mpi_integer, mpi_info_null, mpi_success
          19             :
          20             :   implicit none
          21             :   private
          22             :   save
          23             :
          24             :   ! 8-byte Integer routines
          25             :   public :: LongEdgeVpack, LongEdgeVunpackMIN
          26             :
          27             :   ! 8-byte Real routines
          28             :   public :: zeroEdgeBuffer
          29             :
          30             :   interface initEdgeBuffer
          31             :      module procedure initEdgeBuffer_r8
          32             :      module procedure initEdgeBuffer_i8
          33             :   end interface
          34             :   interface initEdgeSBuffer
          35             :      module procedure initEdgeSbuffer_r8
          36             :   end interface
          37             :   interface freeEdgeBuffer
          38             :      module procedure freeEdgeBuffer_r8
          39             :      module procedure freeEdgeBuffer_i8
          40             :   end interface
          41             :   interface freeGhostBuffer
          42             :      module procedure freeGhostBuffer_r8
          43             :   end interface
          44             :
          45             :   public :: initEdgeBuffer
          46             :   public :: initEdgeSBuffer
          47             :   public :: freeEdgeBuffer
          48             :
          49             :   public :: initGhostBuffer
          50             :   public :: ghostpack, ghostunpack
          51             :   public :: freeGhostBuffer
          52             :
          53             :   !---------------------------------------------------------
          54             :   ! Pack/unpack routines that use the New format Edge buffer
          55             :   !---------------------------------------------------------
          56             :
          57             :   public :: edgeVpack, edgeVunpack
          58             :   public :: edgeVunpackMIN, edgeVunpackMAX
          59             :   public :: edgeDGVpack, edgeDGVunpack
          60             :   public :: edgeVunpackVert
          61             :
          62             :
          63             :   public :: initGhostBuffer3D
          64             :   public :: FreeGhostBuffer3D
          65             :   public :: ghostVpack3D, ghostVunpack3D
          66             :
          67             :   !----------------------------------------------------------------
          68             :   ! Pack/unpack routines that communicate a fixed number values
          69             :   ! per element.  This is used to communicate MIN/MAX values from
          70             :   ! neighboring elemeents
          71             :   !----------------------------------------------------------------
          72             :   interface edgeSpack
          73             :      module procedure edgeSpack_r8
          74             :   end interface
          75             :   public :: edgeSpack
          76             :   public :: edgeSunpackMIN, edgeSunpackMAX
          77             :
          78             :   logical, private :: threadsafe=.true.
          79             :
          80             :   real(kind=r8), parameter, public :: edgeDefaultVal = 1.11e+100_r8
          81             :
          82             : ! NOTE ON ELEMENT ORIENTATION
          83             : !
          84             : ! Element orientation:  index V(i,j)
          85             : !
          86             : !           (1,np) NWEST      (np,np) NEAST
          87             : !
          88             : !           (1,1) SWEST       (np,1) SEAST
          89             : !
          90             : !
          91             : ! for the edge neighbors:
          92             : !    we set the "reverse" flag if two elements who share an edge use a
          93             : !    reverse orientation.  The data is reversed during the *pack* stage
          94             : ! For corner neighbors:
          95             : !    for edge buffers, there is no orientation because two corner neighbors
          96             : !    only share a single point.
          97             : !    For ghost cell data, there is a again two posible orientations. For
          98             : !    this case, we set the "reverse" flag if the corner element is using
          99             : !    the reverse orientation.  In this case, the data is reversed during the
         100             : !    *unpack* stage (not sure why)
         101             : !
         102             : ! The edge orientation is set at startup.  The corner orientation is computed
         103             : ! at run time, via the call to compute_ghost_corner_orientation()
         104             : ! This routine only works for meshes with at most 1 corner element.  It's
         105             : ! not called and the corner orientation flag is not set for unstructured meshes
         106             :
         107             : !
         108             : ! Christoph Erath
         109             : ! pack/unpack partial element of data of size (nx,nx) with user specifed halo size nh
         110             : ! user specifies the sizes when creating the buffer
         111             : ! buffer has 1 extra dimension (as compared to subroutines above) for multiple tracers
         112             : ! input/output arrays are cartesian, and thus assume at most 1 element at each corner
         113             : ! hence currently only supports cube-sphere grids.
         114             : !
         115             : !
         116             : ! routines which including element edge data
         117             : ! (used for FVM arrays where edge data is not shared by neighboring elements)
         118             : ! these routines pack/unpack element data with user specified halo size
         119             :
         120             :   ! Wrap pointer so we can make an array of them.
         121             :   type :: wrap_ptr
         122             :     real (kind=r8), dimension(:,:), pointer :: ptr => null()
         123             :   end type wrap_ptr
         124             :
         125             :   type(wrap_ptr) :: edgebuff_ptrs(0:1)
         126             :
         127             : contains
         128             :
         129           0 :   subroutine initEdgeSBuffer_r8(par,edge,elem,nlyr,bndry_type, nthreads)
         130             :     type (parallel_t),             intent(in)  :: par
         131             :     type (EdgeBuffer_t), target,   intent(out) :: edge
         132             :     type (element_t),              intent(in)  :: elem(:)
         133             :     integer,                       intent(in)  :: nlyr
         134             :     integer ,            optional, intent(in)  :: bndry_type
         135             :     integer,             optional, intent(in)  :: nthreads
         136             :
         137             :
         138             :     call initEdgeBuffer(par,edge,elem,nlyr,bndry_type=bndry_type, &
         139           0 :                          nthreads=nthreads,CardinalLength=1,OrdinalLength=1)
         140             :
         141           0 :   end subroutine initEdgeSBuffer_r8
         142             :
         143       18432 :   subroutine initGhostBuffer(par,edge,elem,nlyr,ndepth, npoints,bndry_type,nthreads)
         144             :
         145             :     type (parallel_t),             intent(in)  :: par
         146             :     type (Edgebuffer_t), target,   intent(out) :: edge
         147             :     type (element_t),              intent(in)  :: elem(:)
         148             :     integer,intent(in)                         :: nlyr,ndepth, npoints
         149             :     integer ,            optional, intent(in)  :: bndry_type
         150             :     integer,             optional, intent(in)  :: nthreads
         151             :
         152             :     call initEdgeBuffer(par,edge,elem,nlyr,bndry_type=bndry_type, &
         153       18432 :                          nthreads=nthreads,CardinalLength=ndepth*npoints,OrdinalLength=ndepth*ndepth)
         154             :     ! set some parameters need to support deep halos
         155       18432 :     edge%ndepth  = ndepth
         156       18432 :     edge%npoints = npoints
         157       18432 :     edge%lb      = 1 - edge%ndepth
         158       18432 :     edge%ub      = edge%npoints + edge%ndepth
         159             :
         160       18432 :   end subroutine initGhostBuffer
         161             :
         162             :
         163             :
         164           0 :   subroutine zeroEdgeBuffer(edge)
         165             :
         166             :     type (EdgeBuffer_t), intent(inout) :: edge
         167             :     integer :: i
         168             :
         169           0 :     do i=1,edge%nbuf
         170           0 :        edge%buf(i)     = 0.0d0
         171           0 :        edge%receive(i) = 0.0d0
         172             :     enddo
         173             :
         174           0 :   end subroutine zeroEdgeBuffer
         175             :
         176             :   ! =========================================
         177             :   ! initEdgeBuffer:
         178             :   !
         179             :   ! create an Real based communication buffer
         180             :   ! =========================================
         181       62976 :   subroutine initEdgeBuffer_r8(par,edge,elem,nlyr, bndry_type,nthreads,CardinalLength,OrdinalLength)
         182             :     use dimensions_mod, only: np, nelemd, max_corner_elem
         183             :     use schedtype_mod,  only: cycle_t, schedule_t, schedule
         184             :     use mpi,            only: MPI_VERSION
         185             :
         186             :     type (parallel_t),             intent(in)  :: par
         187             :     type (EdgeBuffer_t), target,   intent(out) :: edge
         188             :     type (element_t),              intent(in)  :: elem(:)
         189             :     integer,                       intent(in)  :: nlyr
         190             :     integer,             optional, intent(in)  :: bndry_type
         191             :     integer,             optional, intent(in)  :: nthreads
         192             :     integer,             optional, intent(in)  :: CardinalLength
         193             :     integer,             optional, intent(in)  :: OrdinalLength
         194             :
         195             :     ! Notes about the buf_ptr/receive_ptr options:
         196             :     !
         197             :     ! If an EdgeBuffer_t object is initialized from pre-existing storage
         198             :     ! (i.e. buf_ptr is provided and not null), it must *not* be freed,
         199             :     ! and must not be used if the underlying storage has been deallocated.
         200             :     !
         201             :     ! All these restrictions also applied to the old newbuf and newreceive
         202             :     ! options.
         203             :
         204             :     ! Workaround for NAG bug.
         205             :     ! NAG 5.3.1 dies if you use pointer bounds remapping to set
         206             :     ! a pointer that is also a component. So remap to temporary,
         207             :     ! then use that to set component pointer.
         208             :
         209             :     ! Local variables
         210             :     integer                    :: nbuf,ith
         211             :     integer                    :: nSendCycles, nRecvCycles
         212             :     integer                    :: icycle, ierr
         213             :     integer                    :: ie, i
         214             :     integer                    :: edgeid,elemid
         215             :     integer                    :: ptr,llen,moveLength, mLen, tlen
         216             :     type (Cycle_t),    pointer :: pCycle
         217             :     type (Schedule_t), pointer :: pSchedule
         218             :     integer                    :: dest, source, length, tag, iptr
         219             :     integer                    :: nlen, ithr
         220             :
         221             :     integer :: len, lenP,lenS
         222             :     integer :: j,jj,il,mesgid, dst0,src0
         223             :     integer :: moveptr
         224             :     integer :: nbuf2,ilm1,iem1,lenm1
         225       31488 :     integer,allocatable :: putmap2(:,:),getmap2(:,:)
         226             :     integer,allocatable :: scounts(:), rcounts(:)
         227             :     integer,allocatable :: sdispls(:), rdispls(:)
         228             :     integer :: nInter, nIntra
         229             :     integer :: icInter, icIntra
         230             :
         231             :     integer :: maxnsend
         232             :     integer :: tmpnMesg
         233             :     integer :: wintmpnMesg, wintmpDest, wintmpDisp
         234             :     integer(kind=i8) :: winsize
         235             :     integer :: win
         236             :     integer :: sizeofreal
         237             :     integer, allocatable :: tmpDest(:),tmpDisp(:)
         238             :     integer :: nFull
         239             :     integer :: disp, one
         240             :     integer                           :: errorcode,errorlen
         241             :     integer :: CardinalLen, OrdinalLen
         242             :     character(len=80)                 :: errorstring
         243             :     character(len=80), parameter      :: subname='initedgeBuffer'
         244             :
         245       31488 :     if(present(bndry_type)) then
         246        6144 :       if ( MPI_VERSION >= 3 ) then
         247        6144 :         edge%bndry_type = bndry_type
         248             :       else
         249             :         edge%bndry_type = HME_BNDRY_P2P
         250             :       endif
         251             :     else
         252       25344 :        edge%bndry_type = HME_BNDRY_P2P
         253             :     endif
         254             :
         255             :     ! Set the length of the cardinal and ordinal message lengths
         256       31488 :     if(present(CardinalLength)) then
         257       18432 :        CardinalLen = CardinalLength
         258             :     else
         259       13056 :        CardinalLen = np
         260             :     endif
         261       31488 :     if(present(OrdinalLength)) then
         262       18432 :        OrdinalLen = OrdinalLength
         263             :     else
         264       13056 :        OrdinalLen = 1
         265             :     endif
         266             :
         267             : ! DO NOT REMOVE THIS NEXT BARRIER
         268             : ! MT: This initial barrier fixes a long standing issue with Intel compilers on
         269             : ! two different platforms.  Without this barrier, edge buffers initialized from
         270             : ! within the threaded region would not work in a reproducable way with certain
         271             : ! thread combinations.  I cant explain why, but this fixes that issue on Edison
         272             : !$OMP BARRIER
         273             :
         274       31488 :     if (nlyr==0) return  ! tracer code might call initedgebuffer() with zero tracers
         275             :
         276             :
         277             : !$OMP MASTER
         278             :     !
         279             :     ! Keep a counter of how many times initedgebuffer is called.
         280             :     ! This is used to assign a unique message ID for the boundary exchange
         281             :     !
         282       31488 :     initedgebuffer_callid=initedgebuffer_callid+1
         283       31488 :     edge%id  = initedgebuffer_callid
         284       31488 :     edge%tag = BNDRY_TAG_BASE + MODULO(edge%id, MAX_ACTIVE_MSG)
         285             :
         286      125952 :     allocate(edge%putmap(max_neigh_edges,nelemd))
         287       94464 :     allocate(edge%getmap(max_neigh_edges,nelemd))
         288       94464 :     allocate(edge%reverse(max_neigh_edges,nelemd))
         289             :
         290     2024088 :     edge%putmap(:,:)=-1
         291     2024088 :     edge%getmap(:,:)=-1
         292             :
         293      125952 :     allocate(putmap2(max_neigh_edges,nelemd))
         294       94464 :     allocate(getmap2(max_neigh_edges,nelemd))
         295     2024088 :     putmap2(:,:)=-1
         296     2024088 :     getmap2(:,:)=-1
         297      252888 :     do ie=1,nelemd
         298     2024088 :        do i=1,max_neigh_edges
         299     1992600 :           edge%reverse(i,ie) = elem(ie)%desc%reverse(i)
         300             :        enddo
         301             :     enddo
         302             :
         303       31488 :     pSchedule  => Schedule(1)
         304       31488 :     nSendCycles = pSchedule%nSendCycles
         305       31488 :     nRecvCycles = pSchedule%nRecvCycles
         306       31488 :     nInter      = pSchedule%nInter
         307       31488 :     nIntra      = pSchedule%nIntra
         308       31488 :     nFull       = nInter+nIntra
         309             :
         310       31488 :     edge%nInter=nInter
         311       31488 :     edge%nIntra=nIntra
         312             :
         313       31488 :     if(nInter>0) then
         314      125952 :        allocate(edge%rcountsInter(nInter),edge%rdisplsInter(nInter))
         315       94464 :        allocate(edge%scountsInter(nInter),edge%sdisplsInter(nInter))
         316             :     endif
         317       31488 :     if(nIntra>0) then
         318           0 :        allocate(edge%rcountsIntra(nIntra),edge%rdisplsIntra(nIntra))
         319           0 :        allocate(edge%scountsIntra(nIntra),edge%sdisplsIntra(nIntra))
         320             :     endif
         321             :
         322       31488 :     if (nSendCycles>0) then
         323      125952 :        allocate(edge%scountsFull(nSendCycles),edge%sdisplsFull(nSendCycles))
         324       62976 :        allocate(edge%Srequest(nSendCycles))
         325      234274 :        edge%scountsFull(:) = 0
         326             :     endif
         327             :     !
         328             :     ! Setup the data-structures for the sends
         329             :     !
         330       31488 :     j         = 1
         331       31488 :     icycle    = 1
         332       31488 :     dst0      = pSchedule%pIndx(j)%mesgid
         333       31488 :     il        = pSchedule%pIndx(j)%edgeid
         334       31488 :     ie        = pSchedule%pIndx(j)%elemid
         335       31488 :     len       = CalcSegmentLength(pSchedule%pIndx(j),CardinalLen,OrdinalLen,nlyr)
         336       31488 :     edge%putmap(il,ie) = 0
         337       31488 :     if(nSendCycles>0) then
         338       31488 :         edge%sdisplsFull(icycle) = edge%putmap(il,ie)
         339       31488 :         edge%scountsFull(icycle) = len
         340             :     endif
         341       31488 :     ilm1 = il
         342       31488 :     iem1 = ie
         343       31488 :     lenm1 = len
         344             :
         345     1771200 :     do j=2,SIZE(pSchedule%pIndx)
         346     1739712 :        il = pSchedule%pIndx(j)%edgeid
         347     1739712 :        ie = pSchedule%pIndx(j)%elemid
         348     1739712 :        mesgid = pSchedule%pIndx(j)%mesgid
         349     1771200 :        if(il>0 .and. ie >0) then
         350     1738728 :           len     = CalcSegmentLength(pSchedule%pIndx(j),CardinalLen,OrdinalLen,nlyr)
         351     1738728 :           edge%putmap(il,ie) = edge%putmap(ilm1,iem1)+lenm1
         352     1738728 :           if(mesgid .ne. par%rank) then  ! don't enter if this is a move cycle where (mesgid == par%rank)
         353      948740 :              if(mesgid .ne. dst0) then
         354      171298 :                 icycle=icycle+1
         355      171298 :                 if (nSendCycles>0) edge%sdisplsFull(icycle) = edge%putmap(il,ie)
         356             :                 dst0=mesgid
         357             :              endif
         358      948740 :              if (nSendCycles>0) edge%scountsFull(icycle) = edge%scountsFull(icycle)+len
         359             :           endif
         360             :           ilm1=il
         361             :           iem1=ie
         362             :           lenm1=len
         363             :        endif
         364             :     enddo
         365             :
         366             :     icInter=0
         367             :     icIntra=0
         368      234274 :     do icycle=1,nSendCycles
         369      234274 :        if(pSchedule%SendCycle(icycle)%onNode .eqv. .FALSE.) then
         370      202786 :           icInter=icInter+1
         371      202786 :           edge%sdisplsInter(icInter)=edge%sdisplsFull(icycle)
         372      202786 :           edge%scountsInter(icInter)=edge%scountsFull(icycle)
         373             :        else
         374           0 :           icIntra=icIntra+1
         375           0 :           edge%sdisplsIntra(icIntra)=edge%sdisplsFull(icycle)
         376           0 :           edge%scountsIntra(icIntra)=edge%scountsFull(icycle)
         377             :        endif
         378             :     enddo
         379             :
         380       31488 :     if (nRecvCycles>0) then
         381      125952 :        allocate(edge%rcountsFull(nRecvCycles),edge%rdisplsFull(nRecvCycles))
         382       94464 :        allocate(edge%getDisplsFull(nRecvCycles),edge%putDisplsFull(nRecvCycles))
         383      234274 :        edge%rcountsFull(:) = 0
         384             :        ! allocate the MPI Send/Recv request handles
         385       62976 :        allocate(edge%Rrequest(nRecvCycles))
         386       94464 :        allocate(edge%status(HME_status_size,nRecvCycles))
         387             :     endif
         388             :
         389             :     !
         390             :     ! Setup the data-structures for the receives
         391             :     !
         392       31488 :     j      = 1
         393       31488 :     icycle = 1
         394       31488 :     src0    = pSchedule%gIndx(j)%mesgid
         395       31488 :     il      = pSchedule%gIndx(j)%edgeid
         396       31488 :     ie      = pSchedule%gIndx(j)%elemid
         397       31488 :     len     = CalcSegmentLength(pSchedule%gIndx(j),CardinalLen,OrdinalLen,nlyr)
         398       31488 :     edge%getmap(il,ie) = 0
         399       31488 :     if (nRecvCycles>0) then
         400       31488 :        edge%rdisplsFull(icycle) = edge%getmap(il,ie)
         401       31488 :        edge%rcountsFull(icycle) = len
         402             :     endif
         403       31488 :     ilm1=il
         404       31488 :     iem1=ie
         405       31488 :     lenm1=len
         406             :
         407     1771200 :     do j=2,SIZE(pSchedule%gIndx)
         408     1739712 :        il = pSchedule%gIndx(j)%edgeid
         409     1739712 :        ie = pSchedule%gIndx(j)%elemid
         410     1739712 :        mesgid = pSchedule%gIndx(j)%mesgid
         411     1771200 :        if(il>0 .and. ie >0) then
         412     1738728 :           len     = CalcSegmentLength(pSchedule%gIndx(j),CardinalLen,OrdinalLen,nlyr)
         413     1738728 :           edge%getmap(il,ie) = edge%getmap(ilm1,iem1)+lenm1
         414     1738728 :           if(mesgid .ne. par%rank) then ! don't enter if this is a move cycle where (mesgid == par%rank)
         415      948740 :              if(mesgid .ne. src0) then
         416      171298 :                 if (nRecvCycles>0) edge%rdisplsFull(icycle+1) = edge%getmap(il,ie)
         417      171298 :                 icycle=icycle+1
         418      171298 :                 src0=mesgid
         419             :              endif
         420      948740 :              if (nRecvCycles>0) edge%rcountsFull(icycle) = edge%rcountsFull(icycle)+len
         421             :           endif
         422             :           ilm1=il
         423             :           iem1=ie
         424             :           lenm1=len
         425             :        endif
         426             :     enddo
         427             :
         428             :
         429             :     !
         430             :     ! populate the Inter and Intra node communication data-structures
         431             :     !
         432             :     icInter=0
         433             :     icIntra=0
         434      234274 :     do icycle=1,nRecvCycles
         435      234274 :        if(pSchedule%RecvCycle(icycle)%onNode .eqv. .FALSE.) then
         436      202786 :           icInter=icInter+1
         437      202786 :           edge%rdisplsInter(icInter)=edge%rdisplsFull(icycle)
         438      202786 :           edge%rcountsInter(icInter)=edge%rcountsFull(icycle)
         439             :        else
         440           0 :           icIntra=icIntra+1
         441           0 :           edge%rdisplsIntra(icIntra)=edge%rdisplsFull(icycle)
         442           0 :           edge%rcountsIntra(icIntra)=edge%rcountsFull(icycle)
         443             :        endif
         444             :     enddo
         445             :
         446             :
         447             :     ! Setup the data-structures for the on process moves
         448             :     ! Note that this assumes that the data to move is at
         449             :     ! the end of the  message buffer.
         450       31488 :     if(nRecvCycles>0) then
         451       31488 :        moveptr = edge%rdisplsFull(nRecvCycles)+edge%rcountsFull(nRecvCycles)+1
         452             :     else
         453             :        moveptr = 1
         454             :     endif
         455       31488 :     moveLength = 0
         456     1802688 :     do j=1,SIZE(pSchedule%gIndx)
         457     1771200 :        il     = pSchedule%gIndx(j)%edgeid
         458     1771200 :        ie     = pSchedule%gIndx(j)%elemid
         459     1771200 :        mesgid = pSchedule%gIndx(j)%mesgid
         460     1802688 :        if(mesgid == par%rank) then
         461      789988 :           len     = CalcSegmentLength(pSchedule%gIndx(j),CardinalLen,OrdinalLen,nlyr)
         462      789988 :           moveLength = moveLength + len
         463             :        endif
         464             :     enddo
         465             :
         466             :     ! decompose the move data between the available threads
         467       31488 :     if(max_num_threads<=0) then
         468           0 :        nlen = 1
         469             :     else
         470       31488 :        if(present(nthreads)) then
         471       26112 :           if (nthreads > 0) then
         472       26112 :              nlen = nthreads
         473             :           else
         474           0 :              nlen = max_num_threads
         475             :           end if
         476             :        else
         477        5376 :           nlen = 1
         478             :        end if
         479             :     end if
         480       31488 :     call gbarrier_init(edge%gbarrier, nlen)
         481             :
         482       94464 :     allocate(edge%moveLength(nlen))
         483       62976 :     allocate(edge%movePtr(nlen))
         484             :
         485       31488 :     if (nlen > 1) then
         486             :        ! the master thread performs no data movement because it is busy with the
         487             :        ! MPI messaging
         488           0 :        edge%moveLength(1) = -1
         489           0 :        edge%movePtr(1) = 0
         490             :
         491             :        ! Calculate the length of the local copy in bndy_exchange
         492           0 :        llen = ceiling(real(moveLength,kind=r8)/real(nlen-1,kind=r8))
         493           0 :        iptr = moveptr
         494           0 :        mLen = 0
         495           0 :        do i=2,nlen
         496           0 :          if( (mLen+llen) <= moveLength)  then
         497             :             tlen = llen
         498             :          else
         499           0 :             tlen = moveLength - mLen
         500             :          endif
         501           0 :          edge%moveLength(i) = tlen
         502           0 :          edge%movePtr(i)    = iptr
         503           0 :          iptr = iptr + tlen
         504           0 :          mLen = mLen + tLen
         505             :        enddo
         506             :     else
         507       31488 :        edge%moveLength(1) = moveLength
         508       31488 :        edge%movePtr(1) = moveptr
         509             :     endif
         510             :
         511             :     ! Set the maximum length of the message buffer
         512       31488 :     nbuf = movePtr+moveLength
         513             :
         514       31488 :     edge%nlyr=nlyr
         515       31488 :     edge%nbuf=nbuf
         516             :
         517       94464 :     allocate(edge%receive(nbuf))
         518       62976 :     allocate(edge%buf(nbuf))
         519             :
         520             : 21  format('RANK: ',i2, A,8(i6))
         521             :
         522             : !$OMP END MASTER
         523             : ! MT: This next barrier is also needed - threads cannot start using edge()
         524             : ! until MASTER is done initializing it
         525             : !$OMP BARRIER
         526             :
         527       31488 :   end subroutine initEdgeBuffer_r8
         528             :
         529     4330420 :   integer function CalcSegmentLength(pgIndx,CardinalLength,OrdinalLength,nlyr) result(len)
         530             :
         531             :      type(pgindex_t) ::  pgIndx
         532             :      integer, intent(in) :: CardinalLength,OrdinalLength
         533             :      integer, intent(in) :: nlyr
         534             :
         535             :      integer :: rem
         536             :      integer, parameter :: alignment=1  ! align on word boundaries
         537             : !     integer, parameter :: alignment=2  ! align on 2 word boundaries
         538             : !     integer, parameter :: alignment=8  ! align on 8 word boundaries
         539             :
         540     6589110 :      select case(pgIndx%edgeType)
         541             :         CASE(HME_Cardinal)
         542     2258690 :           len = nlyr*CardinalLength
         543             :         CASE(HME_Ordinal)
         544     4330420 :           len = nlyr*OrdinalLength
         545             :      end select
         546             :
         547     4330420 :      rem = MODULO(len,alignment)
         548             :      if(rem .ne. 0) then
         549             :         len = len + (alignment-rem)
         550             :      endif
         551             :
         552     4330420 :   end function calcSegmentLength
         553             :
         554             :   ! =========================================
         555             :   ! initEdgeBuffer:
         556             :   !
         557             :   ! create an Integer based communication buffer
         558             :   ! =========================================
         559        1536 :   subroutine initEdgeBuffer_i8(edge,nlyr)
         560             :     use dimensions_mod, only : np, nelemd, max_corner_elem
         561             :
         562             :     integer,                 intent(in)  :: nlyr
         563             :     type (LongEdgeBuffer_t), intent(out) :: edge
         564             :
         565             :     ! Local variables
         566             :     integer :: nbuf
         567             :
         568             :     ! sanity check for threading
         569        1536 :     if (omp_get_num_threads()>1) then
         570           0 :        call endrun('ERROR: initEdgeBuffer must be called before threaded reagion')
         571             :     endif
         572             :
         573        1536 :     nbuf=4*(np+max_corner_elem)*nelemd
         574        1536 :     edge%nlyr=nlyr
         575        1536 :     edge%nbuf=nbuf
         576        6144 :     allocate(edge%buf(nlyr,nbuf))
         577      433536 :     edge%buf(:,:)=0
         578             :
         579        4608 :     allocate(edge%receive(nlyr,nbuf))
         580      433536 :     edge%receive(:,:)=0
         581             :
         582        1536 :   end subroutine initEdgeBuffer_i8
         583             :   ! =========================================
         584             :   ! edgeDGVpack:
         585             :   !
         586             :   ! Pack edges of v into buf for DG stencil
         587             :   ! =========================================
         588           0 :   subroutine edgeDGVpack(edge,v,vlyr,kptr,ielem)
         589             :     use dimensions_mod, only: np
         590             :
         591             :     type (EdgeBuffer_t)        :: edge
         592             :     integer,        intent(in) :: vlyr
         593             :     real (kind=r8), intent(in) :: v(np,np,vlyr)
         594             :     integer,        intent(in) :: kptr
         595             :     integer,        intent(in) :: ielem
         596             :
         597             :     ! =========================================
         598             :     ! This code is just a wrapper call the
         599             :     !   normal oldedgeVpack
         600             :     ! =========================================
         601           0 :     call edgeVpack(edge,v,vlyr,kptr,ielem)
         602             :
         603           0 :   end subroutine edgeDGVpack
         604             :
         605        4608 :   subroutine FreeGhostBuffer_r8(edge)
         606             :      type (EdgeBuffer_t), intent(inout) :: edge
         607        4608 :      call FreeEdgeBuffer_r8(edge)
         608        4608 :   end subroutine FreeGhostBuffer_r8
         609             :   ! ===========================================
         610             :   !  FreeEdgeBuffer:
         611             :   !
         612             :   !  Freed an edge communication buffer
         613             :   ! =========================================
         614        8448 :   subroutine FreeEdgeBuffer_r8(edge)
         615             :
         616             :     type (EdgeBuffer_t),intent(inout) :: edge
         617             :
         618             : !$OMP BARRIER
         619             : !$OMP MASTER
         620        8448 :     if(allocated(edge%buf))         deallocate(edge%buf)
         621        8448 :     if(allocated(edge%receive))     deallocate(edge%receive)
         622        8448 :     if(associated(edge%putmap))     deallocate(edge%putmap)
         623        8448 :     if(associated(edge%getmap))     deallocate(edge%getmap)
         624        8448 :     if(associated(edge%reverse))    deallocate(edge%reverse)
         625        8448 :     if(associated(edge%moveLength)) deallocate(edge%moveLength)
         626        8448 :     if(associated(edge%movePtr))    deallocate(edge%movePtr)
         627             :
         628             :     ! All MPI communications
         629        8448 :     if(associated(edge%rcountsFull)) deallocate(edge%rcountsFull)
         630        8448 :     if(associated(edge%scountsFull)) deallocate(edge%scountsFull)
         631        8448 :     if(associated(edge%sdisplsFull)) deallocate(edge%sdisplsFull)
         632        8448 :     if(associated(edge%rdisplsFull)) deallocate(edge%rdisplsFull)
         633             :
         634             :     ! Intra-node MPI Communication
         635        8448 :     if(edge%nIntra>0) then
         636           0 :       if(associated(edge%rcountsIntra)) deallocate(edge%rcountsIntra)
         637           0 :       if(associated(edge%scountsIntra)) deallocate(edge%scountsIntra)
         638           0 :       if(associated(edge%sdisplsIntra)) deallocate(edge%sdisplsIntra)
         639           0 :       if(associated(edge%rdisplsIntra)) deallocate(edge%rdisplsIntra)
         640             :     endif
         641             :
         642             :     ! Inter-node MPI Communication
         643        8448 :     if(edge%nInter>0) then
         644        8448 :       if(associated(edge%rcountsInter)) deallocate(edge%rcountsInter)
         645        8448 :       if(associated(edge%scountsInter)) deallocate(edge%scountsInter)
         646        8448 :       if(associated(edge%sdisplsInter)) deallocate(edge%sdisplsInter)
         647        8448 :       if(associated(edge%rdisplsInter)) deallocate(edge%rdisplsInter)
         648             :     endif
         649        8448 :     if(allocated(edge%rRequest)) deallocate(edge%rRequest)
         650        8448 :     if(allocated(edge%sRequest)) deallocate(edge%sRequest)
         651        8448 :     if(allocated(edge%status)) deallocate(edge%status)
         652        8448 :     call gbarrier_delete(edge%gbarrier)
         653             :
         654             : !$OMP END MASTER
         655             :
         656        8448 :   end subroutine FreeEdgeBuffer_r8
         657             :
         658             :   ! ===========================================
         659             :   !  FreeEdgeBuffer:
         660             :   !
         661             :   !  Freed an edge communication buffer
         662             :   ! =========================================
         663        1536 :   subroutine FreeEdgeBuffer_i8(edge)
         664             :
         665             :     type (LongEdgeBuffer_t),intent(inout) :: edge
         666             :
         667        1536 :     edge%nbuf=0
         668        1536 :     edge%nlyr=0
         669        1536 :     deallocate(edge%buf)
         670        1536 :     deallocate(edge%receive)
         671             :
         672        1536 :   end subroutine FreeEdgeBuffer_i8
         673             :
         674             :   ! =========================================
         675             :   !
         676             :   !> @brief Pack edges of v into an edge buffer for boundary exchange.
         677             :   !
         678             :   !> This subroutine packs for one or more vertical layers into an edge
         679             :   !! buffer. If the buffer associated with edge is not large enough to
         680             :   !! hold all vertical layers you intent to pack, the method will
         681             :   !! halt the program with a call to endrum().
         682             :   !! @param[in] edge Edge Buffer into which the data will be packed.
         683             :   !! This buffer must be previously allocated with initEdgeBuffer().
         684             :   !! @param[in] v The data to be packed.
         685             :   !! @param[in] vlyr Number of vertical level coming into the subroutine
         686             :   !! for packing for input v.
         687             :   !! @param[in] kptr Vertical pointer to the place in the edge buffer where
         688             :   !! data will be located.
         689             :   ! =========================================
         690   823548600 :   subroutine edgeVpack(edge,v,vlyr,kptr,ielem)
         691             :     use dimensions_mod, only: np, max_corner_elem
         692             :     use control_mod,    only: north, south, east, west, neast, nwest, seast, swest
         693             :
         694             :     type (EdgeBuffer_t)        :: edge
         695             :     integer,        intent(in) :: vlyr
         696             :     real (kind=r8), intent(in) :: v(np,np,vlyr)
         697             :     integer,        intent(in) :: kptr
         698             :     integer,        intent(in) :: ielem
         699             :
         700             :     ! Local variables
         701             :     integer :: i,k,ir,ll,iptr
         702             :     integer :: is,ie,in,iw,edgeptr
         703             :
         704   823548600 :     is = edge%putmap(south,ielem)
         705   823548600 :     ie = edge%putmap(east,ielem)
         706   823548600 :     in = edge%putmap(north,ielem)
         707   823548600 :     iw = edge%putmap(west,ielem)
         708   823548600 :     if (edge%nlyr < (kptr+vlyr) ) then
         709           0 :        print *,'edge%nlyr = ',edge%nlyr
         710           0 :        print *,'kptr+vlyr = ',kptr+vlyr
         711           0 :        call endrun('edgeVpack: Buffer overflow: size of the vertical dimension must be increased!')
         712             :     endif
         713             :
         714             : !dir$ ivdep
         715 68500355400 :     do k=1,vlyr
         716 67676806800 :        iptr = np*(kptr+k-1)
         717 >33920*10^7 :        do i=1,np
         718 >27070*10^7 :           edge%buf(iptr+ie+i)   = v(np ,i ,k) ! East
         719 >27070*10^7 :           edge%buf(iptr+is+i)   = v(i  ,1 ,k) ! South
         720 >27070*10^7 :           edge%buf(iptr+in+i)   = v(i  ,np,k) ! North
         721 >33838*10^7 :           edge%buf(iptr+iw+i)   = v(1  ,i ,k) ! West
         722             :        enddo
         723             :     enddo
         724             :
         725             :     !  This is really kludgy way to setup the index reversals
         726             :     !  But since it is so a rare event not real need to spend time optimizing
         727             :
         728   823548600 :     if(edge%reverse(south,ielem)) then
         729  1141672590 :        do k=1,vlyr
         730  1127946780 :           iptr = np*(kptr+k-1)+is
         731  5653459710 :           do i=1,np
         732  4511787120 :              ir = np-i+1
         733  5639733900 :              edge%buf(iptr+ir)=v(i,1,k)
         734             :           enddo
         735             :        enddo
         736             :     endif
         737             :
         738   823548600 :     if(edge%reverse(east,ielem)) then
         739   380557530 :        do k=1,vlyr
         740   375982260 :           iptr=np*(kptr+k-1)+ie
         741  1884486570 :           do i=1,np
         742  1503929040 :              ir = np-i+1
         743  1879911300 :              edge%buf(iptr+ir)=v(np,i,k)
         744             :           enddo
         745             :        enddo
         746             :     endif
         747             :
         748   823548600 :     if(edge%reverse(north,ielem)) then
         749  1141672590 :        do k=1,vlyr
         750  1127946780 :           iptr=np*(kptr+k-1)+in
         751  5653459710 :           do i=1,np
         752  4511787120 :              ir = np-i+1
         753  5639733900 :              edge%buf(iptr+ir)=v(i,np,k)
         754             :           enddo
         755             :        enddo
         756             :     endif
         757             :
         758   823548600 :     if(edge%reverse(west,ielem)) then
         759   380557530 :        do k=1,vlyr
         760   375982260 :           iptr=np*(kptr+k-1)+iw
         761  1884486570 :           do i=1,np
         762  1503929040 :              ir = np-i+1
         763  1879911300 :              edge%buf(iptr+ir)=v(1,i,k)
         764             :           enddo
         765             :        enddo
         766             :     endif
         767             :
         768             : ! SWEST
         769  1647097200 :     do ll=swest,swest+max_corner_elem-1
         770  1647097200 :         if (edge%putmap(ll,ielem) /= -1) then
         771   822633546 :             edgeptr = edge%putmap(ll,ielem)+1
         772 68424243894 :             do k=1,vlyr
         773 67601610348 :                 iptr = (kptr+k-1)+edgeptr
         774 67601610348 :                 if (iptr > size(edge%buf)) then
         775           0 :                    write(6, *) 'ERROR SW: ',size(edge%buf),iptr,edge%putmap(ll,ielem)
         776           0 :                    call endrun('pointer bounds ERROR SW')
         777             :                 end if
         778 68424243894 :                 edge%buf(iptr) = v(1, 1, k)
         779             :             end do
         780             :         end if
         781             :     end do
         782             :
         783             : ! SEAST
         784  1647097200 :     do ll=swest+max_corner_elem,swest+2*max_corner_elem-1
         785  1647097200 :         if (edge%putmap(ll,ielem) /= -1) then
         786   822633546 :             edgeptr = edge%putmap(ll,ielem)+1
         787 68424243894 :             do k=1,vlyr
         788 67601610348 :                iptr = (kptr+k-1)+edgeptr
         789 67601610348 :                if (iptr > size(edge%buf)) then
         790           0 :                   write(6, *) 'ERROR SE: ',size(edge%buf),iptr,edge%putmap(ll,ielem)
         791           0 :                   call endrun('pointer bounds ERROR SE')
         792             :                end if
         793 68424243894 :                edge%buf(iptr)=v(np, 1, k)
         794             :             end do
         795             :         end if
         796             :     end do
         797             :
         798             : ! NEAST
         799  1647097200 :     do ll=swest+3*max_corner_elem,swest+4*max_corner_elem-1
         800  1647097200 :         if (edge%putmap(ll,ielem) /= -1) then
         801   822633546 :             edgeptr = edge%putmap(ll,ielem)+1
         802 68424243894 :             do k=1,vlyr
         803 67601610348 :                iptr = (kptr+k-1)+edgeptr
         804 67601610348 :                if (iptr > size(edge%buf)) then
         805           0 :                   write(6, *) 'ERROR NE: ',size(edge%buf),iptr,edge%putmap(ll,ielem)
         806           0 :                   call endrun('pointer bounds ERROR NE')
         807             :                end if
         808 68424243894 :                edge%buf(iptr) = v(np, np, k)
         809             :             end do
         810             :         end if
         811             :     end do
         812             :
         813             : ! NWEST
         814  1647097200 :     do ll=swest+2*max_corner_elem,swest+3*max_corner_elem-1
         815  1647097200 :         if (edge%putmap(ll,ielem) /= -1) then
         816   822633546 :             edgeptr = edge%putmap(ll,ielem)+1
         817 68424243894 :             do k=1,vlyr
         818 67601610348 :                iptr = (kptr+k-1)+edgeptr
         819 67601610348 :                if (iptr > size(edge%buf)) then
         820           0 :                   write(6, *) 'ERROR NW: ',size(edge%buf),iptr,edge%putmap(ll,ielem)
         821           0 :                   call endrun('pointer bounds ERROR NW')
         822             :                end if
         823 68424243894 :                edge%buf(iptr) = v(1, np, k)
         824             :             end do
         825             :         end if
         826             :     end do
         827             :
         828   823548600 :   end subroutine edgeVpack
         829             :
         830           0 :   subroutine edgeSpack_r8(edge,v,vlyr,kptr,ielem)
         831             :     use dimensions_mod, only: np, max_corner_elem
         832             :     use control_mod,    only: north, south, east, west, neast, nwest, seast, swest
         833             :
         834             :     type (EdgeBuffer_t)        :: edge
         835             :     integer,        intent(in) :: vlyr
         836             :     real (kind=r8), intent(in) :: v(vlyr)
         837             :     integer,        intent(in) :: kptr
         838             :     integer,        intent(in) :: ielem
         839             :
         840             :     ! Local variables
         841             :     integer        :: i,k,ir,ll,iptr
         842             :     integer        :: is,ie,in,iw,edgeptr
         843             :     real (kind=r8) :: tmp
         844             :
         845           0 :     is = edge%putmap(south,ielem)
         846           0 :     ie = edge%putmap(east,ielem)
         847           0 :     in = edge%putmap(north,ielem)
         848           0 :     iw = edge%putmap(west,ielem)
         849           0 :     if (edge%nlyr < (kptr+vlyr) ) then
         850           0 :        call endrun('edgeSpack: Buffer overflow: size of the vertical dimension must be increased!')
         851             :     endif
         852             :
         853           0 :     do k=1,vlyr
         854           0 :        iptr = kptr+k-1
         855           0 :        edge%buf(iptr+ie+1)   = v(k) ! East
         856           0 :        edge%buf(iptr+is+1)   = v(k) ! South
         857           0 :        edge%buf(iptr+in+1)   = v(k) ! North
         858           0 :        edge%buf(iptr+iw+1)   = v(k) ! West
         859             :     enddo
         860             :
         861             : ! SWEST
         862           0 :     do ll=swest,swest+max_corner_elem-1
         863           0 :         if (edge%putmap(ll,ielem) /= -1) then
         864           0 :             edgeptr=edge%putmap(ll,ielem)+1
         865           0 :             do k=1,vlyr
         866           0 :                 iptr = (kptr+k-1)+edgeptr
         867           0 :                 edge%buf(iptr)=v(k)
         868             :             end do
         869             :         end if
         870             :     end do
         871             :
         872             : ! SEAST
         873           0 :     do ll=swest+max_corner_elem,swest+2*max_corner_elem-1
         874           0 :         if (edge%putmap(ll,ielem) /= -1) then
         875           0 :             edgeptr=edge%putmap(ll,ielem)+1
         876           0 :             do k=1,vlyr
         877           0 :                 iptr = (kptr+k-1)+edgeptr
         878           0 :                 edge%buf(iptr)=v(k)
         879             :             end do
         880             :         end if
         881             :     end do
         882             :
         883             : ! NEAST
         884           0 :     do ll=swest+3*max_corner_elem,swest+4*max_corner_elem-1
         885           0 :         if (edge%putmap(ll,ielem) /= -1) then
         886           0 :             edgeptr=edge%putmap(ll,ielem)+1
         887           0 :             do k=1,vlyr
         888           0 :                 iptr = (kptr+k-1)+edgeptr
         889           0 :                 edge%buf(iptr)=v(k)
         890             :             end do
         891             :         end if
         892             :     end do
         893             :
         894             : ! NWEST
         895           0 :     do ll=swest+2*max_corner_elem,swest+3*max_corner_elem-1
         896           0 :         if (edge%putmap(ll,ielem) /= -1) then
         897           0 :             edgeptr=edge%putmap(ll,ielem)+1
         898           0 :             do k=1,vlyr
         899           0 :                 iptr = (kptr+k-1)+edgeptr
         900           0 :                 edge%buf(iptr)=v(k)
         901             :             end do
         902             :         end if
         903             :     end do
         904             :
         905           0 :   end subroutine edgeSpack_r8
         906             :
         907             :   ! =========================================
         908             :   ! LongEdgeVpack:
         909             :   !
         910             :   ! Pack edges of v into buf...
         911             :   ! =========================================
         912       10800 :   subroutine LongEdgeVpack(edge,v,vlyr,kptr,desc)
         913             :     use control_mod,    only: north, south, east, west, neast, nwest, seast, swest
         914             :     use dimensions_mod, only: np, max_corner_elem
         915             :
         916             :     type (LongEdgeBuffer_t)             :: edge
         917             :     integer,                 intent(in) :: vlyr
         918             :     integer ,                intent(in) :: v(np,np,vlyr)
         919             :     integer,                 intent(in) :: kptr
         920             :     type (EdgeDescriptor_t), intent(in) :: desc
         921             :
         922             :     ! Local variables
         923             :     logical, parameter :: UseUnroll = .TRUE.
         924             :     integer            :: i,k,ir,l
         925             :     integer            :: is,ie,in,iw
         926             :
         927       10800 :     if(.not. threadsafe) then
         928             : !$OMP BARRIER
         929           0 :        threadsafe=.true.
         930             :     end if
         931             :
         932       10800 :     is = desc%putmapP(south)
         933       10800 :     ie = desc%putmapP(east)
         934       10800 :     in = desc%putmapP(north)
         935       10800 :     iw = desc%putmapP(west)
         936             :
         937             :     if(MODULO(np,2) == 0 .and. UseUnroll) then
         938       21600 :        do k=1,vlyr
         939       10800 :           do i=1,np,2
         940       21600 :              edge%buf(kptr+k,is+i)   = v(i  ,1 ,k)
         941       21600 :              edge%buf(kptr+k,is+i+1) = v(i+1,1 ,k)
         942       21600 :              edge%buf(kptr+k,ie+i)   = v(np ,i ,k)
         943       21600 :              edge%buf(kptr+k,ie+i+1) = v(np ,i+1 ,k)
         944       21600 :              edge%buf(kptr+k,in+i)   = v(i  ,np,k)
         945       21600 :              edge%buf(kptr+k,in+i+1) = v(i+1  ,np,k)
         946       21600 :              edge%buf(kptr+k,iw+i)   = v(1  ,i ,k)
         947       21600 :              edge%buf(kptr+k,iw+i+1) = v(1  ,i+1 ,k)
         948             :
         949             :           enddo
         950             :        end do
         951             :     else
         952             :        do k=1,vlyr
         953             :           do i=1,np
         954             :              edge%buf(kptr+k,is+i)   = v(i  ,1 ,k)
         955             :              edge%buf(kptr+k,ie+i)   = v(np ,i ,k)
         956             :              edge%buf(kptr+k,in+i)   = v(i  ,np,k)
         957             :              edge%buf(kptr+k,iw+i)   = v(1  ,i ,k)
         958             :           enddo
         959             :        end do
         960             :
         961             :     endif
         962             :
         963             :
         964             :     !  This is really kludgy way to setup the index reversals
         965             :     !  But since it is so a rare event not real need to spend time optimizing
         966             :
         967       10800 :     if(desc%reverse(south)) then
         968         180 :        is = desc%putmapP(south)
         969         360 :        do k=1,vlyr
         970        1080 :           do i=1,np
         971         720 :              ir = np-i+1
         972         900 :              edge%buf(kptr+k,is+ir)=v(i,1,k)
         973             :           enddo
         974             :        enddo
         975             :     endif
         976             :
         977       10800 :     if(desc%reverse(east)) then
         978          60 :        ie = desc%putmapP(east)
         979         120 :        do k=1,vlyr
         980         360 :           do i=1,np
         981         240 :              ir = np-i+1
         982         300 :              edge%buf(kptr+k,ie+ir)=v(np,i,k)
         983             :           enddo
         984             :        enddo
         985             :     endif
         986             :
         987       10800 :     if(desc%reverse(north)) then
         988         180 :        in = desc%putmapP(north)
         989         360 :        do k=1,vlyr
         990        1080 :           do i=1,np
         991         720 :              ir = np-i+1
         992         900 :              edge%buf(kptr+k,in+ir)=v(i,np,k)
         993             :           enddo
         994             :        enddo
         995             :     endif
         996             :
         997       10800 :     if(desc%reverse(west)) then
         998          60 :        iw = desc%putmapP(west)
         999         120 :        do k=1,vlyr
        1000         360 :           do i=1,np
        1001         240 :              ir = np-i+1
        1002         300 :              edge%buf(kptr+k,iw+ir)=v(1,i,k)
        1003             :           enddo
        1004             :        enddo
        1005             :     endif
        1006             :
        1007             : ! SWEST
        1008       21600 :     do l=swest,swest+max_corner_elem-1
        1009       21600 :         if (desc%putmapP(l) /= -1) then
        1010       21576 :             do k=1,vlyr
        1011       21576 :                 edge%buf(kptr+k,desc%putmapP(l)+1)=v(1  ,1 ,k)
        1012             :             end do
        1013             :         end if
        1014             :     end do
        1015             :
        1016             : ! SEAST
        1017       21600 :     do l=swest+max_corner_elem,swest+2*max_corner_elem-1
        1018       21600 :         if (desc%putmapP(l) /= -1) then
        1019       21576 :             do k=1,vlyr
        1020       21576 :                 edge%buf(kptr+k,desc%putmapP(l)+1)=v(np ,1 ,k)
        1021             :             end do
        1022             :         end if
        1023             :     end do
        1024             :
        1025             : ! NEAST
        1026       21600 :     do l=swest+3*max_corner_elem,swest+4*max_corner_elem-1
        1027       21600 :         if (desc%putmapP(l) /= -1) then
        1028       21576 :             do k=1,vlyr
        1029       21576 :                 edge%buf(kptr+k,desc%putmapP(l)+1)=v(np ,np,k)
        1030             :             end do
        1031             :         end if
        1032             :     end do
        1033             :
        1034             : ! NWEST
        1035       21600 :     do l=swest+2*max_corner_elem,swest+3*max_corner_elem-1
        1036       21600 :         if (desc%putmapP(l) /= -1) then
        1037       21576 :             do k=1,vlyr
        1038       21576 :                 edge%buf(kptr+k,desc%putmapP(l)+1)=v(1  ,np,k)
        1039             :             end do
        1040             :         end if
        1041             :     end do
        1042             :
        1043       10800 :   end subroutine LongEdgeVpack
        1044             :
        1045   823548600 :   subroutine edgeVunpack(edge,v,vlyr,kptr,ielem,rank)
        1046             :     use dimensions_mod, only: np, max_corner_elem
        1047             :     use control_mod,    only: north, south, east, west, neast, nwest, seast, swest
        1048             :
        1049             :     type (EdgeBuffer_t), intent(in)    :: edge
        1050             :     integer,             intent(in)    :: vlyr
        1051             :     real (kind=r8),      intent(inout) :: v(np,np,vlyr)
        1052             :     integer,             intent(in)    :: kptr
        1053             :     integer,             intent(in)    :: ielem
        1054             :     integer, optional,   intent(in)    :: rank
        1055             :
        1056             :     ! Local
        1057             :     integer :: i,k,ll,iptr
        1058             :     integer :: is,ie,in,iw,edgeptr
        1059             :     integer :: ise,isw,ine,inw
        1060             :     integer :: ks,ke,kblock
        1061             :     logical :: done
        1062             :
        1063   823548600 :     is=edge%getmap(south,ielem)
        1064   823548600 :     ie=edge%getmap(east,ielem)
        1065   823548600 :     in=edge%getmap(north,ielem)
        1066   823548600 :     iw=edge%getmap(west,ielem)
        1067   823548600 :     isw=edge%getmap(swest,ielem)
        1068   823548600 :     ise=edge%getmap(seast,ielem)
        1069   823548600 :     inw=edge%getmap(nwest,ielem)
        1070   823548600 :     ine=edge%getmap(neast,ielem)
        1071             :
        1072             :     !DIR$ IVDEP
        1073 68500355400 :     do k=1,vlyr
        1074 67676806800 :        iptr=np*(kptr+k-1)
        1075 >33920*10^7 :        do i=1,np
        1076 >27070*10^7 :           v(np ,i  ,k) = v(np ,i  ,k)+edge%receive(iptr+i+ie) ! East
        1077 >27070*10^7 :           v(i  ,1  ,k) = v(i  ,1  ,k)+edge%receive(iptr+i+is) ! South
        1078 >27070*10^7 :           v(i  ,np ,k) = v(i  ,np ,k)+edge%receive(iptr+i+in) ! North
        1079 >33838*10^7 :           v(1  ,i  ,k) = v(1  ,i  ,k)+edge%receive(iptr+i+iw) ! West
        1080             :        enddo
        1081             :     enddo
        1082             :
        1083             : ! SWEST
        1084  1647097200 :     do ll=swest,swest+max_corner_elem-1
        1085  1647097200 :         if(edge%getmap(ll,ielem) /= -1) then
        1086   822633546 :             edgeptr=edge%getmap(ll,ielem)+1
        1087 68424243894 :             do k=1,vlyr
        1088 67601610348 :                 iptr = (kptr+k-1)+edgeptr
        1089 68424243894 :                 v(1  ,1 ,k)=v(1 ,1 ,k)+edge%receive(iptr)
        1090             :             enddo
        1091             :         endif
        1092             :     end do
        1093             :
        1094             : ! SEAST
        1095  1647097200 :     do ll=swest+max_corner_elem,swest+2*max_corner_elem-1
        1096  1647097200 :         if(edge%getmap(ll,ielem) /= -1) then
        1097   822633546 :             edgeptr=edge%getmap(ll,ielem)+1
        1098 68424243894 :             do k=1,vlyr
        1099 67601610348 :                 iptr = (kptr+k-1)+edgeptr
        1100 68424243894 :                 v(np ,1 ,k)=v(np,1 ,k)+edge%receive(iptr)
        1101             :             enddo
        1102             :         endif
        1103             :     end do
        1104             :
        1105             : ! NEAST
        1106  1647097200 :     do ll=swest+3*max_corner_elem,swest+4*max_corner_elem-1
        1107  1647097200 :         if(edge%getmap(ll,ielem) /= -1) then
        1108   822633546 :             edgeptr=edge%getmap(ll,ielem)+1
        1109 68424243894 :             do k=1,vlyr
        1110 67601610348 :                 iptr = (kptr+k-1)+edgeptr
        1111 68424243894 :                 v(np ,np,k)=v(np,np,k)+edge%receive(iptr)
        1112             :             enddo
        1113             :         endif
        1114             :     end do
        1115             :
        1116             : ! NWEST
        1117  1647097200 :     do ll=swest+2*max_corner_elem,swest+3*max_corner_elem-1
        1118  1647097200 :         if(edge%getmap(ll,ielem) /= -1) then
        1119   822633546 :             edgeptr=edge%getmap(ll,ielem)+1
        1120 68424243894 :             do k=1,vlyr
        1121 67601610348 :                 iptr = (kptr+k-1)+edgeptr
        1122 68424243894 :                 v(1  ,np,k)=v(1 ,np,k)+edge%receive(iptr)
        1123             :             enddo
        1124             :         endif
        1125             :     end do
        1126             :
        1127             :
        1128   823548600 :   end subroutine edgeVunpack
        1129             : !
        1130           0 :   subroutine edgeVunpackVert(edge,v,ielem)
        1131             :     use control_mod,            only: north, south, east, west, neast, nwest, seast, swest
        1132             :     use dimensions_mod,         only: np, max_corner_elem, ne
        1133             :     use coordinate_systems_mod, only: cartesian3D_t
        1134             :
        1135             :     type (EdgeBuffer_t),  intent(inout) :: edge
        1136             :     type (cartesian3D_t), intent(inout)   :: v(:,:,:)
        1137             :     integer,              intent(in)    :: ielem
        1138             :
        1139             :     ! Local
        1140             :     logical, parameter :: UseUnroll = .TRUE.
        1141             :     integer            :: i,k,l, nce
        1142             :     integer            :: is,ie,in,iw,ine,inw,isw,ise
        1143             :
        1144           0 :     threadsafe=.false.
        1145             :
        1146           0 :     if (max_corner_elem.ne.1 .and. ne==0) then
        1147             :         ! MNL: this is used to construct the dual grid on the cube,
        1148             :         !      currently only supported for the uniform grid. If
        1149             :         !      this is desired on a refined grid, a little bit of
        1150             :         !      work will be required.
        1151           0 :         call endrun("edgeVunpackVert should not be called with unstructured meshes")
        1152             :     end if
        1153             :
        1154           0 :     is=edge%getmap(south,ielem)
        1155           0 :     ie=edge%getmap(east,ielem)
        1156           0 :     in=edge%getmap(north,ielem)
        1157           0 :     iw=edge%getmap(west,ielem)
        1158             :
        1159             :
        1160             :     ! N+S
        1161           0 :     do i=1,np/2
        1162             :        ! North
        1163           0 :        v(3,i ,np)%x = edge%receive(in+i)
        1164           0 :        v(3,i ,np)%y = edge%receive(np+in+i)
        1165           0 :        v(3,i ,np)%z = edge%receive(2*np+in+i)
        1166             :
        1167             :        ! South
        1168           0 :        v(2,i ,1)%x  = edge%receive(is+i)
        1169           0 :        v(2,i ,1)%y  = edge%receive(np+is+i)
        1170           0 :        v(2,i ,1)%z  = edge%receive(2*np+is+i)
        1171             :     enddo
        1172             :
        1173           0 :     do i=np/2+1,np
        1174             :        ! North
        1175           0 :        v(4,i ,np)%x = edge%receive(in+i)
        1176           0 :        v(4,i ,np)%y = edge%receive(np+in+i)
        1177           0 :        v(4,i ,np)%z = edge%receive(2*np+in+i)
        1178             :        ! South
        1179           0 :        v(1,i ,1)%x  = edge%receive(is+i)
        1180           0 :        v(1,i ,1)%y  = edge%receive(np+is+i)
        1181           0 :        v(1,i ,1)%z  = edge%receive(2*np+is+i)
        1182             :     enddo
        1183             :
        1184           0 :     do i=1,np/2
        1185             :        ! East
        1186           0 :        v(3,np,i)%x = edge%receive(ie+i)
        1187           0 :        v(3,np,i)%y = edge%receive(np+ie+i)
        1188           0 :        v(3,np,i)%z = edge%receive(2*np+ie+i)
        1189             :        ! West
        1190           0 :        v(4,1,i)%x  = edge%receive(iw+i)
        1191           0 :        v(4,1,i)%y  = edge%receive(np+iw+i)
        1192           0 :        v(4,1,i)%z  = edge%receive(2*np+iw+i)
        1193             :     end do
        1194             :
        1195           0 :     do i=np/2+1,np
        1196             :        ! East
        1197           0 :        v(2,np,i)%x = edge%receive(ie+i)
        1198           0 :        v(2,np,i)%y = edge%receive(np+ie+i)
        1199           0 :        v(2,np,i)%z = edge%receive(2*np+ie+i)
        1200             :        ! West
        1201           0 :        v(1,1,i)%x  = edge%receive(iw+i)
        1202           0 :        v(1,1,i)%y  = edge%receive(np+iw+i)
        1203           0 :        v(1,1,i)%z  = edge%receive(2*np+iw+i)
        1204             :     end do
        1205             :
        1206             : ! SWEST
        1207           0 :     nce = max_corner_elem
        1208           0 :     do l=swest,swest+max_corner_elem-1
        1209             :        ! find the one active corner, then exist
        1210           0 :         isw=edge%getmap(l,ielem)
        1211           0 :         if(isw /= -1) then
        1212           0 :             v(1,1,1)%x=edge%receive(isw+1)
        1213           0 :             v(1,1,1)%y=edge%receive(nce+isw+1)
        1214           0 :             v(1,1,1)%z=edge%receive(2*nce+isw+1)
        1215           0 :             exit
        1216             :         else
        1217           0 :             v(1,1,1)%x=0_r8
        1218           0 :             v(1,1,1)%y=0_r8
        1219           0 :             v(1,1,1)%z=0_r8
        1220             :         endif
        1221             :     end do
        1222             :
        1223             : ! SEAST
        1224           0 :     do l=swest+max_corner_elem,swest+2*max_corner_elem-1
        1225             :        ! find the one active corner, then exist
        1226           0 :         ise=edge%getmap(l,ielem)
        1227           0 :         if(ise /= -1) then
        1228           0 :             v(2,np,1)%x=edge%receive(ise+1)
        1229           0 :             v(2,np,1)%y=edge%receive(nce+ise+1)
        1230           0 :             v(2,np,1)%z=edge%receive(2*nce+ise+1)
        1231           0 :             exit
        1232             :         else
        1233           0 :             v(2,np,1)%x=0_r8
        1234           0 :             v(2,np,1)%y=0_r8
        1235           0 :             v(2,np,1)%z=0_r8
        1236             :         endif
        1237             :     end do
        1238             :
        1239             : ! NEAST
        1240           0 :     do l=swest+3*max_corner_elem,swest+4*max_corner_elem-1
        1241             :        ! find the one active corner, then exist
        1242           0 :         ine=edge%getmap(l,ielem)
        1243           0 :         if(ine /= -1) then
        1244           0 :             v(3,np,np)%x=edge%receive(ine+1)
        1245           0 :             v(3,np,np)%y=edge%receive(nce+ine+1)
        1246           0 :             v(3,np,np)%z=edge%receive(2*nce+ine+1)
        1247           0 :             exit
        1248             :         else
        1249           0 :             v(3,np,np)%x=0_r8
        1250           0 :             v(3,np,np)%y=0_r8
        1251           0 :             v(3,np,np)%z=0_r8
        1252             :         endif
        1253             :     end do
        1254             :
        1255             : ! NWEST
        1256           0 :     do l=swest+2*max_corner_elem,swest+3*max_corner_elem-1
        1257             :        ! find the one active corner, then exist
        1258           0 :         inw = edge%getmap(l,ielem)
        1259           0 :         if(inw/= -1) then
        1260           0 :             v(4,1,np)%x=edge%receive(inw+1)
        1261           0 :             v(4,1,np)%y=edge%receive(nce+inw+1)
        1262           0 :             v(4,1,np)%z=edge%receive(2*nce+inw+1)
        1263           0 :             exit
        1264             :         else
        1265           0 :             v(4,1,np)%x=0_r8
        1266           0 :             v(4,1,np)%y=0_r8
        1267           0 :             v(4,1,np)%z=0_r8
        1268             :         endif
        1269             :     end do
        1270             :
        1271             :     ! Fill the missing vertex info
        1272             :
        1273           0 :     do i=2,np/2
        1274             :        ! North
        1275           0 :        v(4,i ,np)%x = v(3,i-1 ,np)%x
        1276           0 :        v(4,i ,np)%y = v(3,i-1 ,np)%y
        1277           0 :        v(4,i ,np)%z = v(3,i-1 ,np)%z
        1278             :        ! South
        1279           0 :        v(1,i ,1)%x  = v(2,i-1 ,1)%x
        1280           0 :        v(1,i ,1)%y  = v(2,i-1 ,1)%y
        1281           0 :        v(1,i ,1)%z  = v(2,i-1 ,1)%z
        1282             :     enddo
        1283             :
        1284           0 :     do i=np/2+1,np-1
        1285             :        ! North
        1286           0 :        v(3,i ,np)%x = v(4,i+1 ,np)%x
        1287           0 :        v(3,i ,np)%y = v(4,i+1 ,np)%y
        1288           0 :        v(3,i ,np)%z = v(4,i+1 ,np)%z
        1289             :        ! South
        1290           0 :        v(2,i ,1)%x  = v(1,i+1 ,1)%x
        1291           0 :        v(2,i ,1)%y  = v(1,i+1 ,1)%y
        1292           0 :        v(2,i ,1)%z  = v(1,i+1 ,1)%z
        1293             :     enddo
        1294             :
        1295           0 :     do i=2,np/2
        1296             :        ! East
        1297           0 :        v(2,np,i)%x = v(3,np,i-1)%x
        1298           0 :        v(2,np,i)%y = v(3,np,i-1)%y
        1299           0 :        v(2,np,i)%z = v(3,np,i-1)%z
        1300             :        ! West
        1301           0 :        v(1,1,i)%x  = v(4,1,i-1)%x
        1302           0 :        v(1,1,i)%y  = v(4,1,i-1)%y
        1303           0 :        v(1,1,i)%z  = v(4,1,i-1)%z
        1304             :     end do
        1305             :
        1306           0 :     do i=np/2+1,np-1
        1307             :        ! East
        1308           0 :        v(3,np,i)%x = v(2,np,i+1)%x
        1309           0 :        v(3,np,i)%y = v(2,np,i+1)%y
        1310           0 :        v(3,np,i)%z = v(2,np,i+1)%z
        1311             :        ! West
        1312           0 :        v(4,1,i)%x  = v(1,1,i+1)%x
        1313           0 :        v(4,1,i)%y  = v(1,1,i+1)%y
        1314           0 :        v(4,1,i)%z  = v(1,1,i+1)%z
        1315             :     end do
        1316             :
        1317           0 :   end subroutine edgeVunpackVert
        1318             :   ! ========================================
        1319             :   ! edgeDGVunpack:
        1320             :   !
        1321             :   ! Unpack edges from edge buffer into v...
        1322             :   ! ========================================
        1323             :
        1324   124675200 :   subroutine edgeDGVunpack(edge,v,vlyr,kptr,ielem)
        1325             :     use dimensions_mod, only: np, max_corner_elem
        1326             :     use control_mod,    only: north, south, east, west, neast, nwest, seast, swest
        1327             :
        1328             :     type (EdgeBuffer_t), intent(in)    :: edge
        1329             :     integer,             intent(in)    :: vlyr
        1330             :     real (kind=r8),      intent(inout) :: v(0:np+1,0:np+1,vlyr)
        1331             :     integer,             intent(in)    :: kptr
        1332             :     integer,             intent(in)    :: ielem
        1333             :
        1334             :     ! Local
        1335             :     integer :: i,k,iptr
        1336             :     integer :: is,ie,in,iw
        1337             :
        1338   124675200 :     threadsafe=.false.
        1339             :
        1340   124675200 :     is=edge%getmap(south,ielem)
        1341   124675200 :     ie=edge%getmap(east,ielem)
        1342   124675200 :     in=edge%getmap(north,ielem)
        1343   124675200 :     iw=edge%getmap(west,ielem)
        1344  7558434000 :     do k=1,vlyr
        1345  7433758800 :        iptr=np*(kptr+k-1)
        1346 37293469200 :        do i=1,np
        1347 29735035200 :           v(i   ,0   ,k)=edge%receive(iptr+is+i)
        1348 29735035200 :           v(np+1,i   ,k)=edge%receive(iptr+ie+i)
        1349 29735035200 :           v(i   ,np+1,k)=edge%receive(iptr+in+i)
        1350 37168794000 :           v(0   ,i   ,k)=edge%receive(iptr+iw+i)
        1351             :        end do
        1352             :     end do
        1353             :
        1354   124675200 :     i = swest
        1355   124675200 :     if(edge%getmap(i,ielem) /= -1) then
        1356  7550035740 :       do k=1,vlyr
        1357  7425499068 :         iptr=(kptr+k-1)
        1358  7550035740 :         v(0,0,k) = edge%receive(iptr+edge%getmap(i,ielem)+1)
        1359             :       end do
        1360             :     end if
        1361   124675200 :     i = swest+max_corner_elem
        1362   124675200 :     if(edge%getmap(i,ielem) /= -1) then
        1363  7550035740 :       do k=1,vlyr
        1364  7425499068 :         iptr=(kptr+k-1)
        1365  7550035740 :         v(np+1,0,k) = edge%receive(iptr+edge%getmap(i,ielem)+1)
        1366             :       end do
        1367             :     end if
        1368   124675200 :     i = swest+3*max_corner_elem
        1369   124675200 :     if(edge%getmap(i,ielem) /= -1) then
        1370  7550035740 :       do k=1,vlyr
        1371  7425499068 :         iptr=(kptr+k-1)
        1372  7550035740 :         v(np+1,np+1,k) = edge%receive(iptr+edge%getmap(i,ielem)+1)
        1373             :       end do
        1374             :     end if
        1375   124675200 :     i = swest+2*max_corner_elem
        1376   124675200 :     if(edge%getmap(i,ielem) /= -1) then
        1377  7550035740 :       do k=1,vlyr
        1378  7425499068 :         iptr=(kptr+k-1)
        1379  7550035740 :         v(0,np+1,k) = edge%receive(iptr+edge%getmap(i,ielem)+1)
        1380             :       end do
        1381             :     end if
        1382             :
        1383   124675200 :   end subroutine edgeDGVunpack
        1384             :
        1385             :   ! ========================================
        1386             :   ! edgeVunpackMIN/MAX:
        1387             :   !
        1388             :   ! Finds the Min/Max edges from edge buffer into v...
        1389             :   ! ========================================
        1390           0 :   subroutine edgeVunpackMAX(edge,v,vlyr,kptr,ielem)
        1391             :     use dimensions_mod, only: np, max_corner_elem
        1392             :     use control_mod,    only: north, south, east, west, neast, nwest, seast, swest
        1393             :
        1394             :     type (EdgeBuffer_t), intent(in)    :: edge
        1395             :     integer,             intent(in)    :: vlyr
        1396             :     real (kind=r8),      intent(inout) :: v(np,np,vlyr)
        1397             :     integer,             intent(in)    :: kptr
        1398             :     integer,             intent(in)    :: ielem
        1399             :
        1400             :     ! Local
        1401             :     integer :: i,k,l,iptr
        1402             :     integer :: is,ie,in,iw
        1403             :
        1404           0 :     threadsafe=.false.
        1405             :
        1406           0 :     is=edge%getmap(south,ielem)
        1407           0 :     ie=edge%getmap(east,ielem)
        1408           0 :     in=edge%getmap(north,ielem)
        1409           0 :     iw=edge%getmap(west,ielem)
        1410           0 :     do k=1,vlyr
        1411           0 :        iptr=np*(kptr+k-1)
        1412           0 :        do i=1,np
        1413           0 :           v(np ,i  ,k) = MAX(v(np ,i  ,k),edge%receive(iptr+ie+i  ))
        1414           0 :           v(i  ,1  ,k) = MAX(v(i  ,1  ,k),edge%receive(iptr+is+i  ))
        1415           0 :           v(i  ,np ,k) = MAX(v(i  ,np ,k),edge%receive(iptr+in+i  ))
        1416           0 :           v(1  ,i  ,k) = MAX(v(1  ,i  ,k),edge%receive(iptr+iw+i  ))
        1417             :        end do
        1418             :     end do
        1419             :
        1420             : ! SWEST
        1421           0 :     do l=swest,swest+max_corner_elem-1
        1422           0 :         if(edge%getmap(l,ielem) /= -1) then
        1423           0 :             do k=1,vlyr
        1424           0 :                 v(1  ,1 ,k)=MAX(v(1 ,1 ,k),edge%receive((kptr+k-1)+edge%getmap(l,ielem)+1))
        1425             :             enddo
        1426             :         endif
        1427             :     end do
        1428             :
        1429             : ! SEAST
        1430           0 :     do l=swest+max_corner_elem,swest+2*max_corner_elem-1
        1431           0 :         if(edge%getmap(l,ielem) /= -1) then
        1432           0 :             do k=1,vlyr
        1433           0 :                 v(np ,1 ,k)=MAX(v(np,1 ,k),edge%receive((kptr+k-1)+edge%getmap(l,ielem)+1))
        1434             :             enddo
        1435             :         endif
        1436             :     end do
        1437             :
        1438             : ! NEAST
        1439           0 :     do l=swest+3*max_corner_elem,swest+4*max_corner_elem-1
        1440           0 :         if(edge%getmap(l,ielem) /= -1) then
        1441           0 :             do k=1,vlyr
        1442           0 :                 v(np ,np,k)=MAX(v(np,np,k),edge%receive((kptr+k-1)+edge%getmap(l,ielem)+1))
        1443             :             enddo
        1444             :         endif
        1445             :     end do
        1446             :
        1447             : ! NWEST
        1448           0 :     do l=swest+2*max_corner_elem,swest+3*max_corner_elem-1
        1449           0 :         if(edge%getmap(l,ielem) /= -1) then
        1450           0 :             do k=1,vlyr
        1451           0 :                 v(1  ,np,k)=MAX(v(1 ,np,k),edge%receive((kptr+k-1)+edge%getmap(l,ielem)+1))
        1452             :             enddo
        1453             :         endif
        1454             :     end do
        1455             :
        1456           0 :   end subroutine edgeVunpackMAX
        1457             :
        1458           0 :   subroutine edgeSunpackMAX(edge,v,vlyr,kptr,ielem)
        1459             :     use dimensions_mod, only: np, max_corner_elem
        1460             :     use control_mod,    only: north, south, east, west, neast, nwest, seast, swest
        1461             :
        1462             :     type (EdgeBuffer_t), intent(in)    :: edge
        1463             :     integer,             intent(in)    :: vlyr
        1464             :     real (kind=r8),      intent(inout) :: v(vlyr)
        1465             :     integer,             intent(in)    :: kptr
        1466             :     integer,             intent(in)    :: ielem
        1467             :
        1468             :     ! Local
        1469             :     integer :: i,k,l,iptr
        1470             :     integer :: is,ie,in,iw,edgeptr
        1471             :
        1472           0 :     threadsafe=.false.
        1473             :
        1474           0 :     is=edge%getmap(south,ielem)
        1475           0 :     ie=edge%getmap(east,ielem)
        1476           0 :     in=edge%getmap(north,ielem)
        1477           0 :     iw=edge%getmap(west,ielem)
        1478           0 :     do k=1,vlyr
        1479           0 :        iptr=(kptr+k-1)
        1480           0 :        v(k) = MAX(v(k),edge%receive(iptr+is+1),edge%receive(iptr+ie+1),edge%receive(iptr+in+1),edge%receive(iptr+iw+1))
        1481             :     end do
        1482             :
        1483             : ! SWEST
        1484           0 :     do l=swest,swest+max_corner_elem-1
        1485           0 :         if(edge%getmap(l,ielem) /= -1) then
        1486           0 :             edgeptr = edge%getmap(l,ielem)+1
        1487           0 :             do k=1,vlyr
        1488           0 :                 iptr = (kptr+k-1)+edgeptr
        1489           0 :                 v(k)=MAX(v(k),edge%receive(iptr))
        1490             :             enddo
        1491             :         endif
        1492             :     end do
        1493             :
        1494             : ! SEAST
        1495           0 :     do l=swest+max_corner_elem,swest+2*max_corner_elem-1
        1496           0 :         if(edge%getmap(l,ielem) /= -1) then
        1497           0 :             edgeptr = edge%getmap(l,ielem)+1
        1498           0 :             do k=1,vlyr
        1499           0 :                 iptr = (kptr+k-1)+edgeptr
        1500           0 :                 v(k)=MAX(v(k),edge%receive(iptr))
        1501             :             enddo
        1502             :         endif
        1503             :     end do
        1504             :
        1505             : ! NEAST
        1506           0 :     do l=swest+3*max_corner_elem,swest+4*max_corner_elem-1
        1507           0 :         if(edge%getmap(l,ielem) /= -1) then
        1508           0 :             edgeptr = edge%getmap(l,ielem)+1
        1509           0 :             do k=1,vlyr
        1510           0 :                 iptr = (kptr+k-1)+edgeptr
        1511           0 :                 v(k)=MAX(v(k),edge%receive(iptr))
        1512             :             enddo
        1513             :         endif
        1514             :     end do
        1515             :
        1516             : ! NWEST
        1517           0 :     do l=swest+2*max_corner_elem,swest+3*max_corner_elem-1
        1518           0 :         if(edge%getmap(l,ielem) /= -1) then
        1519           0 :             edgeptr = edge%getmap(l,ielem)+1
        1520           0 :             do k=1,vlyr
        1521           0 :                 iptr = (kptr+k-1)+edgeptr
        1522           0 :                 v(k)=MAX(v(k),edge%receive(iptr))
        1523             :             enddo
        1524             :         endif
        1525             :     end do
        1526             :
        1527           0 :   end subroutine edgeSunpackMAX
        1528             :
        1529           0 :   subroutine edgeSunpackMIN(edge,v,vlyr,kptr,ielem)
        1530             :     use dimensions_mod, only: np, max_corner_elem
        1531             :     use control_mod,    only: north, south, east, west, neast, nwest, seast, swest
        1532             :
        1533             :     type (EdgeBuffer_t), intent(in)    :: edge
        1534             :     integer,             intent(in)    :: vlyr
        1535             :     real (kind=r8),      intent(inout) :: v(vlyr)
        1536             :     integer,             intent(in)    :: kptr
        1537             :     integer,             intent(in)    :: ielem
        1538             :
        1539             :     ! Local
        1540             :     integer :: i,k,l,iptr
        1541             :     integer :: is,ie,in,iw,edgeptr
        1542             :
        1543           0 :     threadsafe=.false.
        1544             :
        1545           0 :     is=edge%getmap(south,ielem)
        1546           0 :     ie=edge%getmap(east,ielem)
        1547           0 :     in=edge%getmap(north,ielem)
        1548           0 :     iw=edge%getmap(west,ielem)
        1549           0 :     do k=1,vlyr
        1550           0 :        iptr=(kptr+k-1)
        1551           0 :        v(k) = MIN(v(k),edge%receive(iptr+is+1),edge%receive(iptr+ie+1),edge%receive(iptr+in+1),edge%receive(iptr+iw+1))
        1552             :     end do
        1553             :
        1554             : ! SWEST
        1555           0 :     do l=swest,swest+max_corner_elem-1
        1556           0 :         if(edge%getmap(l,ielem) /= -1) then
        1557           0 :             edgeptr = edge%getmap(l,ielem)+1
        1558           0 :             do k=1,vlyr
        1559           0 :                 iptr = (kptr+k-1)+edgeptr
        1560           0 :                 v(k)=MiN(v(k),edge%receive(iptr))
        1561             :             enddo
        1562             :         endif
        1563             :     end do
        1564             :
        1565             : ! SEAST
        1566           0 :     do l=swest+max_corner_elem,swest+2*max_corner_elem-1
        1567           0 :         if(edge%getmap(l,ielem) /= -1) then
        1568           0 :             edgeptr = edge%getmap(l,ielem)+1
        1569           0 :             do k=1,vlyr
        1570           0 :                 iptr = (kptr+k-1)+edgeptr
        1571           0 :                 v(k)=MIN(v(k),edge%receive(iptr))
        1572             :             enddo
        1573             :         endif
        1574             :     end do
        1575             :
        1576             : ! NEAST
        1577           0 :     do l=swest+3*max_corner_elem,swest+4*max_corner_elem-1
        1578           0 :         if(edge%getmap(l,ielem) /= -1) then
        1579           0 :             edgeptr = edge%getmap(l,ielem)+1
        1580           0 :             do k=1,vlyr
        1581           0 :                 iptr = (kptr+k-1)+edgeptr
        1582           0 :                 v(k)=MIN(v(k),edge%receive(iptr))
        1583             :             enddo
        1584             :         endif
        1585             :     end do
        1586             :
        1587             : ! NWEST
        1588           0 :     do l=swest+2*max_corner_elem,swest+3*max_corner_elem-1
        1589           0 :         if(edge%getmap(l,ielem) /= -1) then
        1590           0 :             edgeptr = edge%getmap(l,ielem)+1
        1591           0 :             do k=1,vlyr
        1592           0 :                 iptr = (kptr+k-1)+edgeptr
        1593           0 :                 v(k)=MIN(v(k),edge%receive(iptr))
        1594             :             enddo
        1595             :         endif
        1596             :     end do
        1597             :
        1598           0 :   end subroutine edgeSunpackMIN
        1599             :
        1600           0 :   subroutine edgeVunpackMIN(edge,v,vlyr,kptr,ielem)
        1601             :     use dimensions_mod, only: np, max_corner_elem
        1602             :     use control_mod,    only: north, south, east, west, neast, nwest, seast, swest
        1603             :
        1604             :     type (EdgeBuffer_t), intent(in)    :: edge
        1605             :     integer,             intent(in)    :: vlyr
        1606             :     real (kind=r8),      intent(inout) :: v(np,np,vlyr)
        1607             :     integer,             intent(in)    :: kptr
        1608             :     integer,             intent(in)    :: ielem
        1609             :
        1610             :     ! Local
        1611             :     integer :: i,k,l,iptr
        1612             :     integer :: is,ie,in,iw,edgeptr
        1613             :
        1614           0 :     threadsafe=.false.
        1615             :
        1616           0 :     is=edge%getmap(south,ielem)
        1617           0 :     ie=edge%getmap(east,ielem)
        1618           0 :     in=edge%getmap(north,ielem)
        1619           0 :     iw=edge%getmap(west,ielem)
        1620           0 :     do k=1,vlyr
        1621           0 :        iptr = np*(kptr+k-1)
        1622           0 :        do i=1,np
        1623           0 :           v(np ,i  ,k) = MIN(v(np ,i  ,k),edge%receive(iptr+ie+i  ))
        1624           0 :           v(i  ,1  ,k) = MIN(v(i  ,1  ,k),edge%receive(iptr+is+i  ))
        1625           0 :           v(i  ,np ,k) = MIN(v(i  ,np ,k),edge%receive(iptr+in+i  ))
        1626           0 :           v(1  ,i  ,k) = MIN(v(1  ,i  ,k),edge%receive(iptr+iw+i  ))
        1627             :        end do
        1628             :     end do
        1629             :
        1630             : ! SWEST
        1631           0 :     do l=swest,swest+max_corner_elem-1
        1632           0 :         if(edge%getmap(l,ielem) /= -1) then
        1633           0 :             edgeptr=edge%getmap(l,ielem)+1
        1634           0 :             do k=1,vlyr
        1635           0 :                 iptr=(kptr+k-1)+edgeptr
        1636           0 :                 v(1  ,1 ,k)=MIN(v(1 ,1 ,k),edge%receive(iptr))
        1637             :             enddo
        1638             :         endif
        1639             :     end do
        1640             :
        1641             : ! SEAST
        1642           0 :     do l=swest+max_corner_elem,swest+2*max_corner_elem-1
        1643           0 :         if(edge%getmap(l,ielem) /= -1) then
        1644           0 :             edgeptr=edge%getmap(l,ielem)+1
        1645           0 :             do k=1,vlyr
        1646           0 :                 iptr=(kptr+k-1)+edgeptr
        1647           0 :                 v(np ,1 ,k)=MIN(v(np,1 ,k),edge%receive(iptr))
        1648             :             enddo
        1649             :         endif
        1650             :     end do
        1651             :
        1652             : ! NEAST
        1653           0 :     do l=swest+3*max_corner_elem,swest+4*max_corner_elem-1
        1654           0 :         if(edge%getmap(l,ielem) /= -1) then
        1655           0 :             edgeptr=edge%getmap(l,ielem)+1
        1656           0 :             do k=1,vlyr
        1657           0 :                 iptr=(kptr+k-1)+edgeptr
        1658           0 :                 v(np ,np,k)=MIN(v(np,np,k),edge%receive(iptr))
        1659             :             enddo
        1660             :         endif
        1661             :     end do
        1662             :
        1663             : ! NWEST
        1664           0 :     do l=swest+2*max_corner_elem,swest+3*max_corner_elem-1
        1665           0 :         if(edge%getmap(l,ielem) /= -1) then
        1666           0 :             edgeptr=edge%getmap(l,ielem)+1
        1667           0 :             do k=1,vlyr
        1668           0 :                 iptr=(kptr+k-1)+edgeptr
        1669           0 :                 v(1  ,np,k)=MIN(v(1 ,np,k),edge%receive(iptr))
        1670             :             enddo
        1671             :         endif
        1672             :     end do
        1673             :
        1674           0 :   end subroutine edgeVunpackMIN
        1675             :
        1676             :   ! ========================================
        1677             :   ! LongEdgeVunpackMIN:
        1678             :   !
        1679             :   ! Finds the Min edges from edge buffer into v...
        1680             :   ! ========================================
        1681       10800 :   subroutine LongEdgeVunpackMIN(edge,v,vlyr,kptr,desc)
        1682             :     use control_mod,    only: north, south, east, west, neast, nwest, seast, swest
        1683             :     use dimensions_mod, only: np, max_corner_elem
        1684             :
        1685             :     type (LongEdgeBuffer_t), intent(in)    :: edge
        1686             :     integer,                 intent(in)    :: vlyr
        1687             :     integer ,                intent(inout) :: v(np,np,vlyr)
        1688             :     integer,                 intent(in)    :: kptr
        1689             :     type (EdgeDescriptor_t), intent(in)    :: desc
        1690             :
        1691             :     ! Local
        1692             :
        1693             :     integer :: i,k,l
        1694             :     integer :: is,ie,in,iw
        1695             :
        1696       10800 :     threadsafe=.false.
        1697             :
        1698       10800 :     is=desc%getmapP(south)
        1699       10800 :     ie=desc%getmapP(east)
        1700       10800 :     in=desc%getmapP(north)
        1701       10800 :     iw=desc%getmapP(west)
        1702       21600 :     do k=1,vlyr
        1703       64800 :        do i=1,np
        1704       43200 :           v(i  ,1  ,k) = MIN(v(i  ,1  ,k),edge%buf(kptr+k,is+i  ))
        1705       43200 :           v(np ,i  ,k) = MIN(v(np ,i  ,k),edge%buf(kptr+k,ie+i  ))
        1706       43200 :           v(i  ,np ,k) = MIN(v(i  ,np ,k),edge%buf(kptr+k,in+i  ))
        1707       54000 :           v(1  ,i  ,k) = MIN(v(1  ,i  ,k),edge%buf(kptr+k,iw+i  ))
        1708             :        end do
        1709             :     end do
        1710             :
        1711             : ! SWEST
        1712       21600 :     do l=swest,swest+max_corner_elem-1
        1713       21600 :         if(desc%getmapP(l) /= -1) then
        1714       21576 :             do k=1,vlyr
        1715       21576 :                 v(1  ,1 ,k)=MIN(v(1 ,1 ,k),edge%buf(kptr+k,desc%getmapP(l)+1))
        1716             :             enddo
        1717             :         endif
        1718             :     end do
        1719             :
        1720             : ! SEAST
        1721       21600 :     do l=swest+max_corner_elem,swest+2*max_corner_elem-1
        1722       21600 :         if(desc%getmapP(l) /= -1) then
        1723       21576 :             do k=1,vlyr
        1724       21576 :                 v(np ,1 ,k)=MIN(v(np,1 ,k),edge%buf(kptr+k,desc%getmapP(l)+1))
        1725             :             enddo
        1726             :         endif
        1727             :     end do
        1728             :
        1729             : ! NEAST
        1730       21600 :     do l=swest+3*max_corner_elem,swest+4*max_corner_elem-1
        1731       21600 :         if(desc%getmapP(l) /= -1) then
        1732       21576 :             do k=1,vlyr
        1733       21576 :                 v(np ,np,k)=MIN(v(np,np,k),edge%buf(kptr+k,desc%getmapP(l)+1))
        1734             :             enddo
        1735             :         endif
        1736             :     end do
        1737             :
        1738             : ! NWEST
        1739       21600 :     do l=swest+2*max_corner_elem,swest+3*max_corner_elem-1
        1740       21600 :         if(desc%getmapP(l) /= -1) then
        1741       21576 :             do k=1,vlyr
        1742       21576 :                 v(1  ,np,k)=MIN(v(1 ,np,k),edge%buf(kptr+k,desc%getmapP(l)+1))
        1743             :             enddo
        1744             :         endif
        1745             :     end do
        1746             :
        1747       10800 :   end subroutine LongEdgeVunpackMIN
        1748             :
        1749             :
        1750   138018600 : subroutine ghostpack(edge,v,vlyr,kptr,ielem)
        1751             :
        1752             :   use dimensions_mod, only : max_corner_elem
        1753             :   use control_mod, only : north, south, east, west, neast, nwest, seast, swest
        1754             :   use edgetype_mod, only : EdgeDescriptor_t
        1755             :
        1756             :   implicit none
        1757             :
        1758             :   type (Edgebuffer_t)                :: edge
        1759             :   integer,              intent(in)   :: vlyr
        1760             :   integer,              intent(in)   :: kptr
        1761             :
        1762             :   real (kind=r8),intent(in)   :: v(edge%lb:edge%ub,edge%lb:edge%ub,vlyr)
        1763             :   integer,              intent(in)   :: ielem
        1764             :
        1765             :   ! Local variables
        1766             :   integer :: i,j,k,ir,l,itr,ktmp
        1767             :
        1768             :   integer :: is,ie,in,iw,isw,ise,inw,ine
        1769             :   integer :: nhc, npoints
        1770             :   integer :: edgeptr,iptr
        1771             :
        1772   138018600 :   is = edge%putmap(south,ielem)
        1773   138018600 :   ie = edge%putmap(east,ielem)
        1774   138018600 :   in = edge%putmap(north,ielem)
        1775   138018600 :   iw = edge%putmap(west,ielem)
        1776   138018600 :   if (edge%nlyr < (kptr+vlyr) ) then
        1777           0 :        print *,'edge%nlyr = ',edge%nlyr
        1778           0 :        print *,'kptr+vlyr = ',kptr+vlyr
        1779           0 :        call endrun('ghostpack: Buffer overflow: size of the vertical dimension must be increased!')
        1780             :   endif
        1781             :
        1782             :
        1783   138018600 :   nhc     = edge%ndepth
        1784   138018600 :   npoints = edge%npoints
        1785             :
        1786             :     !DIR$ IVDEP
        1787 17289007200 :     do k=1,vlyr
        1788 17150988600 :       ktmp = nhc*(kptr+k-1)
        1789 53282237400 :       do j=1,nhc
        1790 35993230200 :         iptr = npoints*(ktmp + j - 1)
        1791 >16112*10^7 :         do i=1,npoints
        1792 >10797*10^7 :           edge%buf(iptr+is+i)   = v(i  ,j ,k)
        1793 >10797*10^7 :           edge%buf(iptr+ie+i)   = v(npoints-j+1 ,i ,k)
        1794 >10797*10^7 :           edge%buf(iptr+in+i)   = v(i  ,npoints-j+1,k)
        1795 >14397*10^7 :           edge%buf(iptr+iw+i)   = v(j  ,i ,k)
        1796             :         enddo
        1797             :       end do
        1798             :     end do
        1799             :
        1800             :
        1801             :   !  This is really kludgy way to setup the index reversals
        1802             :   !  But since it is so a rare event not real need to spend time optimizing
        1803             :   !  Check if the edge orientation of the recieving element is different
        1804             :   !  if it is, swap the order of data in the edge
        1805   138018600 :   if(edge%reverse(south,ielem)) then
        1806             :      !DIR$ IVDEP
        1807   288150120 :      do k=1,vlyr
        1808   285849810 :        ktmp = nhc*(kptr+k-1)
        1809   888037290 :        do j=1,nhc
        1810   599887170 :          iptr = npoints*(ktmp + j - 1)
        1811  2685398130 :          do i=1,npoints
        1812  1799661150 :            ir = npoints-i+1
        1813  2399548320 :            edge%buf(iptr+is+i)=v(ir,j,k)
        1814             :          enddo
        1815             :        enddo
        1816             :      enddo
        1817             :   endif
        1818             :
        1819   138018600 :   if(edge%reverse(east,ielem)) then
        1820             :      !DIR$ IVDEP
        1821    96050040 :      do k=1,vlyr
        1822    95283270 :        ktmp = nhc*(kptr+k-1)
        1823   296012430 :        do j=1,nhc
        1824   199962390 :          iptr = npoints*(ktmp + j - 1)
        1825   895132710 :          do i=1,npoints
        1826   599887050 :            ir = npoints-i+1
        1827   799849440 :            edge%buf(iptr+ie+i)=v(npoints-j+1,ir,k)
        1828             :           enddo
        1829             :         enddo
        1830             :       enddo
        1831             :   endif
        1832             :
        1833   138018600 :   if(edge%reverse(north,ielem)) then
        1834             :      !DIR$ IVDEP
        1835   288150120 :       do k=1,vlyr
        1836   285849810 :         ktmp = nhc*(kptr+k-1)
        1837   888037290 :         do j=1,nhc
        1838   599887170 :          iptr = npoints*(ktmp + j - 1)
        1839  2685398130 :           do i=1,npoints
        1840  1799661150 :             ir = npoints-i+1
        1841  2399548320 :             edge%buf(iptr+in+i)=v(ir,npoints-j+1,k)
        1842             :           enddo
        1843             :         enddo
        1844             :       enddo
        1845             :   endif
        1846             :
        1847   138018600 :   if(edge%reverse(west,ielem)) then
        1848             :      !DIR$ IVDEP
        1849    96050040 :      do k=1,vlyr
        1850    95283270 :        ktmp = nhc*(kptr+k-1)
        1851   296012430 :        do j=1,nhc
        1852   199962390 :          iptr = npoints*(ktmp + j - 1)
        1853   895132710 :          do i=1,npoints
        1854   599887050 :             ir = npoints-i+1
        1855   799849440 :             edge%buf(iptr+iw+i)=v(j,ir,k)
        1856             :           enddo
        1857             :         enddo
        1858             :       enddo
        1859             :   endif
        1860             :
        1861             :
        1862             :   ! corners.  this is difficult because we dont know the orientaton
        1863             :   ! of the corners, and this which (i,j) dimension maps to which dimension
        1864             : ! SWEST
        1865   276037200 :   do l=swest,swest+max_corner_elem-1
        1866   276037200 :      if (edge%putmap(l,ielem) /= -1) then
        1867             :          isw = edge%putmap(l,ielem)
        1868             :          !DIR$ IVDEP
        1869 17269797192 :          do k=1,vlyr
        1870 17131931946 :           ktmp = nhc*(kptr+k-1)
        1871 53223034914 :           do j=1,nhc
        1872 35953237722 :             iptr = nhc*(ktmp + j - 1)
        1873 >14550*10^7 :             do i=1,nhc
        1874 >12837*10^7 :               edge%buf(iptr+isw+i)=v(i  ,j ,k)
        1875             :             enddo
        1876             :           end do
        1877             :         end do
        1878             :      end if
        1879             :   end do
        1880             :
        1881             : ! SEAST
        1882   276037200 :   do l=swest+max_corner_elem,swest+2*max_corner_elem-1
        1883   276037200 :      if (edge%putmap(l,ielem) /= -1) then
        1884             :          ise = edge%putmap(l,ielem)
        1885             :          !DIR$ IVDEP
        1886 17269797192 :          do k=1,vlyr
        1887 17131931946 :           ktmp = nhc*(kptr+k-1)
        1888 53223034914 :           do j=1,nhc
        1889 35953237722 :             iptr = nhc*(ktmp + j - 1)
        1890 >14550*10^7 :             do i=1,nhc
        1891 >12837*10^7 :               edge%buf(iptr+ise+i)=v(npoints-i+1 ,j ,k)
        1892             :             enddo
        1893             :           end do
        1894             :         end do
        1895             :      end if
        1896             :   end do
        1897             :
        1898             : ! NEAST
        1899   276037200 :   do l=swest+3*max_corner_elem,swest+4*max_corner_elem-1
        1900   276037200 :      if (edge%putmap(l,ielem) /= -1) then
        1901             :         ine = edge%putmap(l,ielem)
        1902             :          !DIR$ IVDEP
        1903 17269797192 :         do k=1,vlyr
        1904 17131931946 :            ktmp = nhc*(kptr+k-1)
        1905 53223034914 :            do j=1,nhc
        1906 35953237722 :               iptr = nhc*(ktmp + j - 1)
        1907 >14550*10^7 :               do i=1,nhc
        1908 >12837*10^7 :                  edge%buf(iptr+ine+i)=v(npoints-i+1,npoints-j+1,k)
        1909             :               enddo
        1910             :             enddo
        1911             :           end do
        1912             :       end if
        1913             :   end do
        1914             :
        1915             : ! NWEST
        1916   276037200 :   do l=swest+2*max_corner_elem,swest+3*max_corner_elem-1
        1917   276037200 :      if (edge%putmap(l,ielem) /= -1) then
        1918             :         inw = edge%putmap(l,ielem)
        1919             :         !DIR$ IVDEP
        1920 17269797192 :         do k=1,vlyr
        1921 17131931946 :           ktmp = nhc*(kptr+k-1)
        1922 53223034914 :           do j=1,nhc
        1923 35953237722 :             iptr = nhc*(ktmp + j - 1)
        1924 >14550*10^7 :             do i=1,nhc
        1925 >12837*10^7 :               edge%buf(iptr+inw+i)=v(i  ,npoints-j+1,k)
        1926             :             enddo
        1927             :           end do
        1928             :         end do
        1929             :      end if
        1930             :   end do
        1931             :
        1932   138018600 : end subroutine ghostpack
        1933             :
        1934   138018600 : subroutine ghostunpack(edge,v,vlyr,kptr,ielem)
        1935             :   use dimensions_mod, only : max_corner_elem
        1936             :   use control_mod, only : north, south, east, west, neast, nwest, seast, swest
        1937             :   type (Edgebuffer_t),         intent(in)  :: edge
        1938             :
        1939             :   integer,               intent(in)  :: vlyr
        1940             :   integer,               intent(in)  :: kptr
        1941             :   integer,               intent(in)  :: ielem
        1942             :
        1943             :   real (kind=r8), intent(inout) :: v(edge%lb:edge%ub,edge%lb:edge%ub,vlyr)
        1944             :
        1945             :
        1946             :   ! Local
        1947             :   logical, parameter :: UseUnroll = .TRUE.
        1948             :   integer :: i,j,k,l,itr, ktmp
        1949             :   integer :: is,ie,in,iw,isw,ise,inw,ine
        1950             :   integer :: nhc,npoints,iptr
        1951             :   logical :: reverse
        1952             :
        1953   138018600 :   threadsafe=.false.
        1954             :
        1955   138018600 :   is=edge%getmap(south,ielem)
        1956   138018600 :   ie=edge%getmap(east,ielem)
        1957   138018600 :   in=edge%getmap(north,ielem)
        1958   138018600 :   iw=edge%getmap(west,ielem)
        1959             :
        1960   138018600 :   nhc     = edge%ndepth
        1961   138018600 :   npoints = edge%npoints
        1962             :
        1963             :   ! example for north buffer
        1964             :   ! first row ('edge') goes in v(:,np+1,k)
        1965             :   ! 2nd   row ('edge') goes in v(:,np+2,k)
        1966             :   ! etc...
        1967             :     !DIR$ IVDEP
        1968 17289007200 :     do k=1,vlyr
        1969 17150988600 :       ktmp = nhc*(kptr+k-1)
        1970 53282237400 :       do j=1,nhc
        1971 35993230200 :         iptr = npoints*(ktmp + j - 1)
        1972 >16112*10^7 :         do i=1,npoints
        1973 >10797*10^7 :           v(i  ,1-j  ,k)      = edge%receive(iptr+is+i)  ! South
        1974 >10797*10^7 :           v(npoints+j ,i  ,k) = edge%receive(iptr+ie+i)  ! East
        1975 >10797*10^7 :           v(i  ,npoints+j ,k) = edge%receive(iptr+in+i)  ! North
        1976 >14397*10^7 :           v(1-j  ,i  ,k)      = edge%receive(iptr+iw+i)  ! West
        1977             :         end do
        1978             :       end do
        1979             :     end do
        1980             :
        1981             :
        1982             : ! SWEST
        1983   276037200 :   do l=swest,swest+max_corner_elem-1
        1984   138018600 :      isw = edge%getmap(l,ielem)
        1985   276037200 :      if(isw /= -1) then
        1986             :         ! note the following is the the correct meaning of reverse in this code.
        1987             :         ! It is  best described as a transponse operation
        1988   137865246 :         if (edge%reverse(l,ielem)) then
        1989   371393024 :            do k=1,vlyr
        1990   368428412 :               ktmp = nhc*(kptr+k-1)
        1991  1144580468 :               do j=1,nhc
        1992   773187444 :                  iptr = nhc*(ktmp + j - 1)
        1993  3129080396 :                  do i=1,nhc
        1994  2760651984 :                     v(1-j,1-i,k)=edge%receive(iptr+isw+i)
        1995             :                  enddo
        1996             :               enddo
        1997             :              enddo
        1998             :         else
        1999 16898404168 :            do k=1,vlyr
        2000 16763503534 :               ktmp = nhc*(kptr+k-1)
        2001 52078454446 :               do i=1,nhc
        2002 35180050278 :                  iptr = nhc*(ktmp + i - 1)
        2003 >14237*10^7 :                  do j=1,nhc
        2004 >12560*10^7 :                     v(1-j,1-i,k)=edge%receive(iptr+isw+j)
        2005             :                  enddo
        2006             :               enddo
        2007             :              enddo
        2008             :         endif
        2009             :      else
        2010    19210008 :          do k=1,vlyr
        2011    59202486 :            do j=1,nhc
        2012   161849070 :              do i=1,nhc
        2013   142792416 :                v(1-i,1-j,k)=edgeDefaultVal
        2014             :              enddo
        2015             :            enddo
        2016             :          enddo
        2017             :      endif
        2018             :   end do
        2019             :
        2020             : ! SEAST
        2021   276037200 :   do l=swest+max_corner_elem,swest+2*max_corner_elem-1
        2022   138018600 :      ise = edge%getmap(l,ielem)
        2023   276037200 :      if(ise /= -1) then
        2024   137865246 :         if (edge%reverse(l,ielem)) then
        2025   371393024 :            do k=1,vlyr
        2026   368428412 :               ktmp = nhc*(kptr+k-1)
        2027  1144580468 :               do i=1,nhc
        2028   773187444 :                  iptr = nhc*(ktmp + i - 1)
        2029  3129080396 :                  do j=1,nhc
        2030  2760651984 :                     v(npoints+i,1-j,k)=edge%receive(iptr+ise+j)
        2031             :                  enddo
        2032             :               enddo
        2033             :             enddo
        2034             :         else
        2035 16898404168 :            do k=1,vlyr
        2036 16763503534 :               ktmp = nhc*(kptr+k-1)
        2037 52078454446 :               do j=1,nhc
        2038 35180050278 :                  iptr = nhc*(ktmp + j - 1)
        2039 >14237*10^7 :                  do i=1,nhc
        2040 >12560*10^7 :                     v(npoints+i ,1-j ,k)=edge%receive(iptr+ise+i)
        2041             :                  enddo
        2042             :               enddo
        2043             :             enddo
        2044             :         endif
        2045             :       else
        2046    19210008 :          do k=1,vlyr
        2047    59202486 :           do j=1,nhc
        2048   161849070 :             do i=1,nhc
        2049   142792416 :               v(npoints+i,1-j,k)=edgeDefaultVal
        2050             :             enddo
        2051             :           enddo
        2052             :          enddo
        2053             :      endif
        2054             :   end do
        2055             :
        2056             : ! NEAST
        2057   276037200 :   do l=swest+3*max_corner_elem,swest+4*max_corner_elem-1
        2058   138018600 :      ine = edge%getmap(l,ielem)
        2059   276037200 :      if(ine /= -1) then
        2060   137865246 :         if (edge%reverse(l,ielem)) then
        2061   371393024 :            do k=1,vlyr
        2062   368428412 :               ktmp = nhc*(kptr+k-1)
        2063  1144580468 :               do j=1,nhc
        2064  3129080396 :                  do i=1,nhc
        2065  1987464540 :                     iptr = nhc*(ktmp + i - 1)
        2066  2760651984 :                     v(npoints+i ,npoints+j,k)=edge%receive(iptr+ine+j)
        2067             :                  enddo
        2068             :               enddo
        2069             :             enddo
        2070             :         else
        2071 16898404168 :            do k=1,vlyr
        2072 16763503534 :               ktmp = nhc*(kptr+k-1)
        2073 52078454446 :               do j=1,nhc
        2074 35180050278 :                  iptr = nhc*(ktmp + j - 1)
        2075 >14237*10^7 :                  do i=1,nhc
        2076 >12560*10^7 :                     v(npoints+i ,npoints+j,k)=edge%receive(iptr+ine+i)
        2077             :                  enddo
        2078             :               enddo
        2079             :             enddo
        2080             :         endif
        2081             :       else
        2082    19210008 :          do k=1,vlyr
        2083    59202486 :           do j=1,nhc
        2084   161849070 :             do i=1,nhc
        2085   142792416 :               v(npoints+i,npoints+j,k)=edgeDefaultVal
        2086             :             enddo
        2087             :           enddo
        2088             :          enddo
        2089             :      endif
        2090             :   end do
        2091             :
        2092             : ! NWEST
        2093   276037200 :   do l=swest+2*max_corner_elem,swest+3*max_corner_elem-1
        2094   138018600 :      inw = edge%getmap(l,ielem)
        2095   276037200 :      if(inw /= -1) then
        2096   137865246 :         if (edge%reverse(l,ielem)) then
        2097   371393024 :            do k=1,vlyr
        2098   368428412 :               ktmp = nhc*(kptr+k-1)
        2099  1144580468 :               do i=1,nhc
        2100   773187444 :                  iptr = nhc*(ktmp + i - 1)
        2101  3129080396 :                  do j=1,nhc
        2102  2760651984 :                     v(1-i ,npoints+j,k)=edge%receive(iptr+inw+j)
        2103             :                  enddo
        2104             :               enddo
        2105             :             enddo
        2106             :         else
        2107 16898404168 :            do k=1,vlyr
        2108 16763503534 :               ktmp = nhc*(kptr+k-1)
        2109 52078454446 :               do j=1,nhc
        2110 35180050278 :                  iptr = nhc*(ktmp + j - 1)
        2111 >14237*10^7 :                  do i=1,nhc
        2112 >12560*10^7 :                     v(1-i ,npoints+j,k)=edge%receive(iptr+inw+i)
        2113             :                  enddo
        2114             :               enddo
        2115             :             enddo
        2116             :         endif
        2117             :       else
        2118    19210008 :          do k=1,vlyr
        2119    59202486 :           do j=1,nhc
        2120   161849070 :             do i=1,nhc
        2121   142792416 :               v(1-i,npoints+j,k)=edgeDefaultVal
        2122             :             enddo
        2123             :           enddo
        2124             :          enddo
        2125             :      endif
        2126             :   end do
        2127             :
        2128   138018600 : end subroutine ghostunpack
        2129             :
        2130             :   ! =========================================
        2131             :   ! initGhostBuffer3d:
        2132             :   ! Author: James Overfelt
        2133             :   ! create an Real based communication buffer
        2134             :   ! npoints is the number of points on one side
        2135             :   ! nhc is the deep of the ghost/halo zone
        2136             :   ! =========================================
        2137           0 :   subroutine initGhostBuffer3d(ghost,nlyr,np,nhc_in)
        2138             :
        2139             :     implicit none
        2140             :     integer,intent(in)                  :: nlyr, np
        2141             :     integer,intent(in),optional         :: nhc_in
        2142             :     type (Ghostbuffer3d_t),intent(out)  :: ghost
        2143             :
        2144             :     ! Local variables
        2145             :
        2146             :     integer :: nbuf,nhc,i
        2147             :
        2148             :     ! sanity check for threading
        2149           0 :     if (omp_get_num_threads()>1) then
        2150           0 :        call endrun('ERROR: initGhostBuffer must be called before threaded region')
        2151             :     endif
        2152             :
        2153           0 :     if (present(nhc_in)) then
        2154           0 :        nhc=nhc_in
        2155             :     else
        2156           0 :        nhc = np-1
        2157             :     endif
        2158             :
        2159           0 :     nbuf=max_neigh_edges*nelemd
        2160             :
        2161           0 :     ghost%nlyr    = nlyr
        2162           0 :     ghost%nhc     = nhc
        2163           0 :     ghost%np      = np
        2164           0 :     ghost%nbuf    = nbuf
        2165           0 :     ghost%elem_size = np*(nhc+1)
        2166           0 :     allocate(ghost%buf    (np,(nhc+1),nlyr,nbuf))
        2167           0 :     allocate(ghost%receive(np,(nhc+1),nlyr,nbuf))
        2168           0 :     ghost%buf=0
        2169           0 :     ghost%receive=0
        2170             :
        2171           0 :   end subroutine initGhostBuffer3d
        2172             :
        2173             :   ! =================================================================================
        2174             :   ! GHOSTVPACK3D
        2175             :   ! AUTHOR: James Overfelt (from a subroutine of Christoph Erath, ghostvpack2D)
        2176             :   ! Pack edges of v into an ghost buffer for boundary exchange.
        2177             :   !
        2178             :   ! This subroutine packs for many vertical layers into an ghost
        2179             :   ! buffer.
        2180             :   ! If the buffer associated with edge is not large enough to
        2181             :   ! hold all vertical layers you intent to pack, the method will
        2182             :   ! halt the program with a call to endrun().
        2183             :   ! INPUT:
        2184             :   ! - ghost Buffer into which the data will be packed.
        2185             :   !   This buffer must be previously allocated with initGhostBuffer().
        2186             :   ! - v The data to be packed.
        2187             :   ! - nhc deep of ghost/halo zone
        2188             :   ! - npoints number of points on on side
        2189             :   ! - kptr Vertical pointer to the place in the edge buffer where
        2190             :   ! data will be located.
        2191             :   ! =================================================================================
        2192           0 :   subroutine ghostVpack3d(ghost, v, vlyr, kptr, desc)
        2193             :     use dimensions_mod, only : max_corner_elem
        2194             :     use control_mod,    only : north, south, east, west, neast, nwest, seast, swest
        2195             :     use edgetype_mod, only : edgedescriptor_t, ghostbuffer3d_t
        2196             :     implicit none
        2197             :
        2198             :     type (Ghostbuffer3d_t)                :: ghost
        2199             :     integer,              intent(in)      :: kptr,vlyr
        2200             :     real (kind=r8),intent(in)      :: v(ghost%np, ghost%np, vlyr)
        2201             :     type (EdgeDescriptor_t),intent(in)    :: desc
        2202             :
        2203             :     integer                               :: nhc, np
        2204             :
        2205             :     ! Local variables
        2206             :     integer :: i,j,k,ir,l,e
        2207             :
        2208             :     integer :: is,ie,in,iw
        2209             :
        2210           0 :     if(.not. threadsafe) then
        2211             : !$OMP BARRIER
        2212           0 :        threadsafe=.true.
        2213             :     end if
        2214             :     ! Example convenction for buffer to the north:
        2215             :     !  buf(:,,:,i,e)
        2216             :     !   each "edge" is a row of data (i=1,np) in the element
        2217             :     !     north most row of data goes into e=1
        2218             :     !     next row of data goes into e=2
        2219             :     !     ....
        2220             :     !     south most row of data goes into e=np
        2221             :     !   We need to pack this way to preserve the orientation
        2222             :     !   so the data can be unpacked correctly
        2223             :
        2224             :     ! note: we think of buf as dimensioned buf(k,is,i,e)
        2225             :     ! but this array is flatted to:   buf(k,is+(i-1)+(e-1)*np)
        2226             :     !
        2227           0 :     nhc = ghost%nhc
        2228           0 :     np = ghost%np
        2229           0 :     is = desc%putmapP_ghost(south)
        2230           0 :     ie = desc%putmapP_ghost(east)
        2231           0 :     in = desc%putmapP_ghost(north)
        2232           0 :     iw = desc%putmapP_ghost(west)
        2233             :
        2234           0 :     do k=1,vlyr
        2235           0 :        do j=1,nhc
        2236           0 :           do i=1,np
        2237           0 :              ghost%buf(i,j,kptr+k,is)   = v(i,    j+1 , k)
        2238           0 :              ghost%buf(i,j,kptr+k,ie)   = v(np-j,   i , k)
        2239           0 :              ghost%buf(i,j,kptr+k,in)   = v(i,   np-j , k)
        2240           0 :              ghost%buf(i,j,kptr+k,iw)   = v(j+1,    i , k)
        2241             :           enddo
        2242             :        end do
        2243             :     end do
        2244             :     !  This is really kludgy way to setup the index reversals
        2245             :     !  But since it is so a rare event not real need to spend time optimizing
        2246             :     !  Check if the edge orientation of the recieving element is different
        2247             :     !  if it is, swap the order of data in the edge
        2248           0 :     if(desc%reverse(south)) then
        2249           0 :        do k=1,vlyr
        2250           0 :           do j=1,nhc
        2251           0 :              do i=1,np
        2252           0 :                 ir = np-i+1
        2253           0 :                 ghost%buf(ir, j, kptr+k, is)=v(i, j+1, k)
        2254             :              enddo
        2255             :           enddo
        2256             :        enddo
        2257             :     endif
        2258             :
        2259           0 :     if(desc%reverse(east)) then
        2260           0 :        do k=1,vlyr
        2261           0 :           do j=1,nhc
        2262           0 :              do i=1,np
        2263           0 :                 ir = np-i+1
        2264           0 :                 ghost%buf(ir, j, kptr+k, ie)=v(np-j, i, k)
        2265             :              enddo
        2266             :           enddo
        2267             :        enddo
        2268             :     endif
        2269             :
        2270           0 :     if(desc%reverse(north)) then
        2271           0 :        do k=1,vlyr
        2272           0 :           do j=1,nhc
        2273           0 :               do i=1,np
        2274           0 :                 ir = np-i+1
        2275           0 :                 ghost%buf(ir, j, kptr+k, in)=v(i, np-j, k)
        2276             :              enddo
        2277             :           enddo
        2278             :        enddo
        2279             :     endif
        2280             :
        2281           0 :     if(desc%reverse(west)) then
        2282           0 :        do k=1,vlyr
        2283           0 :           do j=1,nhc
        2284           0 :              do i=1,np
        2285           0 :                 ir = np-i+1
        2286           0 :                 ghost%buf(ir, j, kptr+k, iw)=v(j+1, i, k)
        2287             :              enddo
        2288             :           enddo
        2289             :        enddo
        2290             :     endif
        2291             :
        2292             :     ! corners.  this is difficult because we dont know the orientaton
        2293             :     ! of the corners, and this which (i,j) dimension maps to which dimension
        2294             : ! SWEST
        2295           0 :     do l=swest, swest+max_corner_elem-1
        2296           0 :        if (desc%putmapP_ghost(l) /= -1) then
        2297           0 :           do k=1,vlyr
        2298           0 :              do j=1,nhc+1
        2299           0 :                 do i=1,nhc+1
        2300           0 :                    ghost%buf(i, j, kptr+k, desc%putmapP_ghost(l))=v(i, j, k)
        2301             :                 enddo
        2302             :              enddo
        2303             :           enddo
        2304             :        end if
        2305             :     end do
        2306             :
        2307             : ! SEAST
        2308           0 :     do l=swest+max_corner_elem,swest+2*max_corner_elem-1
        2309           0 :        if (desc%putmapP_ghost(l) /= -1) then
        2310           0 :           do k=1,vlyr
        2311           0 :              do j=1,nhc+1
        2312           0 :                 do i=1,nhc+1
        2313           0 :                    ghost%buf(i, j, kptr+k, desc%putmapP_ghost(l))=v(np-i+1, j, k)
        2314             :                 enddo
        2315             :              enddo
        2316             :           enddo
        2317             :        end if
        2318             :     end do
        2319             :
        2320             : ! NEAST
        2321           0 :     do l=swest+3*max_corner_elem,swest+4*max_corner_elem-1
        2322           0 :        if (desc%putmapP_ghost(l) /= -1) then
        2323           0 :           do k=1,vlyr
        2324           0 :              do j=1,nhc+1
        2325           0 :                 do i=1,nhc+1
        2326           0 :                    ghost%buf(i, j, kptr+k,desc%putmapP_ghost(l))=v(np-i+1, np-j+1, k)
        2327             :                 enddo
        2328             :              enddo
        2329             :           enddo
        2330             :        end if
        2331             :     end do
        2332             :
        2333             : ! NWEST
        2334           0 :     do l=swest+2*max_corner_elem,swest+3*max_corner_elem-1
        2335           0 :        if (desc%putmapP_ghost(l) /= -1) then
        2336           0 :           do k=1,vlyr
        2337           0 :              do j=1,nhc+1
        2338           0 :                 do i=1,nhc+1
        2339           0 :                    ghost%buf(i, j, kptr+k,desc%putmapP_ghost(l))=v(i, np-j+1, k)
        2340             :                 enddo
        2341             :              enddo
        2342             :           enddo
        2343             :        end if
        2344             :     end do
        2345           0 :   end subroutine ghostVpack3d
        2346             :
        2347             :   ! =================================================================================
        2348             :   ! GHOSTVUNPACK3D
        2349             :   ! AUTHOR: James Overfelt (from a subroutine of Christoph Erath,
        2350             :   ! ghostVunpack2d)
        2351             :   ! Unpack ghost points from ghost buffer into v...
        2352             :   ! It is for cartesian points (v is only two dimensional).
        2353             :   ! INPUT SAME arguments as for GHOSTVPACK
        2354             :   ! =================================================================================
        2355             :
        2356           0 :   subroutine ghostVunpack3d(g, v, vlyr, kptr, desc, sw, se, nw, ne, mult)
        2357             :     use dimensions_mod, only : max_corner_elem
        2358             :     use control_mod, only : north, south, east, west, neast, nwest, seast, swest
        2359             :     use edgetype_mod, only : edgedescriptor_t, ghostbuffer3d_t
        2360             :     implicit none
        2361             :     type (Ghostbuffer3d_t),         intent(in)  :: g
        2362             :
        2363             :     integer,               intent(in)     :: kptr,vlyr
        2364             :     real (kind=r8), intent(inout)  :: v (1-g%nhc : g%np+g%nhc, 1-g%nhc : g%np+g%nhc, vlyr)
        2365             :     integer,               intent(out)    :: mult(5:8)
        2366             :     real (kind=r8), intent(out)    :: sw(1-g%nhc : 1,          1-g%nhc : 1,          vlyr, max_corner_elem-1)
        2367             :     real (kind=r8), intent(out)    :: se(   g%np : g%np+g%nhc, 1-g%nhc : 1,          vlyr, max_corner_elem-1)
        2368             :     real (kind=r8), intent(out)    :: ne(   g%np : g%np+g%nhc,    g%np : g%np+g%nhc, vlyr, max_corner_elem-1)
        2369             :     real (kind=r8), intent(out)    :: nw(1-g%nhc : 1,             g%np : g%np+g%nhc, vlyr, max_corner_elem-1)
        2370             :     type (EdgeDescriptor_t)               :: desc
        2371             :
        2372             :     integer                               :: nhc, np
        2373             :
        2374             :     ! Local
        2375             :     logical, parameter :: UseUnroll = .TRUE.
        2376             :     integer :: i,j,k,l
        2377             :     integer :: is,ie,in,iw,ic
        2378             :     logical :: reverse
        2379             :
        2380           0 :     threadsafe=.false.
        2381             :
        2382           0 :     nhc     = g%nhc
        2383           0 :     np      = g%np
        2384             :
        2385           0 :     is=desc%getmapP_ghost(south)
        2386           0 :     ie=desc%getmapP_ghost(east)
        2387           0 :     in=desc%getmapP_ghost(north)
        2388           0 :     iw=desc%getmapP_ghost(west)
        2389             :
        2390             : ! fill in optional values with edgeDefaultVal
        2391           0 :     do k=1,vlyr
        2392           0 :        do j=1,nhc
        2393           0 :           do i=1,nhc
        2394           0 :              v(1-i,   1-j, k)=edgeDefaultVal
        2395           0 :              v(np+i , 1-j, k)=edgeDefaultVal
        2396           0 :              v(np+i, np+j, k)=edgeDefaultVal
        2397           0 :              v(1-i , np+j, k)=edgeDefaultVal
        2398             :           enddo
        2399             :        enddo
        2400             :     enddo
        2401             :
        2402             :     ! example for north buffer
        2403             :     ! first row ('edge') goes in v(:,np+1)
        2404             :     ! 2nd   row ('edge') goes in v(:,np+2)
        2405             :     ! etc...
        2406             :
        2407           0 :     do k=1,vlyr
        2408           0 :        do j=1,nhc
        2409           0 :           do i=1,np
        2410           0 :              v(i  ,  1-j , k) = g%buf(i,j,kptr+k,is  )
        2411           0 :              v(np+j ,  i , k) = g%buf(i,j,kptr+k,ie  )
        2412           0 :              v(i  , np+j,  k) = g%buf(i,j,kptr+k,in  )
        2413           0 :              v(1-j  ,  i , k) = g%buf(i,j,kptr+k,iw  )
        2414             :           end do
        2415             :        end do
        2416             :     end do
        2417             :
        2418             :     ! four sides are always just one
        2419           0 :     mult(swest) = 0
        2420           0 :     mult(seast) = 0
        2421           0 :     mult(neast) = 0
        2422           0 :     mult(nwest) = 0
        2423             :
        2424             :
        2425             :
        2426             : ! SWEST
        2427           0 :     do l=swest, swest+max_corner_elem-1
        2428           0 :        ic = desc%getmapP_ghost(l)
        2429           0 :        if(ic /= -1) then
        2430           0 :           reverse=desc%reverse(l)
        2431           0 :           if (mult(swest) .eq. 0) then
        2432           0 :             if (reverse) then
        2433           0 :                do k=1,vlyr
        2434           0 :                   do j=1,nhc
        2435           0 :                      do i=1,nhc
        2436           0 :                         v(1-i, 1-j, k)=g%buf(j+1, i+1, kptr+k, ic)
        2437             :                      enddo
        2438             :                   enddo
        2439             :                enddo
        2440             :              else
        2441           0 :                do k=1,vlyr
        2442           0 :                   do j=1,nhc
        2443           0 :                      do i=1,nhc
        2444           0 :                         v(1-i,1-j,k)=g%buf(i+1,j+1,kptr+k,ic)
        2445             :                      enddo
        2446             :                   enddo
        2447             :                enddo
        2448             :             endif
        2449             :           else
        2450           0 :             if (reverse) then
        2451           0 :                do k=1,vlyr
        2452           0 :                   do j=0,nhc
        2453           0 :                      do i=0,nhc
        2454           0 :                         sw(1-i,1-j,k,mult(swest))=g%buf(j+1,i+1,kptr+k,ic)
        2455             :                      enddo
        2456             :                   enddo
        2457             :                enddo
        2458             :              else
        2459           0 :                 do k=1,vlyr
        2460           0 :                    do j=0,nhc
        2461           0 :                       do i=0,nhc
        2462           0 :                          sw(1-i,1-j,k,mult(swest))=g%buf(i+1,j+1,kptr+k,ic)
        2463             :                      enddo
        2464             :                   enddo
        2465             :                enddo
        2466             :             endif
        2467             :           endif
        2468           0 :           mult(swest) = mult(swest) + 1
        2469             :        endif
        2470             :     end do
        2471             :
        2472             : ! SEAST
        2473           0 :     do l=swest+max_corner_elem,swest+2*max_corner_elem-1
        2474           0 :        ic = desc%getmapP_ghost(l)
        2475           0 :        if(ic /= -1) then
        2476           0 :           reverse=desc%reverse(l)
        2477           0 :           if (mult(seast) .eq. 0) then
        2478           0 :             if (reverse) then
        2479           0 :                do k=1,vlyr
        2480           0 :                   do j=1,nhc
        2481           0 :                      do i=1,nhc
        2482           0 :                         v(np+i,1-j,k)=g%buf(j+1,i+1,kptr+k,ic)
        2483             :                      enddo
        2484             :                   enddo
        2485             :                enddo
        2486             :             else
        2487           0 :                do k=1,vlyr
        2488           0 :                   do j=1,nhc
        2489           0 :                      do i=1,nhc
        2490           0 :                         v(np+i ,1-j,k)=g%buf(i+1,j+1,kptr+k,ic)
        2491             :                      enddo
        2492             :                   enddo
        2493             :                enddo
        2494             :             endif
        2495             :           else
        2496           0 :             if (reverse) then
        2497           0 :                do k=1,vlyr
        2498           0 :                   do j=0,nhc
        2499           0 :                      do i=0,nhc
        2500           0 :                         se(np+i,1-j,k,mult(seast))=g%buf(j+1,i+1,kptr+k,ic)
        2501             :                      enddo
        2502             :                   enddo
        2503             :                enddo
        2504             :             else
        2505           0 :                do k=1,vlyr
        2506           0 :                   do j=0,nhc
        2507           0 :                      do i=0,nhc
        2508           0 :                         se(np+i ,1-j,k,mult(seast))=g%buf(i+1,j+1,kptr+k,ic)
        2509             :                      enddo
        2510             :                   enddo
        2511             :                enddo
        2512             :             endif
        2513             :           endif
        2514           0 :           mult(seast) = mult(seast) + 1
        2515             :        endif
        2516             :     end do
        2517             :
        2518             :
        2519             : ! NEAST
        2520           0 :     do l=swest+3*max_corner_elem,swest+4*max_corner_elem-1
        2521           0 :        ic = desc%getmapP_ghost(l)
        2522           0 :        if(ic /= -1) then
        2523           0 :           reverse=desc%reverse(l)
        2524           0 :           if (mult(neast) .eq. 0) then
        2525           0 :             if (reverse) then
        2526           0 :                do k=1,vlyr
        2527           0 :                   do j=1,nhc
        2528           0 :                      do i=1,nhc
        2529           0 :                         v(np+i ,np+j,k)=g%buf(j+1,i+1,kptr+k,ic)
        2530             :                      enddo
        2531             :                   enddo
        2532             :                enddo
        2533             :             else
        2534           0 :                do k=1,vlyr
        2535           0 :                   do j=1,nhc
        2536           0 :                      do i=1,nhc
        2537           0 :                         v(np+i ,np+j,k)=g%buf(i+1,j+1,kptr+k,ic)
        2538             :                      enddo
        2539             :                   enddo
        2540             :                enddo
        2541             :             endif
        2542             :           else
        2543           0 :             if (reverse) then
        2544           0 :                do k=1,vlyr
        2545           0 :                   do j=0,nhc
        2546           0 :                      do i=0,nhc
        2547           0 :                         ne(np+i ,np+j,k,mult(neast))=g%buf(j+1,i+1,kptr+k,ic)
        2548             :                      enddo
        2549             :                   enddo
        2550             :                enddo
        2551             :             else
        2552           0 :                do k=1,vlyr
        2553           0 :                   do j=0,nhc
        2554           0 :                      do i=0,nhc
        2555           0 :                         ne(np+i ,np+j,k,mult(neast))=g%buf(i+1,j+1,kptr+k,ic)
        2556             :                      enddo
        2557             :                   enddo
        2558             :                enddo
        2559             :             endif
        2560             :           endif
        2561           0 :           mult(neast) = mult(neast) + 1
        2562             :        endif
        2563             :     end do
        2564             :
        2565             : ! NWEST
        2566           0 :     do l=swest+2*max_corner_elem,swest+3*max_corner_elem-1
        2567           0 :        ic = desc%getmapP_ghost(l)
        2568           0 :        if(ic /= -1) then
        2569           0 :           reverse=desc%reverse(l)
        2570           0 :           if (mult(nwest) .eq. 0) then
        2571           0 :              if (reverse) then
        2572           0 :                 do k=1,vlyr
        2573           0 :                    do j=1,nhc
        2574           0 :                       do i=1,nhc
        2575           0 :                          v(1-i ,np+j,k)=g%buf(j+1,i+1,kptr+k,ic)
        2576             :                       enddo
        2577             :                    enddo
        2578             :                 enddo
        2579             :              else
        2580           0 :                 do k=1,vlyr
        2581           0 :                    do j=1,nhc
        2582           0 :                       do i=1,nhc
        2583           0 :                          v(1-i ,np+j,k)=g%buf(i+1,j+1,kptr+k,ic)
        2584             :                       enddo
        2585             :                    enddo
        2586             :                 enddo
        2587             :              endif
        2588             :           else
        2589           0 :              if (reverse) then
        2590           0 :                 do k=1,vlyr
        2591           0 :                    do j=0,nhc
        2592           0 :                       do i=0,nhc
        2593           0 :                          nw(1-i ,np+j,k,mult(nwest))=g%buf(j+1,i+1,kptr+k,ic)
        2594             :                       enddo
        2595             :                    enddo
        2596             :                 enddo
        2597             :              else
        2598           0 :                 do k=1,vlyr
        2599           0 :                    do j=0,nhc
        2600           0 :                       do i=0,nhc
        2601           0 :                          nw(1-i ,np+j,k,mult(nwest))=g%buf(i+1,j+1,kptr+k,ic)
        2602             :                       enddo
        2603             :                    enddo
        2604             :                 enddo
        2605             :              endif
        2606             :           endif
        2607           0 :           mult(nwest) = mult(nwest) + 1
        2608             :        endif
        2609             :     end do
        2610             :
        2611           0 :   end subroutine ghostVunpack3d
        2612             :
        2613           0 :   subroutine FreeGhostBuffer3D(buffer)
        2614             :     use edgetype_mod, only : ghostbuffer3d_t
        2615             :     implicit none
        2616             :     type (Ghostbuffer3d_t),intent(inout) :: buffer
        2617             :
        2618             : !$OMP BARRIER
        2619             : !$OMP MASTER
        2620           0 :     buffer%nbuf=0
        2621           0 :     buffer%nlyr=0
        2622           0 :     deallocate(buffer%buf)
        2623           0 :     deallocate(buffer%receive)
        2624             : !$OMP END MASTER
        2625             :
        2626           0 :   end subroutine FreeGhostBuffer3D
        2627             :
        2628             :
        2629           0 : End module edge_mod

Generated by: LCOV version 1.14