Back to home page

MITgcm

 
 

    


File indexing completed on 2018-03-02 18:36:37 UTC

view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
8ab93ea9a3 Jean*0001 #include "PACKAGES_CONFIG.h"
                0002 #include "CPP_OPTIONS.h"
                0003 
842349fcb7 Jean*0004 C--  File diags_rho.F: density & density advection diagnostics
                0005 C--   Contents
                0006 C--   o DIAGS_RHO_L
                0007 C--   o DIAGS_RHO_G
                0008 
                0009 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
8ab93ea9a3 Jean*0010 CBOP
842349fcb7 Jean*0011 C     !ROUTINE: DIAGS_RHO_L
8ab93ea9a3 Jean*0012 C     !INTERFACE:
d497465397 Jean*0013       SUBROUTINE DIAGS_RHO_L(
                0014      I                        diagRho, k, bi, bj,
                0015      I                        rho3d, rhoKm1, wFld,
842349fcb7 Jean*0016      I                        myTime, myIter, myThid )
8ab93ea9a3 Jean*0017 C     !DESCRIPTION: \bv
                0018 C     *==========================================================*
d497465397 Jean*0019 C     | S/R DIAGS_RHO_L
                0020 C     | o Density vertical advective term diagnostics
8ab93ea9a3 Jean*0021 C     *==========================================================*
842349fcb7 Jean*0022 C     | works with local arrays, and called inside k,bi,bj loops
8ab93ea9a3 Jean*0023 C     *==========================================================*
                0024 C     \ev
                0025 
                0026 C     !USES:
                0027       IMPLICIT NONE
                0028 C     == Global variables ==
                0029 #include "SIZE.h"
                0030 #include "EEPARAMS.h"
                0031 #include "PARAMS.h"
                0032 #include "GRID.h"
                0033 
                0034 C     !INPUT/OUTPUT PARAMETERS:
                0035 C     == Routine Arguments ==
d497465397 Jean*0036 C     diagRho  :: select which diags to fill
                0037 C     k        :: level index
                0038 C     bi, bj   :: tile indices
                0039 C     rho3d    :: in-situ density anomaly
                0040 C     rhoKm1   :: density of water @ level above (k-1), evaluated at pressure level k
                0041 C     wFld     :: vertical velocity
                0042 C     myTime   :: Current time
                0043 C     myIter   :: Iteration number
                0044 C     myThid   :: my Thread Id number
                0045       INTEGER diagRho
842349fcb7 Jean*0046       INTEGER k, bi, bj
d497465397 Jean*0047       _RL rho3d (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
8ab93ea9a3 Jean*0048       _RL rhoKm1(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
842349fcb7 Jean*0049       _RL wFld  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
8ab93ea9a3 Jean*0050       _RL myTime
                0051       INTEGER myIter, myThid
                0052 CEOP
                0053 
                0054 #ifdef ALLOW_DIAGNOSTICS
                0055 
                0056 C     !LOCAL VARIABLES:
                0057 C     == Local variables ==
d497465397 Jean*0058 C     i,j      :: Loop counters
                0059 C     tmpFld   :: temporary working array
8ab93ea9a3 Jean*0060       INTEGER i,j
8996cf5a3c Jean*0061       _RL tmpFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
8ab93ea9a3 Jean*0062 
842349fcb7 Jean*0063 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0064 
d497465397 Jean*0065       IF ( k.GE.2 .AND. MOD(diagRho,8).GE.4 ) THEN
                0066 C--   Diagnose Vertical velocity times vertical difference
                0067 C     of potential density reference at level below (i.e. level k)
8ab93ea9a3 Jean*0068         DO j=1,sNy
                0069          DO i=1,sNx
d497465397 Jean*0070            tmpFld(i,j) = wFld(i,j,k,bi,bj)
                0071      &                 *( rho3d(i,j,k) - rhoKm1(i,j) )*rkSign
8ab93ea9a3 Jean*0072          ENDDO
                0073         ENDDO
d497465397 Jean*0074         CALL DIAGNOSTICS_FILL(tmpFld,'WdRHO_P ',k,1,2,bi,bj,myThid)
                0075         IF ( k.EQ.2 ) CALL DIAGNOSTICS_COUNT('WdRHO_P ',bi,bj,myThid)
                0076       ENDIF
                0077 
                0078       IF ( k.GE.2 .AND. diagRho.GE.8 ) THEN
                0079 C--   Diagnose Vertical velocity times vertical difference
                0080 C     of density at fixed Temp & Salt (from level above, i.e. level k-1)
8ab93ea9a3 Jean*0081         DO j=1,sNy
842349fcb7 Jean*0082          DO i=1,sNx
                0083            tmpFld(i,j) = wFld(i,j,k,bi,bj)
d497465397 Jean*0084      &                 *( rhoKm1(i,j) - rho3d(i,j,k-1) )*rkSign
8ab93ea9a3 Jean*0085          ENDDO
                0086         ENDDO
d497465397 Jean*0087         CALL DIAGNOSTICS_FILL(tmpFld,'WdRHOdP ',k,1,2,bi,bj,myThid)
                0088         IF ( k.EQ.2 ) CALL DIAGNOSTICS_COUNT('WdRHOdP ',bi,bj,myThid)
8ab93ea9a3 Jean*0089       ENDIF
                0090 
842349fcb7 Jean*0091 #endif /* ALLOW_DIAGNOSTICS */
                0092 
                0093       RETURN
                0094       END
                0095 
                0096 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0097 CBOP
                0098 C     !ROUTINE: DIAGS_RHO_G
                0099 C     !INTERFACE:
d497465397 Jean*0100       SUBROUTINE DIAGS_RHO_G(
                0101      I                        rho3d, uFld, vFld, wFld,
842349fcb7 Jean*0102      I                        myTime, myIter, myThid )
                0103 C     !DESCRIPTION: \bv
                0104 C     *==========================================================*
d497465397 Jean*0105 C     | S/R DIAGS_RHO_G
                0106 C     | o Density & Density advective Flux diagnostics
842349fcb7 Jean*0107 C     *==========================================================*
                0108 C     | works with global arrays; k,bi,bj loops are done here
                0109 C     *==========================================================*
                0110 C     \ev
                0111 
                0112 C     !USES:
                0113       IMPLICIT NONE
                0114 C     == Global variables ==
                0115 #include "SIZE.h"
                0116 #include "EEPARAMS.h"
                0117 #include "PARAMS.h"
                0118 #include "GRID.h"
                0119 
                0120 C     !INPUT/OUTPUT PARAMETERS:
                0121 C     == Routine Arguments ==
d497465397 Jean*0122 C     rho3d    :: in-situ density anomaly
                0123 C     uFld     :: zonal velocity
                0124 C     vFld     :: meridional velocity
                0125 C     wFld     :: vertical velocity
                0126 C     myTime   :: Current time
                0127 C     myIter   :: Iteration number
                0128 C     myThid   :: my Thread Id number
842349fcb7 Jean*0129       _RL rho3d(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
                0130       _RL uFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
                0131       _RL vFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
d497465397 Jean*0132       _RL wFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
842349fcb7 Jean*0133       _RL myTime
                0134       INTEGER myIter, myThid
                0135 CEOP
                0136 
                0137 #ifdef ALLOW_DIAGNOSTICS
d497465397 Jean*0138 C     !FUNCTIONS:
                0139       LOGICAL  DIAGNOSTICS_IS_ON
                0140       EXTERNAL DIAGNOSTICS_IS_ON
842349fcb7 Jean*0141 
                0142 C     !LOCAL VARIABLES:
                0143 C     == Local variables ==
d497465397 Jean*0144 C     i,j      :: Loop counters
                0145 C     k, bi,bj :: level & tile indices
                0146 C     tmpFld   :: temporary working array
842349fcb7 Jean*0147       INTEGER i,j
                0148       INTEGER k, bi,bj
                0149       _RL tmpFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0150       _RL tmpFac
                0151 
                0152 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0153 
                0154       CALL DIAGNOSTICS_FILL( rho3d, 'RHOAnoma',
                0155      &                               0, Nr, 0, 1, 1, myThid )
                0156       tmpFac = 1. _d 0
                0157       CALL DIAGNOSTICS_SCALE_FILL( rho3d, tmpFac, 2,
                0158      &                   'RHOANOSQ', 0, Nr, 0, 1, 1, myThid )
                0159 
                0160       IF ( DIAGNOSTICS_IS_ON('URHOMASS',myThid) ) THEN
                0161        DO bj=myByLo(myThid),myByHi(myThid)
                0162         DO bi=myBxLo(myThid),myBxHi(myThid)
                0163          DO k=1,Nr
                0164           DO j=1,sNy
                0165            DO i=1,sNx+1
                0166              tmpFld(i,j) = uFld(i,j,k,bi,bj)*_hFacW(i,j,k,bi,bj)
                0167      &                   *(rho3d(i-1,j,k,bi,bj)+rho3d(i,j,k,bi,bj))
                0168      &                   *0.5 _d 0
                0169            ENDDO
                0170           ENDDO
                0171           CALL DIAGNOSTICS_FILL(tmpFld,'URHOMASS',k,1,2,bi,bj,myThid)
8ab93ea9a3 Jean*0172          ENDDO
                0173         ENDDO
842349fcb7 Jean*0174        ENDDO
8ab93ea9a3 Jean*0175       ENDIF
                0176 
842349fcb7 Jean*0177       IF ( DIAGNOSTICS_IS_ON('VRHOMASS',myThid) ) THEN
                0178        DO bj=myByLo(myThid),myByHi(myThid)
                0179         DO bi=myBxLo(myThid),myBxHi(myThid)
                0180          DO k=1,Nr
                0181           DO j=1,sNy+1
                0182            DO i=1,sNx
                0183              tmpFld(i,j) = vFld(i,j,k,bi,bj)*_hFacS(i,j,k,bi,bj)
                0184      &                   *(rho3d(i,j-1,k,bi,bj)+rho3d(i,j,k,bi,bj))
                0185      &                   *0.5 _d 0
                0186            ENDDO
                0187           ENDDO
                0188           CALL DIAGNOSTICS_FILL(tmpFld,'VRHOMASS',k,1,2,bi,bj,myThid)
8ab93ea9a3 Jean*0189          ENDDO
                0190         ENDDO
842349fcb7 Jean*0191        ENDDO
8ab93ea9a3 Jean*0192       ENDIF
                0193 
d497465397 Jean*0194       IF ( DIAGNOSTICS_IS_ON('WRHOMASS',myThid) ) THEN
                0195        DO bj=myByLo(myThid),myByHi(myThid)
                0196         DO bi=myBxLo(myThid),myBxHi(myThid)
                0197          DO k=1,Nr
                0198           IF ( k.EQ.1 ) THEN
                0199            DO j=1,sNy
                0200             DO i=1,sNx
                0201              tmpFld(i,j) = wFld(i,j,k,bi,bj)*rho3d(i,j,k,bi,bj)
                0202             ENDDO
                0203            ENDDO
                0204           ELSE
                0205            DO j=1,sNy
                0206             DO i=1,sNx
                0207              tmpFld(i,j) = wFld(i,j,k,bi,bj)
                0208      &                   *(rho3d(i,j,k-1,bi,bj)+rho3d(i,j,k,bi,bj))
                0209      &                   *0.5 _d 0
                0210             ENDDO
                0211            ENDDO
                0212           ENDIF
                0213           CALL DIAGNOSTICS_FILL(tmpFld,'WRHOMASS',k,1,2,bi,bj,myThid)
                0214          ENDDO
                0215         ENDDO
                0216        ENDDO
                0217       ENDIF
                0218 
8ab93ea9a3 Jean*0219 #endif /* ALLOW_DIAGNOSTICS */
                0220 
                0221       RETURN
                0222       END