/[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.36 by heimbach, Tue Jan 21 19:19:45 2003 UTC revision 1.65 by dimitri, Sun Jan 25 00:31:52 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    #ifdef ALLOW_PTRACERS
7    # include "PTRACERS_OPTIONS.h"
8    #endif
9    
10  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
11  # ifdef ALLOW_GMREDI  # ifdef ALLOW_GMREDI
12  #  include "GMREDI_OPTIONS.h"  #  include "GMREDI_OPTIONS.h"
# Line 9  C $Name$ Line 14  C $Name$
14  # ifdef ALLOW_KPP  # ifdef ALLOW_KPP
15  #  include "KPP_OPTIONS.h"  #  include "KPP_OPTIONS.h"
16  # endif  # endif
 cswdice --- add ----  
 #ifdef ALLOW_THERM_SEAICE  
 #include "ICE.h"  
 #endif  
 cswdice ------  
17  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
18    
19  CBOP  CBOP
# Line 84  C     == Global variables === Line 84  C     == Global variables ===
84  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_PASSIVE_TRACER
85  #include "TR1.h"  #include "TR1.h"
86  #endif  #endif
87    #ifdef ALLOW_PTRACERS
88    #include "PTRACERS.h"
89    #endif
90    #ifdef ALLOW_TIMEAVE
91    #include "TIMEAVE_STATV.h"
92    #endif
93    
94  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
95  # include "tamc.h"  # include "tamc.h"
96  # include "tamc_keys.h"  # include "tamc_keys.h"
# Line 96  C     == Global variables === Line 103  C     == Global variables ===
103  #  include "GMREDI.h"  #  include "GMREDI.h"
104  # endif  # endif
105  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
 #ifdef ALLOW_TIMEAVE  
 #include "TIMEAVE_STATV.h"  
 #endif  
106    
107  C     !INPUT/OUTPUT PARAMETERS:  C     !INPUT/OUTPUT PARAMETERS:
108  C     == Routine arguments ==  C     == Routine arguments ==
# Line 117  C                              transport Line 121  C                              transport
121  C                              o uTrans: Zonal transport  C                              o uTrans: Zonal transport
122  C                              o vTrans: Meridional transport  C                              o vTrans: Meridional transport
123  C                              o rTrans: Vertical transport  C                              o rTrans: Vertical transport
124    C     rTransKp1                o vertical volume transp. at interface k+1
125  C     maskUp                   o maskUp: land/water mask for W points  C     maskUp                   o maskUp: land/water mask for W points
126  C     fVer[STUV]               o fVer: Vertical flux term - note fVer  C     fVer[STUV]               o fVer: Vertical flux term - note fVer
127  C                                      is "pipelined" in the vertical  C                                      is "pipelined" in the vertical
128  C                                      so we need an fVer for each  C                                      so we need an fVer for each
129  C                                      variable.  C                                      variable.
130  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.  
131  C     phiSurfX, - gradient of Surface potentiel (Pressure/rho, ocean)  C     phiSurfX, - gradient of Surface potentiel (Pressure/rho, ocean)
132  C     phiSurfY             or geopotentiel (atmos) in X and Y direction  C     phiSurfY             or geopotentiel (atmos) in X and Y direction
133  C     KappaRT,       - Total diffusion in vertical for T and S.  C     KappaRT,       - Total diffusion in vertical for T and S.
134  C     KappaRS          (background + spatially varying, isopycnal term).  C     KappaRS          (background + spatially varying, isopycnal term).
135    C     useVariableK   = T when vertical diffusion is not constant
136  C     iMin, iMax     - Ranges and sub-block indices on which calculations  C     iMin, iMax     - Ranges and sub-block indices on which calculations
137  C     jMin, jMax       are applied.  C     jMin, jMax       are applied.
138  C     bi, bj  C     bi, bj
# Line 143  C                      index into fVerTe Line 144  C                      index into fVerTe
144        _RL uTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL uTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
145        _RL vTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL vTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
146        _RL rTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
147          _RL rTransKp1(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
148        _RS maskUp  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RS maskUp  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
149        _RL fVerT   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerT   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
150        _RL fVerS   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerS   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
151    #ifdef ALLOW_PASSIVE_TRACER
152        _RL fVerTr1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerTr1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
153        _RL phiHyd  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)  #endif
154    #ifdef ALLOW_PTRACERS
155          _RL fVerP   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2,PTRACERS_num)
156    #endif
157        _RL rhokm1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rhokm1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
158        _RL rhok    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rhok    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
159        _RL phiSurfX(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL phiSurfX(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
# Line 159  C                      index into fVerTe Line 165  C                      index into fVerTe
165        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
166  C     This is currently used by IVDC and Diagnostics  C     This is currently used by IVDC and Diagnostics
167        _RL ConvectCount (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL ConvectCount (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
168          _RL kp1Msk
169          LOGICAL useVariableK
170        INTEGER iMin, iMax        INTEGER iMin, iMax
171        INTEGER jMin, jMax        INTEGER jMin, jMax
172        INTEGER bi, bj        INTEGER bi, bj
173        INTEGER i, j        INTEGER i, j
174        INTEGER k, km1, kup, kDown        INTEGER k, km1, kup, kDown
175          INTEGER iTracer, ip
176    
177  CEOP  CEOP
178    
179    #ifdef ALLOW_DEBUG
180             IF ( debugLevel .GE. debLevB )
181         &    CALL DEBUG_ENTER('THERMODYNAMICS',myThid)
182    #endif
183    
184  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
185  C--   dummy statement to end declaration part  C--   dummy statement to end declaration part
# Line 173  C--   dummy statement to end declaration Line 187  C--   dummy statement to end declaration
187        itdkey = 1        itdkey = 1
188  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
189    
 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  
         rhok   (i,j) = 0. _d 0  
         phiSurfX(i,j) = 0. _d 0  
         phiSurfY(i,j) = 0. _d 0  
        ENDDO  
       ENDDO  
   
   
190  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
191  C--   HPF directive to help TAMC  C--   HPF directive to help TAMC
192  CHPF$ INDEPENDENT  CHPF$ INDEPENDENT
# Line 201  CHPF$ INDEPENDENT Line 197  CHPF$ INDEPENDENT
197  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
198  C--    HPF directive to help TAMC  C--    HPF directive to help TAMC
199  CHPF$  INDEPENDENT, NEW (rTrans,fVerT,fVerS  CHPF$  INDEPENDENT, NEW (rTrans,fVerT,fVerS
200  CHPF$&                  ,phiHyd,utrans,vtrans,xA,yA  CHPF$&                  ,utrans,vtrans,xA,yA
201  CHPF$&                  ,KappaRT,KappaRS  CHPF$&                  ,KappaRT,KappaRS
202  CHPF$&                  )  CHPF$&                  )
203  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
# Line 221  CHPF$&                  ) Line 217  CHPF$&                  )
217       &                      + act4*max1*max2*max3       &                      + act4*max1*max2*max3
218  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
219    
220  C--     Set up work arrays that need valid initial values  C--   Set up work arrays with valid (i.e. not NaN) values
221    C     These inital values do not alter the numerical results. They
222    C     just ensure that all memory references are to valid floating
223    C     point numbers. This prevents spurious hardware signals due to
224    C     uninitialised but inert locations.
225    
226          DO j=1-OLy,sNy+OLy          DO j=1-OLy,sNy+OLy
227           DO i=1-OLx,sNx+OLx           DO i=1-OLx,sNx+OLx
228              xA(i,j)        = 0. _d 0
229              yA(i,j)        = 0. _d 0
230              uTrans(i,j)    = 0. _d 0
231              vTrans(i,j)    = 0. _d 0
232              rhok   (i,j)   = 0. _d 0
233              rhoKM1 (i,j)   = 0. _d 0
234              phiSurfX(i,j)  = 0. _d 0
235              phiSurfY(i,j)  = 0. _d 0
236            rTrans (i,j)   = 0. _d 0            rTrans (i,j)   = 0. _d 0
237              rTransKp1(i,j) = 0. _d 0
238            fVerT  (i,j,1) = 0. _d 0            fVerT  (i,j,1) = 0. _d 0
239            fVerT  (i,j,2) = 0. _d 0            fVerT  (i,j,2) = 0. _d 0
240            fVerS  (i,j,1) = 0. _d 0            fVerS  (i,j,1) = 0. _d 0
241            fVerS  (i,j,2) = 0. _d 0            fVerS  (i,j,2) = 0. _d 0
242    #ifdef ALLOW_PASSIVE_TRACER
243            fVerTr1(i,j,1) = 0. _d 0            fVerTr1(i,j,1) = 0. _d 0
244            fVerTr1(i,j,2) = 0. _d 0            fVerTr1(i,j,2) = 0. _d 0
245            rhoKM1 (i,j)   = 0. _d 0  #endif
246    #ifdef ALLOW_PTRACERS
247              DO ip=1,PTRACERS_num
248               fVerP  (i,j,1,ip) = 0. _d 0
249               fVerP  (i,j,2,ip) = 0. _d 0
250              ENDDO
251    #endif
252           ENDDO           ENDDO
253          ENDDO          ENDDO
254    
# Line 239  C--     Set up work arrays that need val Line 256  C--     Set up work arrays that need val
256           DO j=1-OLy,sNy+OLy           DO j=1-OLy,sNy+OLy
257            DO i=1-OLx,sNx+OLx            DO i=1-OLx,sNx+OLx
258  C This is currently also used by IVDC and Diagnostics  C This is currently also used by IVDC and Diagnostics
            phiHyd(i,j,k) = 0. _d 0  
259             sigmaX(i,j,k) = 0. _d 0             sigmaX(i,j,k) = 0. _d 0
260             sigmaY(i,j,k) = 0. _d 0             sigmaY(i,j,k) = 0. _d 0
261             sigmaR(i,j,k) = 0. _d 0             sigmaR(i,j,k) = 0. _d 0
262             ConvectCount(i,j,k) = 0.             ConvectCount(i,j,k) = 0.
263             KappaRT(i,j,k)    = 0. _d 0             KappaRT(i,j,k)    = 0. _d 0
264             KappaRS(i,j,k)    = 0. _d 0             KappaRS(i,j,k)    = 0. _d 0
265  #ifdef ALLOW_AUTODIFF_TAMC  C- tracer tendency needs to be set to zero (moved here from gad_calc_rhs):
 cph all the following init. are necessary for TAF  
 cph although some of these are re-initialised later.  
266             gT(i,j,k,bi,bj)   = 0. _d 0             gT(i,j,k,bi,bj)   = 0. _d 0
267             gS(i,j,k,bi,bj)   = 0. _d 0             gS(i,j,k,bi,bj)   = 0. _d 0
268  # ifdef ALLOW_PASSIVE_TRACER  # ifdef ALLOW_PASSIVE_TRACER
269    ceh3 needs an IF ( use PASSIVE_TRACER) THEN
270             gTr1(i,j,k,bi,bj) = 0. _d 0             gTr1(i,j,k,bi,bj) = 0. _d 0
271  # endif  # endif
272    # ifdef ALLOW_PTRACERS
273    ceh3 this should have an   IF ( usePTRACERS ) THEN
274               DO iTracer=1,PTRACERS_numInUse
275                gPTr(i,j,k,bi,bj,itracer) = 0. _d 0
276               ENDDO
277    # endif
278    #ifdef ALLOW_AUTODIFF_TAMC
279    cph all the following init. are necessary for TAF
280    cph although some of these are re-initialised later.
281  # ifdef ALLOW_GMREDI  # ifdef ALLOW_GMREDI
282             Kwx(i,j,k,bi,bj)  = 0. _d 0             Kwx(i,j,k,bi,bj)  = 0. _d 0
283             Kwy(i,j,k,bi,bj)  = 0. _d 0             Kwy(i,j,k,bi,bj)  = 0. _d 0
# Line 284  cph although some of these are re-initia Line 308  cph although some of these are re-initia
308          jMin = 1-OLy          jMin = 1-OLy
309          jMax = sNy+OLy          jMax = sNy+OLy
310    
   
311  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
312  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
313  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
314    CADJ STORE totphihyd
315    CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte
316  #ifdef ALLOW_KPP  #ifdef ALLOW_KPP
317  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
318  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
319  #endif  #endif
320  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
321    
322    #ifdef ALLOW_DEBUG
323            IF ( debugLevel .GE. debLevB )
324         &    CALL DEBUG_MSG('ENTERING UPWARD K LOOP',myThid)
325    #endif
326    
327  C--     Start of diagnostic loop  C--     Start of diagnostic loop
328          DO k=Nr,1,-1          DO k=Nr,1,-1
329    
# Line 328  C--       Calculate gradients of potenti Line 358  C--       Calculate gradients of potenti
358  C         slope terms (e.g. GM/Redi tensor or IVDC diffusivity)  C         slope terms (e.g. GM/Redi tensor or IVDC diffusivity)
359  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
360            IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.) ) THEN            IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.) ) THEN
361    #ifdef ALLOW_DEBUG
362                IF ( debugLevel .GE. debLevB )
363         &       CALL DEBUG_CALL('FIND_RHO',myThid)
364    #endif
365  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
366  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
367  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
 CADJ STORE pressure(:,:,k,bi,bj) =  
 CADJ &     comlev1_bibj_k, key=kkey, byte=isbyte  
368  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
369              CALL FIND_RHO(              CALL FIND_RHO(
370       I        bi, bj, iMin, iMax, jMin, jMax, k, k,       I        bi, bj, iMin, iMax, jMin, jMax, k, k,
# Line 344  CADJ &     comlev1_bibj_k, key=kkey, byt Line 376  CADJ &     comlev1_bibj_k, key=kkey, byt
376  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
377  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
378  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
 CADJ STORE pressure(:,:,k-1,bi,bj) =  
 CADJ &     comlev1_bibj_k, key=kkey, byte=isbyte  
379  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
380               CALL FIND_RHO(               CALL FIND_RHO(
381       I        bi, bj, iMin, iMax, jMin, jMax, k-1, k,       I        bi, bj, iMin, iMax, jMin, jMax, k-1, k,
# Line 353  CADJ &     comlev1_bibj_k, key=kkey, byt Line 383  CADJ &     comlev1_bibj_k, key=kkey, byt
383       O        rhoKm1,       O        rhoKm1,
384       I        myThid )       I        myThid )
385              ENDIF              ENDIF
386    #ifdef ALLOW_DEBUG
387                IF ( debugLevel .GE. debLevB )
388         &       CALL DEBUG_CALL('GRAD_SIGMA',myThid)
389    #endif
390              CALL GRAD_SIGMA(              CALL GRAD_SIGMA(
391       I             bi, bj, iMin, iMax, jMin, jMax, k,       I             bi, bj, iMin, iMax, jMin, jMax, k,
392       I             rhoK, rhoKm1, rhoK,       I             rhoK, rhoKm1, rhoK,
# Line 367  CADJ STORE rhokm1 (:,:) = comlev1_bibj_k Line 401  CADJ STORE rhokm1 (:,:) = comlev1_bibj_k
401  C--       Implicit Vertical Diffusion for Convection  C--       Implicit Vertical Diffusion for Convection
402  c ==> should use sigmaR !!!  c ==> should use sigmaR !!!
403            IF (k.GT.1 .AND. ivdc_kappa.NE.0.) THEN            IF (k.GT.1 .AND. ivdc_kappa.NE.0.) THEN
404    #ifdef ALLOW_DEBUG
405                IF ( debugLevel .GE. debLevB )
406         &       CALL DEBUG_CALL('CALC_IVDC',myThid)
407    #endif
408              CALL CALC_IVDC(              CALL CALC_IVDC(
409       I        bi, bj, iMin, iMax, jMin, jMax, k,       I        bi, bj, iMin, iMax, jMin, jMax, k,
410       I        rhoKm1, rhoK,       I        rhoKm1, rhoK,
# Line 382  C--     end of diagnostic k loop (Nr:1) Line 420  C--     end of diagnostic k loop (Nr:1)
420  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
421  cph avoids recomputation of integrate_for_w  cph avoids recomputation of integrate_for_w
422  CADJ STORE wvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE wvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
 CADJ STORE pressure (:,:,:,bi,bj) =  
 CADJ &     comlev1_bibj, key=itdkey, byte=isbyte  
423  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
424    
425  #ifdef  ALLOW_OBCS  #ifdef  ALLOW_OBCS
426  C--     Calculate future values on open boundaries  C--     Calculate future values on open boundaries
427          IF (useOBCS) THEN          IF (useOBCS) THEN
428    #ifdef ALLOW_DEBUG
429              IF ( debugLevel .GE. debLevB )
430         &     CALL DEBUG_CALL('OBCS_CALC',myThid)
431    #endif
432            CALL OBCS_CALC( bi, bj, myTime+deltaT, myIter+1,            CALL OBCS_CALC( bi, bj, myTime+deltaT, myIter+1,
433       I            uVel, vVel, wVel, theta, salt,       I            uVel, vVel, wVel, theta, salt,
434       I            myThid )       I            myThid )
435          ENDIF          ENDIF
436  #endif  /* ALLOW_OBCS */  #endif  /* ALLOW_OBCS */
437    
438            IF ( buoyancyRelation(1:7) .EQ. 'OCEANIC' ) THEN
 c********************************************  
 cswdice --- add ---  
 #ifdef ALLOW_THERM_SEAICE  
 C--     Determines forcing terms based on external fields  
 c--     including effects from ice  
         CALL ICE_FORCING(  
      I             bi, bj, iMin, iMax, jMin, jMax,  
      I             myThid )  
 #else  
   
 cswdice --- end add ---  
   
439  C--     Determines forcing terms based on external fields  C--     Determines forcing terms based on external fields
440  C       relaxation terms, etc.  C       relaxation terms, etc.
441    #ifdef ALLOW_DEBUG
442            IF ( debugLevel .GE. debLevB )
443         &    CALL DEBUG_CALL('EXTERNAL_FORCING_SURF',myThid)
444    #endif
445          CALL EXTERNAL_FORCING_SURF(          CALL EXTERNAL_FORCING_SURF(
446       I             bi, bj, iMin, iMax, jMin, jMax,       I             bi, bj, iMin, iMax, jMin, jMax,
447       I             myThid )       I             myTime, myIter, myThid )
 cswdice --- add ----  
 #endif  
 cswdice --- end add ---  
 c******************************************  
   
   
   
448    
449  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
450  cph needed for KPP  cph needed for KPP
# Line 431  CADJ STORE surfacetendencyS(:,:,bi,bj) Line 456  CADJ STORE surfacetendencyS(:,:,bi,bj)
456  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte
457  CADJ STORE surfacetendencyT(:,:,bi,bj)  CADJ STORE surfacetendencyT(:,:,bi,bj)
458  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte
459    # ifdef ALLOW_SEAICE
460    CADJ STORE surfacetendencyTice(:,:,bi,bj)
461    CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte
462    # endif
463  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
464            ENDIF
465    
466  C--     Attention: by defining "SINGLE_LAYER_MODE" in CPP_OPTIONS.h  C--     Attention: by defining "SINGLE_LAYER_MODE" in CPP_OPTIONS.h
467  C--     MOST of THERMODYNAMICS will be disabled  C--     MOST of THERMODYNAMICS will be disabled
# Line 450  CADJ STORE sigmaR(:,:,:)        = comlev Line 480  CADJ STORE sigmaR(:,:,:)        = comlev
480    
481  C--     Calculate iso-neutral slopes for the GM/Redi parameterisation  C--     Calculate iso-neutral slopes for the GM/Redi parameterisation
482          IF (useGMRedi) THEN          IF (useGMRedi) THEN
483    #ifdef ALLOW_DEBUG
484              IF ( debugLevel .GE. debLevB )
485         &     CALL DEBUG_CALL('GMREDI_CALC_TENSOR',myThid)
486    #endif
487            CALL GMREDI_CALC_TENSOR(            CALL GMREDI_CALC_TENSOR(
488       I             bi, bj, iMin, iMax, jMin, jMax,       I             bi, bj, iMin, iMax, jMin, jMax,
489       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
# Line 474  CADJ STORE Kwz(:,:,:,bi,bj)   = comlev1_ Line 508  CADJ STORE Kwz(:,:,:,bi,bj)   = comlev1_
508  #ifdef  ALLOW_KPP  #ifdef  ALLOW_KPP
509  C--     Compute KPP mixing coefficients  C--     Compute KPP mixing coefficients
510          IF (useKPP) THEN          IF (useKPP) THEN
511    #ifdef ALLOW_DEBUG
512              IF ( debugLevel .GE. debLevB )
513         &     CALL DEBUG_CALL('KPP_CALC',myThid)
514    #endif
515            CALL KPP_CALC(            CALL KPP_CALC(
516       I                  bi, bj, myTime, myThid )       I                  bi, bj, myTime, myThid )
517  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 494  CADJ &                 = comlev1_bibj, k Line 532  CADJ &                 = comlev1_bibj, k
532  #endif  /* ALLOW_KPP */  #endif  /* ALLOW_KPP */
533    
534  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
 CADJ STORE KappaRT(:,:,:)     = comlev1_bibj, key=itdkey, byte=isbyte  
 CADJ STORE KappaRS(:,:,:)     = comlev1_bibj, key=itdkey, byte=isbyte  
535  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
536  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
537  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
# Line 503  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_ Line 539  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_
539  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_PASSIVE_TRACER
540  CADJ STORE tr1  (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE tr1  (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
541  #endif  #endif
542    #ifdef ALLOW_PTRACERS
543    cph-- moved to forward_step to avoid key computation
544    cphCADJ STORE ptracer(:,:,:,bi,bj,itracer) = comlev1_bibj,
545    cphCADJ &                              key=itdkey, byte=isbyte
546    #endif
547  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
548    
549  #ifdef ALLOW_AIM  #ifdef ALLOW_AIM
550  C       AIM - atmospheric intermediate model, physics package code.  C       AIM - atmospheric intermediate model, physics package code.
551          IF ( useAIM ) THEN          IF ( useAIM ) THEN
552    #ifdef ALLOW_DEBUG
553              IF ( debugLevel .GE. debLevB )
554         &     CALL DEBUG_CALL('AIM_DO_PHYSICS',myThid)
555    #endif
556           CALL TIMER_START('AIM_DO_PHYSICS   [THERMODYNAMICS]', myThid)           CALL TIMER_START('AIM_DO_PHYSICS   [THERMODYNAMICS]', myThid)
557           CALL AIM_DO_PHYSICS( bi, bj, myTime, myIter, myThid )           CALL AIM_DO_PHYSICS( bi, bj, myTime, myIter, myThid )
558           CALL TIMER_STOP( 'AIM_DO_PHYSICS   [THERMODYNAMICS]', myThid)           CALL TIMER_STOP( 'AIM_DO_PHYSICS   [THERMODYNAMICS]', myThid)
559          ENDIF          ENDIF
560  #endif /* ALLOW_AIM */  #endif /* ALLOW_AIM */
561    
 #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 */  
   
562  #ifndef DISABLE_MULTIDIM_ADVECTION  #ifndef DISABLE_MULTIDIM_ADVECTION
563  C--     Some advection schemes are better calculated using a multi-dimensional  C--     Some advection schemes are better calculated using a multi-dimensional
564  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 533  C recomputation. It *is* differentiable, Line 571  C recomputation. It *is* differentiable,
571  C Edit GAD_OPTIONS.h and #define DISABLE_MULTIDIM_ADVECTION to  C Edit GAD_OPTIONS.h and #define DISABLE_MULTIDIM_ADVECTION to
572  C disable this section of code.  C disable this section of code.
573          IF (tempMultiDimAdvec) THEN          IF (tempMultiDimAdvec) THEN
574            CALL GAD_ADVECTION(bi,bj,tempAdvScheme,GAD_TEMPERATURE,  #ifdef ALLOW_DEBUG
575       U                      theta,gT,            IF ( debugLevel .GE. debLevB )
576       I                      myTime,myIter,myThid)       &     CALL DEBUG_CALL('GAD_ADVECTION',myThid)
577    #endif
578              CALL GAD_ADVECTION(
579         I             tempImplVertAdv,tempAdvScheme,GAD_TEMPERATURE,
580         I             uVel, vVel, wVel, theta,
581         O             gT,
582         I             bi,bj,myTime,myIter,myThid)
583          ENDIF          ENDIF
584          IF (saltMultiDimAdvec) THEN          IF (saltMultiDimAdvec) THEN
585            CALL GAD_ADVECTION(bi,bj,saltAdvScheme,GAD_SALINITY,  #ifdef ALLOW_DEBUG
586       U                      salt,gS,            IF ( debugLevel .GE. debLevB )
587       I                      myTime,myIter,myThid)       &     CALL DEBUG_CALL('GAD_ADVECTION',myThid)
588    #endif
589              CALL GAD_ADVECTION(
590         I             saltImplVertAdv,saltAdvScheme,GAD_SALINITY,
591         I             uVel, vVel, wVel, salt,
592         O             gS,
593         I             bi,bj,myTime,myIter,myThid)
594          ENDIF          ENDIF
595  C Since passive tracers are configurable separately from T,S we  C Since passive tracers are configurable separately from T,S we
596  C call the multi-dimensional method for PTRACERS regardless  C call the multi-dimensional method for PTRACERS regardless
597  C of whether multiDimAdvection is set or not.  C of whether multiDimAdvection is set or not.
598  #ifdef ALLOW_PTRACERS  #ifdef ALLOW_PTRACERS
599          IF ( usePTRACERS ) THEN          IF ( usePTRACERS ) THEN
600    #ifdef ALLOW_DEBUG
601              IF ( debugLevel .GE. debLevB )
602         &     CALL DEBUG_CALL('PTRACERS_ADVECTION',myThid)
603    #endif
604           CALL PTRACERS_ADVECTION( bi,bj,myIter,myTime,myThid )           CALL PTRACERS_ADVECTION( bi,bj,myIter,myTime,myThid )
605          ENDIF          ENDIF
606  #endif /* ALLOW_PTRACERS */  #endif /* ALLOW_PTRACERS */
607  #endif /* DISABLE_MULTIDIM_ADVECTION */  #endif /* DISABLE_MULTIDIM_ADVECTION */
608    
609    #ifdef ALLOW_DEBUG
610            IF ( debugLevel .GE. debLevB )
611         &    CALL DEBUG_MSG('ENTERING DOWNWARD K LOOP',myThid)
612    #endif
613    
614  C--     Start of thermodynamics loop  C--     Start of thermodynamics loop
615          DO k=Nr,1,-1          DO k=Nr,1,-1
616  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 574  C--       kDown  Cycles through 2,1 to p Line 633  C--       kDown  Cycles through 2,1 to p
633            jMin = 1-OLy            jMin = 1-OLy
634            jMax = sNy+OLy            jMax = sNy+OLy
635    
636              kp1Msk=1.
637              IF (k.EQ.Nr) kp1Msk=0.
638               DO j=1-Oly,sNy+Oly
639                DO i=1-Olx,sNx+Olx
640                 rTransKp1(i,j) = kp1Msk*rTrans(i,j)
641                ENDDO
642               ENDDO
643    
644  C--       Get temporary terms used by tendency routines  C--       Get temporary terms used by tendency routines
645            CALL CALC_COMMON_FACTORS (            CALL CALC_COMMON_FACTORS (
646       I         bi,bj,iMin,iMax,jMin,jMax,k,       I         bi,bj,iMin,iMax,jMin,jMax,k,
647       O         xA,yA,uTrans,vTrans,rTrans,maskUp,       O         xA,yA,uTrans,vTrans,rTrans,maskUp,
648       I         myThid)       I         myThid)
649    
650              IF (k.EQ.1) THEN
651    C- Surface interface :
652               DO j=1-Oly,sNy+Oly
653                DO i=1-Olx,sNx+Olx
654                 rTrans(i,j) = 0.
655                ENDDO
656               ENDDO
657              ELSE
658    C- Interior interface :
659               DO j=1-Oly,sNy+Oly
660                DO i=1-Olx,sNx+Olx
661                 rTrans(i,j) = rTrans(i,j)*maskC(i,j,k-1,bi,bj)
662                ENDDO
663               ENDDO
664              ENDIF
665    
666  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
667    
668  C--   Residual transp = Bolus transp + Eulerian transp  C--   Residual transp = Bolus transp + Eulerian transp
# Line 600  CADJ STORE rTrans(:,:)    = comlev1_bibj Line 683  CADJ STORE rTrans(:,:)    = comlev1_bibj
683    
684  #endif /* ALLOW_GMREDI */  #endif /* ALLOW_GMREDI */
685    
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE KappaRT(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte  
 CADJ STORE KappaRS(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
686  #ifdef  INCLUDE_CALC_DIFFUSIVITY_CALL  #ifdef  INCLUDE_CALC_DIFFUSIVITY_CALL
687  C--      Calculate the total vertical diffusivity  C--      Calculate the total vertical diffusivity
688           CALL CALC_DIFFUSIVITY(           CALL CALC_DIFFUSIVITY(
# Line 612  C--      Calculate the total vertical di Line 690  C--      Calculate the total vertical di
690       I        maskUp,       I        maskUp,
691       O        KappaRT,KappaRS,       O        KappaRT,KappaRS,
692       I        myThid)       I        myThid)
693    # ifdef ALLOW_AUTODIFF_TAMC
694    CADJ STORE KappaRT(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte
695    CADJ STORE KappaRS(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte
696    # endif /* ALLOW_AUTODIFF_TAMC */
697  #endif  #endif
698    
699            iMin = 1-OLx+2            iMin = 1-OLx+2
# Line 624  C        and step forward storing result Line 706  C        and step forward storing result
706           IF ( tempStepping ) THEN           IF ( tempStepping ) THEN
707             CALL CALC_GT(             CALL CALC_GT(
708       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
709       I         xA,yA,uTrans,vTrans,rTrans,maskUp,       I         xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
710       I         KappaRT,       I         KappaRT,
711       U         fVerT,       U         fVerT,
712       I         myTime,myIter,myThid)       I         myTime,myIter,myThid)
# Line 633  C        and step forward storing result Line 715  C        and step forward storing result
715       I         theta, gT,       I         theta, gT,
716       I         myIter, myThid)       I         myIter, myThid)
717           ENDIF           ENDIF
718  cswdice ---- add ---  
 #ifdef ALLOW_THERM_SEAICE  
        if (k.eq.1) then  
         call ICE_FREEZE(bi, bj, iMin, iMax, jMin, jMax, myThid )  
        endif  
 #endif  
 cswdice -- end add ---  
719           IF ( saltStepping ) THEN           IF ( saltStepping ) THEN
720             CALL CALC_GS(             CALL CALC_GS(
721       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
722       I         xA,yA,uTrans,vTrans,rTrans,maskUp,       I         xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
723       I         KappaRS,       I         KappaRS,
724       U         fVerS,       U         fVerS,
725       I         myTime,myIter,myThid)       I         myTime,myIter,myThid)
# Line 653  cswdice -- end add --- Line 729  cswdice -- end add ---
729       I         myIter, myThid)       I         myIter, myThid)
730           ENDIF           ENDIF
731  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_PASSIVE_TRACER
732    ceh3 needs an IF ( usePASSIVE_TRACER ) THEN
733           IF ( tr1Stepping ) THEN           IF ( tr1Stepping ) THEN
734             CALL CALC_GTR1(             CALL CALC_GTR1(
735       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
736       I         xA,yA,uTrans,vTrans,rTrans,maskUp,       I         xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
737       I         KappaRT,       I         KappaRT,
738       U         fVerTr1,       U         fVerTr1,
739       I         myTime,myIter,myThid)       I         myTime,myIter,myThid)
# Line 668  cswdice -- end add --- Line 745  cswdice -- end add ---
745  #endif  #endif
746  #ifdef ALLOW_PTRACERS  #ifdef ALLOW_PTRACERS
747           IF ( usePTRACERS ) THEN           IF ( usePTRACERS ) THEN
748             CALL PTRACERS_INTEGERATE(             CALL PTRACERS_INTEGRATE(
749       I         bi,bj,k,       I         bi,bj,k,
750       I         xA,yA,uTrans,vTrans,rTrans,maskUp,       I         xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
751       X         KappaRS,       X         fVerP, KappaRS,
752       I         myIter,myTime,myThid)       I         myIter,myTime,myThid)
753           ENDIF           ENDIF
754  #endif /* ALLOW_PTRACERS */  #endif /* ALLOW_PTRACERS */
# Line 684  C--      Apply open boundary conditions Line 761  C--      Apply open boundary conditions
761  #endif   /* ALLOW_OBCS */  #endif   /* ALLOW_OBCS */
762    
763  C--      Freeze water  C--      Freeze water
764           IF ( allowFreezing .AND. .NOT. useSEAICE ) THEN  C  this bit of code is left here for backward compatibility.
765    C  freezing at surface level has been moved to FORWARD_STEP
766             IF ( useOldFreezing .AND. .NOT. useSEAICE
767         &       .AND. .NOT.(useThSIce.AND.k.EQ.1) ) THEN
768  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
769  CADJ STORE gT(:,:,k,bi,bj) = comlev1_bibj_k  CADJ STORE gT(:,:,k,bi,bj) = comlev1_bibj_k
770  CADJ &   , key = kkey, byte = isbyte  CADJ &   , key = kkey, byte = isbyte
771  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
772              CALL FREEZE( bi, bj, iMin, iMax, jMin, jMax, k, myThid )              CALL FREEZE( bi, bj, iMin, iMax, jMin, jMax, k, myThid )
773           END IF           ENDIF
774    
775  C--     end of thermodynamic k loop (Nr:1)  C--     end of thermodynamic k loop (Nr:1)
776          ENDDO          ENDDO
777    
 cswdice -- add ---  
 #ifdef ALLOW_THERM_SEAICE  
 c timeaveraging for ice model values  
            CALL ICE_AVE(bi,bj,iMin,iMax,jMin,jMax,myThid )  
 #endif  
 cswdice --- end add ---  
   
   
778    
779    C--     Implicit vertical advection & diffusion
780  C--     Implicit diffusion  #ifdef INCLUDE_IMPLVERTADV_CODE
781          IF (implicitDiffusion) THEN          IF ( tempImplVertAdv ) THEN
782              CALL GAD_IMPLICIT_R(
783           IF (tempStepping) THEN       I         tempImplVertAdv, tempAdvScheme, GAD_TEMPERATURE,
784         I         kappaRT, wVel, theta,
785         U         gT,
786         I         bi, bj, myTime, myIter, myThid )
787            ELSEIF ( tempStepping .AND. implicitDiffusion ) THEN
788    #else /* INCLUDE_IMPLVERTADV_CODE */
789            IF     ( tempStepping .AND. implicitDiffusion ) THEN
790    #endif /* INCLUDE_IMPLVERTADV_CODE */
791  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
792    CADJ STORE KappaRT(:,:,:) = comlev1_bibj , key=itdkey, byte=isbyte
793  CADJ STORE gT(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte  CADJ STORE gT(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte
794  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
795              CALL IMPLDIFF(            CALL IMPLDIFF(
796       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
797       I         deltaTtracer, KappaRT, recip_HFacC,       I         deltaTtracer, KappaRT, recip_HFacC,
798       U         gT,       U         gT,
799       I         myThid )       I         myThid )
800           ENDIF          ENDIF
801    
802           IF (saltStepping) THEN  #ifdef INCLUDE_IMPLVERTADV_CODE
803            IF ( saltImplVertAdv ) THEN
804              CALL GAD_IMPLICIT_R(
805         I         saltImplVertAdv, saltAdvScheme, GAD_SALINITY,
806         I         kappaRS, wVel, salt,
807         U         gS,
808         I         bi, bj, myTime, myIter, myThid )
809            ELSEIF ( saltStepping .AND. implicitDiffusion ) THEN
810    #else /* INCLUDE_IMPLVERTADV_CODE */
811            IF     ( saltStepping .AND. implicitDiffusion ) THEN
812    #endif /* INCLUDE_IMPLVERTADV_CODE */
813  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
814    CADJ STORE KappaRS(:,:,:) = comlev1_bibj , key=itdkey, byte=isbyte
815  CADJ STORE gS(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte  CADJ STORE gS(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte
816  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
817              CALL IMPLDIFF(            CALL IMPLDIFF(
818       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
819       I         deltaTtracer, KappaRS, recip_HFacC,       I         deltaTtracer, KappaRS, recip_HFacC,
820       U         gS,       U         gS,
821       I         myThid )       I         myThid )
822           ENDIF          ENDIF
823    
824  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_PASSIVE_TRACER
825           IF (tr1Stepping) THEN          IF ( tr1Stepping .AND. implicitDiffusion ) THEN
826  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
827  CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte  CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte
828  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
# Line 740  CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_b Line 831  CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_b
831       I      deltaTtracer, KappaRT, recip_HFacC,       I      deltaTtracer, KappaRT, recip_HFacC,
832       U      gTr1,       U      gTr1,
833       I      myThid )       I      myThid )
834           ENDIF          ENDIF
835  #endif  #endif
836    
837  #ifdef ALLOW_PTRACERS  #ifdef ALLOW_PTRACERS
838  C Vertical diffusion (implicit) for passive tracers  c #ifdef INCLUDE_IMPLVERTADV_CODE
839           IF ( usePTRACERS ) THEN  c       IF ( usePTRACERS .AND. ptracerImplVertAdv ) THEN
840    c       ELSEIF ( usePTRACERS .AND. implicitDiffusion ) THEN
841    c #else
842            IF     ( usePTRACERS .AND. implicitDiffusion ) THEN
843    C--     Vertical diffusion (implicit) for passive tracers
844             CALL PTRACERS_IMPLDIFF( bi,bj,KappaRS,myThid )             CALL PTRACERS_IMPLDIFF( bi,bj,KappaRS,myThid )
845           ENDIF          ENDIF
846  #endif /* ALLOW_PTRACERS */  #endif /* ALLOW_PTRACERS */
847    
848  #ifdef   ALLOW_OBCS  #ifdef   ALLOW_OBCS
849  C--      Apply open boundary conditions  C--      Apply open boundary conditions
850           IF (useOBCS) THEN          IF ( ( implicitDiffusion
851         &    .OR. tempImplVertAdv
852         &    .OR. saltImplVertAdv
853         &       ) .AND. useOBCS     ) THEN
854             DO K=1,Nr             DO K=1,Nr
855               CALL OBCS_APPLY_TS( bi, bj, k, gT, gS, myThid )               CALL OBCS_APPLY_TS( bi, bj, k, gT, gS, myThid )
856             ENDDO             ENDDO
857           END IF          ENDIF
858  #endif   /* ALLOW_OBCS */  #endif   /* ALLOW_OBCS */
859    
860  C--     End If implicitDiffusion  #ifdef ALLOW_TIMEAVE
861    #ifndef HRCUBE
862            IF (taveFreq.GT.0. .AND. ivdc_kappa.NE.0.) THEN
863              CALL TIMEAVE_CUMUL_1T(ConvectCountTave, ConvectCount,
864         I                           Nr, deltaTclock, bi, bj, myThid)
865          ENDIF          ENDIF
866            useVariableK = useKPP .OR. useGMredi .OR. ivdc_kappa.NE.0.
867            IF (taveFreq.GT.0. .AND. useVariableK ) THEN
868             IF (implicitDiffusion) THEN
869              CALL TIMEAVE_CUMUL_DIF_1T(TdiffRtave, gT, kappaRT,
870         I                        Nr, 3, deltaTclock, bi, bj, myThid)
871             ELSE
872              CALL TIMEAVE_CUMUL_DIF_1T(TdiffRtave, theta, kappaRT,
873         I                        Nr, 3, deltaTclock, bi, bj, myThid)
874             ENDIF
875            ENDIF
876    #endif /* ndef HRCUBE */
877    #endif /* ALLOW_TIMEAVE */
878    
879  #endif /* SINGLE_LAYER_MODE */  #endif /* SINGLE_LAYER_MODE */
880    
881  Ccs-  C--   end bi,bj loops.
882         ENDDO         ENDDO
883        ENDDO        ENDDO
884    
885  #ifdef ALLOW_AIM  #ifdef ALLOW_DEBUG
 c     IF ( useAIM ) THEN  
 c      CALL AIM_AIM2DYN_EXCHANGES( myTime, myThid )  
 c     ENDIF  
 #endif /* ALLOW_AIM */  
 c     IF ( staggerTimeStep ) THEN  
 c      IF ( useAIM .OR. useCubedSphereExchange ) THEN  
 c        IF (tempStepping) _EXCH_XYZ_R8(gT,myThid)  
 c        IF (saltStepping) _EXCH_XYZ_R8(gS,myThid)  
 c      ELSEIF ( useGMRedi .AND. Oly.LT.4 ) THEN  
 cc       .AND. GM_AdvForm .AND. .NOT.GM_AdvSeparate ) THEN  
 c        IF (tempMultiDimAdvec) _EXCH_XYZ_R8(gT,myThid)  
 c        IF (saltMultiDimAdvec) _EXCH_XYZ_R8(gS,myThid)  
 c      ENDIF    
 c     ENDIF    
   
 #ifndef DISABLE_DEBUGMODE  
886        If (debugMode) THEN        If (debugMode) THEN
887         CALL DEBUG_STATS_RL(Nr,uVel,'Uvel (THERMODYNAMICS)',myThid)         CALL DEBUG_STATS_RL(Nr,uVel,'Uvel (THERMODYNAMICS)',myThid)
888         CALL DEBUG_STATS_RL(Nr,vVel,'Vvel (THERMODYNAMICS)',myThid)         CALL DEBUG_STATS_RL(Nr,vVel,'Vvel (THERMODYNAMICS)',myThid)
# Line 803  c     ENDIF Line 901  c     ENDIF
901        ENDIF        ENDIF
902  #endif  #endif
903    
904    #ifdef ALLOW_DEBUG
905             IF ( debugLevel .GE. debLevB )
906         &    CALL DEBUG_LEAVE('THERMODYNAMICS',myThid)
907    #endif
908    
909        RETURN        RETURN
910        END        END
Legend:
Removed from v.1.36  
changed lines
  Added in v.1.65
  ViewVC Help
Powered by ViewVC 1.1.22