Back to home page

MITgcm

 
 

    


File indexing completed on 2023-06-14 05:10:20 UTC

view on githubraw file Latest commit 5b0716a6 on 2023-06-13 17:16:45 UTC
aef2e94c45 Davi*0001 #include "GGL90_OPTIONS.h"
                0002 
                0003 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
a10c595eb6 Timo*0004 CBOP
aef2e94c45 Davi*0005 C !ROUTINE: GGL90_DIAGNOSTICS_INIT
                0006 C !INTERFACE:
                0007       SUBROUTINE GGL90_DIAGNOSTICS_INIT( myThid )
                0008 
31a3206180 Mart*0009 C     !DESCRIPTION: \bv
                0010 C     *================================================================*
                0011 C     | Initialize list of all available diagnostics
                0012 C     *================================================================*
                0013 C     \ev
aef2e94c45 Davi*0014 
                0015 C     !USES:
                0016       IMPLICIT NONE
31a3206180 Mart*0017 C     === Global variables ===
aef2e94c45 Davi*0018 #include "EEPARAMS.h"
                0019 #include "SIZE.h"
                0020 #include "GGL90.h"
                0021 
                0022 C     !INPUT/OUTPUT PARAMETERS:
31a3206180 Mart*0023 C     == Routine arguments ==
aef2e94c45 Davi*0024 C     myThid ::  my Thread Id number
                0025       INTEGER myThid
                0026 CEOP
                0027 
                0028 #ifdef ALLOW_DIAGNOSTICS
                0029 C     !LOCAL VARIABLES:
                0030 C     === Local variables ===
                0031 C     diagNum   :: diagnostics number in the (long) list of available diag.
                0032 C     diagMate  :: diag. mate number in the (long) list of available diag.
                0033 C     diagName  :: local short name (8c) of a diagnostics
                0034 C     diagCode  :: local parser field with characteristics of the diagnostics
                0035 C              cf head of S/R DIAGNOSTICS_INIT_EARLY or DIAGNOSTICS_MAIN_INIT
                0036 C     diagUnits :: local string (16c): physical units of a diagnostic field
                0037 C     diagTitle :: local string (80c): description of field in diagnostic
                0038       INTEGER       diagNum
24ffd6fe01 Jean*0039 c     INTEGER       diagMate
aef2e94c45 Davi*0040       CHARACTER*8   diagName
                0041       CHARACTER*16  diagCode
                0042       CHARACTER*16  diagUnits
                0043       CHARACTER*(80) diagTitle
31a3206180 Mart*0044 CEOP
aef2e94c45 Davi*0045 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0046 
                0047         diagName  = 'GGL90TKE'
                0048         diagTitle = 'GGL90 sub-grid turbulent kinetic energy'
                0049         diagUnits = 'm^2/s^2         '
                0050         diagCode  = 'SM      LR      '
                0051         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0052      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0053 
5b0716a6b3 Mart*0054         diagName  = 'GGL90Emn'
                0055         diagTitle = 'rate of TKE energy added by applying GGL90TKEmin'
                0056         diagUnits = 'm^2/s^3         '
                0057         diagCode  = 'SM      LR      '
                0058         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0059      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0060 
aef2e94c45 Davi*0061         diagName  = 'GGL90Lmx'
                0062         diagTitle = 'Mixing length scale              '
                0063         diagUnits = 'm               '
                0064         diagCode  = 'SM      LR      '
                0065         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0066      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0067 
                0068         diagName  = 'GGL90Prl'
                0069         diagTitle = 'Prandtl number used in GGL90'
                0070         diagUnits = '1               '
                0071         diagCode  = 'SM      LR      '
                0072         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0073      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0074 
f6b150f7f1 Gael*0075         diagName  = 'GGL90ArU'
                0076         diagTitle = 'GGL90 eddy viscosity at U-point'
                0077         diagUnits = 'm^2/s           '
                0078         diagCode  = 'SM      LR      '
                0079         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0080      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0081 
                0082         diagName  = 'GGL90ArV'
                0083         diagTitle = 'GGL90 eddy viscosity at V-point'
aef2e94c45 Davi*0084         diagUnits = 'm^2/s           '
                0085         diagCode  = 'SM      LR      '
                0086         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0087      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0088 
56d13a40ed Mart*0089         diagName  = 'GGL90Kr '
aef2e94c45 Davi*0090         diagTitle = 'GGL90 diffusion coefficient for temperature'
                0091         diagUnits = 'm^2/s           '
                0092         diagCode  = 'SM      LR      '
                0093         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0094      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0095 
5b0716a6b3 Mart*0096         diagName  = 'GGL90KN2'
                0097         diagTitle = 'GGL90 diffusivity times buoyancy frequency'
                0098         diagUnits = 'm^2/s^3         '
                0099         diagCode  = 'SM      LR      '
                0100         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0101      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0102 
31a3206180 Mart*0103         diagName  = 'GGL90flx'
                0104         diagTitle = 'Surface flux of TKE                       '
                0105         diagUnits = 'm^3/s^3         '
                0106         diagCode  = 'SM      L1      '
                0107         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0108      I       diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0109 
                0110         diagName  = 'GGL90tau'
                0111         diagTitle = 'Work done by the wind                     '
                0112         diagUnits = 'm^3/s^3         '
                0113         diagCode  = 'SM      L1      '
                0114         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0115      I       diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0116 
                0117 #ifdef ALLOW_GGL90_IDEMIX
                0118         diagName  = 'IDEMIX_E'
                0119         diagTitle = 'IDEMIX internal wave energy            '
                0120         diagUnits = 'm^2/s^2         '
                0121         diagCode  = 'SM      LR      '
                0122         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0123      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0124 
5b0716a6b3 Mart*0125         diagName  = 'IDEMgTKE'
                0126         diagTitle = 'IDEMIX tendency: tau_d times IDEMIX_E^2'
                0127         diagUnits = 'm^2/s^3         '
                0128         diagCode  = 'SM      LR      '
                0129         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0130      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0131 
31a3206180 Mart*0132         diagName  = 'IDEMIX_c'
                0133         diagTitle = 'IDEMIX vertical group velocity             '
                0134         diagUnits = 'm/s             '
                0135         diagCode  = 'SM      LR      '
                0136         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0137      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0138 
                0139         diagName  = 'IDEMIX_v'
                0140         diagTitle = 'IDEMIX horizontal group velocity           '
                0141         diagUnits = 'm/s             '
                0142         diagCode  = 'SM      LR      '
                0143         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0144      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0145 
                0146         diagName  = 'IDEMIX_t'  ! m^2/s^3 /(m^4/s^4)
                0147         diagTitle = 'IDEMIX dissipation constant                '
                0148         diagUnits = 's/m^2           '
                0149         diagCode  = 'SM      LR      '
                0150         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0151      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0152 
                0153         diagName  = 'IDEMIX_K'
                0154         diagTitle = 'IDEMIX vertical diffusivity                '
                0155         diagUnits = 'm^2/s           '
                0156         diagCode  = 'SM      LR      '
                0157         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0158      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0159 
                0160         diagName  = 'IDEMIX_F'
                0161         diagTitle = 'IDEMIX Forcing by gm                       '
                0162         diagUnits = 'm^2/s^3         '
                0163         diagCode  = 'SM      LR      '
                0164         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0165      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0166 
                0167         diagName  = 'IDEM_F_b'
                0168         diagTitle = 'Tidal forcing at bottom                    '
                0169         diagUnits = 'm^3/s^3         '
                0170         diagCode  = 'SM      L1      '
                0171         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0172      I       diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0173 
                0174         diagName  = 'IDEM_F_s'
                0175         diagTitle = 'Wind forcing at surface                   '
                0176         diagUnits = 'm^3/s^3         '
                0177         diagCode  = 'SM      L1      '
                0178         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0179      I       diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0180 
                0181         diagName  = 'IDEM_F_g'
                0182         diagTitle = 'Integrated GM forcing                     '
                0183         diagUnits = 'm^3/s^3         '
                0184         diagCode  = 'SM      L1      '
                0185         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0186      I       diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0187 
                0188 #endif
                0189 
a10c595eb6 Timo*0190 #ifdef ALLOW_AUTODIFF
                0191       diagName  = 'ADJtke90'
                0192       diagTitle = 'dJ/dTKE: Sensitivity to GGL90 TKE'
                0193       diagUnits = 'dJ/(m^2/s^2)   '
                0194       diagCode  = 'SMRA    MR      '
                0195       CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0196      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0197 #endif
                0198 
aef2e94c45 Davi*0199 #endif /* ALLOW_DIAGNOSTICS */
                0200 
                0201       RETURN
                0202       END