C $Header: /u/gcmpack/MITgcm/pkg/monitor/mon_calc_advcfl.F,v 1.2 2016/07/22 14:21:20 jmc Exp $
C $Name:  $

#include "MONITOR_OPTIONS.h"

C--  File mon_calc_advcfl.F: Routines to compute Tracer Advective CFL
C--   Contents
C--   o MON_CALC_ADVCFL_TILE
C--   o MON_CALC_ADVCFL_GLOB

C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
CBOP
C     !ROUTINE: MON_CALC_ADVCFL_TILE

C     !INTERFACE:
      SUBROUTINE MON_CALC_ADVCFL_TILE(
     I               myNr, bi, bj,
     I               uFld, vFld, wFld, dT_lev,
     O               maxCFL,
     I               myIter, myThid )

C     Calculate Maximum advective CFL in 3 direction (x,y,z)
C     for current tile bi,bj

C     !USES:
      IMPLICIT NONE
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
c#include "SURFACE.h"
c#include "MONITOR.h"

C     !INPUT/OUTPUT PARAMETERS:
C     myNr    :: number of levels
C     bi, bj  :: tile indices
C     uFld    :: zonal velocity
C     vFld    :: merid velocity
C     wFld    :: vert. velocity
C     dT_lev  :: tracer time-step
C     maxCFL  :: maximum advective CFL in 3 directions
C     myIter  :: my Iteration number
C     myThid  :: my Thread Id number
      INTEGER myNr, bi, bj
      _RL uFld  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,myNr)
      _RL vFld  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,myNr)
      _RL wFld  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,myNr)
      _RL dT_lev(myNr)
      _RL maxCFL(3)
      INTEGER myIter, myThid
CEOP

C     !LOCAL VARIABLES:
      INTEGER i,j,k
      _RL uTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
      _RL vTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
      _RL rTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
      _RL rTransKp1, tmpVal, recVol_dT

C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|

      maxCFL(1) = 0.
      maxCFL(2) = 0.
      maxCFL(3) = 0.

      DO j=1-OLy,sNy+OLy
        DO i=1-OLx,sNx+OLx
          rTrans(i,j) = 0.
        ENDDO
      ENDDO

      DO k=myNr,1,-1

C--   compute horiz. transport and time-step divided by volume
C     (without level thickness "drF" since it cancels)
        DO j=1,sNy+1
         DO i=1,sNx+1
           uTrans(i,j) = uFld(i,j,k)*dyG(i,j,bi,bj)*deepFacC(k)
     &                         *hFacW(i,j,k,bi,bj)
           vTrans(i,j) = vFld(i,j,k)*dxG(i,j,bi,bj)*deepFacC(k)
     &                         *hFacS(i,j,k,bi,bj)
         ENDDO
        ENDDO
        DO j=1,sNy
         DO i=1,sNx
           recVol_dT = dT_lev(k)
     &             *recip_rA(i,j,bi,bj)*recip_deepFac2C(k)
     &                   *recip_hFacC(i,j,k,bi,bj)
           tmpVal = (
     &          + MAX( uTrans(i+1,j), zeroRL )
     &          - MIN( uTrans( i ,j), zeroRL )
     &              )*recVol_dT
           maxCFL(1) = MAX( maxCFL(1), tmpVal )
           tmpVal = (
     &          + MAX( vTrans(i,j+1), zeroRL )
     &          - MIN( vTrans(i, j ), zeroRL )
     &              )*recVol_dT
           maxCFL(2) = MAX( maxCFL(2), tmpVal )
         ENDDO
        ENDDO

C--   compute vert. transport and time-step divided by volume
C     (without grid-cell area "rA" since it cancels)
        DO j=1,sNy
         DO i=1,sNx
           rTransKp1 = rTrans(i,j)
           rTrans(i,j) = wFld(i,j,k)
     &                 *deepFac2F(k)*rhoFacF(k)
           recVol_dT = dT_lev(k)
     &             *recip_deepFac2C(k)*recip_rhoFacC(k)
     &             *recip_drF(k)*recip_hFacC(i,j,k,bi,bj)
           tmpVal = (
     &          + MAX( rTrans(i,j), zeroRL )
     &          - MIN( rTransKp1, zeroRL )
     &              )*recVol_dT
           maxCFL(3) = MAX( maxCFL(3), tmpVal )
         ENDDO
        ENDDO

      ENDDO

      RETURN
      END


C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| CBOP C !ROUTINE: MON_CALC_ADVCFL_GLOB C !INTERFACE: SUBROUTINE MON_CALC_ADVCFL_GLOB( I maxCFL, myIter, myThid ) C !DESCRIPTION: C Calculate Maximum advective CFL in 3 direction (x,y,z) C in global domain (from tile-max value) C !USES: IMPLICIT NONE #include "SIZE.h" #include "EEPARAMS.h" #include "PARAMS.h" #include "MONITOR.h" C !INPUT PARAMETERS: C maxCFL :: maximum advective CFL (per tile) in 3 directions C myIter :: Current iteration number in simulation C myThid :: my Thread Id number _RL maxCFL(3,nSx,nSy) INTEGER myIter INTEGER myThid CEOP C !LOCAL VARIABLES: INTEGER bi, bj _RL uCFL, vCFL, wCFL uCFL = 0. vCFL = 0. wCFL = 0. DO bj=myByLo(myThid),myByHi(myThid) DO bi=myBxLo(myThid),myBxHi(myThid) uCFL = MAX( uCFL, maxCFL(1,bi,bj) ) vCFL = MAX( vCFL, maxCFL(2,bi,bj) ) wCFL = MAX( wCFL, maxCFL(3,bi,bj) ) ENDDO ENDDO _GLOBAL_MAX_RL( uCFL, myThid ) _GLOBAL_MAX_RL( vCFL, myThid ) _GLOBAL_MAX_RL( wCFL, myThid ) C- store values in common bloc _BEGIN_MASTER(myThid) mon_trAdvCFL(1) = uCFL mon_trAdvCFL(2) = vCFL mon_trAdvCFL(3) = wCFL _END_MASTER(myThid) RETURN END