/[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.18 by adcroft, Tue Mar 5 14:15:34 2002 UTC revision 1.79 by jmc, Tue Oct 19 02:39:58 2004 UTC
# Line 1  Line 1 
1  C $Header$  C $Header$
2  C $Name$  C $Name$
3    
4    #include "PACKAGES_CONFIG.h"
5  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
6    
7    #ifdef ALLOW_AUTODIFF_TAMC
8    # ifdef ALLOW_GMREDI
9    #  include "GMREDI_OPTIONS.h"
10    # endif
11    # ifdef ALLOW_KPP
12    #  include "KPP_OPTIONS.h"
13    # endif
14    #endif /* ALLOW_AUTODIFF_TAMC */
15    
16  CBOP  CBOP
17  C     !ROUTINE: THERMODYNAMICS  C     !ROUTINE: THERMODYNAMICS
18  C     !INTERFACE:  C     !INTERFACE:
# Line 68  C     == Global variables === Line 78  C     == Global variables ===
78  #include "DYNVARS.h"  #include "DYNVARS.h"
79  #include "GRID.h"  #include "GRID.h"
80  #include "GAD.h"  #include "GAD.h"
81  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_OFFLINE
82  #include "TR1.h"  #include "OFFLINE.h"
83    #endif
84    #ifdef ALLOW_PTRACERS
85    #include "PTRACERS_SIZE.h"
86    #include "PTRACERS.h"
87    #endif
88    #ifdef ALLOW_TIMEAVE
89    #include "TIMEAVE_STATV.h"
90  #endif  #endif
91    
92  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
93  # include "tamc.h"  # include "tamc.h"
94  # include "tamc_keys.h"  # include "tamc_keys.h"
95  # include "FFIELDS.h"  # include "FFIELDS.h"
96    # include "EOS.h"
97  # ifdef ALLOW_KPP  # ifdef ALLOW_KPP
98  #  include "KPP.h"  #  include "KPP.h"
99  # endif  # endif
100  # ifdef ALLOW_GMREDI  # ifdef ALLOW_GMREDI
101  #  include "GMREDI.h"  #  include "GMREDI.h"
102  # endif  # endif
103    # ifdef ALLOW_EBM
104    #  include "EBM.h"
105    # endif
106  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
107  #ifdef ALLOW_TIMEAVE  
 #include "TIMEAVE_STATV.h"  
 #endif  
108    
109  C     !INPUT/OUTPUT PARAMETERS:  C     !INPUT/OUTPUT PARAMETERS:
110  C     == Routine arguments ==  C     == Routine arguments ==
# Line 103  C                              transport Line 123  C                              transport
123  C                              o uTrans: Zonal transport  C                              o uTrans: Zonal transport
124  C                              o vTrans: Meridional transport  C                              o vTrans: Meridional transport
125  C                              o rTrans: Vertical transport  C                              o rTrans: Vertical transport
126    C     rTransKp1                o vertical volume transp. at interface k+1
127  C     maskUp                   o maskUp: land/water mask for W points  C     maskUp                   o maskUp: land/water mask for W points
128  C     fVer[STUV]               o fVer: Vertical flux term - note fVer  C     fVer[STUV]               o fVer: Vertical flux term - note fVer
129  C                                      is "pipelined" in the vertical  C                                      is "pipelined" in the vertical
130  C                                      so we need an fVer for each  C                                      so we need an fVer for each
131  C                                      variable.  C                                      variable.
 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.  
 C     phiSurfX, - gradient of Surface potentiel (Pressure/rho, ocean)  
 C     phiSurfY             or geopotentiel (atmos) in X and Y direction  
132  C     KappaRT,       - Total diffusion in vertical for T and S.  C     KappaRT,       - Total diffusion in vertical for T and S.
133  C     KappaRS          (background + spatially varying, isopycnal term).  C     KappaRS          (background + spatially varying, isopycnal term).
134    C     useVariableK   = T when vertical diffusion is not constant
135  C     iMin, iMax     - Ranges and sub-block indices on which calculations  C     iMin, iMax     - Ranges and sub-block indices on which calculations
136  C     jMin, jMax       are applied.  C     jMin, jMax       are applied.
137  C     bi, bj  C     bi, bj
# Line 129  C                      index into fVerTe Line 143  C                      index into fVerTe
143        _RL uTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL uTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
144        _RL vTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL vTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
145        _RL rTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
146          _RL rTransKp1(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
147        _RS maskUp  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RS maskUp  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
148        _RL fVerT   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerT   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
149        _RL fVerS   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerS   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
150        _RL fVerTr1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)  #ifdef ALLOW_PTRACERS
151        _RL phiHyd  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL fVerP   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2,PTRACERS_num)
152        _RL rhokm1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  #endif
       _RL rhok    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
       _RL phiSurfX(1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
       _RL phiSurfY(1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
153        _RL KappaRT (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL KappaRT (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
154        _RL KappaRS (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL KappaRS (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
155        _RL sigmaX  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaX  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
156        _RL sigmaY  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaY  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
157        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
158  C     This is currently used by IVDC and Diagnostics        _RL kp1Msk
159        _RL ConvectCount (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        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, ip
166    
167  CEOP  CEOP
168    
169    #ifdef ALLOW_DEBUG
170             IF ( debugLevel .GE. debLevB )
171         &    CALL DEBUG_ENTER('THERMODYNAMICS',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            rTrans (i,j)   = 0. _d 0            rTrans (i,j)   = 0. _d 0
223              rTransKp1(i,j) = 0. _d 0
224            fVerT  (i,j,1) = 0. _d 0            fVerT  (i,j,1) = 0. _d 0
225            fVerT  (i,j,2) = 0. _d 0            fVerT  (i,j,2) = 0. _d 0
226            fVerS  (i,j,1) = 0. _d 0            fVerS  (i,j,1) = 0. _d 0
227            fVerS  (i,j,2) = 0. _d 0            fVerS  (i,j,2) = 0. _d 0
228            fVerTr1(i,j,1) = 0. _d 0  #ifdef ALLOW_PTRACERS
229            fVerTr1(i,j,2) = 0. _d 0            DO ip=1,PTRACERS_num
230               fVerP  (i,j,1,ip) = 0. _d 0
231               fVerP  (i,j,2,ip) = 0. _d 0
232              ENDDO
233    #endif
234           ENDDO           ENDDO
235          ENDDO          ENDDO
236    
# Line 230  C--     Set up work arrays that need val Line 238  C--     Set up work arrays that need val
238           DO j=1-OLy,sNy+OLy           DO j=1-OLy,sNy+OLy
239            DO i=1-OLx,sNx+OLx            DO i=1-OLx,sNx+OLx
240  C This is currently also used by IVDC and Diagnostics  C This is currently also used by IVDC and Diagnostics
241             ConvectCount(i,j,k) = 0.             KappaRT(i,j,k)    = 0. _d 0
242             KappaRT(i,j,k) = 0. _d 0             KappaRS(i,j,k)    = 0. _d 0
243             KappaRS(i,j,k) = 0. _d 0  C- tracer tendency needs to be set to zero (moved here from gad_calc_rhs):
244  #ifdef ALLOW_AUTODIFF_TAMC             gT(i,j,k,bi,bj)   = 0. _d 0
245             gT(i,j,k,bi,bj) = 0. _d 0             gS(i,j,k,bi,bj)   = 0. _d 0
246             gS(i,j,k,bi,bj) = 0. _d 0  # ifdef ALLOW_PTRACERS
247  #ifdef ALLOW_PASSIVE_TRACER  ceh3 this should have an   IF ( usePTRACERS ) THEN
248             gTr1(i,j,k,bi,bj) = 0. _d 0             DO iTracer=1,PTRACERS_numInUse
249  #endif              gPTr(i,j,k,bi,bj,itracer) = 0. _d 0
250  #endif             ENDDO
251    # endif
252            ENDDO            ENDDO
253           ENDDO           ENDDO
254          ENDDO          ENDDO
255    
256          iMin = 1-OLx+1  c       iMin = 1-OLx
257          iMax = sNx+OLx  c       iMax = sNx+OLx
258          jMin = 1-OLy+1  c       jMin = 1-OLy
259          jMax = sNy+OLy  c       jMax = sNy+OLy
   
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  
 CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  
 #ifdef ALLOW_KPP  
 CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  
 CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  
 #endif  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
 C--     Start of diagnostic loop  
         DO k=Nr,1,-1  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 C? Patrick, is this formula correct now that we change the loop range?  
 C? Do we still need this?  
 cph kkey formula corrected.  
 cph Needed for rhok, rhokm1, in the case useGMREDI.  
          kkey = (ikey-1)*Nr + k  
 CADJ STORE rhokm1(:,:) = comlev1_bibj_k ,       key=kkey, byte=isbyte  
 CADJ STORE rhok  (:,:) = comlev1_bibj_k ,       key=kkey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
 C--       Integrate continuity vertically for vertical velocity  
           CALL INTEGRATE_FOR_W(  
      I                         bi, bj, k, uVel, vVel,  
      O                         wVel,  
      I                         myThid )  
   
 #ifdef    ALLOW_OBCS  
 #ifdef    ALLOW_NONHYDROSTATIC  
 C--       Apply OBC to W if in N-H mode  
           IF (useOBCS.AND.nonHydrostatic) THEN  
             CALL OBCS_APPLY_W( bi, bj, k, wVel, myThid )  
           ENDIF  
 #endif    /* ALLOW_NONHYDROSTATIC */  
 #endif    /* ALLOW_OBCS */  
   
 C--       Calculate gradients of potential density for isoneutral  
 C         slope terms (e.g. GM/Redi tensor or IVDC diffusivity)  
 c         IF ( k.GT.1 .AND. (useGMRedi.OR.ivdc_kappa.NE.0.) ) THEN  
           IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.) ) THEN  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  
 CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
             CALL FIND_RHO(  
      I        bi, bj, iMin, iMax, jMin, jMax, k, k, eosType,  
      I        theta, salt,  
      O        rhoK,  
      I        myThid )  
             IF (k.GT.1) THEN  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE theta(:,:,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  
 #endif /* ALLOW_AUTODIFF_TAMC */  
              CALL FIND_RHO(  
      I        bi, bj, iMin, iMax, jMin, jMax, k-1, k, eosType,  
      I        theta, salt,  
      O        rhoKm1,  
      I        myThid )  
             ENDIF  
             CALL GRAD_SIGMA(  
      I             bi, bj, iMin, iMax, jMin, jMax, k,  
      I             rhoK, rhoKm1, rhoK,  
      O             sigmaX, sigmaY, sigmaR,  
      I             myThid )  
           ENDIF  
   
 C--       Implicit Vertical Diffusion for Convection  
 c ==> should use sigmaR !!!  
           IF (k.GT.1 .AND. ivdc_kappa.NE.0.) THEN  
             CALL CALC_IVDC(  
      I        bi, bj, iMin, iMax, jMin, jMax, k,  
      I        rhoKm1, rhoK,  
      U        ConvectCount, KappaRT, KappaRS,  
      I        myTime, myIter, myThid)  
           ENDIF  
   
 C--     end of diagnostic k loop (Nr:1)  
         ENDDO  
260    
261  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
262  cph avoids recomputation of integrate_for_w  cph avoids recomputation of integrate_for_w
263  CADJ STORE wvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE wvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
264  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
265    
266  #ifdef  ALLOW_OBCS  C--     Attention: by defining "SINGLE_LAYER_MODE" in CPP_OPTIONS.h
267  C--     Calculate future values on open boundaries  C--     MOST of THERMODYNAMICS will be disabled
268          IF (useOBCS) THEN  #ifndef SINGLE_LAYER_MODE
           CALL OBCS_CALC( bi, bj, myTime+deltaT, myIter+1,  
      I            uVel, vVel, wVel, theta, salt,  
      I            myThid )  
         ENDIF  
 #endif  /* ALLOW_OBCS */  
   
 C--     Determines forcing terms based on external fields  
 C       relaxation terms, etc.  
         CALL EXTERNAL_FORCING_SURF(  
      I             bi, bj, iMin, iMax, jMin, jMax,  
      I             myThid )  
 #ifdef ALLOW_AUTODIFF_TAMC  
 cph needed for KPP  
 CADJ STORE surfacetendencyU(:,:,bi,bj)  
 CADJ &     = comlev1_bibj, key=ikey, byte=isbyte  
 CADJ STORE surfacetendencyV(:,:,bi,bj)  
 CADJ &     = comlev1_bibj, key=ikey, byte=isbyte  
 CADJ STORE surfacetendencyS(:,:,bi,bj)  
 CADJ &     = comlev1_bibj, key=ikey, byte=isbyte  
 CADJ STORE surfacetendencyT(:,:,bi,bj)  
 CADJ &     = comlev1_bibj, key=ikey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
 #ifdef  ALLOW_GMREDI  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE sigmaX(:,:,:) = comlev1, key=ikey, byte=isbyte  
 CADJ STORE sigmaY(:,:,:) = comlev1, key=ikey, byte=isbyte  
 CADJ STORE sigmaR(:,:,:) = comlev1, key=ikey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
 C--     Calculate iso-neutral slopes for the GM/Redi parameterisation  
         IF (useGMRedi) THEN  
           CALL GMREDI_CALC_TENSOR(  
      I             bi, bj, iMin, iMax, jMin, jMax,  
      I             sigmaX, sigmaY, sigmaR,  
      I             myThid )  
 #ifdef ALLOW_AUTODIFF_TAMC  
         ELSE  
           CALL GMREDI_CALC_TENSOR_DUMMY(  
      I             bi, bj, iMin, iMax, jMin, jMax,  
      I             sigmaX, sigmaY, sigmaR,  
      I             myThid )  
 #endif /* ALLOW_AUTODIFF_TAMC */  
         ENDIF  
269    
270  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
271  CADJ STORE Kwx(:,:,:,bi,bj)   = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
272  CADJ STORE Kwy(:,:,:,bi,bj)   = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
273  CADJ STORE Kwz(:,:,:,bi,bj)   = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
274  #endif /* ALLOW_AUTODIFF_TAMC */  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
275    #ifdef ALLOW_PTRACERS
276  #endif  /* ALLOW_GMREDI */  cph-- moved to forward_step to avoid key computation
277    cphCADJ STORE ptracer(:,:,:,bi,bj,itracer) = comlev1_bibj,
278  #ifdef  ALLOW_KPP  cphCADJ &                              key=itdkey, byte=isbyte
 C--     Compute KPP mixing coefficients  
         IF (useKPP) THEN  
           CALL KPP_CALC(  
      I                  bi, bj, myTime, myThid )  
 #ifdef ALLOW_AUTODIFF_TAMC  
         ELSE  
           CALL KPP_CALC_DUMMY(  
      I                  bi, bj, myTime, myThid )  
 #endif /* ALLOW_AUTODIFF_TAMC */  
         ENDIF  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE KPPghat   (:,:,:,bi,bj)  
 CADJ &   , KPPdiffKzT(:,:,:,bi,bj)  
 CADJ &   , KPPdiffKzS(:,:,:,bi,bj)  
 CADJ &   , KPPfrac   (:,:  ,bi,bj)  
 CADJ &                 = comlev1_bibj, key=ikey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
 #endif  /* ALLOW_KPP */  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE KappaRT(:,:,:)     = comlev1_bibj, key=ikey, byte=isbyte  
 CADJ STORE KappaRS(:,:,:)     = comlev1_bibj, key=ikey, byte=isbyte  
 CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  
 CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  
 CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  
 CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  
 #ifdef ALLOW_PASSIVE_TRACER  
 CADJ STORE tr1  (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte  
279  #endif  #endif
280  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
281    
 #ifdef ALLOW_AIM  
 C       AIM - atmospheric intermediate model, physics package code.  
 C note(jmc) : phiHyd=0 at this point but is not really used in Molteni Physics  
         IF ( useAIM ) THEN  
          CALL TIMER_START('AIM_DO_ATMOS_PHYS      [DYNAMICS]', myThid)  
          CALL AIM_DO_ATMOS_PHYSICS( phiHyd, bi, bj, myTime, myThid )  
          CALL TIMER_STOP ('AIM_DO_ATMOS_PHYS      [DYNAMICS]', myThid)  
         ENDIF  
 #endif /* ALLOW_AIM */  
   
 #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 */  
   
282  #ifndef DISABLE_MULTIDIM_ADVECTION  #ifndef DISABLE_MULTIDIM_ADVECTION
283  C--     Some advection schemes are better calculated using a multi-dimensional  C--     Some advection schemes are better calculated using a multi-dimensional
284  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 454  C to be able to exclude this scheme to a Line 290  C to be able to exclude this scheme to a
290  C recomputation. It *is* differentiable, if you need it.  C recomputation. It *is* differentiable, if you need it.
291  C Edit GAD_OPTIONS.h and #define DISABLE_MULTIDIM_ADVECTION to  C Edit GAD_OPTIONS.h and #define DISABLE_MULTIDIM_ADVECTION to
292  C disable this section of code.  C disable this section of code.
293          IF (multiDimAdvection) THEN  #ifndef ALLOW_OFFLINE
294           IF (tempStepping .AND.          IF (tempMultiDimAdvec) THEN
295       &       tempAdvScheme.NE.ENUM_CENTERED_2ND .AND.  #ifdef ALLOW_DEBUG
296       &       tempAdvScheme.NE.ENUM_UPWIND_3RD .AND.            IF ( debugLevel .GE. debLevB )
297       &       tempAdvScheme.NE.ENUM_CENTERED_4TH ) THEN       &     CALL DEBUG_CALL('GAD_ADVECTION',myThid)
298            CALL GAD_ADVECTION(bi,bj,tempAdvScheme,GAD_TEMPERATURE,  #endif
299       U                      theta,gT,            CALL GAD_ADVECTION(
300       I                      myTime,myIter,myThid)       I             tempImplVertAdv, tempAdvScheme, tempVertAdvScheme,
301           ENDIF       I             GAD_TEMPERATURE,
302           IF (saltStepping .AND.       I             uVel, vVel, wVel, theta,
303       &       saltAdvScheme.NE.ENUM_CENTERED_2ND .AND.       O             gT,
304       &       saltAdvScheme.NE.ENUM_UPWIND_3RD .AND.       I             bi,bj,myTime,myIter,myThid)
305       &       saltAdvScheme.NE.ENUM_CENTERED_4TH ) THEN          ENDIF
306            CALL GAD_ADVECTION(bi,bj,saltAdvScheme,GAD_SALINITY,  #endif
307       U                      salt,gS,  #ifndef ALLOW_OFFLINE
308       I                      myTime,myIter,myThid)          IF (saltMultiDimAdvec) THEN
309           ENDIF  #ifdef ALLOW_DEBUG
310              IF ( debugLevel .GE. debLevB )
311         &     CALL DEBUG_CALL('GAD_ADVECTION',myThid)
312    #endif
313              CALL GAD_ADVECTION(
314         I             saltImplVertAdv, saltAdvScheme, saltVertAdvScheme,
315         I             GAD_SALINITY,
316         I             uVel, vVel, wVel, salt,
317         O             gS,
318         I             bi,bj,myTime,myIter,myThid)
319          ENDIF          ENDIF
320    #endif
321  C Since passive tracers are configurable separately from T,S we  C Since passive tracers are configurable separately from T,S we
322  C call the multi-dimensional method for PTRACERS regardless  C call the multi-dimensional method for PTRACERS regardless
323  C of whether multiDimAdvection is set or not.  C of whether multiDimAdvection is set or not.
324  #ifdef ALLOW_PTRACERS  #ifdef ALLOW_PTRACERS
325          IF ( usePTRACERS ) THEN          IF ( usePTRACERS ) THEN
326    #ifdef ALLOW_DEBUG
327              IF ( debugLevel .GE. debLevB )
328         &     CALL DEBUG_CALL('PTRACERS_ADVECTION',myThid)
329    #endif
330           CALL PTRACERS_ADVECTION( bi,bj,myIter,myTime,myThid )           CALL PTRACERS_ADVECTION( bi,bj,myIter,myTime,myThid )
331          ENDIF          ENDIF
332  #endif /* ALLOW_PTRACERS */  #endif /* ALLOW_PTRACERS */
333  #endif /* DISABLE_MULTIDIM_ADVECTION */  #endif /* DISABLE_MULTIDIM_ADVECTION */
334    
335    #ifdef ALLOW_DEBUG
336            IF ( debugLevel .GE. debLevB )
337         &    CALL DEBUG_MSG('ENTERING DOWNWARD K LOOP',myThid)
338    #endif
339    
340  C--     Start of thermodynamics loop  C--     Start of thermodynamics loop
341          DO k=Nr,1,-1          DO k=Nr,1,-1
342  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
343  C? Patrick Is this formula correct?  C? Patrick Is this formula correct?
344  cph Yes, but I rewrote it.  cph Yes, but I rewrote it.
345  cph Also, the KappaR? need the index and subscript k!  cph Also, the KappaR? need the index and subscript k!
346           kkey = (ikey-1)*Nr + k           kkey = (itdkey-1)*Nr + k
347  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
348    
349  C--       km1    Points to level above k (=k-1)  C--       km1    Points to level above k (=k-1)
# Line 504  C--       kDown  Cycles through 2,1 to p Line 359  C--       kDown  Cycles through 2,1 to p
359            jMin = 1-OLy            jMin = 1-OLy
360            jMax = sNy+OLy            jMax = sNy+OLy
361    
362              kp1Msk=1.
363              IF (k.EQ.Nr) kp1Msk=0.
364               DO j=1-Oly,sNy+Oly
365                DO i=1-Olx,sNx+Olx
366                 rTransKp1(i,j) = kp1Msk*rTrans(i,j)
367                ENDDO
368               ENDDO
369    #ifdef ALLOW_AUTODIFF_TAMC
370    CADJ STORE rTransKp1(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
371    #endif
372    
373  C--       Get temporary terms used by tendency routines  C--       Get temporary terms used by tendency routines
374            CALL CALC_COMMON_FACTORS (            CALL CALC_COMMON_FACTORS (
375       I         bi,bj,iMin,iMax,jMin,jMax,k,       I         bi,bj,iMin,iMax,jMin,jMax,k,
376       O         xA,yA,uTrans,vTrans,rTrans,maskUp,       O         xA,yA,uTrans,vTrans,rTrans,maskUp,
377       I         myThid)       I         myThid)
378    
379  #ifdef ALLOW_AUTODIFF_TAMC            IF (k.EQ.1) THEN
380  CADJ STORE KappaRT(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte  C- Surface interface :
381  CADJ STORE KappaRS(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte             DO j=1-Oly,sNy+Oly
382                DO i=1-Olx,sNx+Olx
383                 rTrans(i,j) = 0.
384                ENDDO
385               ENDDO
386              ELSE
387    C- Interior interface :
388               DO j=1-Oly,sNy+Oly
389                DO i=1-Olx,sNx+Olx
390                 rTrans(i,j) = rTrans(i,j)*maskC(i,j,k-1,bi,bj)
391                ENDDO
392               ENDDO
393              ENDIF
394    
395    #ifdef ALLOW_GMREDI
396    
397    C--   Residual transp = Bolus transp + Eulerian transp
398              IF (useGMRedi) THEN
399                CALL GMREDI_CALC_UVFLOW(
400         &            uTrans, vTrans, bi, bj, k, myThid)
401                IF (K.GE.2) CALL GMREDI_CALC_WFLOW(
402         &                    rTrans, bi, bj, k, myThid)
403              ENDIF
404    
405    #ifdef ALLOW_AUTODIFF_TAMC
406    CADJ STORE rTrans(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
407    #ifdef GM_BOLUS_ADVEC
408    CADJ STORE uTrans(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
409    CADJ STORE vTrans(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
410    #endif
411  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
412    
413    #endif /* ALLOW_GMREDI */
414    
415  #ifdef  INCLUDE_CALC_DIFFUSIVITY_CALL  #ifdef  INCLUDE_CALC_DIFFUSIVITY_CALL
416  C--      Calculate the total vertical diffusivity  C--      Calculate the total vertical diffusivity
417           CALL CALC_DIFFUSIVITY(           CALL CALC_DIFFUSIVITY(
# Line 522  C--      Calculate the total vertical di Line 419  C--      Calculate the total vertical di
419       I        maskUp,       I        maskUp,
420       O        KappaRT,KappaRS,       O        KappaRT,KappaRS,
421       I        myThid)       I        myThid)
422    # ifdef ALLOW_AUTODIFF_TAMC
423    CADJ STORE KappaRT(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte
424    CADJ STORE KappaRS(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte
425    # endif /* ALLOW_AUTODIFF_TAMC */
426  #endif  #endif
427    
428            iMin = 1-OLx+2            iMin = 1-OLx+2
# Line 531  C--      Calculate the total vertical di Line 432  C--      Calculate the total vertical di
432    
433  C--      Calculate active tracer tendencies (gT,gS,...)  C--      Calculate active tracer tendencies (gT,gS,...)
434  C        and step forward storing result in gTnm1, gSnm1, etc.  C        and step forward storing result in gTnm1, gSnm1, etc.
435    #ifndef ALLOW_OFFLINE
436           IF ( tempStepping ) THEN           IF ( tempStepping ) THEN
437             CALL CALC_GT(             CALL CALC_GT(
438       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
439       I         xA,yA,uTrans,vTrans,rTrans,maskUp,       I         xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
440       I         KappaRT,       I         KappaRT,
441       U         fVerT,       U         fVerT,
442       I         myTime,myIter,myThid)       I         myTime,myIter,myThid)
# Line 543  C        and step forward storing result Line 445  C        and step forward storing result
445       I         theta, gT,       I         theta, gT,
446       I         myIter, myThid)       I         myIter, myThid)
447           ENDIF           ENDIF
448    #endif
449    
450    #ifndef ALLOW_OFFLINE
451           IF ( saltStepping ) THEN           IF ( saltStepping ) THEN
452             CALL CALC_GS(             CALL CALC_GS(
453       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
454       I         xA,yA,uTrans,vTrans,rTrans,maskUp,       I         xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
455       I         KappaRS,       I         KappaRS,
456       U         fVerS,       U         fVerS,
457       I         myTime,myIter,myThid)       I         myTime,myIter,myThid)
# Line 555  C        and step forward storing result Line 460  C        and step forward storing result
460       I         salt, gS,       I         salt, gS,
461       I         myIter, myThid)       I         myIter, myThid)
462           ENDIF           ENDIF
 #ifdef ALLOW_PASSIVE_TRACER  
          IF ( tr1Stepping ) THEN  
            CALL CALC_GTR1(  
      I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,  
      I         xA,yA,uTrans,vTrans,rTrans,maskUp,  
      I         KappaRT,  
      U         fVerTr1,  
      I         myTime,myIter,myThid)  
            CALL TIMESTEP_TRACER(  
      I         bi,bj,iMin,iMax,jMin,jMax,k,tracerAdvScheme,  
      I         Tr1, gTr1,  
      I         myIter,myThid)  
          ENDIF  
463  #endif  #endif
464  #ifdef ALLOW_PTRACERS  #ifdef ALLOW_PTRACERS
465           IF ( usePTRACERS ) THEN           IF ( usePTRACERS ) THEN
466             CALL PTRACERS_INTEGERATE(             CALL PTRACERS_INTEGRATE(
467       I         bi,bj,k,       I         bi,bj,k,
468       I         xA,yA,uTrans,vTrans,rTrans,maskUp,       I         xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
469       X         KappaRS,       X         fVerP, KappaRS,
470       I         myIter,myTime,myThid)       I         myIter,myTime,myThid)
471           ENDIF           ENDIF
472  #endif /* ALLOW_PTRACERS */  #endif /* ALLOW_PTRACERS */
# Line 587  C--      Apply open boundary conditions Line 479  C--      Apply open boundary conditions
479  #endif   /* ALLOW_OBCS */  #endif   /* ALLOW_OBCS */
480    
481  C--      Freeze water  C--      Freeze water
482           IF (allowFreezing) THEN  C  this bit of code is left here for backward compatibility.
483    C  freezing at surface level has been moved to FORWARD_STEP
484    #ifndef ALLOW_OFFLINE
485             IF ( useOldFreezing .AND. .NOT. useSEAICE
486         &       .AND. .NOT.(useThSIce.AND.k.EQ.1) ) THEN
487  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
488  CADJ STORE gT(:,:,k,bi,bj) = comlev1_bibj_k  CADJ STORE gT(:,:,k,bi,bj) = comlev1_bibj_k
489  CADJ &   , key = kkey, byte = isbyte  CADJ &   , key = kkey, byte = isbyte
490  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
491              CALL FREEZE( bi, bj, iMin, iMax, jMin, jMax, k, myThid )              CALL FREEZE( bi, bj, iMin, iMax, jMin, jMax, k, myThid )
492           END IF           ENDIF
493    #endif
494    
495  C--     end of thermodynamic k loop (Nr:1)  C--     end of thermodynamic k loop (Nr:1)
496          ENDDO          ENDDO
497    
498    
499    C--     Implicit vertical advection & diffusion
500    #ifndef ALLOW_OFFLINE
501    #ifdef INCLUDE_IMPLVERTADV_CODE
502            IF ( tempImplVertAdv ) THEN
503              CALL GAD_IMPLICIT_R(
504         I         tempImplVertAdv, tempAdvScheme, GAD_TEMPERATURE,
505         I         kappaRT, wVel, theta,
506         U         gT,
507         I         bi, bj, myTime, myIter, myThid )
508            ELSEIF ( tempStepping .AND. implicitDiffusion ) THEN
509    #else /* INCLUDE_IMPLVERTADV_CODE */
510            IF     ( tempStepping .AND. implicitDiffusion ) THEN
511    #endif /* INCLUDE_IMPLVERTADV_CODE */
512  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
513  C? Patrick? What about this one?  CADJ STORE KappaRT(:,:,:) = comlev1_bibj , key=itdkey, byte=isbyte
514  cph Keys iikey and idkey dont seem to be needed  CADJ STORE gT(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte
 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 */  
   
 C--     Implicit diffusion  
         IF (implicitDiffusion) THEN  
   
          IF (tempStepping) THEN  
 #ifdef ALLOW_AUTODIFF_TAMC  
             idkey = iikey + 1  
 CADJ STORE gT(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte  
515  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
516              CALL IMPLDIFF(            CALL IMPLDIFF(
517       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
518       I         deltaTtracer, KappaRT, recip_HFacC,       I         deltaTtracer, KappaRT, recip_HFacC,
519       U         gT,       U         gT,
520       I         myThid )       I         myThid )
521           ENDIF          ENDIF
522    #endif
523    
524           IF (saltStepping) THEN  #ifndef ALLOW_OFFLINE
525    #ifdef INCLUDE_IMPLVERTADV_CODE
526            IF ( saltImplVertAdv ) THEN
527              CALL GAD_IMPLICIT_R(
528         I         saltImplVertAdv, saltAdvScheme, GAD_SALINITY,
529         I         kappaRS, wVel, salt,
530         U         gS,
531         I         bi, bj, myTime, myIter, myThid )
532            ELSEIF ( saltStepping .AND. implicitDiffusion ) THEN
533    #else /* INCLUDE_IMPLVERTADV_CODE */
534            IF     ( saltStepping .AND. implicitDiffusion ) THEN
535    #endif /* INCLUDE_IMPLVERTADV_CODE */
536  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
537           idkey = iikey + 2  CADJ STORE KappaRS(:,:,:) = comlev1_bibj , key=itdkey, byte=isbyte
538  CADJ STORE gS(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte  CADJ STORE gS(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte
539  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
540              CALL IMPLDIFF(            CALL IMPLDIFF(
541       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
542       I         deltaTtracer, KappaRS, recip_HFacC,       I         deltaTtracer, KappaRS, recip_HFacC,
543       U         gS,       U         gS,
544       I         myThid )       I         myThid )
545           ENDIF          ENDIF
   
 #ifdef ALLOW_PASSIVE_TRACER  
          IF (tr1Stepping) THEN  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
           CALL IMPLDIFF(  
      I      bi, bj, iMin, iMax, jMin, jMax,  
      I      deltaTtracer, KappaRT, recip_HFacC,  
      U      gTr1,  
      I      myThid )  
          ENDIF  
546  #endif  #endif
547    
548  #ifdef ALLOW_PTRACERS  #ifdef ALLOW_PTRACERS
549  C Vertical diffusion (implicit) for passive tracers  c #ifdef INCLUDE_IMPLVERTADV_CODE
550           IF ( usePTRACERS ) THEN  c       IF ( usePTRACERS .AND. ptracerImplVertAdv ) THEN
551    c       ELSEIF ( usePTRACERS .AND. implicitDiffusion ) THEN
552    c #else
553            IF     ( usePTRACERS .AND. implicitDiffusion ) THEN
554    C--     Vertical diffusion (implicit) for passive tracers
555             CALL PTRACERS_IMPLDIFF( bi,bj,KappaRS,myThid )             CALL PTRACERS_IMPLDIFF( bi,bj,KappaRS,myThid )
556           ENDIF          ENDIF
557  #endif /* ALLOW_PTRACERS */  #endif /* ALLOW_PTRACERS */
558    
559  #ifdef   ALLOW_OBCS  #ifdef   ALLOW_OBCS
560  C--      Apply open boundary conditions  C--      Apply open boundary conditions
561           IF (useOBCS) THEN          IF ( ( implicitDiffusion
562         &    .OR. tempImplVertAdv
563         &    .OR. saltImplVertAdv
564         &       ) .AND. useOBCS     ) THEN
565             DO K=1,Nr             DO K=1,Nr
566               CALL OBCS_APPLY_TS( bi, bj, k, gT, gS, myThid )               CALL OBCS_APPLY_TS( bi, bj, k, gT, gS, myThid )
567             ENDDO             ENDDO
568           END IF          ENDIF
569  #endif   /* ALLOW_OBCS */  #endif   /* ALLOW_OBCS */
570    
571  C--     End If implicitDiffusion  #ifdef ALLOW_TIMEAVE
572            IF ( taveFreq.GT. 0. _d 0 .AND. fluidIsWater ) THEN
573              CALL TIMEAVE_SURF_FLUX( bi, bj, myTime, myIter, myThid)
574            ENDIF
575    #ifndef HRCUBE
576            IF (taveFreq.GT.0. .AND. ivdc_kappa.NE.0.) THEN
577              CALL TIMEAVE_CUMULATE(ConvectCountTave, IVDConvCount,
578         I                           Nr, deltaTclock, bi, bj, myThid)
579          ENDIF          ENDIF
580            useVariableK = useKPP .OR. usePP81 .OR. useMY82 .OR. useGGL90
581         &       .OR. useGMredi .OR. ivdc_kappa.NE.0.
582            IF (taveFreq.GT.0. .AND. useVariableK ) THEN
583             IF (implicitDiffusion) THEN
584              CALL TIMEAVE_CUMUL_DIF_1T(TdiffRtave, gT, kappaRT,
585         I                        Nr, 3, deltaTclock, bi, bj, myThid)
586             ELSE
587              CALL TIMEAVE_CUMUL_DIF_1T(TdiffRtave, theta, kappaRT,
588         I                        Nr, 3, deltaTclock, bi, bj, myThid)
589             ENDIF
590            ENDIF
591    #endif /* ndef HRCUBE */
592    #endif /* ALLOW_TIMEAVE */
593    
594  Ccs-  #endif /* SINGLE_LAYER_MODE */
595    
596    C--   end bi,bj loops.
597         ENDDO         ENDDO
598        ENDDO        ENDDO
599    
600  #ifdef ALLOW_AIM  #ifdef ALLOW_DEBUG
       IF ( useAIM ) THEN  
        CALL AIM_AIM2DYN_EXCHANGES( myTime, myThid )  
       ENDIF  
        _EXCH_XYZ_R8(gT,myThid)  
        _EXCH_XYZ_R8(gS,myThid)  
 #else  
       IF (staggerTimeStep.AND.useCubedSphereExchange) THEN  
        _EXCH_XYZ_R8(gT,myThid)  
        _EXCH_XYZ_R8(gS,myThid)  
       ENDIF  
 #endif /* ALLOW_AIM */  
   
 #ifndef DISABLE_DEBUGMODE  
601        If (debugMode) THEN        If (debugMode) THEN
602         CALL DEBUG_STATS_RL(Nr,uVel,'Uvel (THERMODYNAMICS)',myThid)         CALL DEBUG_STATS_RL(Nr,uVel,'Uvel (THERMODYNAMICS)',myThid)
603         CALL DEBUG_STATS_RL(Nr,vVel,'Vvel (THERMODYNAMICS)',myThid)         CALL DEBUG_STATS_RL(Nr,vVel,'Vvel (THERMODYNAMICS)',myThid)
# Line 704  Ccs- Line 616  Ccs-
616        ENDIF        ENDIF
617  #endif  #endif
618    
619    #ifdef ALLOW_DEBUG
620             IF ( debugLevel .GE. debLevB )
621         &    CALL DEBUG_LEAVE('THERMODYNAMICS',myThid)
622    #endif
623    
624        RETURN        RETURN
625        END        END

Legend:
Removed from v.1.18  
changed lines
  Added in v.1.79

  ViewVC Help
Powered by ViewVC 1.1.22