/[MITgcm]/MITgcm/model/src/thermodynamics.F
ViewVC logotype

Diff of /MITgcm/model/src/thermodynamics.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.14 by jmc, Fri Nov 16 03:25:41 2001 UTC revision 1.44 by jmc, Fri Jul 11 16:00:08 2003 UTC
# Line 2  C $Header$ Line 2  C $Header$
2  C $Name$  C $Name$
3    
4  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
5    #ifdef ALLOW_AUTODIFF_TAMC
6    # ifdef ALLOW_GMREDI
7    #  include "GMREDI_OPTIONS.h"
8    # endif
9    # ifdef ALLOW_KPP
10    #  include "KPP_OPTIONS.h"
11    # endif
12    # ifdef ALLOW_PTRACERS
13    #  include "PTRACERS_OPTIONS.h"
14    # endif
15    #endif /* ALLOW_AUTODIFF_TAMC */
16    
17  CBOP  CBOP
18  C     !ROUTINE: THERMODYNAMICS  C     !ROUTINE: THERMODYNAMICS
# Line 71  C     == Global variables === Line 82  C     == Global variables ===
82  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_PASSIVE_TRACER
83  #include "TR1.h"  #include "TR1.h"
84  #endif  #endif
85    #ifdef ALLOW_TIMEAVE
86    #include "TIMEAVE_STATV.h"
87    #endif
88    
89  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
90  # include "tamc.h"  # include "tamc.h"
91  # include "tamc_keys.h"  # include "tamc_keys.h"
92  # include "FFIELDS.h"  # include "FFIELDS.h"
93    # include "EOS.h"
94  # ifdef ALLOW_KPP  # ifdef ALLOW_KPP
95  #  include "KPP.h"  #  include "KPP.h"
96  # endif  # endif
97  # ifdef ALLOW_GMREDI  # ifdef ALLOW_GMREDI
98  #  include "GMREDI.h"  #  include "GMREDI.h"
99  # endif  # endif
100    # ifdef ALLOW_PTRACERS
101    #  include "PTRACERS.h"
102    # endif
103  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
 #ifdef ALLOW_TIMEAVE  
 #include "TIMEAVE_STATV.h"  
 #endif  
104    
105  C     !INPUT/OUTPUT PARAMETERS:  C     !INPUT/OUTPUT PARAMETERS:
106  C     == Routine arguments ==  C     == Routine arguments ==
# Line 109  C                                      i Line 125  C                                      i
125  C                                      so we need an fVer for each  C                                      so we need an fVer for each
126  C                                      variable.  C                                      variable.
127  C     rhoK, rhoKM1   - Density at current level, and level above  C     rhoK, rhoKM1   - Density at current level, and level above
 C     phiHyd         - Hydrostatic part of the potential phiHydi.  
 C                      In z coords phiHydiHyd is the hydrostatic  
 C                      Potential (=pressure/rho0) anomaly  
 C                      In p coords phiHydiHyd is the geopotential  
 C                      surface height anomaly.  
128  C     phiSurfX, - gradient of Surface potentiel (Pressure/rho, ocean)  C     phiSurfX, - gradient of Surface potentiel (Pressure/rho, ocean)
129  C     phiSurfY             or geopotentiel (atmos) in X and Y direction  C     phiSurfY             or geopotentiel (atmos) in X and Y direction
130  C     KappaRT,       - Total diffusion in vertical for T and S.  C     KappaRT,       - Total diffusion in vertical for T and S.
131  C     KappaRS          (background + spatially varying, isopycnal term).  C     KappaRS          (background + spatially varying, isopycnal term).
132    C     useVariableK   = T when vertical diffusion is not constant
133  C     iMin, iMax     - Ranges and sub-block indices on which calculations  C     iMin, iMax     - Ranges and sub-block indices on which calculations
134  C     jMin, jMax       are applied.  C     jMin, jMax       are applied.
135  C     bi, bj  C     bi, bj
# Line 133  C                      index into fVerTe Line 145  C                      index into fVerTe
145        _RL fVerT   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerT   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
146        _RL fVerS   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerS   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
147        _RL fVerTr1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerTr1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
       _RL phiHyd  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)  
148        _RL rhokm1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rhokm1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
149        _RL rhok    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rhok    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
150        _RL phiSurfX(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL phiSurfX(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
# Line 145  C                      index into fVerTe Line 156  C                      index into fVerTe
156        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
157  C     This is currently used by IVDC and Diagnostics  C     This is currently used by IVDC and Diagnostics
158        _RL ConvectCount (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL ConvectCount (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
159          LOGICAL useVariableK
160        INTEGER iMin, iMax        INTEGER iMin, iMax
161        INTEGER jMin, jMax        INTEGER jMin, jMax
162        INTEGER bi, bj        INTEGER bi, bj
163        INTEGER i, j        INTEGER i, j
164        INTEGER k, km1, kup, kDown        INTEGER k, km1, kup, kDown
165          INTEGER iTracer
166    
167  CEOP  CEOP
168    
169    #ifndef DISABLE_DEBUGMODE
170             IF ( debugLevel .GE. debLevB )
171         &    CALL DEBUG_ENTER('FORWARD_STEP',myThid)
172    #endif
173    
174  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
175  C--   dummy statement to end declaration part  C--   dummy statement to end declaration part
176        ikey = 1        ikey = 1
177          itdkey = 1
178  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
179    
 C--   Set up work arrays with valid (i.e. not NaN) values  
 C     These inital values do not alter the numerical results. They  
 C     just ensure that all memory references are to valid floating  
 C     point numbers. This prevents spurious hardware signals due to  
 C     uninitialised but inert locations.  
       DO j=1-OLy,sNy+OLy  
        DO i=1-OLx,sNx+OLx  
         xA(i,j)      = 0. _d 0  
         yA(i,j)      = 0. _d 0  
         uTrans(i,j)  = 0. _d 0  
         vTrans(i,j)  = 0. _d 0  
         DO k=1,Nr  
          phiHyd(i,j,k)  = 0. _d 0  
          sigmaX(i,j,k) = 0. _d 0  
          sigmaY(i,j,k) = 0. _d 0  
          sigmaR(i,j,k) = 0. _d 0  
         ENDDO  
         rhoKM1 (i,j) = 0. _d 0  
         rhok   (i,j) = 0. _d 0  
         phiSurfX(i,j) = 0. _d 0  
         phiSurfY(i,j) = 0. _d 0  
        ENDDO  
       ENDDO  
   
   
180  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
181  C--   HPF directive to help TAMC  C--   HPF directive to help TAMC
182  CHPF$ INDEPENDENT  CHPF$ INDEPENDENT
# Line 193  CHPF$ INDEPENDENT Line 187  CHPF$ INDEPENDENT
187  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
188  C--    HPF directive to help TAMC  C--    HPF directive to help TAMC
189  CHPF$  INDEPENDENT, NEW (rTrans,fVerT,fVerS  CHPF$  INDEPENDENT, NEW (rTrans,fVerT,fVerS
190  CHPF$&                  ,phiHyd,utrans,vtrans,xA,yA  CHPF$&                  ,utrans,vtrans,xA,yA
191  CHPF$&                  ,KappaRT,KappaRS  CHPF$&                  ,KappaRT,KappaRS
192  CHPF$&                  )  CHPF$&                  )
193  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
# Line 208  CHPF$&                  ) Line 202  CHPF$&                  )
202            act3 = myThid - 1            act3 = myThid - 1
203            max3 = nTx*nTy            max3 = nTx*nTy
204            act4 = ikey_dynamics - 1            act4 = ikey_dynamics - 1
205            ikey = (act1 + 1) + act2*max1            itdkey = (act1 + 1) + act2*max1
206       &                      + act3*max1*max2       &                      + act3*max1*max2
207       &                      + act4*max1*max2*max3       &                      + act4*max1*max2*max3
208  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
209    
210  C--     Set up work arrays that need valid initial values  C--   Set up work arrays with valid (i.e. not NaN) values
211    C     These inital values do not alter the numerical results. They
212    C     just ensure that all memory references are to valid floating
213    C     point numbers. This prevents spurious hardware signals due to
214    C     uninitialised but inert locations.
215    
216          DO j=1-OLy,sNy+OLy          DO j=1-OLy,sNy+OLy
217           DO i=1-OLx,sNx+OLx           DO i=1-OLx,sNx+OLx
218              xA(i,j)        = 0. _d 0
219              yA(i,j)        = 0. _d 0
220              uTrans(i,j)    = 0. _d 0
221              vTrans(i,j)    = 0. _d 0
222              rhok   (i,j)   = 0. _d 0
223              rhoKM1 (i,j)   = 0. _d 0
224              phiSurfX(i,j)  = 0. _d 0
225              phiSurfY(i,j)  = 0. _d 0
226            rTrans (i,j)   = 0. _d 0            rTrans (i,j)   = 0. _d 0
227            fVerT  (i,j,1) = 0. _d 0            fVerT  (i,j,1) = 0. _d 0
228            fVerT  (i,j,2) = 0. _d 0            fVerT  (i,j,2) = 0. _d 0
# Line 230  C--     Set up work arrays that need val Line 237  C--     Set up work arrays that need val
237           DO j=1-OLy,sNy+OLy           DO j=1-OLy,sNy+OLy
238            DO i=1-OLx,sNx+OLx            DO i=1-OLx,sNx+OLx
239  C This is currently also used by IVDC and Diagnostics  C This is currently also used by IVDC and Diagnostics
240               sigmaX(i,j,k) = 0. _d 0
241               sigmaY(i,j,k) = 0. _d 0
242               sigmaR(i,j,k) = 0. _d 0
243             ConvectCount(i,j,k) = 0.             ConvectCount(i,j,k) = 0.
244             KappaRT(i,j,k) = 0. _d 0             KappaRT(i,j,k)    = 0. _d 0
245             KappaRS(i,j,k) = 0. _d 0             KappaRS(i,j,k)    = 0. _d 0
246  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
247             gT(i,j,k,bi,bj) = 0. _d 0  cph all the following init. are necessary for TAF
248             gS(i,j,k,bi,bj) = 0. _d 0  cph although some of these are re-initialised later.
249  #ifdef ALLOW_PASSIVE_TRACER             gT(i,j,k,bi,bj)   = 0. _d 0
250               gS(i,j,k,bi,bj)   = 0. _d 0
251    # ifdef ALLOW_PASSIVE_TRACER
252             gTr1(i,j,k,bi,bj) = 0. _d 0             gTr1(i,j,k,bi,bj) = 0. _d 0
253  #endif  # endif
254  #endif  # ifdef ALLOW_PTRACERS
255               DO iTracer=1,PTRACERS_numInUse
256                gPTr(i,j,k,bi,bj,itracer) = 0. _d 0
257               ENDDO
258    # endif
259    # ifdef ALLOW_GMREDI
260               Kwx(i,j,k,bi,bj)  = 0. _d 0
261               Kwy(i,j,k,bi,bj)  = 0. _d 0
262               Kwz(i,j,k,bi,bj)  = 0. _d 0
263    #  ifdef GM_NON_UNITY_DIAGONAL
264               Kux(i,j,k,bi,bj)  = 0. _d 0
265               Kvy(i,j,k,bi,bj)  = 0. _d 0
266    #  endif
267    #  ifdef GM_EXTRA_DIAGONAL
268               Kuz(i,j,k,bi,bj)  = 0. _d 0
269               Kvz(i,j,k,bi,bj)  = 0. _d 0
270    #  endif
271    #  ifdef GM_BOLUS_ADVEC
272               GM_PsiX(i,j,k,bi,bj)  = 0. _d 0
273               GM_PsiY(i,j,k,bi,bj)  = 0. _d 0
274    #  endif
275    #  ifdef GM_VISBECK_VARIABLE_K
276               VisbeckK(i,j,bi,bj)   = 0. _d 0
277    #  endif
278    # endif /* ALLOW_GMREDI */
279    #endif /* ALLOW_AUTODIFF_TAMC */
280            ENDDO            ENDDO
281           ENDDO           ENDDO
282          ENDDO          ENDDO
283    
284          iMin = 1-OLx+1          iMin = 1-OLx
285          iMax = sNx+OLx          iMax = sNx+OLx
286          jMin = 1-OLy+1          jMin = 1-OLy
287          jMax = sNy+OLy          jMax = sNy+OLy
288    
   
289  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
290  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
291  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
292    CADJ STORE totphihyd
293    CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte
294  #ifdef ALLOW_KPP  #ifdef ALLOW_KPP
295  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
296  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
297  #endif  #endif
298  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
299    
300    #ifndef DISABLE_DEBUGMODE
301            IF ( debugLevel .GE. debLevB )
302         &    CALL DEBUG_MSG('ENTERING UPWARD K LOOP',myThid)
303    #endif
304    
305  C--     Start of diagnostic loop  C--     Start of diagnostic loop
306          DO k=Nr,1,-1          DO k=Nr,1,-1
307    
# Line 267  C? Patrick, is this formula correct now Line 310  C? Patrick, is this formula correct now
310  C? Do we still need this?  C? Do we still need this?
311  cph kkey formula corrected.  cph kkey formula corrected.
312  cph Needed for rhok, rhokm1, in the case useGMREDI.  cph Needed for rhok, rhokm1, in the case useGMREDI.
313           kkey = (ikey-1)*Nr + k           kkey = (itdkey-1)*Nr + k
 CADJ STORE rhokm1(:,:) = comlev1_bibj_k ,       key=kkey, byte=isbyte  
 CADJ STORE rhok  (:,:) = comlev1_bibj_k ,       key=kkey, byte=isbyte  
314  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
315    
316  C--       Integrate continuity vertically for vertical velocity  C--       Integrate continuity vertically for vertical velocity
317            CALL INTEGRATE_FOR_W(  c         CALL INTEGRATE_FOR_W(
318       I                         bi, bj, k, uVel, vVel,  c    I                         bi, bj, k, uVel, vVel,
319       O                         wVel,  c    O                         wVel,
320       I                         myThid )  c    I                         myThid )
321    
322  #ifdef    ALLOW_OBCS  #ifdef    ALLOW_OBCS
323  #ifdef    ALLOW_NONHYDROSTATIC  #ifdef    ALLOW_NONHYDROSTATIC
324  C--       Apply OBC to W if in N-H mode  C--       Apply OBC to W if in N-H mode
325            IF (useOBCS.AND.nonHydrostatic) THEN  c         IF (useOBCS.AND.nonHydrostatic) THEN
326              CALL OBCS_APPLY_W( bi, bj, k, wVel, myThid )  c           CALL OBCS_APPLY_W( bi, bj, k, wVel, myThid )
327            ENDIF  c         ENDIF
328  #endif    /* ALLOW_NONHYDROSTATIC */  #endif    /* ALLOW_NONHYDROSTATIC */
329  #endif    /* ALLOW_OBCS */  #endif    /* ALLOW_OBCS */
330    
331    C--       Attention: by defining "SINGLE_LAYER_MODE" in CPP_OPTIONS.h
332    C--       MOST of THERMODYNAMICS will be disabled
333    #ifndef SINGLE_LAYER_MODE
334    
335  C--       Calculate gradients of potential density for isoneutral  C--       Calculate gradients of potential density for isoneutral
336  C         slope terms (e.g. GM/Redi tensor or IVDC diffusivity)  C         slope terms (e.g. GM/Redi tensor or IVDC diffusivity)
337  c         IF ( k.GT.1 .AND. (useGMRedi.OR.ivdc_kappa.NE.0.) ) THEN  c         IF ( k.GT.1 .AND. (useGMRedi.OR.ivdc_kappa.NE.0.) ) THEN
338            IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.) ) THEN            IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.) ) THEN
339    #ifndef DISABLE_DEBUGMODE
340                IF ( debugLevel .GE. debLevB )
341         &       CALL DEBUG_CALL('FIND_RHO',myThid)
342    #endif
343  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
344  CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
345  CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
346  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
347              CALL FIND_RHO(              CALL FIND_RHO(
348       I        bi, bj, iMin, iMax, jMin, jMax, k, k, eosType,       I        bi, bj, iMin, iMax, jMin, jMax, k, k,
349       I        theta, salt,       I        theta, salt,
350       O        rhoK,       O        rhoK,
351       I        myThid )       I        myThid )
352    
353              IF (k.GT.1) THEN              IF (k.GT.1) THEN
354  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
355  CADJ STORE theta(:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE theta(:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
356  CADJ STORE salt (:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE salt (:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
357  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
358               CALL FIND_RHO(               CALL FIND_RHO(
359       I        bi, bj, iMin, iMax, jMin, jMax, k-1, k, eosType,       I        bi, bj, iMin, iMax, jMin, jMax, k-1, k,
360       I        theta, salt,       I        theta, salt,
361       O        rhoKm1,       O        rhoKm1,
362       I        myThid )       I        myThid )
363              ENDIF              ENDIF
364    #ifndef DISABLE_DEBUGMODE
365                IF ( debugLevel .GE. debLevB )
366         &       CALL DEBUG_CALL('GRAD_SIGMA',myThid)
367    #endif
368              CALL GRAD_SIGMA(              CALL GRAD_SIGMA(
369       I             bi, bj, iMin, iMax, jMin, jMax, k,       I             bi, bj, iMin, iMax, jMin, jMax, k,
370       I             rhoK, rhoKm1, rhoK,       I             rhoK, rhoKm1, rhoK,
# Line 318  CADJ STORE salt (:,:,k-1,bi,bj) = comlev Line 372  CADJ STORE salt (:,:,k-1,bi,bj) = comlev
372       I             myThid )       I             myThid )
373            ENDIF            ENDIF
374    
375    #ifdef ALLOW_AUTODIFF_TAMC
376    CADJ STORE rhok   (:,:) = comlev1_bibj_k ,       key=kkey, byte=isbyte
377    CADJ STORE rhokm1 (:,:) = comlev1_bibj_k ,       key=kkey, byte=isbyte
378    #endif /* ALLOW_AUTODIFF_TAMC */
379  C--       Implicit Vertical Diffusion for Convection  C--       Implicit Vertical Diffusion for Convection
380  c ==> should use sigmaR !!!  c ==> should use sigmaR !!!
381            IF (k.GT.1 .AND. ivdc_kappa.NE.0.) THEN            IF (k.GT.1 .AND. ivdc_kappa.NE.0.) THEN
382    #ifndef DISABLE_DEBUGMODE
383                IF ( debugLevel .GE. debLevB )
384         &       CALL DEBUG_CALL('CALC_IVDC',myThid)
385    #endif
386              CALL CALC_IVDC(              CALL CALC_IVDC(
387       I        bi, bj, iMin, iMax, jMin, jMax, k,       I        bi, bj, iMin, iMax, jMin, jMax, k,
388       I        rhoKm1, rhoK,       I        rhoKm1, rhoK,
# Line 328  c ==> should use sigmaR !!! Line 390  c ==> should use sigmaR !!!
390       I        myTime, myIter, myThid)       I        myTime, myIter, myThid)
391            ENDIF            ENDIF
392    
393    #endif /* SINGLE_LAYER_MODE */
394    
395  C--     end of diagnostic k loop (Nr:1)  C--     end of diagnostic k loop (Nr:1)
396          ENDDO          ENDDO
397    
398  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
399  cph avoids recomputation of integrate_for_w  cph avoids recomputation of integrate_for_w
400  CADJ STORE wvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE wvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
401  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
402    
403  #ifdef  ALLOW_OBCS  #ifdef  ALLOW_OBCS
404  C--     Calculate future values on open boundaries  C--     Calculate future values on open boundaries
405          IF (useOBCS) THEN          IF (useOBCS) THEN
406            CALL OBCS_CALC( bi, bj, myTime+deltaT,  #ifndef DISABLE_DEBUGMODE
407              IF ( debugLevel .GE. debLevB )
408         &     CALL DEBUG_CALL('OBCS_CALC',myThid)
409    #endif
410              CALL OBCS_CALC( bi, bj, myTime+deltaT, myIter+1,
411       I            uVel, vVel, wVel, theta, salt,       I            uVel, vVel, wVel, theta, salt,
412       I            myThid )       I            myThid )
413          ENDIF          ENDIF
414  #endif  /* ALLOW_OBCS */  #endif  /* ALLOW_OBCS */
415    
416    
417    #ifdef ALLOW_THERM_SEAICE
418           IF (useThermSeaIce) THEN
419    #ifndef DISABLE_DEBUGMODE
420            IF ( debugLevel .GE. debLevB )
421         &    CALL DEBUG_CALL('ICE_FORCING',myThid)
422    #endif
423    C--     Determines forcing terms based on external fields
424    C       including effects from ice
425            CALL ICE_FORCING(
426         I             bi, bj, iMin, iMax, jMin, jMax,
427         I             myThid )
428           ELSE
429    #else  /* ALLOW_THERM_SEAICE */
430           IF (.TRUE.) THEN
431    #endif /* ALLOW_THERM_SEAICE */
432    
433  C--     Determines forcing terms based on external fields  C--     Determines forcing terms based on external fields
434  C       relaxation terms, etc.  C       relaxation terms, etc.
435          CALL EXTERNAL_FORCING_SURF(  #ifndef DISABLE_DEBUGMODE
436            IF ( debugLevel .GE. debLevB )
437         &    CALL DEBUG_CALL('EXTERNAL_FORCING_SURF',myThid)
438    #endif
439            CALL EXTERNAL_FORCING_SURF(
440       I             bi, bj, iMin, iMax, jMin, jMax,       I             bi, bj, iMin, iMax, jMin, jMax,
441       I             myThid )       I             myThid )
442    
443    C--    end of if/else block useThermSeaIce --
444           ENDIF
445    
446    
447  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
448  cph needed for KPP  cph needed for KPP
449  CADJ STORE surfacetendencyU(:,:,bi,bj)  CADJ STORE surfacetendencyU(:,:,bi,bj)
450  CADJ &     = comlev1_bibj, key=ikey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte
451  CADJ STORE surfacetendencyV(:,:,bi,bj)  CADJ STORE surfacetendencyV(:,:,bi,bj)
452  CADJ &     = comlev1_bibj, key=ikey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte
453  CADJ STORE surfacetendencyS(:,:,bi,bj)  CADJ STORE surfacetendencyS(:,:,bi,bj)
454  CADJ &     = comlev1_bibj, key=ikey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte
455  CADJ STORE surfacetendencyT(:,:,bi,bj)  CADJ STORE surfacetendencyT(:,:,bi,bj)
456  CADJ &     = comlev1_bibj, key=ikey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte
457  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
458    
459    C--     Attention: by defining "SINGLE_LAYER_MODE" in CPP_OPTIONS.h
460    C--     MOST of THERMODYNAMICS will be disabled
461    #ifndef SINGLE_LAYER_MODE
462    
463  #ifdef  ALLOW_GMREDI  #ifdef  ALLOW_GMREDI
464    
465  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
466  CADJ STORE sigmaX(:,:,:) = comlev1, key=ikey, byte=isbyte  cph storing here is needed only for one GMREDI_OPTIONS:
467  CADJ STORE sigmaY(:,:,:) = comlev1, key=ikey, byte=isbyte  cph define GM_BOLUS_ADVEC
468  CADJ STORE sigmaR(:,:,:) = comlev1, key=ikey, byte=isbyte  cph but I've avoided the #ifdef for now, in case more things change
469    CADJ STORE sigmaX(:,:,:)        = comlev1_bibj, key=itdkey, byte=isbyte
470    CADJ STORE sigmaY(:,:,:)        = comlev1_bibj, key=itdkey, byte=isbyte
471    CADJ STORE sigmaR(:,:,:)        = comlev1_bibj, key=itdkey, byte=isbyte
472  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
473    
474  C--     Calculate iso-neutral slopes for the GM/Redi parameterisation  C--     Calculate iso-neutral slopes for the GM/Redi parameterisation
475          IF (useGMRedi) THEN          IF (useGMRedi) THEN
476            DO k=1,Nr  #ifndef DISABLE_DEBUGMODE
477              CALL GMREDI_CALC_TENSOR(            IF ( debugLevel .GE. debLevB )
478       I             bi, bj, iMin, iMax, jMin, jMax, k,       &     CALL DEBUG_CALL('GMREDI_CALC_TENSOR',myThid)
479    #endif
480              CALL GMREDI_CALC_TENSOR(
481         I             bi, bj, iMin, iMax, jMin, jMax,
482       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
483       I             myThid )       I             myThid )
           ENDDO  
484  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
485          ELSE          ELSE
486            DO k=1, Nr            CALL GMREDI_CALC_TENSOR_DUMMY(
487              CALL GMREDI_CALC_TENSOR_DUMMY(       I             bi, bj, iMin, iMax, jMin, jMax,
      I             bi, bj, iMin, iMax, jMin, jMax, k,  
488       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
489       I             myThid )       I             myThid )
           ENDDO  
490  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
491          ENDIF          ENDIF
492    
493  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
494  CADJ STORE Kwx(:,:,:,bi,bj)   = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE Kwx(:,:,:,bi,bj)   = comlev1_bibj, key=itdkey, byte=isbyte
495  CADJ STORE Kwy(:,:,:,bi,bj)   = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE Kwy(:,:,:,bi,bj)   = comlev1_bibj, key=itdkey, byte=isbyte
496  CADJ STORE Kwz(:,:,:,bi,bj)   = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE Kwz(:,:,:,bi,bj)   = comlev1_bibj, key=itdkey, byte=isbyte
497  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
498    
499  #endif  /* ALLOW_GMREDI */  #endif  /* ALLOW_GMREDI */
# Line 399  CADJ STORE Kwz(:,:,:,bi,bj)   = comlev1_ Line 501  CADJ STORE Kwz(:,:,:,bi,bj)   = comlev1_
501  #ifdef  ALLOW_KPP  #ifdef  ALLOW_KPP
502  C--     Compute KPP mixing coefficients  C--     Compute KPP mixing coefficients
503          IF (useKPP) THEN          IF (useKPP) THEN
504    #ifndef DISABLE_DEBUGMODE
505              IF ( debugLevel .GE. debLevB )
506         &     CALL DEBUG_CALL('KPP_CALC',myThid)
507    #endif
508            CALL KPP_CALC(            CALL KPP_CALC(
509       I                  bi, bj, myTime, myThid )       I                  bi, bj, myTime, myThid )
510  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 413  CADJ STORE KPPghat   (:,:,:,bi,bj) Line 519  CADJ STORE KPPghat   (:,:,:,bi,bj)
519  CADJ &   , KPPdiffKzT(:,:,:,bi,bj)  CADJ &   , KPPdiffKzT(:,:,:,bi,bj)
520  CADJ &   , KPPdiffKzS(:,:,:,bi,bj)  CADJ &   , KPPdiffKzS(:,:,:,bi,bj)
521  CADJ &   , KPPfrac   (:,:  ,bi,bj)  CADJ &   , KPPfrac   (:,:  ,bi,bj)
522  CADJ &                 = comlev1_bibj, key=ikey, byte=isbyte  CADJ &                 = comlev1_bibj, key=itdkey, byte=isbyte
523  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
524    
525  #endif  /* ALLOW_KPP */  #endif  /* ALLOW_KPP */
526    
527  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
528  CADJ STORE KappaRT(:,:,:)     = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE KappaRT(:,:,:)     = comlev1_bibj, key=itdkey, byte=isbyte
529  CADJ STORE KappaRS(:,:,:)     = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE KappaRS(:,:,:)     = comlev1_bibj, key=itdkey, byte=isbyte
530  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
531  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
532  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
533  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
534  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_PASSIVE_TRACER
535  CADJ STORE tr1  (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE tr1  (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
536    #endif
537    #ifdef ALLOW_PTRACERS
538    cph-- moved to forward_step to avoid key computation
539    cphCADJ STORE ptracer(:,:,:,bi,bj,itracer) = comlev1_bibj,
540    cphCADJ &                              key=itdkey, byte=isbyte
541  #endif  #endif
542  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
543    
544  #ifdef ALLOW_AIM  #ifdef ALLOW_AIM
545  C       AIM - atmospheric intermediate model, physics package code.  C       AIM - atmospheric intermediate model, physics package code.
 C note(jmc) : phiHyd=0 at this point but is not really used in Molteni Physics  
546          IF ( useAIM ) THEN          IF ( useAIM ) THEN
547           CALL TIMER_START('AIM_DO_ATMOS_PHYS      [DYNAMICS]', myThid)  #ifndef DISABLE_DEBUGMODE
548           CALL AIM_DO_ATMOS_PHYSICS( phiHyd, bi, bj, myTime, myThid )            IF ( debugLevel .GE. debLevB )
549           CALL TIMER_STOP ('AIM_DO_ATMOS_PHYS      [DYNAMICS]', myThid)       &     CALL DEBUG_CALL('AIM_DO_PHYSICS',myThid)
550    #endif
551             CALL TIMER_START('AIM_DO_PHYSICS   [THERMODYNAMICS]', myThid)
552             CALL AIM_DO_PHYSICS( bi, bj, myTime, myIter, myThid )
553             CALL TIMER_STOP( 'AIM_DO_PHYSICS   [THERMODYNAMICS]', myThid)
554          ENDIF          ENDIF
555  #endif /* ALLOW_AIM */  #endif /* ALLOW_AIM */
556    
 #ifdef ALLOW_TIMEAVE  
         IF (taveFreq.GT.0. .AND. ivdc_kappa.NE.0.) THEN  
           CALL TIMEAVE_CUMULATE(ConvectCountTave, ConvectCount, Nr,  
      I                               deltaTclock, bi, bj, myThid)  
         ENDIF  
 #endif /* ALLOW_TIMEAVE */  
   
557  #ifndef DISABLE_MULTIDIM_ADVECTION  #ifndef DISABLE_MULTIDIM_ADVECTION
558  C--     Some advection schemes are better calculated using a multi-dimensional  C--     Some advection schemes are better calculated using a multi-dimensional
559  C       method in the absence of any other terms and, if used, is done here.  C       method in the absence of any other terms and, if used, is done here.
# Line 458  C to be able to exclude this scheme to a Line 565  C to be able to exclude this scheme to a
565  C recomputation. It *is* differentiable, if you need it.  C recomputation. It *is* differentiable, if you need it.
566  C Edit GAD_OPTIONS.h and #define DISABLE_MULTIDIM_ADVECTION to  C Edit GAD_OPTIONS.h and #define DISABLE_MULTIDIM_ADVECTION to
567  C disable this section of code.  C disable this section of code.
568          IF (multiDimAdvection) THEN          IF (tempMultiDimAdvec) THEN
569           IF (tempStepping .AND.  #ifndef DISABLE_DEBUGMODE
570       &       tempAdvScheme.NE.ENUM_CENTERED_2ND .AND.            IF ( debugLevel .GE. debLevB )
571       &       tempAdvScheme.NE.ENUM_UPWIND_3RD .AND.       &     CALL DEBUG_CALL('GAD_ADVECTION',myThid)
572       &       tempAdvScheme.NE.ENUM_CENTERED_4TH ) THEN  #endif
573            CALL GAD_ADVECTION(bi,bj,tempAdvScheme,GAD_TEMPERATURE,            CALL GAD_ADVECTION(bi,bj,tempAdvScheme,GAD_TEMPERATURE,
574       U                      theta,gT,       U                      theta,gT,
575       I                      myTime,myIter,myThid)       I                      myTime,myIter,myThid)
576           ENDIF          ENDIF
577           IF (saltStepping .AND.          IF (saltMultiDimAdvec) THEN
578       &       saltAdvScheme.NE.ENUM_CENTERED_2ND .AND.  #ifndef DISABLE_DEBUGMODE
579       &       saltAdvScheme.NE.ENUM_UPWIND_3RD .AND.            IF ( debugLevel .GE. debLevB )
580       &       saltAdvScheme.NE.ENUM_CENTERED_4TH ) THEN       &     CALL DEBUG_CALL('GAD_ADVECTION',myThid)
581    #endif
582            CALL GAD_ADVECTION(bi,bj,saltAdvScheme,GAD_SALINITY,            CALL GAD_ADVECTION(bi,bj,saltAdvScheme,GAD_SALINITY,
583       U                      salt,gS,       U                      salt,gS,
584       I                      myTime,myIter,myThid)       I                      myTime,myIter,myThid)
          ENDIF  
585          ENDIF          ENDIF
586    C Since passive tracers are configurable separately from T,S we
587    C call the multi-dimensional method for PTRACERS regardless
588    C of whether multiDimAdvection is set or not.
589    #ifdef ALLOW_PTRACERS
590            IF ( usePTRACERS ) THEN
591    #ifndef DISABLE_DEBUGMODE
592              IF ( debugLevel .GE. debLevB )
593         &     CALL DEBUG_CALL('PTRACERS_ADVECTION',myThid)
594    #endif
595             CALL PTRACERS_ADVECTION( bi,bj,myIter,myTime,myThid )
596            ENDIF
597    #endif /* ALLOW_PTRACERS */
598  #endif /* DISABLE_MULTIDIM_ADVECTION */  #endif /* DISABLE_MULTIDIM_ADVECTION */
599    
600    #ifndef DISABLE_DEBUGMODE
601           IF ( debugLevel .GE. debLevB )
602         &    CALL DEBUG_MSG('ENTERING DOWNWARD K LOOP',myThid)
603    #endif
604    
605  C--     Start of thermodynamics loop  C--     Start of thermodynamics loop
606          DO k=Nr,1,-1          DO k=Nr,1,-1
607  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
608  C? Patrick Is this formula correct?  C? Patrick Is this formula correct?
609  cph Yes, but I rewrote it.  cph Yes, but I rewrote it.
610  cph Also, the KappaR? need the index and subscript k!  cph Also, the KappaR? need the index and subscript k!
611           kkey = (ikey-1)*Nr + k           kkey = (itdkey-1)*Nr + k
612  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
613    
614  C--       km1    Points to level above k (=k-1)  C--       km1    Points to level above k (=k-1)
# Line 506  C--       Get temporary terms used by te Line 630  C--       Get temporary terms used by te
630       O         xA,yA,uTrans,vTrans,rTrans,maskUp,       O         xA,yA,uTrans,vTrans,rTrans,maskUp,
631       I         myThid)       I         myThid)
632    
633    #ifdef ALLOW_GMREDI
634    
635    C--   Residual transp = Bolus transp + Eulerian transp
636              IF (useGMRedi) THEN
637                CALL GMREDI_CALC_UVFLOW(
638         &            uTrans, vTrans, bi, bj, k, myThid)
639                IF (K.GE.2) CALL GMREDI_CALC_WFLOW(
640         &                    rTrans, bi, bj, k, myThid)
641              ENDIF
642    
643    #ifdef ALLOW_AUTODIFF_TAMC
644    #ifdef GM_BOLUS_ADVEC
645    CADJ STORE uTrans(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
646    CADJ STORE vTrans(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
647    CADJ STORE rTrans(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
648    #endif
649    #endif /* ALLOW_AUTODIFF_TAMC */
650    
651    #endif /* ALLOW_GMREDI */
652    
653  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
654  CADJ STORE KappaRT(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE KappaRT(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte
655  CADJ STORE KappaRS(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE KappaRS(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte
# Line 539  C        and step forward storing result Line 683  C        and step forward storing result
683       I         theta, gT,       I         theta, gT,
684       I         myIter, myThid)       I         myIter, myThid)
685           ENDIF           ENDIF
686    
687    #ifdef ALLOW_THERM_SEAICE
688             IF (useThermSeaIce .AND. k.EQ.1) THEN
689               CALL ICE_FREEZE( bi,bj, iMin,iMax,jMin,jMax, myThid )
690             ENDIF
691    #endif
692    
693           IF ( saltStepping ) THEN           IF ( saltStepping ) THEN
694             CALL CALC_GS(             CALL CALC_GS(
695       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
# Line 565  C        and step forward storing result Line 716  C        and step forward storing result
716       I         myIter,myThid)       I         myIter,myThid)
717           ENDIF           ENDIF
718  #endif  #endif
719    #ifdef ALLOW_PTRACERS
720             IF ( usePTRACERS ) THEN
721               CALL PTRACERS_INTEGRATE(
722         I         bi,bj,k,
723         I         xA,yA,uTrans,vTrans,rTrans,maskUp,
724         X         KappaRS,
725         I         myIter,myTime,myThid)
726             ENDIF
727    #endif /* ALLOW_PTRACERS */
728    
729  #ifdef   ALLOW_OBCS  #ifdef   ALLOW_OBCS
730  C--      Apply open boundary conditions  C--      Apply open boundary conditions
# Line 574  C--      Apply open boundary conditions Line 734  C--      Apply open boundary conditions
734  #endif   /* ALLOW_OBCS */  #endif   /* ALLOW_OBCS */
735    
736  C--      Freeze water  C--      Freeze water
737           IF (allowFreezing) THEN           IF ( allowFreezing .AND. .NOT. useSEAICE
738         &       .AND. .NOT.(useThermSeaIce.AND.k.EQ.1) ) THEN
739  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
740  CADJ STORE gT(:,:,k,bi,bj) = comlev1_bibj_k  CADJ STORE gT(:,:,k,bi,bj) = comlev1_bibj_k
741  CADJ &   , key = kkey, byte = isbyte  CADJ &   , key = kkey, byte = isbyte
742  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
743              CALL FREEZE( bi, bj, iMin, iMax, jMin, jMax, k, myThid )              CALL FREEZE( bi, bj, iMin, iMax, jMin, jMax, k, myThid )
744           END IF           ENDIF
745    
746  C--     end of thermodynamic k loop (Nr:1)  C--     end of thermodynamic k loop (Nr:1)
747          ENDDO          ENDDO
748    
749    cswdice -- add ---
750    #ifdef ALLOW_THERM_SEAICE
751    c timeaveraging for ice model values
752               CALL ICE_AVE(bi,bj,iMin,iMax,jMin,jMax,myThid )
753    #endif
754    cswdice --- end add ---
755    
756    
757    
 #ifdef ALLOW_AUTODIFF_TAMC  
 C? Patrick? What about this one?  
 cph Keys iikey and idkey dont seem to be needed  
 cph since storing occurs on different tape for each  
 cph impldiff call anyways.  
 cph Thus, common block comlev1_impl isnt needed either.  
 cph Storing below needed in the case useGMREDI.  
         iikey = (ikey-1)*maximpl  
 #endif /* ALLOW_AUTODIFF_TAMC */  
758    
759  C--     Implicit diffusion  C--     Implicit diffusion
760          IF (implicitDiffusion) THEN          IF (implicitDiffusion) THEN
761    
762           IF (tempStepping) THEN           IF (tempStepping) THEN
763  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
764              idkey = iikey + 1  CADJ STORE gT(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte
 CADJ STORE gT(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte  
765  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
766              CALL IMPLDIFF(              CALL IMPLDIFF(
767       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
# Line 613  CADJ STORE gT(:,:,:,bi,bj) = comlev1_bib Line 772  CADJ STORE gT(:,:,:,bi,bj) = comlev1_bib
772    
773           IF (saltStepping) THEN           IF (saltStepping) THEN
774  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
775           idkey = iikey + 2  CADJ STORE gS(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte
 CADJ STORE gS(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte  
776  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
777              CALL IMPLDIFF(              CALL IMPLDIFF(
778       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
# Line 626  CADJ STORE gS(:,:,:,bi,bj) = comlev1_bib Line 784  CADJ STORE gS(:,:,:,bi,bj) = comlev1_bib
784  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_PASSIVE_TRACER
785           IF (tr1Stepping) THEN           IF (tr1Stepping) THEN
786  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
787  CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte  CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte
788  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
789            CALL IMPLDIFF(            CALL IMPLDIFF(
790       I      bi, bj, iMin, iMax, jMin, jMax,       I      bi, bj, iMin, iMax, jMin, jMax,
# Line 636  CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_b Line 794  CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_b
794           ENDIF           ENDIF
795  #endif  #endif
796    
797    #ifdef ALLOW_PTRACERS
798    C Vertical diffusion (implicit) for passive tracers
799             IF ( usePTRACERS ) THEN
800               CALL PTRACERS_IMPLDIFF( bi,bj,KappaRS,myThid )
801             ENDIF
802    #endif /* ALLOW_PTRACERS */
803    
804  #ifdef   ALLOW_OBCS  #ifdef   ALLOW_OBCS
805  C--      Apply open boundary conditions  C--      Apply open boundary conditions
806           IF (useOBCS) THEN           IF (useOBCS) THEN
# Line 648  C--      Apply open boundary conditions Line 813  C--      Apply open boundary conditions
813  C--     End If implicitDiffusion  C--     End If implicitDiffusion
814          ENDIF          ENDIF
815    
816  Ccs-  #ifdef ALLOW_TIMEAVE
817            IF (taveFreq.GT.0. .AND. ivdc_kappa.NE.0.) THEN
818              CALL TIMEAVE_CUMUL_1T(ConvectCountTave, ConvectCount,
819         I                           Nr, deltaTclock, bi, bj, myThid)
820            ENDIF
821            useVariableK = useKPP .OR. useGMredi .OR. ivdc_kappa.NE.0.
822            IF (taveFreq.GT.0. .AND. useVariableK ) THEN
823             IF (implicitDiffusion) THEN
824              CALL TIMEAVE_CUMUL_DIF_1T(TdiffRtave, gT, kappaRT,
825         I                        Nr, 3, deltaTclock, bi, bj, myThid)
826             ELSE
827              CALL TIMEAVE_CUMUL_DIF_1T(TdiffRtave, theta, kappaRT,
828         I                        Nr, 3, deltaTclock, bi, bj, myThid)
829             ENDIF
830            ENDIF
831    #endif /* ALLOW_TIMEAVE */
832    
833    #endif /* SINGLE_LAYER_MODE */
834    
835    C--   end bi,bj loops.
836         ENDDO         ENDDO
837        ENDDO        ENDDO
838    
839  #ifdef ALLOW_AIM  #ifndef DISABLE_DEBUGMODE
840        IF ( useAIM ) THEN        If (debugMode) THEN
841         CALL AIM_AIM2DYN_EXCHANGES( myTime, myThid )         CALL DEBUG_STATS_RL(Nr,uVel,'Uvel (THERMODYNAMICS)',myThid)
842        ENDIF         CALL DEBUG_STATS_RL(Nr,vVel,'Vvel (THERMODYNAMICS)',myThid)
843         _EXCH_XYZ_R8(gT,myThid)         CALL DEBUG_STATS_RL(Nr,wVel,'Wvel (THERMODYNAMICS)',myThid)
844         _EXCH_XYZ_R8(gS,myThid)         CALL DEBUG_STATS_RL(Nr,theta,'Theta (THERMODYNAMICS)',myThid)
845  #else         CALL DEBUG_STATS_RL(Nr,salt,'Salt (THERMODYNAMICS)',myThid)
846        IF (staggerTimeStep.AND.useCubedSphereExchange) THEN         CALL DEBUG_STATS_RL(Nr,Gt,'Gt (THERMODYNAMICS)',myThid)
847         _EXCH_XYZ_R8(gT,myThid)         CALL DEBUG_STATS_RL(Nr,Gs,'Gs (THERMODYNAMICS)',myThid)
848         _EXCH_XYZ_R8(gS,myThid)         CALL DEBUG_STATS_RL(Nr,GtNm1,'GtNm1 (THERMODYNAMICS)',myThid)
849           CALL DEBUG_STATS_RL(Nr,GsNm1,'GsNm1 (THERMODYNAMICS)',myThid)
850    #ifdef ALLOW_PTRACERS
851           IF ( usePTRACERS ) THEN
852             CALL PTRACERS_DEBUG(myThid)
853           ENDIF
854    #endif /* ALLOW_PTRACERS */
855        ENDIF        ENDIF
856  #endif /* ALLOW_AIM */  #endif
857    
858    #ifndef DISABLE_DEBUGMODE
859             IF ( debugLevel .GE. debLevB )
860         &    CALL DEBUG_LEAVE('FORWARD_STEP',myThid)
861    #endif
862    
863        RETURN        RETURN
864        END        END

Legend:
Removed from v.1.14  
changed lines
  Added in v.1.44

  ViewVC Help
Powered by ViewVC 1.1.22