C $Header: /u/gcmpack/MITgcm/pkg/generic_advdiff/gad_ppm_flx_y.F,v 1.1 2016/03/13 01:44:03 jmc Exp $
C $Name:  $

#     include "GAD_OPTIONS.h"

      SUBROUTINE GAD_PPM_FLX_Y(bi,bj,kk,ix,
     &           calc_CFL,delT,vvel,
     &           vfac,fhat,flux,myThid )
C     |================================================================|
C     | PPM_FLX_Y: evaluate PPM flux on grid-cell edges.               |
C     |================================================================|

          implicit none

C     =============================================== global variables
#         include "SIZE.h"
#         include "GRID.h"
#         include "GAD.h"

C     ================================================================
C       bi,bj     :: tile indexing.
C       kk        :: r-index.
C       ix        :: x-index.
C       calc_CFL  :: TRUE to calc. CFL from vel.
C       delT      :: time-step.
C       vvel      :: vel.-comp in y-direction.
C       vfac      :: vel.-flux in y-direction.
C       fhat      :: row of poly. coeff.
C       flux      :: adv.-flux in y-direction.
C       myThid    :: thread number.
C     ================================================================
          integer bi,bj,kk,ix
          logical calc_CFL
          _RL delT
          _RL vvel(1-OLx:sNx+OLx,
     &             1-OLy:sNy+OLy)
          _RL vfac(1-OLx:sNx+OLx,
     &             1-OLy:sNy+OLy)
          _RL fhat(1:3,
     &             1-OLy:sNy+OLy)
          _RL flux(1-OLx:sNx+OLx,
     &             1-OLy:sNy+OLy)
          integer myThid

C     ================================================================
C       iy        :: y-indexing.
C       vCFL      :: CFL number.
C       intF      :: upwind tracer edge-value.
C       ss11,ss22 :: int. endpoints.
C       ivec      :: int. basis vec.
C     ================================================================
          integer iy
          _RL vCFL,intF
          _RL ss11,ss22
          _RL ivec(1:3)

C     ================================================================
C       (1): calc. "departure-points" for each grid-cell edge by int-
C            egrating edge position backward in time over one single
C            time-step. This is a "single-cell" implementation: requ-
C            ires CFL <= 1.0.
C       (2): calc. flux as the integral of the upwind grid-cell poly-
C            nomial over the deformation interval found in (1).
C     ================================================================

          do  iy = 1-OLy+3, sNy+OLy-2

              if (vvel(ix,iy) .eq. 0. _d 0) then

                  flux(ix,iy)    = 0. _d 0

              else

              if (vvel(ix,iy) .gt. 0. _d 0) then

C     ==================== integrate PPM profile over upwind cell IY-1
              if ( calc_CFL ) then
              vCFL = vvel(ix,iy) * delT
     &             * recip_dyF(ix,iy-1,bi,bj)
     &             * recip_deepFacC(kk)
              else
              vCFL = vvel(ix,iy)
              end


if ss11 = +1. _d 0 - 2. _d 0 * vCFL ss22 = +1. _d 0 C ==================== integrate profile over region swept by face ivec(1) = ss22 - ss11 ivec(2) =(ss22 ** 2 & - ss11 ** 2)*(1. _d 0 / 2. _d 0) ivec(3) =(ss22 ** 3 & - ss11 ** 3)*(1. _d 0 / 3. _d 0) intF = ivec(1) * fhat(1,iy-1) & + ivec(2) * fhat(2,iy-1) & + ivec(3) * fhat(3,iy-1) intF = intF / (ss22 - ss11) else C ==================== integrate PPM profile over upwind cell IY+0 if ( calc_CFL ) then vCFL = vvel(ix,iy) * delT & * recip_dyF(ix,iy-0,bi,bj) & * recip_deepFacC(kk) else vCFL = vvel(ix,iy) end


if ss11 = -1. _d 0 - 2. _d 0 * vCFL ss22 = -1. _d 0 C ==================== integrate profile over region swept by face ivec(1) = ss22 - ss11 ivec(2) =(ss22 ** 2 & - ss11 ** 2)*(1. _d 0 / 2. _d 0) ivec(3) =(ss22 ** 3 & - ss11 ** 3)*(1. _d 0 / 3. _d 0) intF = ivec(1) * fhat(1,iy-0) & + ivec(2) * fhat(2,iy-0) & + ivec(3) * fhat(3,iy-0) intF = intF / (ss22 - ss11) end


if C ==================== calc. flux = upwind tracer * face-transport flux(ix,iy) = + vfac(ix,iy) * intF end


if end


do return c end subroutine GAD_PPM_FLX_Y end