C $Header: /u/gcmpack/MITgcm/model/src/diags_rho.F,v 1.5 2011/11/10 21:03:35 jmc Exp $
C $Name: $
#include "PACKAGES_CONFIG.h"
#include "CPP_OPTIONS.h"
C-- File diags_rho.F: density & density advection diagnostics
C-- Contents
C-- o DIAGS_RHO_L
C-- o DIAGS_RHO_G
C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
CBOP
C !ROUTINE: DIAGS_RHO_L
C !INTERFACE:
SUBROUTINE DIAGS_RHO_L(
I diagRho, k, bi, bj,
I rho3d, rhoKm1, wFld,
I myTime, myIter, myThid )
C !DESCRIPTION: \bv
C *==========================================================*
C | S/R DIAGS_RHO_L
C | o Density vertical advective term diagnostics
C *==========================================================*
C | works with local arrays, and called inside k,bi,bj loops
C *==========================================================*
C \ev
C !USES:
IMPLICIT NONE
C == Global variables ==
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
C !INPUT/OUTPUT PARAMETERS:
C == Routine Arguments ==
C diagRho :: select which diags to fill
C k :: level index
C bi, bj :: tile indices
C rho3d :: in-situ density anomaly
C rhoKm1 :: density of water @ level above (k-1), evaluated at pressure level k
C wFld :: vertical velocity
C myTime :: Current time
C myIter :: Iteration number
C myThid :: my Thread Id number
INTEGER diagRho
INTEGER k, bi, bj
_RL rho3d (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
_RL rhoKm1(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
_RL wFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
_RL myTime
INTEGER myIter, myThid
CEOP
#ifdef ALLOW_DIAGNOSTICS
C !LOCAL VARIABLES:
C == Local variables ==
C i,j :: Loop counters
C tmpFld :: temporary working array
INTEGER i,j
_RL tmpFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
IF ( k.GE.2 .AND. MOD(diagRho,8).GE.4 ) THEN
C-- Diagnose Vertical velocity times vertical difference
C of potential density reference at level below (i.e. level k)
DO j=1,sNy
DO i=1,sNx
tmpFld(i,j) = wFld(i,j,k,bi,bj)
& *( rho3d(i,j,k) - rhoKm1(i,j) )*rkSign
ENDDO
ENDDO
CALL DIAGNOSTICS_FILL(tmpFld,'WdRHO_P ',k,1,2,bi,bj,myThid)
IF ( k.EQ.2 ) CALL DIAGNOSTICS_COUNT('WdRHO_P ',bi,bj,myThid)
ENDIF
IF ( k.GE.2 .AND. diagRho.GE.8 ) THEN
C-- Diagnose Vertical velocity times vertical difference
C of density at fixed Temp & Salt (from level above, i.e. level k-1)
DO j=1,sNy
DO i=1,sNx
tmpFld(i,j) = wFld(i,j,k,bi,bj)
& *( rhoKm1(i,j) - rho3d(i,j,k-1) )*rkSign
ENDDO
ENDDO
CALL DIAGNOSTICS_FILL(tmpFld,'WdRHOdP ',k,1,2,bi,bj,myThid)
IF ( k.EQ.2 ) CALL DIAGNOSTICS_COUNT('WdRHOdP ',bi,bj,myThid)
ENDIF
#endif /* ALLOW_DIAGNOSTICS */
RETURN
END
C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
CBOP
C !ROUTINE: DIAGS_RHO_G
C !INTERFACE:
SUBROUTINE DIAGS_RHO_G(
I rho3d, uFld, vFld, wFld,
I myTime, myIter, myThid )
C !DESCRIPTION: \bv
C *==========================================================*
C | S/R DIAGS_RHO_G
C | o Density & Density advective Flux diagnostics
C *==========================================================*
C | works with global arrays; k,bi,bj loops are done here
C *==========================================================*
C \ev
C !USES:
IMPLICIT NONE
C == Global variables ==
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
C !INPUT/OUTPUT PARAMETERS:
C == Routine Arguments ==
C rho3d :: in-situ density anomaly
C uFld :: zonal velocity
C vFld :: meridional velocity
C wFld :: vertical velocity
C myTime :: Current time
C myIter :: Iteration number
C myThid :: my Thread Id number
_RL rho3d(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
_RL uFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
_RL vFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
_RL wFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
_RL myTime
INTEGER myIter, myThid
CEOP
#ifdef ALLOW_DIAGNOSTICS
C !FUNCTIONS:
LOGICAL DIAGNOSTICS_IS_ON
EXTERNAL
C !LOCAL VARIABLES:
C == Local variables ==
C i,j :: Loop counters
C k, bi,bj :: level & tile indices
C tmpFld :: temporary working array
INTEGER i,j
INTEGER k, bi,bj
_RL tmpFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
_RL tmpFac
C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
CALL DIAGNOSTICS_FILL( rho3d, 'RHOAnoma',
& 0, Nr, 0, 1, 1, myThid )
tmpFac = 1. _d 0
CALL DIAGNOSTICS_SCALE_FILL( rho3d, tmpFac, 2,
& 'RHOANOSQ', 0, Nr, 0, 1, 1, myThid )
IF ( DIAGNOSTICS_IS_ON('URHOMASS',myThid) ) THEN
DO bj=myByLo(myThid),myByHi(myThid)
DO bi=myBxLo(myThid),myBxHi(myThid)
DO k=1,Nr
DO j=1,sNy
DO i=1,sNx+1
tmpFld(i,j) = uFld(i,j,k,bi,bj)*_hFacW(i,j,k,bi,bj)
& *(rho3d(i-1,j,k,bi,bj)+rho3d(i,j,k,bi,bj))
& *0.5 _d 0
ENDDO
ENDDO
CALL DIAGNOSTICS_FILL(tmpFld,'URHOMASS',k,1,2,bi,bj,myThid)
ENDDO
ENDDO
ENDDO
ENDIF
IF ( DIAGNOSTICS_IS_ON('VRHOMASS',myThid) ) THEN
DO bj=myByLo(myThid),myByHi(myThid)
DO bi=myBxLo(myThid),myBxHi(myThid)
DO k=1,Nr
DO j=1,sNy+1
DO i=1,sNx
tmpFld(i,j) = vFld(i,j,k,bi,bj)*_hFacS(i,j,k,bi,bj)
& *(rho3d(i,j-1,k,bi,bj)+rho3d(i,j,k,bi,bj))
& *0.5 _d 0
ENDDO
ENDDO
CALL DIAGNOSTICS_FILL(tmpFld,'VRHOMASS',k,1,2,bi,bj,myThid)
ENDDO
ENDDO
ENDDO
ENDIF
IF ( DIAGNOSTICS_IS_ON('WRHOMASS',myThid) ) THEN
DO bj=myByLo(myThid),myByHi(myThid)
DO bi=myBxLo(myThid),myBxHi(myThid)
DO k=1,Nr
IF ( k.EQ.1 ) THEN
DO j=1,sNy
DO i=1,sNx
tmpFld(i,j) = wFld(i,j,k,bi,bj)*rho3d(i,j,k,bi,bj)
ENDDO
ENDDO
ELSE
DO j=1,sNy
DO i=1,sNx
tmpFld(i,j) = wFld(i,j,k,bi,bj)
& *(rho3d(i,j,k-1,bi,bj)+rho3d(i,j,k,bi,bj))
& *0.5 _d 0
ENDDO
ENDDO
ENDIF
CALL DIAGNOSTICS_FILL(tmpFld,'WRHOMASS',k,1,2,bi,bj,myThid)
ENDDO
ENDDO
ENDDO
ENDIF
#endif /* ALLOW_DIAGNOSTICS */
RETURN
END