Back to home page

MITgcm

 
 

    


File indexing completed on 2025-06-05 05:08:28 UTC

view on githubraw file Latest commit 6a6c83f9 on 2025-06-04 22:00:11 UTC
4e66ab0b67 Oliv*0001 #include "LONGSTEP_OPTIONS.h"
                0002 
                0003 CBOP
                0004 C     !ROUTINE: LONGSTEP_DIAGNOSTICS_INIT
                0005 C     !INTERFACE:
                0006       SUBROUTINE LONGSTEP_DIAGNOSTICS_INIT( myThid )
                0007 C     !DESCRIPTION:
                0008 C     Routine to initialize longstep diagnostics
                0009 
                0010 C     !USES:
                0011       IMPLICIT NONE
                0012 C     === Global variables ===
                0013 #include "SIZE.h"
                0014 #include "EEPARAMS.h"
                0015 #include "PARAMS.h"
                0016 #include "LONGSTEP_PARAMS.h"
                0017 
                0018 C     !INPUT/OUTPUT PARAMETERS:
                0019 C     === Routine arguments ===
                0020 C     myThid -  Number of this instance of LONGSTEP_DIAGNOSTICS_INIT
                0021       INTEGER myThid
                0022 CEOP
                0023 
                0024 #ifdef ALLOW_DIAGNOSTICS
                0025 C     !LOCAL VARIABLES:
                0026 C     === Local variables ===
8e7e785cad Jean*0027 C     msgBuf      - Informational/error message buffer
4e66ab0b67 Oliv*0028 c     CHARACTER*(MAX_LEN_MBUF) msgBuf
                0029 
                0030       INTEGER       diagNum
                0031       INTEGER       diagMate
                0032       CHARACTER*8   diagName
                0033       CHARACTER*16  diagCode
                0034       CHARACTER*16  diagUnits
                0035       CHARACTER*(80) diagTitle
                0036 
                0037 C--   Add diagnostics to the (long) list
                0038 
                0039       diagName  = 'LSwVel  '
                0040       diagTitle = 'Vertical Component of Velocity (m/s)'
                0041       diagUnits = 'm/s             '
                0042       diagCode  = 'WM      LR      '
                0043       CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0044      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0045 
                0046       diagName  = 'LSuVel  '
                0047       diagTitle = 'Zonal Component of Velocity (m/s)'
                0048       diagUnits = 'm/s             '
8e7e785cad Jean*0049       diagCode  = 'UUR     MR      '
4e66ab0b67 Oliv*0050       diagMate  = diagNum + 2
                0051       CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0052      I   diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
                0053 
                0054       diagName  = 'LSvVel  '
                0055       diagTitle = 'Meridional Component of Velocity (m/s)'
                0056       diagUnits = 'm/s             '
8e7e785cad Jean*0057       diagCode  = 'VVR     MR      '
4e66ab0b67 Oliv*0058       diagMate  = diagNum
                0059       CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0060      I   diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
                0061 
                0062       diagName  = 'LStheta '
                0063       diagTitle = 'Potential Temperature'
                0064       diagUnits = 'degC            '
8e7e785cad Jean*0065       diagCode  = 'SMR     MR      '
4e66ab0b67 Oliv*0066       CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0067      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0068 
                0069       diagName  = 'LSsalt  '
                0070       diagTitle = 'Salinity'
ba0b047096 Mart*0071       diagUnits = 'g/kg            '
8e7e785cad Jean*0072       diagCode  = 'SMR     MR      '
4e66ab0b67 Oliv*0073       CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0074      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0075 
                0076       IF ( ivdc_kappa .NE. 0. _d 0 ) THEN
                0077         diagName  = 'LScnvAdj'
                0078         diagTitle = 'Convective Adjustment Index [0-1] '
                0079         diagUnits = 'fraction        '
                0080         diagCode  = 'SM      LR      '
                0081         CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0082      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0083       ENDIF
                0084 
                0085 #ifdef SHORTWAVE_HEATING
                0086       diagName  = 'LSqsw   '
                0087       diagTitle = 'net Short-Wave radiation (+=up)'
                0088       diagUnits = 'W/m^2           '
                0089       diagCode  = 'SM      U1      '
                0090       CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0091      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0092 #endif
                0093 
fda3710353 Oliv*0094       diagName  = 'LSfwFlux'
                0095       diagTitle = 'net surface Fresh-Water flux into the ocean'
                0096       diagUnits = 'kg/m^2/s        '
                0097       diagCode  = 'SM      U1      '
                0098       CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0099      I          diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0100 
4e66ab0b67 Oliv*0101 #ifdef ALLOW_GMREDI
                0102       IF ( useGMRedi ) THEN
                0103        diagName  = 'LSkwx   '
                0104        diagTitle = 'K_31 element (W.point, X.dir) of GM-Redi tensor'
                0105        diagUnits = 'm^2/s           '
8e7e785cad Jean*0106        diagCode  = 'UM      LR      '
                0107        diagMate  = diagNum + 2
                0108        CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0109      I    diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
4e66ab0b67 Oliv*0110 
                0111        diagName  = 'LSkwy   '
                0112        diagTitle = 'K_32 element (W.point, Y.dir) of GM-Redi tensor'
                0113        diagUnits = 'm^2/s           '
8e7e785cad Jean*0114        diagCode  = 'VM      LR      '
                0115        diagMate  = diagNum
                0116        CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0117      I    diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
4e66ab0b67 Oliv*0118 
                0119        diagName  = 'LSkwz   '
                0120        diagTitle = 'K_33 element (W.point, Z.dir) of GM-Redi tensor'
                0121        diagUnits = 'm^2/s           '
                0122        diagCode  = 'WM P    LR      '
8e7e785cad Jean*0123        CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0124      I       diagName, diagCode, diagUnits, diagTitle, 0, myThid )
6a6c83f9ac Hajo*0125 
                0126        diagName  = 'LSkux   '
                0127        diagTitle = 'K_11 element (U.point, X.dir) of GM-Redi tensor'
                0128        diagUnits = 'm^2/s           '
                0129        diagCode  = 'UU P    MR      '
                0130        CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0131      I       diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0132 
                0133        diagName  = 'LSkvy   '
                0134        diagTitle = 'K_22 element (V.point, Y.dir) of GM-Redi tensor'
                0135        diagUnits = 'm^2/s           '
                0136        diagCode  = 'VV P    MR      '
                0137        CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0138      I       diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0139 
                0140        diagName  = 'LSkuz   '
                0141        diagTitle = 'K_13 element (U.point, Z.dir) of GM-Redi tensor'
                0142        diagUnits = 'm^2/s           '
                0143        diagCode  = 'UU      MR      '
                0144        CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0145      I       diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0146 
                0147        diagName  = 'LSkvz   '
                0148        diagTitle = 'K_23 element (V.point, Z.dir) of GM-Redi tensor'
                0149        diagUnits = 'm^2/s           '
                0150        diagCode  = 'VV      MR      '
                0151        CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0152      I       diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0153 
                0154        diagName  = 'LS_PsiX '
                0155        diagTitle = 'GM Bolus transport Psi in longstep : U component'
                0156        diagUnits = 'm^2/s           '
                0157        diagCode  = 'UU      LR      '
                0158        CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0159      I       diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0160 
                0161        diagName  = 'LS_PsiY '
                0162        diagTitle = 'GM Bolus transport Psi in longstep : V component'
                0163        diagUnits = 'm^2/s           '
                0164        diagCode  = 'VV      LR      '
                0165        CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0166      I       diagName, diagCode, diagUnits, diagTitle, 0, myThid )
                0167 
4e66ab0b67 Oliv*0168       ENDIF
                0169 #endif
                0170 
                0171 #ifdef ALLOW_KPP
                0172       IF ( useKPP ) THEN
                0173        diagName  = 'LSKPPdfS'
                0174        diagTitle = 'Vertical diffusion coefficient for salt & tracers'
                0175        diagUnits = 'm^2/s           '
8e7e785cad Jean*0176        diagCode  = 'SM P    LR      '
                0177        CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0178      I       diagName, diagCode, diagUnits, diagTitle, 0, myThid )
4e66ab0b67 Oliv*0179 
                0180        diagName  = 'LSKPPght'
                0181        diagTitle = 'Nonlocal transport coefficient'
                0182        diagUnits = 's/m^2           '
                0183        diagCode  = 'SM P    LR      '
8e7e785cad Jean*0184        CALL DIAGNOSTICS_ADDTOLIST( diagNum,
                0185      I       diagName, diagCode, diagUnits, diagTitle, 0, myThid )
4e66ab0b67 Oliv*0186       ENDIF
                0187 #endif /* ALLOW_KPP */
                0188 
                0189 #endif /* ALLOW_DIAGNOSTICS */
                0190 
                0191       RETURN
                0192       END