pkg/seaice/seaice_model.F

809c36b928 Patr*0001 #include "SEAICE_OPTIONS.h"
8377b8ee87 Mart*0002 #ifdef ALLOW_EXF
                0003 # include "EXF_OPTIONS.h"
                0004 #endif
772b2ed80e Gael*0005 #ifdef ALLOW_AUTODIFF
                0006 # include "AUTODIFF_OPTIONS.h"
                0007 #endif
7e1af563d8 Jean*0008
ef080e1d37 Dimi*0009 CBOP
                0010 C !ROUTINE: SEAICE_MODEL
                0011
                0012 C !INTERFACE: ==========================================================
b92ebad606 Jean*0013       SUBROUTINE SEAICE_MODEL( myTime, myIter, myThid )
ef080e1d37 Dimi*0014
                0015 C !DESCRIPTION: \bv
ff8e41478f Jean*0016 C     *===========================================================*
09510da3bb Dimi*0017 C     | SUBROUTINE SEAICE_MODEL                                   |
                0018 C     | o Time stepping of a dynamic/thermodynamic sea ice model. |
                0019 C     |  Dynamics solver: Zhang/Hibler, JGR, 102, 8691-8702, 1997 |
                0020 C     |  Thermodynamics:        Hibler, MWR, 108, 1943-1973, 1980 |
                0021 C     |  Rheology:              Hibler, JPO,   9,  815- 846, 1979 |
                0022 C     |  Snow:          Zhang et al.  , JPO,  28,  191- 217, 1998 |
                0023 C     |  Parallel forward ice model written by Jinlun Zhang PSC/UW|
                0024 C     |  & coupled into MITgcm by Dimitris Menemenlis (JPL) 2/2001|
                0025 C     |  zhang@apl.washington.edu / menemenlis@jpl.nasa.gov       |
45315406aa Mart*0026 C     | o The code has been rewritten substantially to use the    |
                0027 C     |   MITgcm C-grid, see Losch et al. OM, 33,  129- 144, 2010 |
ff8e41478f Jean*0028 C     *===========================================================*
809c36b928 Patr*0029       IMPLICIT NONE
7e1af563d8 Jean*0030 C \ev
                0031
ef080e1d37 Dimi*0032 C !USES: ===============================================================
809c36b928 Patr*0033 #include "SIZE.h"
                0034 #include "EEPARAMS.h"
                0035 #include "DYNVARS.h"
                0036 #include "PARAMS.h"
4634ed4578 Dimi*0037 #include "GRID.h"
809c36b928 Patr*0038 #include "FFIELDS.h"
ccaa3c61f4 Patr*0039 #include "SEAICE_SIZE.h"
809c36b928 Patr*0040 #include "SEAICE_PARAMS.h"
ccaa3c61f4 Patr*0041 #include "SEAICE.h"
                0042 #include "SEAICE_TRACER.h"
ae1fb66b64 Dimi*0043 #ifdef ALLOW_EXF
                0044 # include "EXF_FIELDS.h"
                0045 #endif
7109a141b2 Patr*0046 #ifdef ALLOW_AUTODIFF_TAMC
                0047 # include "tamc.h"
6060ec2938 Dimi*0048 #endif
baa476eeba Dimi*0049
ef080e1d37 Dimi*0050 C !INPUT PARAMETERS: ===================================================
f003822524 Jean*0051 C     myTime :: Current time in simulation
                0052 C     myIter :: Current iteration number in simulation
                0053 C     myThid :: my Thread Id number
809c36b928 Patr*0054       _RL     myTime
                0055       INTEGER myIter
                0056       INTEGER myThid
                0057
ef080e1d37 Dimi*0058 C !LOCAL VARIABLES: ====================================================
a1ab12d5e7 Dimi*0059 C     i,j,bi,bj :: Loop counters
e8a5a3baa4 Jean*0060       INTEGER i, j
                0061       INTEGER bi, bj
40d541aac0 Jean*0062 #ifdef ALLOW_EXF
                0063       INTEGER grpDiag
                0064 #endif
f50f58ec54 Gael*0065 #ifdef ALLOW_SITRACER
                0066       INTEGER iTr
                0067 #endif
45315406aa Mart*0068 #ifdef SEAICE_BGRID_DYNAMICS
f003822524 Jean*0069       _RL uLoc(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0070       _RL vLoc(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0071 #endif
ef080e1d37 Dimi*0072 CEOP
5fe78992ba Mart*0073 #ifdef ALLOW_DIAGNOSTICS
                0074       LOGICAL  DIAGNOSTICS_IS_ON
                0075       EXTERNAL DIAGNOSTICS_IS_ON
                0076       LOGICAL  diag_SIenph_isOn
                0077       _RL hTmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0078       _RL recip_deltaTtherm
                0079 #endif
b92ebad606 Jean*0080 #ifdef ALLOW_DEBUG
a8b75793c4 Jean*0081       IF (debugMode) CALL DEBUG_ENTER( 'SEAICE_MODEL', myThid )
b92ebad606 Jean*0082 #endif
                0083
0320e25227 Mart*0084       IF ( usingPCoords ) THEN
                0085 C     In z-coordinates, phiHydLow is just a diagnostics and hence its
                0086 C     overlaps are never filled properly
                0087 C     in p-coordinates, phiHydLow is the sea level elevation that is required
                0088 C     for the sea level tilt term in the momentum equation, so we need
                0089 C     to fill the overlaps properly here
                0090 C     (or elsewhere, not sure where it would be best)
                0091        _EXCH_XY_RL( phiHydLow, myThid )
                0092       ENDIF
47852c9c0c Mart*0093 #ifdef ALLOW_THSICE
                0094       IF ( useThSice ) THEN
20e951135e Mart*0095 C--   Map thSice-variables to HEFF, HSNOW, and AREA, because they are
                0096 C     needed in S/R SEAICE_DYNSOLVER
47852c9c0c Mart*0097        CALL SEAICE_MAP_THSICE( myTime, myIter, myThid )
                0098       ENDIF
                0099 #endif /* ALLOW_THSICE */
d5254b4e3d Mart*0100 #ifdef ALLOW_DIAGNOSTICS
                0101       IF ( useDiagnostics ) THEN
                0102        diag_SIenph_isOn = DIAGNOSTICS_IS_ON('SIenph  ',myThid)
                0103        IF ( diag_SIenph_isOn ) THEN
                0104         DO bj=myByLo(myThid),myByHi(myThid)
                0105          DO bi=myBxLo(myThid),myBxHi(myThid)
                0106           DO j=1-OLy,sNy+OLy
                0107            DO i=1-OLx,sNx+OLx
                0108             hTmp(i,j,bi,bj)=HEFF(i,j,bi,bj)
                0109            ENDDO
                0110           ENDDO
                0111          ENDDO
                0112         ENDDO
                0113        ENDIF
                0114       ENDIF
                0115 #endif /* ALLOW_DIAGNOSTICS */
47852c9c0c Mart*0116
40d541aac0 Jean*0117 #ifdef ALLOW_EXF
                0118       IF ( useEXF ) THEN
                0119 C--   Winds are from pkg/exf, which does not update edges.
                0120        CALL EXCH_UV_AGRID_3D_RL( uwind, vwind, .TRUE., 1, myThid )
                0121        IF ( useDiagnostics ) THEN
                0122 C-    Fill-in EXF wind-stess diags, weighted by open-ocean fraction
                0123         grpDiag = -1
                0124         IF ( SEAICEuseDYNAMICS ) grpDiag = 1
                0125         CALL EXF_WEIGHT_SFX_DIAGS(
                0126      I                  AREA, grpDiag, myTime, myIter, myThid )
                0127        ENDIF
                0128       ENDIF
                0129 #endif /* ALLOW_EXF */
                0130
8377b8ee87 Mart*0131 #ifdef ALLOW_AUTODIFF
c20cddf271 Patr*0132       DO bj=myByLo(myThid),myByHi(myThid)
                0133        DO bi=myBxLo(myThid),myBxHi(myThid)
d24acf3cfc Jean*0134         DO j=1-OLy,sNy+OLy
                0135          DO i=1-OLx,sNx+OLx
c20cddf271 Patr*0136           uIceNm1(i,j,bi,bj) = 0. _d 0
                0137           vIceNm1(i,j,bi,bj) = 0. _d 0
aef1621d33 Patr*0138 # ifdef ALLOW_SITRACER
                0139           DO iTr = 1, SItrMaxNum
                0140            SItrBucket(i,j,bi,bj,iTr) = 0. _d 0
                0141           ENDDO
                0142 # endif
c20cddf271 Patr*0143          ENDDO
                0144         ENDDO
                0145        ENDDO
                0146       ENDDO
8377b8ee87 Mart*0147 #endif /* ALLOW_AUTODIFF */
                0148 #ifdef ALLOW_AUTODIFF_TAMC
95c72ef3a1 Patr*0149 CADJ STORE heff  = comlev1, key=ikey_dynamics, kind=isbyte
                0150 CADJ STORE area  = comlev1, key=ikey_dynamics, kind=isbyte
3cb8290a3c Mart*0151 CADJ STORE hsnow = comlev1, key=ikey_dynamics, kind=isbyte
45315406aa Mart*0152 # ifdef SEAICE_CGRID
95c72ef3a1 Patr*0153 CADJ STORE fu    = comlev1, key=ikey_dynamics, kind=isbyte
                0154 CADJ STORE fv    = comlev1, key=ikey_dynamics, kind=isbyte
d5254b4e3d Mart*0155 #  ifdef SEAICE_ALLOW_EVP
b81c343615 Patr*0156 CADJ STORE seaice_sigma1  = comlev1, key=ikey_dynamics, kind=isbyte
                0157 CADJ STORE seaice_sigma2  = comlev1, key=ikey_dynamics, kind=isbyte
                0158 CADJ STORE seaice_sigma12 = comlev1, key=ikey_dynamics, kind=isbyte
d5254b4e3d Mart*0159 #  endif /* SEAICE_ALLOW_EVP */
45315406aa Mart*0160 # endif /* SEAICE_CGRID */
aef1621d33 Patr*0161 # ifdef ALLOW_SITRACER
d5254b4e3d Mart*0162 CADJ STORE siceload = comlev1, key=ikey_dynamics, kind=isbyte
                0163 CADJ STORE sitracer = comlev1, key=ikey_dynamics, kind=isbyte
aef1621d33 Patr*0164 # endif
7109a141b2 Patr*0165 #endif /* ALLOW_AUTODIFF_TAMC */
                0166
45315406aa Mart*0167 C--   Solve ice momentum equations and calculate ocean surface stress.
                0168 C     The surface stress always needs to be updated, even if neither B-
                0169 C     or C-grid dynamics are compiled, and SEAICEuseDYNAMICS = F.
b92ebad606 Jean*0170 #ifdef ALLOW_DEBUG
a8b75793c4 Jean*0171       IF (debugMode) CALL DEBUG_CALL( 'SEAICE_DYNSOLVER', myThid )
b92ebad606 Jean*0172 #endif
45315406aa Mart*0173 #ifdef SEAICE_BGRID_DYNAMICS
02aabb8fe5 Dimi*0174       CALL TIMER_START('DYNSOLVER          [SEAICE_MODEL]',myThid)
                0175       CALL DYNSOLVER ( myTime, myIter, myThid )
                0176       CALL TIMER_STOP ('DYNSOLVER          [SEAICE_MODEL]',myThid)
45315406aa Mart*0177 #else /* use default C-grid solver */
                0178       CALL TIMER_START('SEAICE_DYNSOLVER   [SEAICE_MODEL]',myThid)
                0179       CALL SEAICE_DYNSOLVER ( myTime, myIter, myThid )
                0180       CALL TIMER_STOP ('SEAICE_DYNSOLVER   [SEAICE_MODEL]',myThid)
                0181 #endif
809c36b928 Patr*0182
0f64a7d6fb Dimi*0183 C--   Apply ice velocity open boundary conditions
5459643feb Dimi*0184 #ifdef ALLOW_OBCS
0f64a7d6fb Dimi*0185 # ifndef DISABLE_SEAICE_OBCS
d5254b4e3d Mart*0186       IF ( useOBCS ) CALL OBCS_ADJUST_UVICE( uice, vice, myThid )
0f64a7d6fb Dimi*0187 # endif /* DISABLE_SEAICE_OBCS */
5459643feb Dimi*0188 #endif /* ALLOW_OBCS */
                0189
45315406aa Mart*0190 #if ( defined ALLOW_AUTODIFF_TAMC &&  defined SEAICE_CGRID )
d5254b4e3d Mart*0191 CADJ STORE uice  = comlev1, key=ikey_dynamics, kind=isbyte
                0192 CADJ STORE vice  = comlev1, key=ikey_dynamics, kind=isbyte
                0193 C     Note: Storing u/vice **after** seaice_dynsolver (and obcs_adjust_uvice)
                0194 C     has the effect that seaice_dynsolver is not called from seaice_model_ad
                0195 C     anymore. This is important because with the numerous tricks to avoid
                0196 C     complicated code in the backward integration (see pkg/autodiff),
                0197 C     the extra call would update the ice velocities with the wrong
                0198 C     set of flag values.
45315406aa Mart*0199 #endif
d5254b4e3d Mart*0200
4a98994c14 Jean*0201 #ifdef ALLOW_THSICE
de80e32bba Jean*0202       IF ( useThSice ) THEN
                0203 #ifdef ALLOW_DEBUG
d5254b4e3d Mart*0204        IF (debugMode) CALL DEBUG_CALL( 'THSICE_DO_ADVECT', myThid )
de80e32bba Jean*0205 #endif
d5254b4e3d Mart*0206        CALL THSICE_DO_ADVECT( 0, 0, myTime, myIter, myThid )
de80e32bba Jean*0207       ELSE
bfae2ee370 jm-c 0208 #endif /* ALLOW_THSICE */
ff8e41478f Jean*0209 C--   Only call advection of heff, area, snow, and salt and
09e3b53265 Mart*0210 C--   growth for the generic 0-layer thermodynamics of seaice
                0211 C--   if useThSice=.false., otherwise the 3-layer Winton thermodynamics
                0212 C--   (called from DO_OCEANIC_PHYSICS) take care of this
                0213 C NOW DO ADVECTION and DIFFUSION
d5254b4e3d Mart*0214        IF ( SEAICEadvHeff .OR. SEAICEadvArea .OR. SEAICEadvSnow
e54fe3e1f9 Gael*0215      &        .OR. SEAICEadvSalt ) THEN
b92ebad606 Jean*0216 #ifdef ALLOW_DEBUG
d5254b4e3d Mart*0217         IF (debugMode) CALL DEBUG_CALL( 'SEAICE_ADVDIFF', myThid )
b92ebad606 Jean*0218 #endif
45315406aa Mart*0219 C--   There always needs to be advection, even if neither B- or C-grid
                0220 C     dynamics are compiled.
                0221 #ifdef SEAICE_BGRID_DYNAMICS
d5254b4e3d Mart*0222         CALL SEAICE_ADVDIFF( uLoc, vLoc, myTime, myIter, myThid )
45315406aa Mart*0223 #else /* default C-grid advection */
                0224         CALL SEAICE_ADVDIFF( uIce, vIce, myTime, myIter, myThid )
                0225 #endif
d5254b4e3d Mart*0226        ENDIF
                0227 #ifdef ALLOW_AUTODIFF_TAMC
                0228 CADJ STORE heff  = comlev1, key=ikey_dynamics, kind=isbyte
                0229 CADJ STORE area  = comlev1, key=ikey_dynamics, kind=isbyte
                0230 CADJ STORE hsnow = comlev1, key=ikey_dynamics, kind=isbyte
                0231 # ifdef SEAICE_VARIABLE_SALINITY
                0232 CADJ STORE hsalt = comlev1, key=ikey_dynamics, kind=isbyte
                0233 # endif
                0234 C     Note: This store has the effect that seaice_advdiff is not called
                0235 C     from seaice_model_ad anymore. Instead, the stored values are used.
                0236 #endif /* ALLOW_AUTODIFF_TAMC */
809c36b928 Patr*0237
1cf549c217 Mart*0238 C     After advection, the sea ice variables may have unphysical values
20970ee22f Jean*0239 C     e.g., < 0, that are regularized here. Concentration as a special case
1cf549c217 Mart*0240 C     may be > 1 in convergent motion and a ridging algorithm redistributes
                0241 C     the ice to limit the concentration to 1.
d5254b4e3d Mart*0242        CALL SEAICE_REG_RIDGE( myTime, myIter, myThid )
1cf549c217 Mart*0243
40d541aac0 Jean*0244 #ifdef ALLOW_EXF
d5254b4e3d Mart*0245        IF ( useEXF .AND. useDiagnostics ) THEN
40d541aac0 Jean*0246 C-    Fill-in EXF surface flux diags, weighted by open-ocean fraction
                0247         grpDiag = -2
                0248         IF ( usePW79thermodynamics ) grpDiag = 2
                0249         CALL EXF_WEIGHT_SFX_DIAGS(
                0250      I                  AREA, grpDiag, myTime, myIter, myThid )
d5254b4e3d Mart*0251        ENDIF
40d541aac0 Jean*0252 #endif /* ALLOW_EXF */
                0253
bfae2ee370 jm-c 0254 #ifdef DISABLE_SEAICE_GROWTH
                0255        IF ( .TRUE. ) THEN
                0256 #else /* DISABLE_SEAICE_GROWTH */
09e3b53265 Mart*0257 C     thermodynamics growth
ff8e41478f Jean*0258 C     must call growth after calling advection
09e3b53265 Mart*0259 C     because of ugly time level business
                0260        IF ( usePW79thermodynamics ) THEN
45315406aa Mart*0261 # ifdef SEAICE_USE_GROWTH_ADX
                0262 #  ifdef ALLOW_DEBUG
dc54d31829 Ian *0263         IF (debugMode) CALL DEBUG_CALL( 'SEAICE_GROWTH_ADX', myThid )
45315406aa Mart*0264 #  endif
dc54d31829 Ian *0265         CALL SEAICE_GROWTH_ADX( myTime, myIter, myThid )
382462ccb5 Mart*0266 # else /* SEAICE_USE_GROWTH_ADX */
45315406aa Mart*0267 #  ifdef ALLOW_DEBUG
a8b75793c4 Jean*0268         IF (debugMode) CALL DEBUG_CALL( 'SEAICE_GROWTH', myThid )
45315406aa Mart*0269 #  endif
09e3b53265 Mart*0270         CALL SEAICE_GROWTH( myTime, myIter, myThid )
382462ccb5 Mart*0271 # endif /* SEAICE_USE_GROWTH_ADX */
bfae2ee370 jm-c 0272        ELSE
4ce2104710 Patr*0273 #endif /* DISABLE_SEAICE_GROWTH */
bfae2ee370 jm-c 0274         DO bj=myByLo(myThid),myByHi(myThid)
                0275          DO bi=myBxLo(myThid),myBxHi(myThid)
                0276           DO j=1,sNy
                0277            DO i=1,sNx
                0278             sIceLoad(i,j,bi,bj) = HEFF(i,j,bi,bj)*SEAICE_rhoIce
                0279      &                         + HSNOW(i,j,bi,bj)*SEAICE_rhoSnow
                0280            ENDDO
                0281           ENDDO
                0282 c#ifdef SEAICE_CAP_ICELOAD
                0283 c         sIceTooHeavy = rhoConst*drF(1) / 5. _d 0
                0284 c         DO j=1,sNy
                0285 c          DO i=1,sNx
                0286 c           sIceLoad(i,j,bi,bj) = MIN( sIceLoad(i,j,bi,bj),
                0287 c    &                                 sIceTooHeavy )
                0288 c          ENDDO
                0289 c         ENDDO
                0290 c#endif
                0291          ENDDO
                0292         ENDDO
                0293        ENDIF
809c36b928 Patr*0294
f50f58ec54 Gael*0295 #ifdef ALLOW_SITRACER
aef1621d33 Patr*0296 # ifdef ALLOW_AUTODIFF_TAMC
                0297 CADJ STORE sitracer  = comlev1, key=ikey_dynamics, kind=isbyte
                0298 # endif
a746f35887 Jean*0299        CALL SEAICE_TRACER_PHYS ( myTime, myIter, myThid )
f50f58ec54 Gael*0300 #endif
                0301
0f64a7d6fb Dimi*0302 C--   Apply ice tracer open boundary conditions
a1ab12d5e7 Dimi*0303 #ifdef ALLOW_OBCS
0f64a7d6fb Dimi*0304 # ifndef DISABLE_SEAICE_OBCS
                0305        IF ( useOBCS ) CALL OBCS_APPLY_SEAICE( myThid )
                0306 # endif /* DISABLE_SEAICE_OBCS */
a1ab12d5e7 Dimi*0307 #endif /* ALLOW_OBCS */
                0308
809c36b928 Patr*0309 C--   Update overlap regions for a bunch of stuff
772590b63c Mart*0310        _EXCH_XY_RL( HEFF,  myThid )
                0311        _EXCH_XY_RL( AREA,  myThid )
                0312        _EXCH_XY_RL( HSNOW, myThid )
963c7c7f95 Torg*0313 #ifdef SEAICE_ITD
                0314        CALL EXCH_3D_RL( HEFFITD,  nITD, myThid )
                0315        CALL EXCH_3D_RL( AREAITD,  nITD, myThid )
                0316        CALL EXCH_3D_RL( HSNOWITD, nITD, myThid )
                0317 #endif
a98c4b8072 Ian *0318 #ifdef SEAICE_VARIABLE_SALINITY
772590b63c Mart*0319        _EXCH_XY_RL( HSALT, myThid )
fdfa8e151f Dimi*0320 #endif
f50f58ec54 Gael*0321 #ifdef ALLOW_SITRACER
be02c52974 Gael*0322        DO iTr = 1, SItrNumInUse
d24acf3cfc Jean*0323         _EXCH_XY_RL( SItracer(1-OLx,1-OLy,1,1,iTr),myThid )
f50f58ec54 Gael*0324        ENDDO
                0325 #endif
14cee76775 Jean*0326        _EXCH_XY_RS(EmPmR, myThid )
                0327        _EXCH_XY_RS(saltFlux, myThid )
                0328        _EXCH_XY_RS(Qnet , myThid )
809c36b928 Patr*0329 #ifdef SHORTWAVE_HEATING
14cee76775 Jean*0330        _EXCH_XY_RS(Qsw  , myThid )
e4775240e5 Dimi*0331 #endif /* SHORTWAVE_HEATING */
f3ce416a61 Jean*0332 #ifdef ATMOSPHERIC_LOADING
7e1af563d8 Jean*0333        IF ( useRealFreshWaterFlux )
09e3b53265 Mart*0334      &      _EXCH_XY_RS( sIceLoad, myThid )
f3ce416a61 Jean*0335 #endif
809c36b928 Patr*0336
4f6e464f38 Jean*0337 #ifdef ALLOW_OBCS
3352db77f4 Jean*0338 C--   In case we use scheme with a large stencil that extends into overlap:
                0339 C     no longer needed with the right masking in advection & diffusion S/R.
                0340 c      IF ( useOBCS ) THEN
                0341 c       DO bj=myByLo(myThid),myByHi(myThid)
                0342 c        DO bi=myBxLo(myThid),myBxHi(myThid)
d24acf3cfc Jean*0343 c          CALL OBCS_COPY_TRACER( HEFF(1-OLx,1-OLy,bi,bj),
3352db77f4 Jean*0344 c    I                            1, bi, bj, myThid )
d24acf3cfc Jean*0345 c          CALL OBCS_COPY_TRACER( AREA(1-OLx,1-OLy,bi,bj),
3352db77f4 Jean*0346 c    I                            1, bi, bj, myThid )
d24acf3cfc Jean*0347 c          CALL OBCS_COPY_TRACER( HSNOW(1-OLx,1-OLy,bi,bj),
3352db77f4 Jean*0348 c    I                            1, bi, bj, myThid )
a98c4b8072 Ian *0349 #ifdef SEAICE_VARIABLE_SALINITY
d24acf3cfc Jean*0350 c          CALL OBCS_COPY_TRACER( HSALT(1-OLx,1-OLy,bi,bj),
3352db77f4 Jean*0351 c    I                            1, bi, bj, myThid )
4f6e464f38 Jean*0352 #endif
3352db77f4 Jean*0353 c        ENDDO
                0354 c       ENDDO
                0355 c      ENDIF
4f6e464f38 Jean*0356 #endif /* ALLOW_OBCS */
                0357
d74b8b60ab Patr*0358 #ifdef ALLOW_DIAGNOSTICS
                0359        IF ( useDiagnostics ) THEN
9ef71a93d2 Mart*0360 C     diagnostics for "non-state variables" that are modified by
d5254b4e3d Mart*0361 C     the seaice model ...
ead4e7560e Jean*0362         CALL DIAGNOSTICS_FILL_RS(EmPmR,'SIempmr ',0,1 ,0,1,1,myThid)
                0363         CALL DIAGNOSTICS_FILL_RS(Qnet ,'SIqnet  ',0,1 ,0,1,1,myThid)
                0364         CALL DIAGNOSTICS_FILL_RS(Qsw  ,'SIqsw   ',0,1 ,0,1,1,myThid)
d5254b4e3d Mart*0365 C     ... and energy diagnostic
5fe78992ba Mart*0366         IF ( diag_SIenph_isOn ) THEN
                0367          recip_deltaTtherm = 1. _d 0 / SEAICE_deltaTtherm
                0368          DO bj=myByLo(myThid),myByHi(myThid)
                0369           DO bi=myBxLo(myThid),myBxHi(myThid)
                0370            DO j=1,sNy
                0371             DO i=1,sNx
                0372              hTmp(i,j,bi,bj) = 0.5 _d 0 * SEAICE_rhoIce
                0373      &            * 0.5 _d 0 * (
                0374      &              uIce (i,j,bi,bj)**2 + uIce (i+1,j,  bi,bj)**2
                0375      &            + vIce (i,j,bi,bj)**2 + vIce (i,  j+1,bi,bj)**2 )
                0376      &            * (HEFF(i,j,bi,bj) - hTmp(i,j,bi,bj))
                0377      &             *recip_deltaTtherm
                0378             ENDDO
                0379            ENDDO
                0380           ENDDO
                0381          ENDDO
                0382          CALL DIAGNOSTICS_FILL(hTmp,'SIenph  ',0,1,0,1,1,myThid)
                0383         ENDIF
d74b8b60ab Patr*0384        ENDIF
                0385 #endif /* ALLOW_DIAGNOSTICS */
                0386
09e3b53265 Mart*0387 #ifdef ALLOW_THSICE
                0388 C     endif .not.useThSice
                0389       ENDIF
                0390 #endif /* ALLOW_THSICE */
                0391 CML   This has already been done in seaice_ocean_stress/ostres, so why repeat?
                0392 CML   CALL EXCH_UV_XY_RS(fu,fv,.TRUE.,myThid)
                0393
fb1912d055 Patr*0394 #ifdef ALLOW_EXF
5d96c49326 Jean*0395 # ifdef ALLOW_AUTODIFF
                0396 #  ifdef ALLOW_AUTODIFF_MONITOR
d5254b4e3d Mart*0397       CALL EXF_ADJOINT_SNAPSHOTS( 3, myTime, myIter, myThid )
fb1912d055 Patr*0398 #  endif
                0399 # endif
                0400 #endif
                0401
b92ebad606 Jean*0402 #ifdef ALLOW_DEBUG
a8b75793c4 Jean*0403       IF (debugMode) CALL DEBUG_LEAVE( 'SEAICE_MODEL', myThid )
b92ebad606 Jean*0404 #endif
                0405
809c36b928 Patr*0406       RETURN
                0407       END