/[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.42 by heimbach, Fri Jun 27 01:51:10 2003 UTC revision 1.76 by mlosch, Thu Sep 2 09:13:49 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  #ifdef ALLOW_AUTODIFF_TAMC
8  # ifdef ALLOW_GMREDI  # ifdef ALLOW_GMREDI
9  #  include "GMREDI_OPTIONS.h"  #  include "GMREDI_OPTIONS.h"
# Line 9  C $Name$ Line 11  C $Name$
11  # ifdef ALLOW_KPP  # ifdef ALLOW_KPP
12  #  include "KPP_OPTIONS.h"  #  include "KPP_OPTIONS.h"
13  # endif  # endif
 # ifdef ALLOW_PTRACERS  
 #  include "PTRACERS_OPTIONS.h"  
 # endif  
14  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
15    
16  CBOP  CBOP
# Line 82  C     == Global variables === Line 81  C     == Global variables ===
81  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_PASSIVE_TRACER
82  #include "TR1.h"  #include "TR1.h"
83  #endif  #endif
84    #ifdef ALLOW_OFFLINE
85    #include "OFFLINE.h"
86    #endif
87    #ifdef ALLOW_PTRACERS
88    #include "PTRACERS_SIZE.h"
89    #include "PTRACERS.h"
90    #endif
91  #ifdef ALLOW_TIMEAVE  #ifdef ALLOW_TIMEAVE
92  #include "TIMEAVE_STATV.h"  #include "TIMEAVE_STATV.h"
93  #endif  #endif
# Line 97  C     == Global variables === Line 103  C     == Global variables ===
103  # ifdef ALLOW_GMREDI  # ifdef ALLOW_GMREDI
104  #  include "GMREDI.h"  #  include "GMREDI.h"
105  # endif  # endif
106  # ifdef ALLOW_PTRACERS  # ifdef ALLOW_EBM
107  #  include "PTRACERS.h"  #  include "EBM.h"
 # endif  
 cswdice --- add ----  
 # ifdef ALLOW_THERM_SEAICE  
 #  include "ICE.h"  
108  # endif  # endif
 cswdice ------  
109  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
110    
111    
112  C     !INPUT/OUTPUT PARAMETERS:  C     !INPUT/OUTPUT PARAMETERS:
113  C     == Routine arguments ==  C     == Routine arguments ==
114  C     myTime - Current time in simulation  C     myTime - Current time in simulation
# Line 124  C                              transport Line 126  C                              transport
126  C                              o uTrans: Zonal transport  C                              o uTrans: Zonal transport
127  C                              o vTrans: Meridional transport  C                              o vTrans: Meridional transport
128  C                              o rTrans: Vertical transport  C                              o rTrans: Vertical transport
129    C     rTransKp1                o vertical volume transp. at interface k+1
130  C     maskUp                   o maskUp: land/water mask for W points  C     maskUp                   o maskUp: land/water mask for W points
131  C     fVer[STUV]               o fVer: Vertical flux term - note fVer  C     fVer[STUV]               o fVer: Vertical flux term - note fVer
132  C                                      is "pipelined" in the vertical  C                                      is "pipelined" in the vertical
133  C                                      so we need an fVer for each  C                                      so we need an fVer for each
134  C                                      variable.  C                                      variable.
 C     rhoK, rhoKM1   - Density at current level, and level above  
 C     phiSurfX, - gradient of Surface potentiel (Pressure/rho, ocean)  
 C     phiSurfY             or geopotentiel (atmos) in X and Y direction  
135  C     KappaRT,       - Total diffusion in vertical for T and S.  C     KappaRT,       - Total diffusion in vertical for T and S.
136  C     KappaRS          (background + spatially varying, isopycnal term).  C     KappaRS          (background + spatially varying, isopycnal term).
137  C     useVariableK   = T when vertical diffusion is not constant  C     useVariableK   = T when vertical diffusion is not constant
# Line 146  C                      index into fVerTe Line 146  C                      index into fVerTe
146        _RL uTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL uTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
147        _RL vTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL vTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
148        _RL rTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
149          _RL rTransKp1(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
150        _RS maskUp  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RS maskUp  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
151        _RL fVerT   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerT   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
152        _RL fVerS   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerS   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
153    #ifdef ALLOW_PASSIVE_TRACER
154        _RL fVerTr1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerTr1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
155        _RL rhokm1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  #endif
156        _RL rhok    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  #ifdef ALLOW_PTRACERS
157        _RL phiSurfX(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL fVerP   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2,PTRACERS_num)
158        _RL phiSurfY(1-OLx:sNx+OLx,1-OLy:sNy+OLy)  #endif
159        _RL KappaRT (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL KappaRT (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
160        _RL KappaRS (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL KappaRS (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
161        _RL sigmaX  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaX  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
162        _RL sigmaY  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaY  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
163        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
164  C     This is currently used by IVDC and Diagnostics        _RL kp1Msk
       _RL ConvectCount (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)  
165        LOGICAL useVariableK        LOGICAL useVariableK
166        INTEGER iMin, iMax        INTEGER iMin, iMax
167        INTEGER jMin, jMax        INTEGER jMin, jMax
168        INTEGER bi, bj        INTEGER bi, bj
169        INTEGER i, j        INTEGER i, j
170        INTEGER k, km1, kup, kDown        INTEGER k, km1, kup, kDown
171        INTEGER iTracer        INTEGER iTracer, ip
172    
173  CEOP  CEOP
174    
175  #ifndef DISABLE_DEBUGMODE  #ifdef ALLOW_DEBUG
176           IF (debugMode) CALL DEBUG_ENTER('FORWARD_STEP',myThid)           IF ( debugLevel .GE. debLevB )
177         &    CALL DEBUG_ENTER('THERMODYNAMICS',myThid)
178  #endif  #endif
179    
180  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 223  C     uninitialised but inert locations. Line 225  C     uninitialised but inert locations.
225            yA(i,j)        = 0. _d 0            yA(i,j)        = 0. _d 0
226            uTrans(i,j)    = 0. _d 0            uTrans(i,j)    = 0. _d 0
227            vTrans(i,j)    = 0. _d 0            vTrans(i,j)    = 0. _d 0
           rhok   (i,j)   = 0. _d 0  
           rhoKM1 (i,j)   = 0. _d 0  
           phiSurfX(i,j)  = 0. _d 0  
           phiSurfY(i,j)  = 0. _d 0  
228            rTrans (i,j)   = 0. _d 0            rTrans (i,j)   = 0. _d 0
229              rTransKp1(i,j) = 0. _d 0
230            fVerT  (i,j,1) = 0. _d 0            fVerT  (i,j,1) = 0. _d 0
231            fVerT  (i,j,2) = 0. _d 0            fVerT  (i,j,2) = 0. _d 0
232            fVerS  (i,j,1) = 0. _d 0            fVerS  (i,j,1) = 0. _d 0
233            fVerS  (i,j,2) = 0. _d 0            fVerS  (i,j,2) = 0. _d 0
234    #ifdef ALLOW_PASSIVE_TRACER
235            fVerTr1(i,j,1) = 0. _d 0            fVerTr1(i,j,1) = 0. _d 0
236            fVerTr1(i,j,2) = 0. _d 0            fVerTr1(i,j,2) = 0. _d 0
237    #endif
238    #ifdef ALLOW_PTRACERS
239              DO ip=1,PTRACERS_num
240               fVerP  (i,j,1,ip) = 0. _d 0
241               fVerP  (i,j,2,ip) = 0. _d 0
242              ENDDO
243    #endif
244           ENDDO           ENDDO
245          ENDDO          ENDDO
246    
# Line 241  C     uninitialised but inert locations. Line 248  C     uninitialised but inert locations.
248           DO j=1-OLy,sNy+OLy           DO j=1-OLy,sNy+OLy
249            DO i=1-OLx,sNx+OLx            DO i=1-OLx,sNx+OLx
250  C This is currently also used by IVDC and Diagnostics  C This is currently also used by IVDC and Diagnostics
            sigmaX(i,j,k) = 0. _d 0  
            sigmaY(i,j,k) = 0. _d 0  
            sigmaR(i,j,k) = 0. _d 0  
            ConvectCount(i,j,k) = 0.  
251             KappaRT(i,j,k)    = 0. _d 0             KappaRT(i,j,k)    = 0. _d 0
252             KappaRS(i,j,k)    = 0. _d 0             KappaRS(i,j,k)    = 0. _d 0
253  #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.  
254             gT(i,j,k,bi,bj)   = 0. _d 0             gT(i,j,k,bi,bj)   = 0. _d 0
255             gS(i,j,k,bi,bj)   = 0. _d 0             gS(i,j,k,bi,bj)   = 0. _d 0
256  # ifdef ALLOW_PASSIVE_TRACER  # ifdef ALLOW_PASSIVE_TRACER
257    ceh3 needs an IF ( use PASSIVE_TRACER) THEN
258             gTr1(i,j,k,bi,bj) = 0. _d 0             gTr1(i,j,k,bi,bj) = 0. _d 0
259  # endif  # endif
260  # ifdef ALLOW_PTRACERS  # ifdef ALLOW_PTRACERS
261    ceh3 this should have an   IF ( usePTRACERS ) THEN
262             DO iTracer=1,PTRACERS_numInUse             DO iTracer=1,PTRACERS_numInUse
263              gPTr(i,j,k,bi,bj,itracer) = 0. _d 0              gPTr(i,j,k,bi,bj,itracer) = 0. _d 0
264             ENDDO             ENDDO
265  # endif  # endif
 # ifdef ALLOW_GMREDI  
            Kwx(i,j,k,bi,bj)  = 0. _d 0  
            Kwy(i,j,k,bi,bj)  = 0. _d 0  
            Kwz(i,j,k,bi,bj)  = 0. _d 0  
 #  ifdef GM_NON_UNITY_DIAGONAL  
            Kux(i,j,k,bi,bj)  = 0. _d 0  
            Kvy(i,j,k,bi,bj)  = 0. _d 0  
 #  endif  
 #  ifdef GM_EXTRA_DIAGONAL  
            Kuz(i,j,k,bi,bj)  = 0. _d 0  
            Kvz(i,j,k,bi,bj)  = 0. _d 0  
 #  endif  
 #  ifdef GM_BOLUS_ADVEC  
            GM_PsiX(i,j,k,bi,bj)  = 0. _d 0  
            GM_PsiY(i,j,k,bi,bj)  = 0. _d 0  
 #  endif  
 #  ifdef GM_VISBECK_VARIABLE_K  
            VisbeckK(i,j,bi,bj)   = 0. _d 0  
 #  endif  
 # endif /* ALLOW_GMREDI */  
 #endif /* ALLOW_AUTODIFF_TAMC */  
266            ENDDO            ENDDO
267           ENDDO           ENDDO
268          ENDDO          ENDDO
269    
270          iMin = 1-OLx  c       iMin = 1-OLx
271          iMax = sNx+OLx  c       iMax = sNx+OLx
272          jMin = 1-OLy  c       jMin = 1-OLy
273          jMax = sNy+OLy  c       jMax = sNy+OLy
274    
275  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
276  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
# Line 301  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_ Line 283  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_
283  #endif  #endif
284  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
285    
 #ifndef DISABLE_DEBUGMODE  
         IF (debugMode) CALL DEBUG_MSG('ENTERING UPWARD K LOOP',myThid)  
 #endif  
   
 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 = (itdkey-1)*Nr + k  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
 C--       Integrate continuity vertically for vertical velocity  
 c         CALL INTEGRATE_FOR_W(  
 c    I                         bi, bj, k, uVel, vVel,  
 c    O                         wVel,  
 c    I                         myThid )  
   
 #ifdef    ALLOW_OBCS  
 #ifdef    ALLOW_NONHYDROSTATIC  
 C--       Apply OBC to W if in N-H mode  
 c         IF (useOBCS.AND.nonHydrostatic) THEN  
 c           CALL OBCS_APPLY_W( bi, bj, k, wVel, myThid )  
 c         ENDIF  
 #endif    /* ALLOW_NONHYDROSTATIC */  
 #endif    /* ALLOW_OBCS */  
   
 C--       Attention: by defining "SINGLE_LAYER_MODE" in CPP_OPTIONS.h  
 C--       MOST of THERMODYNAMICS will be disabled  
 #ifndef SINGLE_LAYER_MODE  
   
 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  
 #ifndef DISABLE_DEBUGMODE  
             IF (debugMode) CALL DEBUG_CALL('FIND_RHO',myThid)  
 #endif  
 #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,  
      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,  
      I        theta, salt,  
      O        rhoKm1,  
      I        myThid )  
             ENDIF  
 #ifndef DISABLE_DEBUGMODE  
             IF (debugMode) CALL DEBUG_CALL('GRAD_SIGMA',myThid)  
 #endif  
             CALL GRAD_SIGMA(  
      I             bi, bj, iMin, iMax, jMin, jMax, k,  
      I             rhoK, rhoKm1, rhoK,  
      O             sigmaX, sigmaY, sigmaR,  
      I             myThid )  
           ENDIF  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE rhok   (:,:) = comlev1_bibj_k ,       key=kkey, byte=isbyte  
 CADJ STORE rhokm1 (:,:) = comlev1_bibj_k ,       key=kkey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
 C--       Implicit Vertical Diffusion for Convection  
 c ==> should use sigmaR !!!  
           IF (k.GT.1 .AND. ivdc_kappa.NE.0.) THEN  
 #ifndef DISABLE_DEBUGMODE  
             IF (debugMode) CALL DEBUG_CALL('CALC_IVDC',myThid)  
 #endif  
             CALL CALC_IVDC(  
      I        bi, bj, iMin, iMax, jMin, jMax, k,  
      I        rhoKm1, rhoK,  
      U        ConvectCount, KappaRT, KappaRS,  
      I        myTime, myIter, myThid)  
           ENDIF  
   
 #endif /* SINGLE_LAYER_MODE */  
   
 C--     end of diagnostic k loop (Nr:1)  
         ENDDO  
   
286  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
287  cph avoids recomputation of integrate_for_w  cph avoids recomputation of integrate_for_w
288  CADJ STORE wvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE wvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
289  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
290    
 #ifdef  ALLOW_OBCS  
 C--     Calculate future values on open boundaries  
         IF (useOBCS) THEN  
 #ifndef DISABLE_DEBUGMODE  
           IF (debugMode) CALL DEBUG_CALL('OBCS_CALC',myThid)  
 #endif  
           CALL OBCS_CALC( bi, bj, myTime+deltaT, myIter+1,  
      I            uVel, vVel, wVel, theta, salt,  
      I            myThid )  
         ENDIF  
 #endif  /* ALLOW_OBCS */  
   
   
 c********************************************  
 cswdice --- add ---  
 #ifdef ALLOW_THERM_SEAICE  
 #ifndef DISABLE_DEBUGMODE  
         IF (debugMode) CALL DEBUG_CALL('ICE_FORCING',myThid)  
 #endif  
 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 ---  
   
 C--     Determines forcing terms based on external fields  
 C       relaxation terms, etc.  
 #ifndef DISABLE_DEBUGMODE  
         IF (debugMode) CALL DEBUG_CALL('EXTERNAL_FORCING_SURF',myThid)  
 #endif  
         CALL EXTERNAL_FORCING_SURF(  
      I             bi, bj, iMin, iMax, jMin, jMax,  
      I             myThid )  
 cswdice --- add ----  
 #endif  
 cswdice --- end add ---  
 c******************************************  
   
   
   
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 cph needed for KPP  
 CADJ STORE surfacetendencyU(:,:,bi,bj)  
 CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  
 CADJ STORE surfacetendencyV(:,:,bi,bj)  
 CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  
 CADJ STORE surfacetendencyS(:,:,bi,bj)  
 CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  
 CADJ STORE surfacetendencyT(:,:,bi,bj)  
 CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
291  C--     Attention: by defining "SINGLE_LAYER_MODE" in CPP_OPTIONS.h  C--     Attention: by defining "SINGLE_LAYER_MODE" in CPP_OPTIONS.h
292  C--     MOST of THERMODYNAMICS will be disabled  C--     MOST of THERMODYNAMICS will be disabled
293  #ifndef SINGLE_LAYER_MODE  #ifndef SINGLE_LAYER_MODE
294    
 #ifdef  ALLOW_GMREDI  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 cph storing here is needed only for one GMREDI_OPTIONS:  
 cph define GM_BOLUS_ADVEC  
 cph but I've avoided the #ifdef for now, in case more things change  
 CADJ STORE sigmaX(:,:,:)        = comlev1_bibj, key=itdkey, byte=isbyte  
 CADJ STORE sigmaY(:,:,:)        = comlev1_bibj, key=itdkey, byte=isbyte  
 CADJ STORE sigmaR(:,:,:)        = comlev1_bibj, key=itdkey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
 C--     Calculate iso-neutral slopes for the GM/Redi parameterisation  
         IF (useGMRedi) THEN  
 #ifndef DISABLE_DEBUGMODE  
           IF (debugMode) CALL DEBUG_CALL('GMREDI_CALC_TENSOR',myThid)  
 #endif  
           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  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE Kwx(:,:,:,bi,bj)   = comlev1_bibj, key=itdkey, byte=isbyte  
 CADJ STORE Kwy(:,:,:,bi,bj)   = comlev1_bibj, key=itdkey, byte=isbyte  
 CADJ STORE Kwz(:,:,:,bi,bj)   = comlev1_bibj, key=itdkey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
 #endif  /* ALLOW_GMREDI */  
   
 #ifdef  ALLOW_KPP  
 C--     Compute KPP mixing coefficients  
         IF (useKPP) THEN  
 #ifndef DISABLE_DEBUGMODE  
           IF (debugMode) CALL DEBUG_CALL('KPP_CALC',myThid)  
 #endif  
           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  
   
295  #ifdef ALLOW_AUTODIFF_TAMC  #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=itdkey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
 #endif  /* ALLOW_KPP */  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE KappaRT(:,:,:)     = comlev1_bibj, key=itdkey, byte=isbyte  
 CADJ STORE KappaRS(:,:,:)     = comlev1_bibj, key=itdkey, byte=isbyte  
296  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
297  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
298  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
# Line 539  cphCADJ &                              k Line 307  cphCADJ &                              k
307  #endif  #endif
308  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
309    
 #ifdef ALLOW_AIM  
 C       AIM - atmospheric intermediate model, physics package code.  
         IF ( useAIM ) THEN  
 #ifndef DISABLE_DEBUGMODE  
           IF (debugMode) CALL DEBUG_CALL('AIM_DO_PHYSICS',myThid)  
 #endif  
          CALL TIMER_START('AIM_DO_PHYSICS   [THERMODYNAMICS]', myThid)  
          CALL AIM_DO_PHYSICS( bi, bj, myTime, myIter, myThid )  
          CALL TIMER_STOP( 'AIM_DO_PHYSICS   [THERMODYNAMICS]', myThid)  
         ENDIF  
 #endif /* ALLOW_AIM */  
   
310  #ifndef DISABLE_MULTIDIM_ADVECTION  #ifndef DISABLE_MULTIDIM_ADVECTION
311  C--     Some advection schemes are better calculated using a multi-dimensional  C--     Some advection schemes are better calculated using a multi-dimensional
312  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 562  C to be able to exclude this scheme to a Line 318  C to be able to exclude this scheme to a
318  C recomputation. It *is* differentiable, if you need it.  C recomputation. It *is* differentiable, if you need it.
319  C Edit GAD_OPTIONS.h and #define DISABLE_MULTIDIM_ADVECTION to  C Edit GAD_OPTIONS.h and #define DISABLE_MULTIDIM_ADVECTION to
320  C disable this section of code.  C disable this section of code.
321    #ifndef ALLOW_OFFLINE
322          IF (tempMultiDimAdvec) THEN          IF (tempMultiDimAdvec) THEN
323  #ifndef DISABLE_DEBUGMODE  #ifdef ALLOW_DEBUG
324            IF (debugMode) CALL DEBUG_CALL('GAD_ADVECTION',myThid)            IF ( debugLevel .GE. debLevB )
325  #endif       &     CALL DEBUG_CALL('GAD_ADVECTION',myThid)
326            CALL GAD_ADVECTION(bi,bj,tempAdvScheme,GAD_TEMPERATURE,  #endif
327       U                      theta,gT,            CALL GAD_ADVECTION(
328       I                      myTime,myIter,myThid)       I             tempImplVertAdv, tempAdvScheme, tempVertAdvScheme,
329         I             GAD_TEMPERATURE,
330         I             uVel, vVel, wVel, theta,
331         O             gT,
332         I             bi,bj,myTime,myIter,myThid)
333          ENDIF          ENDIF
         IF (saltMultiDimAdvec) THEN  
 #ifndef DISABLE_DEBUGMODE  
           IF (debugMode) CALL DEBUG_CALL('GAD_ADVECTION',myThid)  
334  #endif  #endif
335            CALL GAD_ADVECTION(bi,bj,saltAdvScheme,GAD_SALINITY,  #ifndef ALLOW_OFFLINE
336       U                      salt,gS,          IF (saltMultiDimAdvec) THEN
337       I                      myTime,myIter,myThid)  #ifdef ALLOW_DEBUG
338              IF ( debugLevel .GE. debLevB )
339         &     CALL DEBUG_CALL('GAD_ADVECTION',myThid)
340    #endif
341              CALL GAD_ADVECTION(
342         I             saltImplVertAdv, saltAdvScheme, saltVertAdvScheme,
343         I             GAD_SALINITY,
344         I             uVel, vVel, wVel, salt,
345         O             gS,
346         I             bi,bj,myTime,myIter,myThid)
347          ENDIF          ENDIF
348    #endif
349  C Since passive tracers are configurable separately from T,S we  C Since passive tracers are configurable separately from T,S we
350  C call the multi-dimensional method for PTRACERS regardless  C call the multi-dimensional method for PTRACERS regardless
351  C of whether multiDimAdvection is set or not.  C of whether multiDimAdvection is set or not.
352  #ifdef ALLOW_PTRACERS  #ifdef ALLOW_PTRACERS
353          IF ( usePTRACERS ) THEN          IF ( usePTRACERS ) THEN
354  #ifndef DISABLE_DEBUGMODE  #ifdef ALLOW_DEBUG
355            IF (debugMode) CALL DEBUG_CALL('PTRACERS_ADVECTION',myThid)            IF ( debugLevel .GE. debLevB )
356         &     CALL DEBUG_CALL('PTRACERS_ADVECTION',myThid)
357  #endif  #endif
358           CALL PTRACERS_ADVECTION( bi,bj,myIter,myTime,myThid )           CALL PTRACERS_ADVECTION( bi,bj,myIter,myTime,myThid )
359          ENDIF          ENDIF
360  #endif /* ALLOW_PTRACERS */  #endif /* ALLOW_PTRACERS */
361  #endif /* DISABLE_MULTIDIM_ADVECTION */  #endif /* DISABLE_MULTIDIM_ADVECTION */
362    
363  #ifndef DISABLE_DEBUGMODE  #ifdef ALLOW_DEBUG
364         IF (debugMode) CALL DEBUG_MSG('ENTERING DOWNWARD K LOOP',myThid)          IF ( debugLevel .GE. debLevB )
365         &    CALL DEBUG_MSG('ENTERING DOWNWARD K LOOP',myThid)
366  #endif  #endif
367    
368  C--     Start of thermodynamics loop  C--     Start of thermodynamics loop
# Line 617  C--       kDown  Cycles through 2,1 to p Line 387  C--       kDown  Cycles through 2,1 to p
387            jMin = 1-OLy            jMin = 1-OLy
388            jMax = sNy+OLy            jMax = sNy+OLy
389    
390              kp1Msk=1.
391              IF (k.EQ.Nr) kp1Msk=0.
392               DO j=1-Oly,sNy+Oly
393                DO i=1-Olx,sNx+Olx
394                 rTransKp1(i,j) = kp1Msk*rTrans(i,j)
395                ENDDO
396               ENDDO
397    #ifdef ALLOW_AUTODIFF_TAMC
398    CADJ STORE rTransKp1(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
399    #endif
400    
401  C--       Get temporary terms used by tendency routines  C--       Get temporary terms used by tendency routines
402            CALL CALC_COMMON_FACTORS (            CALL CALC_COMMON_FACTORS (
403       I         bi,bj,iMin,iMax,jMin,jMax,k,       I         bi,bj,iMin,iMax,jMin,jMax,k,
404       O         xA,yA,uTrans,vTrans,rTrans,maskUp,       O         xA,yA,uTrans,vTrans,rTrans,maskUp,
405       I         myThid)       I         myThid)
406    
407              IF (k.EQ.1) THEN
408    C- Surface interface :
409               DO j=1-Oly,sNy+Oly
410                DO i=1-Olx,sNx+Olx
411                 rTrans(i,j) = 0.
412                ENDDO
413               ENDDO
414              ELSE
415    C- Interior interface :
416               DO j=1-Oly,sNy+Oly
417                DO i=1-Olx,sNx+Olx
418                 rTrans(i,j) = rTrans(i,j)*maskC(i,j,k-1,bi,bj)
419                ENDDO
420               ENDDO
421              ENDIF
422    
423  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
424    
425  C--   Residual transp = Bolus transp + Eulerian transp  C--   Residual transp = Bolus transp + Eulerian transp
# Line 633  C--   Residual transp = Bolus transp + E Line 430  C--   Residual transp = Bolus transp + E
430       &                    rTrans, bi, bj, k, myThid)       &                    rTrans, bi, bj, k, myThid)
431            ENDIF            ENDIF
432    
433  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
434    CADJ STORE rTrans(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
435  #ifdef GM_BOLUS_ADVEC  #ifdef GM_BOLUS_ADVEC
436  CADJ STORE uTrans(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE uTrans(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
437  CADJ STORE vTrans(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE vTrans(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
 CADJ STORE rTrans(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte  
438  #endif  #endif
439  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
440    
441  #endif /* ALLOW_GMREDI */  #endif /* ALLOW_GMREDI */
442    
 #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 */  
   
443  #ifdef  INCLUDE_CALC_DIFFUSIVITY_CALL  #ifdef  INCLUDE_CALC_DIFFUSIVITY_CALL
444  C--      Calculate the total vertical diffusivity  C--      Calculate the total vertical diffusivity
445           CALL CALC_DIFFUSIVITY(           CALL CALC_DIFFUSIVITY(
# Line 655  C--      Calculate the total vertical di Line 447  C--      Calculate the total vertical di
447       I        maskUp,       I        maskUp,
448       O        KappaRT,KappaRS,       O        KappaRT,KappaRS,
449       I        myThid)       I        myThid)
450    # ifdef ALLOW_AUTODIFF_TAMC
451    CADJ STORE KappaRT(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte
452    CADJ STORE KappaRS(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte
453    # endif /* ALLOW_AUTODIFF_TAMC */
454  #endif  #endif
455    
456            iMin = 1-OLx+2            iMin = 1-OLx+2
# Line 664  C--      Calculate the total vertical di Line 460  C--      Calculate the total vertical di
460    
461  C--      Calculate active tracer tendencies (gT,gS,...)  C--      Calculate active tracer tendencies (gT,gS,...)
462  C        and step forward storing result in gTnm1, gSnm1, etc.  C        and step forward storing result in gTnm1, gSnm1, etc.
463    #ifndef ALLOW_OFFLINE
464           IF ( tempStepping ) THEN           IF ( tempStepping ) THEN
465             CALL CALC_GT(             CALL CALC_GT(
466       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
467       I         xA,yA,uTrans,vTrans,rTrans,maskUp,       I         xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
468       I         KappaRT,       I         KappaRT,
469       U         fVerT,       U         fVerT,
470       I         myTime,myIter,myThid)       I         myTime,myIter,myThid)
# Line 676  C        and step forward storing result Line 473  C        and step forward storing result
473       I         theta, gT,       I         theta, gT,
474       I         myIter, myThid)       I         myIter, myThid)
475           ENDIF           ENDIF
 cswdice ---- add ---  
 #ifdef ALLOW_THERM_SEAICE  
        if (k.eq.1) then  
         call ICE_FREEZE(bi, bj, iMin, iMax, jMin, jMax, myThid )  
        endif  
476  #endif  #endif
477  cswdice -- end add ---  
478    #ifndef ALLOW_OFFLINE
479           IF ( saltStepping ) THEN           IF ( saltStepping ) THEN
480             CALL CALC_GS(             CALL CALC_GS(
481       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
482       I         xA,yA,uTrans,vTrans,rTrans,maskUp,       I         xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
483       I         KappaRS,       I         KappaRS,
484       U         fVerS,       U         fVerS,
485       I         myTime,myIter,myThid)       I         myTime,myIter,myThid)
# Line 695  cswdice -- end add --- Line 488  cswdice -- end add ---
488       I         salt, gS,       I         salt, gS,
489       I         myIter, myThid)       I         myIter, myThid)
490           ENDIF           ENDIF
491    #endif
492  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_PASSIVE_TRACER
493    ceh3 needs an IF ( usePASSIVE_TRACER ) THEN
494           IF ( tr1Stepping ) THEN           IF ( tr1Stepping ) THEN
495             CALL CALC_GTR1(             CALL CALC_GTR1(
496       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
497       I         xA,yA,uTrans,vTrans,rTrans,maskUp,       I         xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
498       I         KappaRT,       I         KappaRT,
499       U         fVerTr1,       U         fVerTr1,
500       I         myTime,myIter,myThid)       I         myTime,myIter,myThid)
# Line 713  cswdice -- end add --- Line 508  cswdice -- end add ---
508           IF ( usePTRACERS ) THEN           IF ( usePTRACERS ) THEN
509             CALL PTRACERS_INTEGRATE(             CALL PTRACERS_INTEGRATE(
510       I         bi,bj,k,       I         bi,bj,k,
511       I         xA,yA,uTrans,vTrans,rTrans,maskUp,       I         xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
512       X         KappaRS,       X         fVerP, KappaRS,
513       I         myIter,myTime,myThid)       I         myIter,myTime,myThid)
514           ENDIF           ENDIF
515  #endif /* ALLOW_PTRACERS */  #endif /* ALLOW_PTRACERS */
# Line 727  C--      Apply open boundary conditions Line 522  C--      Apply open boundary conditions
522  #endif   /* ALLOW_OBCS */  #endif   /* ALLOW_OBCS */
523    
524  C--      Freeze water  C--      Freeze water
525           IF ( allowFreezing .AND. .NOT. useSEAICE ) THEN  C  this bit of code is left here for backward compatibility.
526    C  freezing at surface level has been moved to FORWARD_STEP
527    #ifndef ALLOW_OFFLINE
528             IF ( useOldFreezing .AND. .NOT. useSEAICE
529         &       .AND. .NOT.(useThSIce.AND.k.EQ.1) ) THEN
530  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
531  CADJ STORE gT(:,:,k,bi,bj) = comlev1_bibj_k  CADJ STORE gT(:,:,k,bi,bj) = comlev1_bibj_k
532  CADJ &   , key = kkey, byte = isbyte  CADJ &   , key = kkey, byte = isbyte
533  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
534              CALL FREEZE( bi, bj, iMin, iMax, jMin, jMax, k, myThid )              CALL FREEZE( bi, bj, iMin, iMax, jMin, jMax, k, myThid )
535           END IF           ENDIF
536    #endif
537    
538  C--     end of thermodynamic k loop (Nr:1)  C--     end of thermodynamic k loop (Nr:1)
539          ENDDO          ENDDO
540    
 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 ---  
541    
542    C--     Implicit vertical advection & diffusion
543    #ifndef ALLOW_OFFLINE
544    #ifdef INCLUDE_IMPLVERTADV_CODE
545  C--     Implicit diffusion          IF ( tempImplVertAdv ) THEN
546          IF (implicitDiffusion) THEN            CALL GAD_IMPLICIT_R(
547         I         tempImplVertAdv, tempAdvScheme, GAD_TEMPERATURE,
548           IF (tempStepping) THEN       I         kappaRT, wVel, theta,
549         U         gT,
550         I         bi, bj, myTime, myIter, myThid )
551            ELSEIF ( tempStepping .AND. implicitDiffusion ) THEN
552    #else /* INCLUDE_IMPLVERTADV_CODE */
553            IF     ( tempStepping .AND. implicitDiffusion ) THEN
554    #endif /* INCLUDE_IMPLVERTADV_CODE */
555  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
556    CADJ STORE KappaRT(:,:,:) = comlev1_bibj , key=itdkey, byte=isbyte
557  CADJ STORE gT(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte  CADJ STORE gT(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte
558  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
559              CALL IMPLDIFF(            CALL IMPLDIFF(
560       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
561       I         deltaTtracer, KappaRT, recip_HFacC,       I         deltaTtracer, KappaRT, recip_HFacC,
562       U         gT,       U         gT,
563       I         myThid )       I         myThid )
564           ENDIF          ENDIF
565    #endif
566    
567           IF (saltStepping) THEN  #ifndef ALLOW_OFFLINE
568    #ifdef INCLUDE_IMPLVERTADV_CODE
569            IF ( saltImplVertAdv ) THEN
570              CALL GAD_IMPLICIT_R(
571         I         saltImplVertAdv, saltAdvScheme, GAD_SALINITY,
572         I         kappaRS, wVel, salt,
573         U         gS,
574         I         bi, bj, myTime, myIter, myThid )
575            ELSEIF ( saltStepping .AND. implicitDiffusion ) THEN
576    #else /* INCLUDE_IMPLVERTADV_CODE */
577            IF     ( saltStepping .AND. implicitDiffusion ) THEN
578    #endif /* INCLUDE_IMPLVERTADV_CODE */
579  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
580    CADJ STORE KappaRS(:,:,:) = comlev1_bibj , key=itdkey, byte=isbyte
581  CADJ STORE gS(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte  CADJ STORE gS(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte
582  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
583              CALL IMPLDIFF(            CALL IMPLDIFF(
584       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
585       I         deltaTtracer, KappaRS, recip_HFacC,       I         deltaTtracer, KappaRS, recip_HFacC,
586       U         gS,       U         gS,
587       I         myThid )       I         myThid )
588           ENDIF          ENDIF
589    #endif
590    
591  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_PASSIVE_TRACER
592           IF (tr1Stepping) THEN          IF ( tr1Stepping .AND. implicitDiffusion ) THEN
593  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
594  CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte  CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte
595  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
# Line 783  CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_b Line 598  CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_b
598       I      deltaTtracer, KappaRT, recip_HFacC,       I      deltaTtracer, KappaRT, recip_HFacC,
599       U      gTr1,       U      gTr1,
600       I      myThid )       I      myThid )
601           ENDIF          ENDIF
602  #endif  #endif
603    
604  #ifdef ALLOW_PTRACERS  #ifdef ALLOW_PTRACERS
605  C Vertical diffusion (implicit) for passive tracers  c #ifdef INCLUDE_IMPLVERTADV_CODE
606           IF ( usePTRACERS ) THEN  c       IF ( usePTRACERS .AND. ptracerImplVertAdv ) THEN
607    c       ELSEIF ( usePTRACERS .AND. implicitDiffusion ) THEN
608    c #else
609            IF     ( usePTRACERS .AND. implicitDiffusion ) THEN
610    C--     Vertical diffusion (implicit) for passive tracers
611             CALL PTRACERS_IMPLDIFF( bi,bj,KappaRS,myThid )             CALL PTRACERS_IMPLDIFF( bi,bj,KappaRS,myThid )
612           ENDIF          ENDIF
613  #endif /* ALLOW_PTRACERS */  #endif /* ALLOW_PTRACERS */
614    
615  #ifdef   ALLOW_OBCS  #ifdef   ALLOW_OBCS
616  C--      Apply open boundary conditions  C--      Apply open boundary conditions
617           IF (useOBCS) THEN          IF ( ( implicitDiffusion
618         &    .OR. tempImplVertAdv
619         &    .OR. saltImplVertAdv
620         &       ) .AND. useOBCS     ) THEN
621             DO K=1,Nr             DO K=1,Nr
622               CALL OBCS_APPLY_TS( bi, bj, k, gT, gS, myThid )               CALL OBCS_APPLY_TS( bi, bj, k, gT, gS, myThid )
623             ENDDO             ENDDO
          END IF  
 #endif   /* ALLOW_OBCS */  
   
 C--     End If implicitDiffusion  
624          ENDIF          ENDIF
625    #endif   /* ALLOW_OBCS */
626    
627  #ifdef ALLOW_TIMEAVE  #ifdef ALLOW_TIMEAVE
628            IF ( taveFreq.GT. 0. _d 0 .AND.
629         &       buoyancyRelation(1:7) .EQ. 'OCEANIC' ) THEN
630              CALL TIMEAVE_SURF_FLUX( bi, bj, myTime, myIter, myThid)
631            ENDIF
632    #ifndef HRCUBE
633          IF (taveFreq.GT.0. .AND. ivdc_kappa.NE.0.) THEN          IF (taveFreq.GT.0. .AND. ivdc_kappa.NE.0.) THEN
634            CALL TIMEAVE_CUMUL_1T(ConvectCountTave, ConvectCount,            CALL TIMEAVE_CUMULATE(ConvectCountTave, IVDConvCount,
635       I                           Nr, deltaTclock, bi, bj, myThid)       I                           Nr, deltaTclock, bi, bj, myThid)
636          ENDIF          ENDIF
637          useVariableK = useKPP .OR. useGMredi .OR. ivdc_kappa.NE.0.          useVariableK = useKPP .OR. usePP81 .OR. useMY82
638         &       .OR. useGMredi .OR. ivdc_kappa.NE.0.
639          IF (taveFreq.GT.0. .AND. useVariableK ) THEN          IF (taveFreq.GT.0. .AND. useVariableK ) THEN
640           IF (implicitDiffusion) THEN           IF (implicitDiffusion) THEN
641            CALL TIMEAVE_CUMUL_DIF_1T(TdiffRtave, gT, kappaRT,            CALL TIMEAVE_CUMUL_DIF_1T(TdiffRtave, gT, kappaRT,
# Line 820  C--     End If implicitDiffusion Line 645  C--     End If implicitDiffusion
645       I                        Nr, 3, deltaTclock, bi, bj, myThid)       I                        Nr, 3, deltaTclock, bi, bj, myThid)
646           ENDIF           ENDIF
647          ENDIF          ENDIF
648    #endif /* ndef HRCUBE */
649  #endif /* ALLOW_TIMEAVE */  #endif /* ALLOW_TIMEAVE */
650    
651  #endif /* SINGLE_LAYER_MODE */  #endif /* SINGLE_LAYER_MODE */
# Line 828  C--   end bi,bj loops. Line 654  C--   end bi,bj loops.
654         ENDDO         ENDDO
655        ENDDO        ENDDO
656    
657  #ifndef DISABLE_DEBUGMODE  #ifdef ALLOW_DEBUG
658        If (debugMode) THEN        If (debugMode) THEN
659         CALL DEBUG_STATS_RL(Nr,uVel,'Uvel (THERMODYNAMICS)',myThid)         CALL DEBUG_STATS_RL(Nr,uVel,'Uvel (THERMODYNAMICS)',myThid)
660         CALL DEBUG_STATS_RL(Nr,vVel,'Vvel (THERMODYNAMICS)',myThid)         CALL DEBUG_STATS_RL(Nr,vVel,'Vvel (THERMODYNAMICS)',myThid)
# Line 847  C--   end bi,bj loops. Line 673  C--   end bi,bj loops.
673        ENDIF        ENDIF
674  #endif  #endif
675    
676  #ifndef DISABLE_DEBUGMODE  #ifdef ALLOW_DEBUG
677           IF (debugMode) CALL DEBUG_LEAVE('FORWARD_STEP',myThid)           IF ( debugLevel .GE. debLevB )
678         &    CALL DEBUG_LEAVE('THERMODYNAMICS',myThid)
679  #endif  #endif
680    
681        RETURN        RETURN
Legend:
Removed from v.1.42  
changed lines
  Added in v.1.76
  ViewVC Help
Powered by ViewVC 1.1.22