Back to home page

MITgcm

 
 

    


File indexing completed on 2023-03-10 06:09:59 UTC

view on githubraw file Latest commit 884e6cdd on 2023-03-09 18:21:20 UTC
785a077159 Alis*0001 #include "PTRACERS_OPTIONS.h"
d12f0d0164 Jean*0002 #ifdef ALLOW_AUTODIFF
                0003 # include "AUTODIFF_OPTIONS.h"
                0004 #endif
785a077159 Alis*0005 
                0006 CBOP
0dcd950fdd Patr*0007 C !ROUTINE: PTRACERS_INTEGRATE
785a077159 Alis*0008 
                0009 C !INTERFACE: ==========================================================
0dcd950fdd Patr*0010       SUBROUTINE PTRACERS_INTEGRATE(
9428e03c42 Jean*0011      I                    bi, bj, recip_hFac,
e4acd7d7df Jean*0012      I                    uFld, vFld, wFld,
                0013      U                    KappaRk,
                0014      I                    myTime, myIter, myThid )
785a077159 Alis*0015 
                0016 C !DESCRIPTION:
5c40528187 Jean*0017 C     Calculates tendency for passive tracers and integrates forward in
                0018 C     time. The tracer array is updated here while adjustments (filters,
                0019 C     conv.adjustment) are applied later, in S/R TRACERS_CORRECTION_STEP
                0020 
785a077159 Alis*0021 C !USES: ===============================================================
5101983465 Jean*0022 #include "PTRACERS_MOD.h"
785a077159 Alis*0023       IMPLICIT NONE
                0024 #include "SIZE.h"
                0025 #include "EEPARAMS.h"
5d45aaff6a Alis*0026 #include "PARAMS.h"
9428e03c42 Jean*0027 #include "GRID.h"
4e66ab0b67 Oliv*0028 #ifdef ALLOW_LONGSTEP
                0029 #include "LONGSTEP_PARAMS.h"
                0030 #endif
636477d15b Jean*0031 #include "PTRACERS_SIZE.h"
0a278985fd Jean*0032 #include "PTRACERS_PARAMS.h"
a2791aaebc Jean*0033 #include "PTRACERS_START.h"
0a278985fd Jean*0034 #include "PTRACERS_FIELDS.h"
785a077159 Alis*0035 #include "GAD.h"
9829fd7e0c Patr*0036 #ifdef ALLOW_AUTODIFF_TAMC
                0037 # include "tamc.h"
                0038 #endif
785a077159 Alis*0039 
                0040 C !INPUT PARAMETERS: ===================================================
9428e03c42 Jean*0041 C  bi, bj           :: tile indices
                0042 C  recip_hFac       :: reciprocal of cell open-depth factor (@ next iter)
                0043 C  uFld, vFld, wFld :: Local copy of velocity field (3 components)
                0044 C  KappaRk          :: vertical diffusion used for one passive tracer
                0045 C  myTime           :: model time
                0046 C  myIter           :: time-step number
                0047 C  myThid           :: thread number
e4acd7d7df Jean*0048       INTEGER bi, bj
9428e03c42 Jean*0049       _RS recip_hFac(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0050       _RL uFld      (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0051       _RL vFld      (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0052       _RL wFld      (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0053       _RL KappaRk   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
785a077159 Alis*0054       _RL myTime
0845f1a203 Jean*0055       INTEGER myIter
785a077159 Alis*0056       INTEGER myThid
                0057 
                0058 C !OUTPUT PARAMETERS: ==================================================
                0059 C  none
                0060 
                0061 #ifdef ALLOW_PTRACERS
ecf3fdef88 Jean*0062 #ifdef ALLOW_DIAGNOSTICS
                0063 C     !FUNCTIONS:
                0064       LOGICAL  DIAGNOSTICS_IS_ON
                0065       EXTERNAL DIAGNOSTICS_IS_ON
                0066       CHARACTER*4 GAD_DIAG_SUFX
                0067       EXTERNAL    GAD_DIAG_SUFX
                0068 #endif /* ALLOW_DIAGNOSTICS */
785a077159 Alis*0069 
                0070 C !LOCAL VARIABLES: ====================================================
9428e03c42 Jean*0071 C  iTracer          :: tracer index
                0072 C  iMin, iMax       :: 1rst index loop range
                0073 C  jMin, jMax       :: 2nd  index loop range
                0074 C  k                :: vertical level number
                0075 C  kUp,kDown        :: toggle indices for even/odd level fluxes
                0076 C  kM1              :: =min(1,k-1)
                0077 C  GAD_TR           :: passive tracer id (GAD_TR1+iTracer-1)
                0078 C  xA               :: face area at U points in level k
                0079 C  yA               :: face area at V points in level k
                0080 C  maskUp           :: mask for vertical transport
                0081 C  uTrans           :: zonal transport in level k
                0082 C  vTrans           :: meridional transport in level k
                0083 C  rTrans           :: vertical volume transport at interface k
                0084 C  rTransKp         :: vertical volume transport at interface k+1
                0085 C  fZon             :: passive tracer zonal flux
                0086 C  fMer             :: passive tracer meridional flux
                0087 C  fVer             :: passive tracer vertical flux
1ece1facde Jean*0088 C  gTracer          :: passive tracer tendency
ecf3fdef88 Jean*0089 C  gTrForc          :: passive tracer forcing tendency
e07d0b88b1 Jean*0090 C  gTr_AB           :: Adams-Bashforth tracer tendency increment
b1c7b8e45b Jean*0091 C  gTr_trp          :: passive tracer total transport tendency before gchem
ef794e16d2 Jean*0092       INTEGER iTracer
785a077159 Alis*0093       INTEGER iMin,iMax,jMin,jMax
e4acd7d7df Jean*0094       INTEGER i, j, k
                0095       INTEGER kUp, kDown, kM1
981f32514d Dimi*0096       INTEGER GAD_TR
e4acd7d7df Jean*0097       _RS xA      (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0098       _RS yA      (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0099       _RS maskUp  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0100       _RL uTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0101       _RL vTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0102       _RL rTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0103       _RL rTransKp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
9428e03c42 Jean*0104       _RL fZon    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0105       _RL fMer    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0106       _RL fVer    (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
1ece1facde Jean*0107       _RL gTracer (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
ecf3fdef88 Jean*0108       _RL gTrForc (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
e4acd7d7df Jean*0109       _RL gTr_AB  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
f8b016827f Jean*0110       LOGICAL calcAdvection
d31276c95b Jean*0111       INTEGER iterNb
e14953c89a Jean*0112       _RL dummy(Nr)
7c50f07931 Mart*0113 #ifdef ALLOW_AUTODIFF_TAMC
d9efc637ba Mart*0114 C     tkey :: tape key (depends on tiles)
                0115 C     ikey :: tape key (depends on tracer and tiles)
edb6656069 Mart*0116 C     kkey :: tape key (depends on levels, tracer, and tiles)
d9efc637ba Mart*0117       INTEGER  tkey, ikey, kkey
7c50f07931 Mart*0118 #endif
2dd0aa05a0 Jean*0119 #ifdef ALLOW_DIAGNOSTICS
b1c7b8e45b Jean*0120       _RL gTr_trp (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0121       _RL recip_dt
2dd0aa05a0 Jean*0122       CHARACTER*8 diagName
                0123       CHARACTER*4 diagSufx
1cbab8a41c Oliv*0124       LOGICAL diagForcing, diagAB_tend, diagTrp_tend
2dd0aa05a0 Jean*0125 #endif /* ALLOW_DIAGNOSTICS */
785a077159 Alis*0126 CEOP
                0127 
5101983465 Jean*0128 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0129 
e4acd7d7df Jean*0130 C-    Compute iter at beginning of ptracer time step
                0131 #ifdef ALLOW_LONGSTEP
                0132       iterNb = myIter - LS_nIter + 1
                0133       IF (LS_whenToSample.GE.2) iterNb = myIter - LS_nIter
                0134 #else
                0135       iterNb = myIter
                0136       IF (staggerTimeStep) iterNb = myIter - 1
                0137 #endif
9829fd7e0c Patr*0138 
d3548146a8 Jean*0139 C-    Loop ranges for daughter routines
                0140 c     iMin = 1
                0141 c     iMax = sNx
                0142 c     jMin = 1
                0143 c     jMax = sNy
                0144 C     Regarding model dynamics, only needs to get correct tracer tendency
                0145 C     (gTracer) in tile interior (1:sNx,1:sNy);
                0146 C     However, for some diagnostics, we may want to get valid tendency
                0147 C     extended over 1 point in tile halo region (0:sNx+1,0:sNy=1).
                0148       iMin = 0
                0149       iMax = sNx+1
                0150       jMin = 0
                0151       jMax = sNy+1
                0152 
e4acd7d7df Jean*0153 C--   Loop over tracers
785a077159 Alis*0154       DO iTracer=1,PTRACERS_numInUse
1ece1facde Jean*0155        IF ( PTRACERS_StepFwd(iTracer) ) THEN
                0156         GAD_TR = GAD_TR1 + iTracer - 1
5101983465 Jean*0157 
                0158 C-    Initialise tracer tendency to zero
1ece1facde Jean*0159         DO k=1,Nr
                0160          DO j=1-OLy,sNy+OLy
                0161           DO i=1-OLx,sNx+OLx
                0162            gTracer(i,j,k) = 0. _d 0
                0163           ENDDO
5101983465 Jean*0164          ENDDO
                0165         ENDDO
785a077159 Alis*0166 
ecf3fdef88 Jean*0167 #ifdef ALLOW_DIAGNOSTICS
                0168         diagForcing = .FALSE.
                0169         diagAB_tend = .FALSE.
1cbab8a41c Oliv*0170         diagTrp_tend = .FALSE.
ecf3fdef88 Jean*0171         IF ( useDiagnostics ) THEN
                0172           diagSufx = GAD_DIAG_SUFX( GAD_TR, myThid )
                0173           diagName = 'Forc'//diagSufx
                0174           diagForcing = DIAGNOSTICS_IS_ON( diagName, myThid )
                0175           diagName = 'AB_g'//diagSufx
                0176           IF ( PTRACERS_AdamsBashGtr(iTracer) )
                0177      &    diagAB_tend = DIAGNOSTICS_IS_ON( diagName, myThid )
1cbab8a41c Oliv*0178           diagName = 'Tp_g'//diagSufx
                0179           diagTrp_tend = DIAGNOSTICS_IS_ON( diagName, myThid )
ecf3fdef88 Jean*0180         ENDIF
                0181 #endif
                0182 
9829fd7e0c Patr*0183 #ifdef ALLOW_AUTODIFF_TAMC
884e6cddfb mjlo*0184         tkey = bi + (bj-1)*nSx + (ikey_dynamics-1)*nSx*nSy
d9efc637ba Mart*0185         ikey = iTracer + (tkey-1)*PTRACERS_num
e4acd7d7df Jean*0186 #endif /* ALLOW_AUTODIFF_TAMC */
                0187 
fc10d43a89 Jean*0188 C-    Apply AB on Tracer :
                0189         IF ( PTRACERS_AdamsBash_Tr(iTracer) ) THEN
                0190 C     compute pTr^n+1/2 (stored in gpTrNm1) extrapolating pTr forward in time
                0191           CALL ADAMS_BASHFORTH2(
                0192      I                      bi, bj, 0, Nr,
                0193      I                      pTracer(1-OLx,1-OLy,1,bi,bj,iTracer),
                0194      U                      gpTrNm1(1-OLx,1-OLy,1,bi,bj,iTracer),
                0195      O                      gTr_AB,
                0196      I                      PTRACERS_startAB(iTracer), iterNb, myThid )
                0197         ENDIF
                0198 
e4acd7d7df Jean*0199         DO j=1-OLy,sNy+OLy
                0200          DO i=1-OLx,sNx+OLx
9428e03c42 Jean*0201            fVer(i,j,1) = 0. _d 0
                0202            fVer(i,j,2) = 0. _d 0
e4acd7d7df Jean*0203          ENDDO
                0204         ENDDO
d8354445a3 Jean*0205 #ifdef ALLOW_AUTODIFF
e4acd7d7df Jean*0206         DO k=1,Nr
                0207          DO j=1-OLy,sNy+OLy
                0208           DO i=1-OLx,sNx+OLx
                0209            kappaRk(i,j,k) = 0. _d 0
                0210           ENDDO
                0211          ENDDO
                0212         ENDDO
d8354445a3 Jean*0213 #endif /* ALLOW_AUTODIFF */
e4acd7d7df Jean*0214 
d9efc637ba Mart*0215 C     salt_integrate.F uses INCLUDE_CALC_DIFFUSIVITY_CALL
                0216 c#ifdef INCLUDE_CALC_DIFFUSIVITY_CALL
9428e03c42 Jean*0217         CALL CALC_3D_DIFFUSIVITY(
e4acd7d7df Jean*0218      I         bi, bj, iMin,iMax,jMin,jMax,
                0219      I         GAD_TR,
                0220      I         PTRACERS_useGMRedi(iTracer), PTRACERS_useKPP(iTracer),
                0221      O         kappaRk,
                0222      I         myThid)
d9efc637ba Mart*0223 # ifdef ALLOW_AUTODIFF_TAMC
                0224 CADJ STORE kappaRk = comlev1_bibj_ptracers, key = ikey, byte = isbyte
                0225 # endif
                0226 c#endif /* INCLUDE_CALC_DIFFUSIVITY_CALL */
e4acd7d7df Jean*0227 
5101983465 Jean*0228 #ifndef DISABLE_MULTIDIM_ADVECTION
                0229 #ifdef ALLOW_AUTODIFF_TAMC
                0230 CADJ STORE pTracer(:,:,:,bi,bj,iTracer)
d9efc637ba Mart*0231 CADJ &      = comlev1_bibj_ptracers, key=ikey, byte=isbyte
5101983465 Jean*0232 #endif /* ALLOW_AUTODIFF_TAMC */
                0233 
                0234 #ifdef PTRACERS_ALLOW_DYN_STATE
                0235         IF ( PTRACERS_SOM_Advection(iTracer) ) THEN
                0236 # ifdef ALLOW_DEBUG
                0237          IF (debugMode) CALL DEBUG_CALL('GAD_SOM_ADVECT',myThid)
                0238 # endif
                0239           CALL GAD_SOM_ADVECT(
1ece1facde Jean*0240      I                        PTRACERS_ImplVertAdv(iTracer),
                0241      I                        PTRACERS_advScheme(iTracer),
                0242      I                        PTRACERS_advScheme(iTracer),
                0243      I                        GAD_TR,
                0244      I                        PTRACERS_dTLev, uFld, vFld, wFld,
                0245      I                        pTracer(1-OLx,1-OLy,1,1,1,iTracer),
                0246      U                        _Ptracers_som(:,:,:,:,:,:,iTracer),
                0247      O                        gTracer,
                0248      I                        bi, bj, myTime, myIter, myThid )
5101983465 Jean*0249         ELSEIF ( PTRACERS_MultiDimAdv(iTracer) ) THEN
                0250 #else /* PTRACERS_ALLOW_DYN_STATE */
                0251         IF ( PTRACERS_MultiDimAdv(iTracer) ) THEN
                0252 #endif /* PTRACERS_ALLOW_DYN_STATE */
                0253 # ifdef ALLOW_DEBUG
                0254           IF (debugMode) CALL DEBUG_CALL('GAD_ADVECTION',myThid)
                0255 # endif
                0256           CALL GAD_ADVECTION(
                0257      I                        PTRACERS_ImplVertAdv(iTracer),
                0258      I                        PTRACERS_advScheme(iTracer),
                0259      I                        PTRACERS_advScheme(iTracer),
e07d0b88b1 Jean*0260      I                        GAD_TR,
                0261      I                        PTRACERS_dTLev, uFld, vFld, wFld,
5101983465 Jean*0262      I                        pTracer(1-OLx,1-OLy,1,1,1,iTracer),
1ece1facde Jean*0263      O                        gTracer,
5101983465 Jean*0264      I                        bi, bj, myTime, myIter, myThid )
                0265         ENDIF
                0266 #endif /* DISABLE_MULTIDIM_ADVECTION */
                0267 
                0268 C-    Start vertical index (k) loop (Nr:1)
                0269         calcAdvection = .NOT.PTRACERS_MultiDimAdv(iTracer)
bd50fe7d15 Gael*0270      &                  .AND. (PTRACERS_advScheme(iTracer).NE.0)
e4acd7d7df Jean*0271         DO k=Nr,1,-1
                0272 #ifdef ALLOW_AUTODIFF_TAMC
d9efc637ba Mart*0273           kkey = k + (ikey-1)*Nr
9829fd7e0c Patr*0274 #endif /* ALLOW_AUTODIFF_TAMC */
785a077159 Alis*0275 
e4acd7d7df Jean*0276           kM1  = MAX(1,k-1)
                0277           kUp  = 1+MOD(k+1,2)
                0278           kDown= 1+MOD(k,2)
                0279 
9829fd7e0c Patr*0280 #ifdef ALLOW_AUTODIFF_TAMC
d9efc637ba Mart*0281 CADJ STORE fVer(:,:,:)                  = comlev1_bibj_k_ptracers,
                0282 CADJ &     key = kkey, kind = isbyte
                0283 CADJ STORE gTracer(:,:,k)               = comlev1_bibj_k_ptracers,
                0284 CADJ &     key = kkey, kind = isbyte
e4acd7d7df Jean*0285 CADJ STORE gpTrNm1(:,:,k,bi,bj,iTracer) = comlev1_bibj_k_ptracers,
d9efc637ba Mart*0286 CADJ &     key = kkey, kind = isbyte
9829fd7e0c Patr*0287 #endif /* ALLOW_AUTODIFF_TAMC */
e4acd7d7df Jean*0288           CALL CALC_ADV_FLOW(
                0289      I                uFld, vFld, wFld,
                0290      U                rTrans,
                0291      O                uTrans, vTrans, rTransKp,
                0292      O                maskUp, xA, yA,
                0293      I                k, bi, bj, myThid )
0845f1a203 Jean*0294 
ecf3fdef88 Jean*0295 C--   Collect forcing term in local array gTrForc:
                0296           DO j=1-OLy,sNy+OLy
                0297            DO i=1-OLx,sNx+OLx
                0298             gTrForc(i,j) = 0. _d 0
                0299            ENDDO
                0300           ENDDO
                0301           CALL PTRACERS_APPLY_FORCING(
                0302      U                    gTrForc,
                0303      I                    surfaceForcingPTr(1-OLx,1-OLy,bi,bj,iTracer),
                0304      I                    iMin,iMax,jMin,jMax, k, bi, bj,
                0305      I                    iTracer, myTime, myIter, myThid )
                0306 #ifdef ALLOW_DIAGNOSTICS
                0307           IF ( diagForcing ) THEN
                0308             diagName = 'Forc'//diagSufx
                0309             CALL DIAGNOSTICS_FILL(gTrForc,diagName,k,1,2,bi,bj,myThid)
                0310           ENDIF
                0311 #endif /* ALLOW_DIAGNOSTICS */
                0312 
1ece1facde Jean*0313 C-    Calculate active tracer tendencies (gTracer) due to internal processes
e4acd7d7df Jean*0314 C      (advection, [explicit] diffusion, parameterizations,...)
                0315           CALL GAD_CALC_RHS(
                0316      I             bi,bj, iMin,iMax,jMin,jMax, k,kM1, kUp,kDown,
                0317      I             xA, yA, maskUp, uFld(1-OLx,1-OLy,k),
                0318      I             vFld(1-OLx,1-OLy,k), wFld(1-OLx,1-OLy,k),
                0319      I             uTrans, vTrans, rTrans, rTransKp,
                0320      I             PTRACERS_diffKh(iTracer),
                0321      I             PTRACERS_diffK4(iTracer),
                0322      I             KappaRk(1-OLx,1-OLy,k), dummy,
cb30c1718d Jean*0323      I             pTracer(1-OLx,1-OLy,1,bi,bj,iTracer),
                0324      I             gpTrNm1(1-OLx,1-OLy,1,bi,bj,iTracer),
e4acd7d7df Jean*0325      I             PTRACERS_dTLev, GAD_TR,
                0326      I             PTRACERS_advScheme(iTracer),
                0327      I             PTRACERS_advScheme(iTracer),
                0328      I             calcAdvection, PTRACERS_ImplVertAdv(iTracer),
fc10d43a89 Jean*0329      I             PTRACERS_AdamsBash_Tr(iTracer), .FALSE.,
e4acd7d7df Jean*0330      I             PTRACERS_useGMRedi(iTracer),
                0331      I             PTRACERS_useKPP(iTracer),
46918f1b26 Jean*0332      I             PTRACERS_stayPositive(iTracer),
9428e03c42 Jean*0333      O             fZon, fMer,
1ece1facde Jean*0334      U             fVer, gTracer,
e4acd7d7df Jean*0335      I             myTime, myIter, myThid )
                0336 
ecf3fdef88 Jean*0337 C-    External forcing term(s) inside Adams-Bashforth:
                0338           IF ( tracForcingOutAB.NE.1 ) THEN
                0339             DO j=1-OLy,sNy+OLy
                0340              DO i=1-OLx,sNx+OLx
1ece1facde Jean*0341               gTracer(i,j,k) = gTracer(i,j,k) + gTrForc(i,j)
ecf3fdef88 Jean*0342              ENDDO
                0343             ENDDO
                0344           ENDIF
e4acd7d7df Jean*0345 
                0346 C-    If using Adams-Bashforth II, then extrapolate tendencies
1ece1facde Jean*0347 C     gTracer is now the tracer tendency for explicit advection/diffusion
e4acd7d7df Jean*0348 
                0349 C     If matrix is being computed, skip call to S/R ADAMS_BASHFORTH2 to
1ece1facde Jean*0350 C     prevent gTracer from being replaced by the average of gTracer and gpTrNm1.
e4acd7d7df Jean*0351           IF ( .NOT.useMATRIX .AND.
                0352      &         PTRACERS_AdamsBashGtr(iTracer) ) THEN
                0353            CALL ADAMS_BASHFORTH2(
94e3f14b29 Jean*0354      I                      bi, bj, k, Nr,
1ece1facde Jean*0355      U                      gTracer,
                0356      U                      gpTrNm1(1-OLx,1-OLy,1,bi,bj,iTracer),
                0357      O                      gTr_AB,
d31276c95b Jean*0358      I                      PTRACERS_startAB(iTracer), iterNb, myThid )
2dd0aa05a0 Jean*0359 #ifdef ALLOW_DIAGNOSTICS
ecf3fdef88 Jean*0360            IF ( diagAB_tend ) THEN
e4acd7d7df Jean*0361              diagName = 'AB_g'//diagSufx
                0362              CALL DIAGNOSTICS_FILL(gTr_AB,diagName,k,1,2,bi,bj,myThid)
                0363            ENDIF
2dd0aa05a0 Jean*0364 #endif /* ALLOW_DIAGNOSTICS */
e4acd7d7df Jean*0365           ENDIF
785a077159 Alis*0366 
ecf3fdef88 Jean*0367 C-    External forcing term(s) outside Adams-Bashforth:
                0368           IF ( tracForcingOutAB.EQ.1 ) THEN
                0369             DO j=1-OLy,sNy+OLy
                0370              DO i=1-OLx,sNx+OLx
1ece1facde Jean*0371               gTracer(i,j,k) = gTracer(i,j,k) + gTrForc(i,j)
ecf3fdef88 Jean*0372              ENDDO
                0373             ENDDO
                0374           ENDIF
f8b016827f Jean*0375 
785a077159 Alis*0376 #ifdef NONLIN_FRSURF
d9efc637ba Mart*0377 # ifdef ALLOW_AUTODIFF_TAMC
                0378 CADJ STORE gTracer(:,:,k) = comlev1_bibj_k, key = kkey, kind = isbyte
                0379 # endif
e4acd7d7df Jean*0380 C-    Account for change in level thickness
                0381           IF (nonlinFreeSurf.GT.0) THEN
                0382            CALL FREESURF_RESCALE_G(
1ece1facde Jean*0383      I                           bi, bj, k,
                0384      U                           gTracer,
                0385      I                           myThid )
d9efc637ba Mart*0386            IF ( PTRACERS_AdamsBashGtr(iTracer) ) THEN
                0387 # ifdef ALLOW_AUTODIFF_TAMC
                0388 CADJ STORE gpTrNm1(:,:,k,bi,bj,iTracer) = comlev1_bibj_k_ptracers,
                0389 CADJ &     key = kkey, kind = isbyte
                0390 # endif
                0391             CALL FREESURF_RESCALE_G(
1ece1facde Jean*0392      I                           bi, bj, k,
                0393      U                           gpTrNm1(1-OLx,1-OLy,1,bi,bj,iTracer),
                0394      I                           myThid )
d9efc637ba Mart*0395            ENDIF
e4acd7d7df Jean*0396           ENDIF
785a077159 Alis*0397 #endif /* NONLIN_FRSURF */
                0398 
e4acd7d7df Jean*0399 C-    end of vertical index (k) loop (Nr:1)
                0400         ENDDO
785a077159 Alis*0401 
9428e03c42 Jean*0402 #ifdef ALLOW_DOWN_SLOPE
                0403         IF ( PTRACERS_useDWNSLP(iTracer) ) THEN
d9efc637ba Mart*0404 # ifdef ALLOW_AUTODIFF_TAMC
                0405 C     Here we can store on comlev1_bibj to save memory, because
                0406 C     recip_hfac is the same for all tracers
                0407 CADJ STORE recip_hfac = comlev1_bibj, key = tkey, kind = isbyte
                0408 # endif
9428e03c42 Jean*0409           IF ( usingPCoords ) THEN
                0410             CALL DWNSLP_APPLY(
                0411      I                  GAD_TR, bi, bj, kSurfC,
6fca64f237 Jean*0412      I                  pTracer(1-OLx,1-OLy,1,bi,bj,iTracer),
1ece1facde Jean*0413      U                  gTracer,
6fca64f237 Jean*0414      I                  recip_hFac, recip_rA, recip_drF,
                0415      I                  PTRACERS_dTLev, myTime, myIter, myThid )
9428e03c42 Jean*0416           ELSE
                0417             CALL DWNSLP_APPLY(
                0418      I                  GAD_TR, bi, bj, kLowC,
6fca64f237 Jean*0419      I                  pTracer(1-OLx,1-OLy,1,bi,bj,iTracer),
1ece1facde Jean*0420      U                  gTracer,
6fca64f237 Jean*0421      I                  recip_hFac, recip_rA, recip_drF,
                0422      I                  PTRACERS_dTLev, myTime, myIter, myThid )
9428e03c42 Jean*0423           ENDIF
                0424         ENDIF
                0425 #endif /* ALLOW_DOWN_SLOPE */
                0426 
1ece1facde Jean*0427 C-    Integrate forward in time, storing in gTracer:  gTr <= pTr + dt*gTr
6fca64f237 Jean*0428         CALL TIMESTEP_TRACER(
                0429      I                bi, bj, PTRACERS_dTLev,
                0430      I                pTracer(1-OLx,1-OLy,1,bi,bj,iTracer),
1ece1facde Jean*0431      U                gTracer,
6fca64f237 Jean*0432      I                myTime, myIter, myThid )
                0433 
cb30c1718d Jean*0434 C     All explicit advection/diffusion/sources should now be
                0435 C     done. The updated tracer field is in gTracer.
9428e03c42 Jean*0436 #ifdef ALLOW_MATRIX
1ece1facde Jean*0437 C       Accumalate explicit tendency and also reset gTracer to initial
9428e03c42 Jean*0438 C       tracer field for implicit matrix calculation
3581a69c62 Jean*0439         IF ( useMATRIX ) THEN
                0440           CALL MATRIX_STORE_TENDENCY_EXP(
                0441      I                      iTracer, bi, bj,
1ece1facde Jean*0442      U                      gTracer,
3581a69c62 Jean*0443      I                      myTime, myIter, myThid )
                0444         ENDIF
9428e03c42 Jean*0445 #endif /* ALLOW_MATRIX */
                0446 
cb30c1718d Jean*0447 C--   Vertical advection & diffusion (implicit) for passive tracers
9428e03c42 Jean*0448 
                0449 #ifdef ALLOW_AUTODIFF_TAMC
d9efc637ba Mart*0450 CADJ STORE gTracer(:,:,:) = comlev1_bibj_ptracers, key=ikey, byte=isbyte
9428e03c42 Jean*0451 #endif /* ALLOW_AUTODIFF_TAMC */
                0452 
                0453 #ifdef INCLUDE_IMPLVERTADV_CODE
37ba176f24 Jean*0454         IF ( PTRACERS_ImplVertAdv(iTracer) .OR. implicitDiffusion ) THEN
                0455 C      to recover older (prior to 2016-10-05) results:
                0456 c       IF ( PTRACERS_ImplVertAdv(iTracer) ) THEN
9428e03c42 Jean*0457           CALL GAD_IMPLICIT_R(
                0458      I         PTRACERS_ImplVertAdv(iTracer),
                0459      I         PTRACERS_advScheme(iTracer), GAD_TR,
                0460      I         PTRACERS_dTLev, kappaRk, recip_hFac, wFld,
1ece1facde Jean*0461      I         pTracer(1-OLx,1-OLy,1,bi,bj,iTracer),
                0462      U         gTracer,
9428e03c42 Jean*0463      I         bi, bj, myTime, myIter, myThid )
                0464 
                0465         ELSEIF ( implicitDiffusion ) THEN
                0466 #else /* INCLUDE_IMPLVERTADV_CODE */
                0467         IF ( implicitDiffusion ) THEN
                0468 #endif /* INCLUDE_IMPLVERTADV_CODE */
                0469           CALL IMPLDIFF(
                0470      I         bi, bj, iMin, iMax, jMin, jMax,
                0471      I         GAD_TR, kappaRk, recip_hFac,
1ece1facde Jean*0472      U         gTracer,
9428e03c42 Jean*0473      I         myThid )
                0474         ENDIF
                0475 
1cbab8a41c Oliv*0476 #ifdef ALLOW_DIAGNOSTICS
                0477         IF ( diagTrp_tend ) THEN
                0478          diagName = 'Tp_g'//diagSufx
                0479          DO k=1,Nr
b1c7b8e45b Jean*0480           IF ( PTRACERS_dTLev(k).NE.zeroRL ) THEN
                0481             recip_dt = oneRL/PTRACERS_dTLev(k)
                0482           ELSE
                0483             recip_dt = 0. _d 0
                0484           ENDIF
1cbab8a41c Oliv*0485           DO j=1,sNy
                0486            DO i=1,sNx
b1c7b8e45b Jean*0487             gTr_trp(i,j) = ( gTracer(i,j,k)
                0488      &                     - pTracer(i,j,k,bi,bj,iTracer) )*recip_dt
1cbab8a41c Oliv*0489            ENDDO
                0490           ENDDO
b1c7b8e45b Jean*0491           CALL DIAGNOSTICS_FILL( gTr_trp,diagName,k,1,2,bi,bj,myThid )
1cbab8a41c Oliv*0492          ENDDO
                0493         ENDIF
                0494 #endif /* ALLOW_DIAGNOSTICS */
                0495 
fc10d43a89 Jean*0496         IF ( PTRACERS_AdamsBash_Tr(iTracer) ) THEN
                0497 C--   Save current tracer field (for AB on tracer) and then update tracer
                0498           CALL CYCLE_AB_TRACER(
                0499      I               bi, bj, gTracer,
                0500      U               pTracer(1-OLx,1-OLy,1,bi,bj,iTracer),
                0501      O               gpTrNm1(1-OLx,1-OLy,1,bi,bj,iTracer),
                0502      I               myTime, myIter, myThid )
                0503         ELSE
                0504 C--   Update tracer fields:  pTr(n) = pTr**
                0505           CALL CYCLE_TRACER(
                0506      I               bi, bj,
                0507      O               pTracer(1-OLx,1-OLy,1,bi,bj,iTracer),
                0508      I               gTracer, myTime, myIter, myThid )
                0509         ENDIF
5c40528187 Jean*0510 
9428e03c42 Jean*0511 #ifdef ALLOW_OBCS
cb30c1718d Jean*0512 C--   Apply open boundary conditions for each passive tracer
9428e03c42 Jean*0513         IF ( useOBCS ) THEN
                0514           CALL OBCS_APPLY_PTRACER(
                0515      I         bi, bj, 0, iTracer,
5c40528187 Jean*0516      U         pTracer(1-OLx,1-OLy,1,bi,bj,iTracer),
9428e03c42 Jean*0517      I         myThid )
                0518         ENDIF
                0519 #endif /* ALLOW_OBCS */
                0520 
e4acd7d7df Jean*0521 C--   end of tracer loop
                0522        ENDIF
785a077159 Alis*0523       ENDDO
                0524 
                0525 #endif /* ALLOW_PTRACERS */
                0526 
                0527       RETURN
                0528       END