model/src/temp_integrate.F

4794324d0b Jean*0001 #include "PACKAGES_CONFIG.h"
                0002 #include "CPP_OPTIONS.h"
517dbdc414 Jean*0003 #ifdef ALLOW_AUTODIFF
                0004 # include "AUTODIFF_OPTIONS.h"
                0005 #endif
fb8688a588 Jean*0006 #ifdef ALLOW_GENERIC_ADVDIFF
                0007 # include "GAD_OPTIONS.h"
                0008 #endif
4794324d0b Jean*0009
                0010 CBOP
                0011 C     !ROUTINE: TEMP_INTEGRATE
                0012 C     !INTERFACE:
                0013       SUBROUTINE TEMP_INTEGRATE(
0b191c5f5a Jean*0014      I           bi, bj, recip_hFac,
fb8688a588 Jean*0015      I           uFld, vFld, wFld,
                0016      U           KappaRk,
4794324d0b Jean*0017      I           myTime, myIter, myThid )
                0018 C     !DESCRIPTION: \bv
                0019 C     *==========================================================*
                0020 C     | SUBROUTINE TEMP_INTEGRATE
bd27360393 Jean*0021 C     | o Calculate tendency for temperature and integrates
                0022 C     |   forward in time. The temperature array is updated here
                0023 C     |   while adjustments (filters, conv.adjustment) are applied
                0024 C     |   later, in S/R TRACERS_CORRECTION_STEP.
4794324d0b Jean*0025 C     *==========================================================*
e874fa47e5 Jean*0026 C     | A procedure called APPLY_FORCING_T is called from
4794324d0b Jean*0027 C     | here. These procedures can be used to add per problem
                0028 C     | heat flux source terms.
                0029 C     | Note: Although it is slightly counter-intuitive the
                0030 C     |       EXTERNAL_FORCING routine is not the place to put
                0031 C     |       file I/O. Instead files that are required to
                0032 C     |       calculate the external source terms are generally
                0033 C     |       read during the model main loop. This makes the
                0034 C     |       logistics of multi-processing simpler and also
                0035 C     |       makes the adjoint generation simpler. It also
                0036 C     |       allows for I/O to overlap computation where that
                0037 C     |       is supported by hardware.
                0038 C     | Aside from the problem specific term the code here
                0039 C     | forms the tendency terms due to advection and mixing
                0040 C     | The baseline implementation here uses a centered
                0041 C     | difference form for the advection term and a tensorial
                0042 C     | divergence of a flux form for the diffusive term. The
fb8688a588 Jean*0043 C     | diffusive term is formulated so that isopycnal mixing
                0044 C     | and GM-style subgrid-scale terms can be incorporated by
                0045 C     | simply setting the diffusion tensor terms appropriately.
4794324d0b Jean*0046 C     *==========================================================*
                0047 C     \ev
                0048
                0049 C     !USES:
                0050       IMPLICIT NONE
                0051 C     == GLobal variables ==
                0052 #include "SIZE.h"
                0053 #include "EEPARAMS.h"
                0054 #include "PARAMS.h"
0b191c5f5a Jean*0055 #include "GRID.h"
                0056 #include "DYNVARS.h"
4794324d0b Jean*0057 #include "RESTART.h"
                0058 #ifdef ALLOW_GENERIC_ADVDIFF
fb8688a588 Jean*0059 # include "GAD.h"
                0060 # include "GAD_SOM_VARS.h"
4794324d0b Jean*0061 #endif
7c50f07931 Mart*0062 #ifdef ALLOW_AUTODIFF_TAMC
4794324d0b Jean*0063 # include "tamc.h"
                0064 #endif
                0065
                0066 C     !INPUT/OUTPUT PARAMETERS:
                0067 C     == Routine arguments ==
0b191c5f5a Jean*0068 C     bi, bj,    :: tile indices
                0069 C     recip_hFac :: reciprocal of cell open-depth factor (@ next iter)
                0070 C     uFld,vFld  :: Local copy of horizontal velocity field
                0071 C     wFld       :: Local copy of vertical velocity field
                0072 C     KappaRk    :: Vertical diffusion for Tempertature
                0073 C     myTime     :: current time
                0074 C     myIter     :: current iteration number
                0075 C     myThid     :: my Thread Id. number
fb8688a588 Jean*0076       INTEGER bi, bj
0b191c5f5a Jean*0077       _RS recip_hFac(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0078       _RL uFld      (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0079       _RL vFld      (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0080       _RL wFld      (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0081       _RL KappaRk   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
4794324d0b Jean*0082       _RL     myTime
                0083       INTEGER myIter
                0084       INTEGER myThid
                0085 CEOP
                0086
                0087 #ifdef ALLOW_GENERIC_ADVDIFF
31067d68ad Jean*0088 #ifdef ALLOW_DIAGNOSTICS
                0089 C     !FUNCTIONS:
                0090       LOGICAL  DIAGNOSTICS_IS_ON
                0091       EXTERNAL DIAGNOSTICS_IS_ON
                0092 #endif /* ALLOW_DIAGNOSTICS */
                0093
4794324d0b Jean*0094 C     !LOCAL VARIABLES:
0b191c5f5a Jean*0095 C     iMin, iMax :: 1rst index loop range
                0096 C     jMin, jMax :: 2nd  index loop range
                0097 C     k          :: vertical index
                0098 C     kM1        :: =k-1 for k>1, =1 for k=1
                0099 C     kUp        :: index into 2 1/2D array, toggles between 1|2
                0100 C     kDown      :: index into 2 1/2D array, toggles between 2|1
                0101 C     xA         :: Tracer cell face area normal to X
                0102 C     yA         :: Tracer cell face area normal to X
                0103 C     maskUp     :: Land/water mask for Wvel points (interface k)
                0104 C     uTrans     :: Zonal volume transport through cell face
                0105 C     vTrans     :: Meridional volume transport through cell face
                0106 C     rTrans     ::   Vertical volume transport at interface k
                0107 C     rTransKp   :: Vertical volume transport at inteface k+1
                0108 C     fZon       :: Flux of temperature (T) in the zonal direction
                0109 C     fMer       :: Flux of temperature (T) in the meridional direction
                0110 C     fVer       :: Flux of temperature (T) in the vertical direction
                0111 C                   at the upper(U) and lower(D) faces of a cell.
1ece1facde Jean*0112 C     gT_loc     :: Temperature tendency (local to this S/R)
31067d68ad Jean*0113 C     gtForc     :: Temperature forcing tendency
e874fa47e5 Jean*0114 C     gt_AB      :: Adams-Bashforth temperature tendency increment
0b191c5f5a Jean*0115 C   useVariableK :: T when vertical diffusion is not constant
fb8688a588 Jean*0116       INTEGER iMin, iMax, jMin, jMax
4794324d0b Jean*0117       INTEGER i, j, k
                0118       INTEGER kUp, kDown, kM1
                0119       _RS xA      (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0120       _RS yA      (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0121       _RS maskUp  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0122       _RL uTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0123       _RL vTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0124       _RL rTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0125       _RL rTransKp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0b191c5f5a Jean*0126       _RL fZon    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0127       _RL fMer    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0128       _RL fVer    (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
1ece1facde Jean*0129       _RL gT_loc  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
31067d68ad Jean*0130       _RL gtForc  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
4794324d0b Jean*0131       _RL gt_AB   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
31067d68ad Jean*0132 #ifdef ALLOW_DIAGNOSTICS
                0133       LOGICAL diagForcing, diagAB_tend
                0134 #endif
4794324d0b Jean*0135       LOGICAL calcAdvection
                0136       INTEGER iterNb
                0137 #ifdef ALLOW_ADAMSBASHFORTH_3
972c0130ec Jean*0138       INTEGER m2
4794324d0b Jean*0139 #endif
7c50f07931 Mart*0140 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0141 C     tkey :: tape key (tile dependent)
                0142 C     kkey :: tape key (level and tile dependent)
                0143       INTEGER tkey, kkey
7c50f07931 Mart*0144 #endif
0b191c5f5a Jean*0145
fb8688a588 Jean*0146 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0147
                0148       iterNb = myIter
0b191c5f5a Jean*0149       IF (staggerTimeStep) iterNb = myIter - 1
4794324d0b Jean*0150
d3548146a8 Jean*0151 C-    Loop ranges for daughter routines:
                0152 c     iMin = 1
                0153 c     iMax = sNx
                0154 c     jMin = 1
                0155 c     jMax = sNy
                0156 C     Regarding model dynamics, only needs to get correct tracer tendency
                0157 C     (gT_loc) in tile interior (1:sNx,1:sNy);
                0158 C     However, for some diagnostics, we may want to get valid tendency
                0159 C     extended over 1 point in tile halo region (0:sNx+1,0:sNy=1).
                0160       iMin = 0
                0161       iMax = sNx+1
                0162       jMin = 0
                0163       jMax = sNy+1
                0164
31067d68ad Jean*0165 #ifdef ALLOW_DIAGNOSTICS
                0166       diagForcing = .FALSE.
                0167       diagAB_tend = .FALSE.
                0168       IF ( useDiagnostics .AND. tempForcing )
                0169      &     diagForcing = DIAGNOSTICS_IS_ON( 'gT_Forc ', myThid )
                0170       IF ( useDiagnostics .AND. AdamsBashforthGt )
                0171      &     diagAB_tend = DIAGNOSTICS_IS_ON( 'AB_gT   ', myThid )
                0172 #endif
                0173
4794324d0b Jean*0174 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0175       tkey = bi + (bj-1)*nSx + (ikey_dynamics-1)*nSx*nSy
4794324d0b Jean*0176 #endif /* ALLOW_AUTODIFF_TAMC */
                0177
ff1bd1eafd Jean*0178 C-    Apply AB on T :
                0179       IF ( AdamsBashforth_T ) THEN
                0180 C     compute T^n+1/2 (stored in gtNm) extrapolating T forward in time
                0181 #ifdef ALLOW_ADAMSBASHFORTH_3
                0182 c         m1 = 1 + MOD(iterNb+1,2)
                0183 c         m2 = 1 + MOD( iterNb ,2)
                0184           CALL ADAMS_BASHFORTH3(
                0185      I                           bi, bj, 0, Nr,
                0186      I                           theta(1-OLx,1-OLy,1,bi,bj),
                0187      U                           gtNm, gt_AB,
                0188      I                           tempStartAB, iterNb, myThid )
                0189 #else /* ALLOW_ADAMSBASHFORTH_3 */
                0190           CALL ADAMS_BASHFORTH2(
                0191      I                           bi, bj, 0, Nr,
                0192      I                           theta(1-OLx,1-OLy,1,bi,bj),
                0193      U                           gtNm1(1-OLx,1-OLy,1,bi,bj), gt_AB,
                0194      I                           tempStartAB, iterNb, myThid )
                0195 #endif /* ALLOW_ADAMSBASHFORTH_3 */
                0196       ENDIF
                0197
fb8688a588 Jean*0198 C-    Tracer tendency needs to be set to zero (moved here from gad_calc_rhs):
                0199       DO k=1,Nr
                0200        DO j=1-OLy,sNy+OLy
                0201         DO i=1-OLx,sNx+OLx
1ece1facde Jean*0202          gT_loc(i,j,k) = 0. _d 0
fb8688a588 Jean*0203         ENDDO
                0204        ENDDO
                0205       ENDDO
4794324d0b Jean*0206       DO j=1-OLy,sNy+OLy
                0207        DO i=1-OLx,sNx+OLx
0b191c5f5a Jean*0208          fVer(i,j,1) = 0. _d 0
                0209          fVer(i,j,2) = 0. _d 0
4794324d0b Jean*0210        ENDDO
                0211       ENDDO
fb8688a588 Jean*0212 #ifdef ALLOW_AUTODIFF
                0213       DO k=1,Nr
                0214        DO j=1-OLy,sNy+OLy
                0215         DO i=1-OLx,sNx+OLx
                0216          kappaRk(i,j,k) = 0. _d 0
                0217         ENDDO
                0218        ENDDO
                0219       ENDDO
4e4ad91a39 Jean*0220 #endif /* ALLOW_AUTODIFF */
                0221 #ifdef ALLOW_AUTODIFF_TAMC
598aebfcee Mart*0222 CADJ STORE wFld(:,:,:)         = comlev1_bibj, key=tkey, kind=isbyte
                0223 CADJ STORE theta(:,:,:,bi,bj)  = comlev1_bibj, key=tkey, kind=isbyte
ff1bd1eafd Jean*0224 # ifdef ALLOW_ADAMSBASHFORTH_3
598aebfcee Mart*0225 CADJ STORE gtNm(:,:,:,bi,bj,1) = comlev1_bibj, key=tkey, kind=isbyte
                0226 CADJ STORE gtNm(:,:,:,bi,bj,2) = comlev1_bibj, key=tkey, kind=isbyte
ff1bd1eafd Jean*0227 # else
598aebfcee Mart*0228 CADJ STORE gtNm1(:,:,:,bi,bj)  = comlev1_bibj, key=tkey, kind=isbyte
ff1bd1eafd Jean*0229 # endif
4e4ad91a39 Jean*0230 #endif /* ALLOW_AUTODIFF_TAMC */
4794324d0b Jean*0231
fb8688a588 Jean*0232 #ifdef INCLUDE_CALC_DIFFUSIVITY_CALL
0b191c5f5a Jean*0233       CALL CALC_3D_DIFFUSIVITY(
fb8688a588 Jean*0234      I         bi, bj, iMin, iMax, jMin, jMax,
                0235      I         GAD_TEMPERATURE, useGMredi, useKPP,
                0236      O         kappaRk,
                0237      I         myThid )
d9efc637ba Mart*0238 # ifdef ALLOW_AUTODIFF_TAMC
                0239 CADJ STORE kappaRk = comlev1_bibj, key = tkey, kind = isbyte
                0240 # endif
fb8688a588 Jean*0241 #endif /* INCLUDE_CALC_DIFFUSIVITY_CALL */
                0242
                0243 #ifndef DISABLE_MULTIDIM_ADVECTION
                0244 C--     Some advection schemes are better calculated using a multi-dimensional
                0245 C       method in the absence of any other terms and, if used, is done here.
                0246 C
                0247 C The CPP flag DISABLE_MULTIDIM_ADVECTION is currently unset in GAD_OPTIONS.h
                0248 C The default is to use multi-dimensinal advection for non-linear advection
                0249 C schemes. However, for the sake of efficiency of the adjoint it is necessary
                0250 C to be able to exclude this scheme to avoid excessive storage and
                0251 C recomputation. It *is* differentiable, if you need it.
                0252 C Edit GAD_OPTIONS.h and #define DISABLE_MULTIDIM_ADVECTION to
                0253 C disable this section of code.
                0254 #ifdef GAD_ALLOW_TS_SOM_ADV
                0255 # ifdef ALLOW_AUTODIFF_TAMC
598aebfcee Mart*0256 C This looks awkward from an F77 point of view, but TAF does the right thing.
                0257 CADJ STORE som_T(:,:,:,bi,bj,:) = comlev1_bibj, key=tkey, kind=isbyte
fb8688a588 Jean*0258 # endif
                0259       IF ( tempSOM_Advection ) THEN
                0260 # ifdef ALLOW_DEBUG
                0261         IF (debugMode) CALL DEBUG_CALL('GAD_SOM_ADVECT',myThid)
                0262 # endif
                0263         CALL GAD_SOM_ADVECT(
e874fa47e5 Jean*0264      I             tempImplVertAdv,
                0265      I             tempAdvScheme, tempVertAdvScheme, GAD_TEMPERATURE,
                0266      I             dTtracerLev, uFld, vFld, wFld, theta,
fb8688a588 Jean*0267      U             som_T,
1ece1facde Jean*0268      O             gT_loc,
fb8688a588 Jean*0269      I             bi, bj, myTime, myIter, myThid )
                0270       ELSEIF (tempMultiDimAdvec) THEN
                0271 #else /* GAD_ALLOW_TS_SOM_ADV */
                0272       IF (tempMultiDimAdvec) THEN
                0273 #endif /* GAD_ALLOW_TS_SOM_ADV */
                0274 # ifdef ALLOW_DEBUG
                0275         IF (debugMode) CALL DEBUG_CALL('GAD_ADVECTION',myThid)
                0276 # endif
                0277         CALL GAD_ADVECTION(
e874fa47e5 Jean*0278      I             tempImplVertAdv,
                0279      I             tempAdvScheme, tempVertAdvScheme, GAD_TEMPERATURE,
                0280      I             dTtracerLev, uFld, vFld, wFld, theta,
1ece1facde Jean*0281      O             gT_loc,
fb8688a588 Jean*0282      I             bi, bj, myTime, myIter, myThid )
                0283       ENDIF
                0284 #endif /* DISABLE_MULTIDIM_ADVECTION */
                0285
                0286 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0287
                0288 C-    Start vertical index (k) loop (Nr:1)
                0289       calcAdvection = tempAdvection .AND. .NOT.tempMultiDimAdvec
4794324d0b Jean*0290       DO k=Nr,1,-1
                0291 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0292         kkey = k + (tkey-1)*Nr
4794324d0b Jean*0293 #endif /* ALLOW_AUTODIFF_TAMC */
                0294         kM1  = MAX(1,k-1)
                0295         kUp  = 1+MOD(k+1,2)
                0296         kDown= 1+MOD(k,2)
                0297
                0298 #ifdef ALLOW_AUTODIFF_TAMC
d9efc637ba Mart*0299 CADJ STORE fVer(:,:,:)         = comlev1_bibj_k, key=kkey, kind = isbyte
                0300 CADJ STORE gT_loc(:,:,k)       = comlev1_bibj_k, key=kkey, kind = isbyte
4794324d0b Jean*0301 # ifdef ALLOW_ADAMSBASHFORTH_3
d9efc637ba Mart*0302 CADJ STORE gtNm(:,:,k,bi,bj,1) = comlev1_bibj_k, key=kkey, kind = isbyte
                0303 CADJ STORE gtNm(:,:,k,bi,bj,2) = comlev1_bibj_k, key=kkey, kind = isbyte
4794324d0b Jean*0304 # else
d9efc637ba Mart*0305 CADJ STORE gtNm1(:,:,k,bi,bj)  = comlev1_bibj_k, key=kkey, kind = isbyte
4794324d0b Jean*0306 # endif
                0307 #endif /* ALLOW_AUTODIFF_TAMC */
                0308         CALL CALC_ADV_FLOW(
                0309      I                uFld, vFld, wFld,
                0310      U                rTrans,
                0311      O                uTrans, vTrans, rTransKp,
                0312      O                maskUp, xA, yA,
                0313      I                k, bi, bj, myThid )
                0314
31067d68ad Jean*0315 C--   Collect forcing term in local array gtForc:
                0316         DO j=1-OLy,sNy+OLy
                0317          DO i=1-OLx,sNx+OLx
                0318           gtForc(i,j) = 0. _d 0
                0319          ENDDO
                0320         ENDDO
                0321         IF ( tempForcing ) THEN
                0322           CALL APPLY_FORCING_T(
                0323      U                        gtForc,
                0324      I                        iMin,iMax,jMin,jMax, k, bi,bj,
                0325      I                        myTime, myIter, myThid )
                0326 #ifdef ALLOW_DIAGNOSTICS
                0327           IF ( diagForcing ) THEN
                0328             CALL DIAGNOSTICS_FILL(gtForc,'gT_Forc ',k,1,2,bi,bj,myThid)
                0329           ENDIF
                0330 #endif /* ALLOW_DIAGNOSTICS */
                0331         ENDIF
                0332
4794324d0b Jean*0333 #ifdef ALLOW_ADAMSBASHFORTH_3
972c0130ec Jean*0334 c       m1 = 1 + MOD(iterNb+1,2)
4794324d0b Jean*0335         m2 = 1 + MOD( iterNb ,2)
                0336         CALL GAD_CALC_RHS(
                0337      I           bi, bj, iMin,iMax,jMin,jMax, k, kM1, kUp, kDown,
                0338      I           xA, yA, maskUp, uFld(1-OLx,1-OLy,k),
                0339      I           vFld(1-OLx,1-OLy,k), wFld(1-OLx,1-OLy,k),
                0340      I           uTrans, vTrans, rTrans, rTransKp,
                0341      I           diffKhT, diffK4T, KappaRk(1-OLx,1-OLy,k), diffKr4T,
cb30c1718d Jean*0342      I           theta(1-OLx,1-OLy,1,bi,bj),
                0343      I           gtNm(1-OLx,1-OLy,1,bi,bj,m2), dTtracerLev,
4794324d0b Jean*0344      I           GAD_TEMPERATURE, tempAdvScheme, tempVertAdvScheme,
                0345      I           calcAdvection, tempImplVertAdv, AdamsBashforth_T,
46918f1b26 Jean*0346      I           tempVertDiff4, useGMRedi, useKPP, temp_stayPositive,
0b191c5f5a Jean*0347      O           fZon, fMer,
1ece1facde Jean*0348      U           fVer, gT_loc,
4794324d0b Jean*0349      I           myTime, myIter, myThid )
                0350 #else /* ALLOW_ADAMSBASHFORTH_3 */
                0351         CALL GAD_CALC_RHS(
                0352      I           bi, bj, iMin,iMax,jMin,jMax, k, kM1, kUp, kDown,
                0353      I           xA, yA, maskUp, uFld(1-OLx,1-OLy,k),
                0354      I           vFld(1-OLx,1-OLy,k), wFld(1-OLx,1-OLy,k),
                0355      I           uTrans, vTrans, rTrans, rTransKp,
                0356      I           diffKhT, diffK4T, KappaRk(1-OLx,1-OLy,k), diffKr4T,
cb30c1718d Jean*0357      I           theta(1-OLx,1-OLy,1,bi,bj),
                0358      I           gtNm1(1-OLx,1-OLy,1,bi,bj), dTtracerLev,
4794324d0b Jean*0359      I           GAD_TEMPERATURE, tempAdvScheme, tempVertAdvScheme,
                0360      I           calcAdvection, tempImplVertAdv, AdamsBashforth_T,
46918f1b26 Jean*0361      I           tempVertDiff4, useGMRedi, useKPP, temp_stayPositive,
0b191c5f5a Jean*0362      O           fZon, fMer,
1ece1facde Jean*0363      U           fVer, gT_loc,
4794324d0b Jean*0364      I           myTime, myIter, myThid )
46918f1b26 Jean*0365 #endif /* ALLOW_ADAMSBASHFORTH_3 */
4794324d0b Jean*0366
                0367 C--   External thermal forcing term(s) inside Adams-Bashforth:
31067d68ad Jean*0368         IF ( tempForcing .AND. tracForcingOutAB.NE.1 ) THEN
                0369           DO j=1-OLy,sNy+OLy
                0370            DO i=1-OLx,sNx+OLx
1ece1facde Jean*0371             gT_loc(i,j,k) = gT_loc(i,j,k) + gtForc(i,j)
31067d68ad Jean*0372            ENDDO
                0373           ENDDO
                0374         ENDIF
4794324d0b Jean*0375
                0376         IF ( AdamsBashforthGt ) THEN
                0377 #ifdef ALLOW_ADAMSBASHFORTH_3
                0378           CALL ADAMS_BASHFORTH3(
                0379      I                          bi, bj, k, Nr,
1ece1facde Jean*0380      U                          gT_loc, gtNm,
                0381      O                          gt_AB,
4794324d0b Jean*0382      I                          tempStartAB, iterNb, myThid )
                0383 #else
                0384           CALL ADAMS_BASHFORTH2(
                0385      I                          bi, bj, k, Nr,
1ece1facde Jean*0386      U                          gT_loc, gtNm1(1-OLx,1-OLy,1,bi,bj),
                0387      O                          gt_AB,
4794324d0b Jean*0388      I                          tempStartAB, iterNb, myThid )
                0389 #endif
                0390 #ifdef ALLOW_DIAGNOSTICS
31067d68ad Jean*0391           IF ( diagAB_tend ) THEN
4794324d0b Jean*0392             CALL DIAGNOSTICS_FILL(gt_AB,'AB_gT   ',k,1,2,bi,bj,myThid)
                0393           ENDIF
                0394 #endif /* ALLOW_DIAGNOSTICS */
                0395         ENDIF
                0396
                0397 C--   External thermal forcing term(s) outside Adams-Bashforth:
31067d68ad Jean*0398         IF ( tempForcing .AND. tracForcingOutAB.EQ.1 ) THEN
                0399           DO j=1-OLy,sNy+OLy
                0400            DO i=1-OLx,sNx+OLx
1ece1facde Jean*0401             gT_loc(i,j,k) = gT_loc(i,j,k) + gtForc(i,j)
31067d68ad Jean*0402            ENDDO
                0403           ENDDO
                0404         ENDIF
4794324d0b Jean*0405
                0406 #ifdef NONLIN_FRSURF
d9efc637ba Mart*0407 # ifdef ALLOW_AUTODIFF_TAMC
                0408 CADJ STORE gT_loc(:,:,k) = comlev1_bibj_k, key = kkey, kind = isbyte
                0409 # endif
4794324d0b Jean*0410         IF (nonlinFreeSurf.GT.0) THEN
                0411           CALL FREESURF_RESCALE_G(
                0412      I                            bi, bj, k,
1ece1facde Jean*0413      U                            gT_loc,
4794324d0b Jean*0414      I                            myThid )
                0415          IF ( AdamsBashforthGt ) THEN
                0416 #ifdef ALLOW_ADAMSBASHFORTH_3
                0417 # ifdef ALLOW_AUTODIFF_TAMC
d9efc637ba Mart*0418 CADJ STORE gtNm(:,:,k,bi,bj,1) = comlev1_bibj_k, key=kkey, kind=isbyte
                0419 CADJ STORE gtNm(:,:,k,bi,bj,2) = comlev1_bibj_k, key=kkey, kind=isbyte
4794324d0b Jean*0420 # endif
                0421           CALL FREESURF_RESCALE_G(
                0422      I                            bi, bj, k,
1ece1facde Jean*0423      U                            gtNm(1-OLx,1-OLy,1,bi,bj,1),
4794324d0b Jean*0424      I                            myThid )
                0425           CALL FREESURF_RESCALE_G(
                0426      I                            bi, bj, k,
1ece1facde Jean*0427      U                            gtNm(1-OLx,1-OLy,1,bi,bj,2),
4794324d0b Jean*0428      I                            myThid )
                0429 #else
d9efc637ba Mart*0430 # ifdef ALLOW_AUTODIFF_TAMC
                0431 CADJ STORE gtNm1(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, kind = isbyte
                0432 # endif
4794324d0b Jean*0433           CALL FREESURF_RESCALE_G(
                0434      I                            bi, bj, k,
1ece1facde Jean*0435      U                            gtNm1(1-OLx,1-OLy,1,bi,bj),
4794324d0b Jean*0436      I                            myThid )
                0437 #endif
                0438          ENDIF
                0439         ENDIF
                0440 #endif /* NONLIN_FRSURF */
                0441
                0442 C-    end of vertical index (k) loop (Nr:1)
                0443       ENDDO
                0444
0b191c5f5a Jean*0445 #ifdef ALLOW_DOWN_SLOPE
                0446       IF ( useDOWN_SLOPE ) THEN
d9efc637ba Mart*0447 # ifdef ALLOW_AUTODIFF_TAMC
                0448 CADJ STORE recip_hFac = comlev1_bibj, key = tkey, kind = isbyte
                0449 # endif
0b191c5f5a Jean*0450         IF ( usingPCoords ) THEN
                0451           CALL DWNSLP_APPLY(
                0452      I                  GAD_TEMPERATURE, bi, bj, kSurfC,
6fca64f237 Jean*0453      I                  theta(1-OLx,1-OLy,1,bi,bj),
1ece1facde Jean*0454      U                  gT_loc,
6fca64f237 Jean*0455      I                  recip_hFac, recip_rA, recip_drF,
                0456      I                  dTtracerLev, myTime, myIter, myThid )
0b191c5f5a Jean*0457         ELSE
                0458           CALL DWNSLP_APPLY(
                0459      I                  GAD_TEMPERATURE, bi, bj, kLowC,
6fca64f237 Jean*0460      I                  theta(1-OLx,1-OLy,1,bi,bj),
1ece1facde Jean*0461      U                  gT_loc,
6fca64f237 Jean*0462      I                  recip_hFac, recip_rA, recip_drF,
                0463      I                  dTtracerLev, myTime, myIter, myThid )
0b191c5f5a Jean*0464         ENDIF
                0465       ENDIF
                0466 #endif /* ALLOW_DOWN_SLOPE */
                0467
1ece1facde Jean*0468 C-    Integrate forward in time, storing in gT_loc:  gT <= T + dt*gT
6fca64f237 Jean*0469       CALL TIMESTEP_TRACER(
                0470      I                  bi, bj, dTtracerLev,
                0471      I                  theta(1-OLx,1-OLy,1,bi,bj),
1ece1facde Jean*0472      U                  gT_loc,
6fca64f237 Jean*0473      I                  myTime, myIter, myThid )
                0474
0b191c5f5a Jean*0475 C--   Implicit vertical advection & diffusion
                0476
d9efc637ba Mart*0477 #ifdef ALLOW_AUTODIFF_TAMC
598aebfcee Mart*0478 CADJ STORE gT_loc (:,:,:)     = comlev1_bibj, key = tkey, kind = isbyte
d9efc637ba Mart*0479 #endif /* ALLOW_AUTODIFF_TAMC */
0b191c5f5a Jean*0480 #ifdef INCLUDE_IMPLVERTADV_CODE
5685a5c914 Jean*0481       IF ( tempImplVertAdv .OR. implicitDiffusion ) THEN
                0482 C     to recover older (prior to 2016-10-05) results:
                0483 c     IF ( tempImplVertAdv ) THEN
0b191c5f5a Jean*0484 #ifdef ALLOW_AUTODIFF_TAMC
598aebfcee Mart*0485 CADJ STORE wFld(:,:,:)        = comlev1_bibj, key = tkey, kind = isbyte
                0486 CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key = tkey, kind = isbyte
                0487 CADJ STORE recip_hFac(:,:,:)  = comlev1_bibj, key = tkey, kind = isbyte
0b191c5f5a Jean*0488 #endif /* ALLOW_AUTODIFF_TAMC */
                0489         CALL GAD_IMPLICIT_R(
                0490      I         tempImplVertAdv, tempVertAdvScheme, GAD_TEMPERATURE,
1ece1facde Jean*0491      I         dTtracerLev, kappaRk, recip_hFac, wFld,
                0492      I         theta(1-OLx,1-OLy,1,bi,bj),
                0493      U         gT_loc,
0b191c5f5a Jean*0494      I         bi, bj, myTime, myIter, myThid )
                0495       ELSEIF ( implicitDiffusion ) THEN
                0496 #else /* INCLUDE_IMPLVERTADV_CODE */
                0497       IF     ( implicitDiffusion ) THEN
                0498 #endif /* INCLUDE_IMPLVERTADV_CODE */
                0499         CALL IMPLDIFF(
                0500      I         bi, bj, iMin, iMax, jMin, jMax,
                0501      I         GAD_TEMPERATURE, kappaRk, recip_hFac,
1ece1facde Jean*0502      U         gT_loc,
0b191c5f5a Jean*0503      I         myThid )
                0504       ENDIF
                0505
bd27360393 Jean*0506       IF ( AdamsBashforth_T ) THEN
                0507 C-    Save current tracer field (for AB on tracer) and then update tracer
ff1bd1eafd Jean*0508 #ifdef ALLOW_ADAMSBASHFORTH_3
bd27360393 Jean*0509         CALL CYCLE_AB_TRACER(
cb30c1718d Jean*0510      I             bi, bj, gT_loc,
                0511      U             theta(1-OLx,1-OLy,1,bi,bj),
                0512      O             gtNm(1-OLx,1-OLy,1,bi,bj,m2),
bd27360393 Jean*0513      I             myTime, myIter, myThid )
                0514 #else /* ALLOW_ADAMSBASHFORTH_3 */
ff1bd1eafd Jean*0515         CALL CYCLE_AB_TRACER(
                0516      I             bi, bj, gT_loc,
                0517      U             theta(1-OLx,1-OLy,1,bi,bj),
                0518      O             gtNm1(1-OLx,1-OLy,1,bi,bj),
                0519      I             myTime, myIter, myThid )
bd27360393 Jean*0520 #endif /* ALLOW_ADAMSBASHFORTH_3 */
ff1bd1eafd Jean*0521       ELSE
bd27360393 Jean*0522 C-    Update tracer fields:  T(n) = T**
                0523         CALL CYCLE_TRACER(
                0524      I             bi, bj,
cb30c1718d Jean*0525      O             theta(1-OLx,1-OLy,1,bi,bj),
                0526      I             gT_loc, myTime, myIter, myThid )
bd27360393 Jean*0527       ENDIF
                0528
4794324d0b Jean*0529 #endif /* ALLOW_GENERIC_ADVDIFF */
                0530
                0531       RETURN
                0532       END