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

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

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

revision 1.57 by heimbach, Thu Feb 1 19:32:02 2001 UTC revision 1.74 by heimbach, Mon Jul 30 20:37:45 2001 UTC
# Line 1  Line 1 
1  C $Header$  C $Header$
2    C $Name$
3    
4  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
5    
# Line 25  C     \================================= Line 26  C     \=================================
26  C     == Global variables ===  C     == Global variables ===
27  #include "SIZE.h"  #include "SIZE.h"
28  #include "EEPARAMS.h"  #include "EEPARAMS.h"
 #include "CG2D.h"  
29  #include "PARAMS.h"  #include "PARAMS.h"
30  #include "DYNVARS.h"  #include "DYNVARS.h"
31  #include "GRID.h"  #include "GRID.h"
32    #ifdef ALLOW_PASSIVE_TRACER
33    #include "TR1.h"
34    #endif
35    
36  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
37  # include "tamc.h"  # include "tamc.h"
# Line 42  C     == Global variables === Line 45  C     == Global variables ===
45  # endif  # endif
46  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
47    
48    #ifdef ALLOW_TIMEAVE
49    #include "TIMEAVE_STATV.h"
50    #endif
51    
52  C     == Routine arguments ==  C     == Routine arguments ==
53  C     myTime - Current time in simulation  C     myTime - Current time in simulation
54  C     myIter - Current iteration number in simulation  C     myIter - Current iteration number in simulation
# Line 54  C     == Local variables Line 61  C     == Local variables
61  C     xA, yA                 - Per block temporaries holding face areas  C     xA, yA                 - Per block temporaries holding face areas
62  C     uTrans, vTrans, rTrans - Per block temporaries holding flow  C     uTrans, vTrans, rTrans - Per block temporaries holding flow
63  C                              transport  C                              transport
64  C     rVel                     o uTrans: Zonal transport  C                              o uTrans: Zonal transport
65  C                              o vTrans: Meridional transport  C                              o vTrans: Meridional transport
66  C                              o rTrans: Vertical transport  C                              o rTrans: Vertical transport
67  C                              o rVel:   Vertical velocity at upper and  C     maskUp                   o maskUp: land/water mask for W points
68  C                                        lower cell faces.  C     fVer[STUV]               o fVer: Vertical flux term - note fVer
 C     maskC,maskUp             o maskC: land/water mask for tracer cells  
 C                              o maskUp: land/water mask for W points  
 C     aTerm, xTerm, cTerm    - Work arrays for holding separate terms in  
 C     mTerm, pTerm,            tendency equations.  
 C     fZon, fMer, fVer[STUV]   o aTerm: Advection term  
 C                              o xTerm: Mixing term  
 C                              o cTerm: Coriolis term  
 C                              o mTerm: Metric term  
 C                              o pTerm: Pressure term  
 C                              o fZon: Zonal flux term  
 C                              o fMer: Meridional flux term  
 C                              o fVer: Vertical flux term - note fVer  
69  C                                      is "pipelined" in the vertical  C                                      is "pipelined" in the vertical
70  C                                      so we need an fVer for each  C                                      so we need an fVer for each
71  C                                      variable.  C                                      variable.
72  C     rhoK, rhoKM1   - Density at current level, level above and level  C     rhoK, rhoKM1   - Density at current level, and level above
 C                      below.  
 C     rhoKP1                                                                    
 C     buoyK, buoyKM1 - Buoyancy at current level and level above.  
73  C     phiHyd         - Hydrostatic part of the potential phiHydi.  C     phiHyd         - Hydrostatic part of the potential phiHydi.
74  C                      In z coords phiHydiHyd is the hydrostatic  C                      In z coords phiHydiHyd is the hydrostatic
75  C                      pressure anomaly  C                      Potential (=pressure/rho0) anomaly
76  C                      In p coords phiHydiHyd is the geopotential  C                      In p coords phiHydiHyd is the geopotential
77  C                      surface height  C                      surface height anomaly.
78  C                      anomaly.  C     phiSurfX, - gradient of Surface potentiel (Pressure/rho, ocean)
79  C     etaSurfX,      - Holds surface elevation gradient in X and Y.  C     phiSurfY             or geopotentiel (atmos) in X and Y direction
 C     etaSurfY  
80  C     KappaRT,       - Total diffusion in vertical for T and S.  C     KappaRT,       - Total diffusion in vertical for T and S.
81  C     KappaRS          (background + spatially varying, isopycnal term).  C     KappaRS          (background + spatially varying, isopycnal term).
82  C     iMin, iMax     - Ranges and sub-block indices on which calculations  C     iMin, iMax     - Ranges and sub-block indices on which calculations
# Line 94  C     bi, bj Line 85  C     bi, bj
85  C     k, kup,        - Index for layer above and below. kup and kDown  C     k, kup,        - Index for layer above and below. kup and kDown
86  C     kDown, km1       are switched with layer to be the appropriate  C     kDown, km1       are switched with layer to be the appropriate
87  C                      index into fVerTerm.  C                      index into fVerTerm.
88    C     tauAB - Adams-Bashforth timestepping weight: 0=forward ; 1/2=Adams-Bashf.
89        _RS xA      (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RS xA      (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
90        _RS yA      (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RS yA      (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
91        _RL uTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL uTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
92        _RL vTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL vTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
93        _RL rTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
       _RL rVel    (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)  
       _RS maskC   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
94        _RS maskUp  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RS maskUp  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
       _RL aTerm   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
       _RL xTerm   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
       _RL cTerm   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
       _RL mTerm   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
       _RL pTerm   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
       _RL fZon    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
       _RL fMer    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
95        _RL fVerT   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerT   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
96        _RL fVerS   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerS   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
97          _RL fVerTr1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
98        _RL fVerU   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerU   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
99        _RL fVerV   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerV   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
100        _RL phiHyd  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL phiHyd  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
101        _RL rhokm1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rhokm1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
       _RL rhokp1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
102        _RL rhok    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rhok    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
103        _RL buoyKM1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL phiSurfX(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
104        _RL buoyK   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL phiSurfY(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
       _RL rhotmp  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
       _RL etaSurfX(1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
       _RL etaSurfY(1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
105        _RL KappaRT (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL KappaRT (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
106        _RL KappaRS (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL KappaRS (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
107        _RL KappaRU (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL KappaRU (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
# Line 129  C                      index into fVerTe Line 109  C                      index into fVerTe
109        _RL sigmaX  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaX  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
110        _RL sigmaY  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaY  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
111        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
112          _RL tauAB
113    
114  C This is currently also used by IVDC and Diagnostics  C This is currently used by IVDC and Diagnostics
 C #ifdef INCLUDE_CONVECT_CALL  
115        _RL ConvectCount (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL ConvectCount (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
 C #endif  
116    
117        INTEGER iMin, iMax        INTEGER iMin, iMax
118        INTEGER jMin, jMax        INTEGER jMin, jMax
119        INTEGER bi, bj        INTEGER bi, bj
120        INTEGER i, j        INTEGER i, j
121        INTEGER k, km1, kup, kDown        INTEGER k, km1, kup, kDown
       LOGICAL BOTTOM_LAYER  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
       INTEGER    isbyte  
       PARAMETER( isbyte = 4 )  
   
       INTEGER act1, act2, act3, act4  
       INTEGER max1, max2, max3  
       INTEGER iikey, kkey  
       INTEGER maximpl  
 #endif /* ALLOW_AUTODIFF_TAMC */  
122    
123    Cjmc : add for phiHyd output <- but not working if multi tile per CPU
124    c     CHARACTER*(MAX_LEN_MBUF) suff
125    c     LOGICAL  DIFFERENT_MULTIPLE
126    c     EXTERNAL DIFFERENT_MULTIPLE
127    Cjmc(end)
128    
129  C---    The algorithm...  C---    The algorithm...
130  C  C
131  C       "Correction Step"  C       "Correction Step"
# Line 166  C       "Calculation of Gs" Line 140  C       "Calculation of Gs"
140  C       ===================  C       ===================
141  C       This is where all the accelerations and tendencies (ie.  C       This is where all the accelerations and tendencies (ie.
142  C       physics, parameterizations etc...) are calculated  C       physics, parameterizations etc...) are calculated
 C         rVel = sum_r ( div. u[n] )  
143  C         rho = rho ( theta[n], salt[n] )  C         rho = rho ( theta[n], salt[n] )
144  C         b   = b(rho, theta)  C         b   = b(rho, theta)
145  C         K31 = K31 ( rho )  C         K31 = K31 ( rho )
146  C         Gu[n] = Gu( u[n], v[n], rVel, b, ... )  C         Gu[n] = Gu( u[n], v[n], wVel, b, ... )
147  C         Gv[n] = Gv( u[n], v[n], rVel, b, ... )  C         Gv[n] = Gv( u[n], v[n], wVel, b, ... )
148  C         Gt[n] = Gt( theta[n], u[n], v[n], rVel, K31, ... )  C         Gt[n] = Gt( theta[n], u[n], v[n], wVel, K31, ... )
149  C         Gs[n] = Gs( salt[n], u[n], v[n], rVel, K31, ... )  C         Gs[n] = Gs( salt[n], u[n], v[n], wVel, K31, ... )
150  C  C
151  C       "Time-stepping" or "Prediction"  C       "Time-stepping" or "Prediction"
152  C       ================================  C       ================================
# Line 202  C--   dummy statement to end declaration Line 175  C--   dummy statement to end declaration
175        ikey = 1        ikey = 1
176  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
177    
   
178  C--   Set up work arrays with valid (i.e. not NaN) values  C--   Set up work arrays with valid (i.e. not NaN) values
179  C     These inital values do not alter the numerical results. They  C     These inital values do not alter the numerical results. They
180  C     just ensure that all memory references are to valid floating  C     just ensure that all memory references are to valid floating
# Line 214  C     uninitialised but inert locations. Line 186  C     uninitialised but inert locations.
186          yA(i,j)      = 0. _d 0          yA(i,j)      = 0. _d 0
187          uTrans(i,j)  = 0. _d 0          uTrans(i,j)  = 0. _d 0
188          vTrans(i,j)  = 0. _d 0          vTrans(i,j)  = 0. _d 0
         aTerm(i,j)   = 0. _d 0  
         xTerm(i,j)   = 0. _d 0  
         cTerm(i,j)   = 0. _d 0  
         mTerm(i,j)   = 0. _d 0  
         pTerm(i,j)   = 0. _d 0  
         fZon(i,j)    = 0. _d 0  
         fMer(i,j)    = 0. _d 0  
189          DO k=1,Nr          DO k=1,Nr
190           phiHyd (i,j,k)  = 0. _d 0           phiHyd(i,j,k)  = 0. _d 0
191           KappaRU(i,j,k) = 0. _d 0           KappaRU(i,j,k) = 0. _d 0
192           KappaRV(i,j,k) = 0. _d 0           KappaRV(i,j,k) = 0. _d 0
193           sigmaX(i,j,k) = 0. _d 0           sigmaX(i,j,k) = 0. _d 0
# Line 231  C     uninitialised but inert locations. Line 196  C     uninitialised but inert locations.
196          ENDDO          ENDDO
197          rhoKM1 (i,j) = 0. _d 0          rhoKM1 (i,j) = 0. _d 0
198          rhok   (i,j) = 0. _d 0          rhok   (i,j) = 0. _d 0
199          rhoKP1 (i,j) = 0. _d 0          phiSurfX(i,j) = 0. _d 0
200          rhoTMP (i,j) = 0. _d 0          phiSurfY(i,j) = 0. _d 0
         buoyKM1(i,j) = 0. _d 0  
         buoyK  (i,j) = 0. _d 0  
         maskC  (i,j) = 0. _d 0  
201         ENDDO         ENDDO
202        ENDDO        ENDDO
203    
# Line 249  CHPF$ INDEPENDENT Line 211  CHPF$ INDEPENDENT
211    
212  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
213  C--    HPF directive to help TAMC  C--    HPF directive to help TAMC
214  CHPF$  INDEPENDENT, NEW (rTrans,rVel,fVerT,fVerS,fVerU,fVerV  CHPF$  INDEPENDENT, NEW (rTrans,fVerT,fVerS,fVerU,fVerV
215  CHPF$&                  ,phiHyd,utrans,vtrans,maskc,xA,yA  CHPF$&                  ,phiHyd,utrans,vtrans,xA,yA
216  CHPF$&                  ,KappaRT,KappaRS,KappaRU,KappaRV  CHPF$&                  ,KappaRT,KappaRS,KappaRU,KappaRV
217  CHPF$&                  )  CHPF$&                  )
218  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
# Line 277  CHPF$&                  ) Line 239  CHPF$&                  )
239  C--     Set up work arrays that need valid initial values  C--     Set up work arrays that need valid initial values
240          DO j=1-OLy,sNy+OLy          DO j=1-OLy,sNy+OLy
241           DO i=1-OLx,sNx+OLx           DO i=1-OLx,sNx+OLx
242            rTrans(i,j)   = 0. _d 0            rTrans (i,j)   = 0. _d 0
243            rVel  (i,j,1) = 0. _d 0            fVerT  (i,j,1) = 0. _d 0
244            rVel  (i,j,2) = 0. _d 0            fVerT  (i,j,2) = 0. _d 0
245            fVerT (i,j,1) = 0. _d 0            fVerS  (i,j,1) = 0. _d 0
246            fVerT (i,j,2) = 0. _d 0            fVerS  (i,j,2) = 0. _d 0
247            fVerS (i,j,1) = 0. _d 0            fVerTr1(i,j,1) = 0. _d 0
248            fVerS (i,j,2) = 0. _d 0            fVerTr1(i,j,2) = 0. _d 0
249            fVerU (i,j,1) = 0. _d 0            fVerU  (i,j,1) = 0. _d 0
250            fVerU (i,j,2) = 0. _d 0            fVerU  (i,j,2) = 0. _d 0
251            fVerV (i,j,1) = 0. _d 0            fVerV  (i,j,1) = 0. _d 0
252            fVerV (i,j,2) = 0. _d 0            fVerV  (i,j,2) = 0. _d 0
           phiHyd(i,j,1) = 0. _d 0  
253           ENDDO           ENDDO
254          ENDDO          ENDDO
255    
256          DO k=1,Nr          DO k=1,Nr
257           DO j=1-OLy,sNy+OLy           DO j=1-OLy,sNy+OLy
258            DO i=1-OLx,sNx+OLx            DO i=1-OLx,sNx+OLx
259  #ifdef INCLUDE_CONVECT_CALL  C This is currently also used by IVDC and Diagnostics
260             ConvectCount(i,j,k) = 0.             ConvectCount(i,j,k) = 0.
 #endif  
261             KappaRT(i,j,k) = 0. _d 0             KappaRT(i,j,k) = 0. _d 0
262             KappaRS(i,j,k) = 0. _d 0             KappaRS(i,j,k) = 0. _d 0
263            ENDDO            ENDDO
# Line 309  C--     Set up work arrays that need val Line 269  C--     Set up work arrays that need val
269          jMin = 1-OLy+1          jMin = 1-OLy+1
270          jMax = sNy+OLy          jMax = sNy+OLy
271    
         k = 1  
         BOTTOM_LAYER = k .EQ. Nr  
   
 #ifdef DO_PIPELINED_CORRECTION_STEP  
 C--     Calculate gradient of surface pressure  
         CALL CALC_GRAD_ETA_SURF(  
      I       bi,bj,iMin,iMax,jMin,jMax,  
      O       etaSurfX,etaSurfY,  
      I       myThid)  
 C--     Update fields in top level according to tendency terms  
         CALL CORRECTION_STEP(  
      I       bi,bj,iMin,iMax,jMin,jMax,k,  
      I       etaSurfX,etaSurfY,myTime,myThid)  
   
 #ifdef ALLOW_OBCS  
         IF (openBoundaries) THEN  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE uvel (:,:,k,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  
 CADJ STORE vvel (:,:,k,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  
 CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  
 CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
            CALL APPLY_OBCS1( bi, bj, k, myThid )  
         END IF  
 #endif  
   
         IF ( .NOT. BOTTOM_LAYER ) THEN  
 C--      Update fields in layer below according to tendency terms  
          CALL CORRECTION_STEP(  
      I        bi,bj,iMin,iMax,jMin,jMax,k+1,  
      I        etaSurfX,etaSurfY,myTime,myThid)  
 #ifdef ALLOW_OBCS  
          IF (openBoundaries) THEN  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE uvel (:,:,k,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  
 CADJ STORE vvel (:,:,k,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  
 CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  
 CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
             CALL APPLY_OBCS1( bi, bj, k+1, myThid )  
          END IF  
 #endif  
         ENDIF  
 #endif  
272    
 C--     Density of 1st level (below W(1)) reference to level 1  
 #ifdef  INCLUDE_FIND_RHO_CALL  
273  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
274  CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
275  CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
276  #endif /* ALLOW_AUTODIFF_TAMC */  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
277          CALL FIND_RHO(  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
278       I     bi, bj, iMin, iMax, jMin, jMax, k, k, eosType,  #ifdef ALLOW_PASSIVE_TRACER
279       O     rhoKm1,  CADJ STORE tr1  (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
      I     myThid )  
 #endif  
   
         IF (.NOT. BOTTOM_LAYER) THEN  
   
 C--      Check static stability with layer below  
 C--      and mix as needed.  
 #ifdef  INCLUDE_FIND_RHO_CALL  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE theta(:,:,k+1,bi,bj) = comlev1_bibj  
 CADJ &   , key = ikey, byte = isbyte  
 CADJ STORE salt (:,:,k+1,bi,bj) = comlev1_bibj  
 CADJ &   , key = ikey, byte = isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
          CALL FIND_RHO(  
      I      bi, bj, iMin, iMax, jMin, jMax, k+1, k, eosType,  
      O      rhoKp1,  
      I      myThid )  
280  #endif  #endif
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE rhoKm1(:,:)  = comlev1_bibj, key = ikey, byte = isbyte  
 CADJ STORE rhoKp1(:,:)  = comlev1_bibj, key = ikey, byte = isbyte  
281  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
282    
283  #ifdef  INCLUDE_CONVECT_CALL  C--     Start of diagnostic loop
284            DO k=Nr,1,-1
285           CALL CONVECT(  
286       I       bi,bj,iMin,iMax,jMin,jMax,k+1,rhoKm1,rhoKp1,  #ifdef ALLOW_AUTODIFF_TAMC
287       U       ConvectCount,  C? Patrick, is this formula correct now that we change the loop range?
288       I       myTime,myIter,myThid)  C? Do we still need this?
289    cph kkey formula corrected.
290  #ifdef ALLOW_AUTODIFF_TAMC  cph Needed for rhok, rhokm1, in the case useGMREDI.
291  CADJ STORE theta(:,:,k+1,bi,bj),theta(:,:,k,bi,bj)           kkey = (ikey-1)*Nr + k
292  CADJ &     = comlev1_bibj, key = ikey, byte = isbyte  CADJ STORE rhokm1(:,:) = comlev1_bibj_k ,       key=kkey, byte=isbyte
293  CADJ STORE salt (:,:,k+1,bi,bj),salt (:,:,k,bi,bj)  CADJ STORE rhok  (:,:) = comlev1_bibj_k ,       key=kkey, byte=isbyte
294  CADJ &     = comlev1_bibj, key = ikey, byte = isbyte  #endif /* ALLOW_AUTODIFF_TAMC */
295  #endif /* ALLOW_AUTODIFF_TAMC */  
296    C--       Integrate continuity vertically for vertical velocity
297  #endif            CALL INTEGRATE_FOR_W(
298         I                         bi, bj, k, uVel, vVel,
299  C--      Implicit Vertical Diffusion for Convection       O                         wVel,
300           IF (ivdc_kappa.NE.0.) THEN       I                         myThid )
301              CALL CALC_IVDC(  
302       I       bi,bj,iMin,iMax,jMin,jMax,k+1,rhoKm1,rhoKp1,  #ifdef    ALLOW_OBCS
303       U       ConvectCount, KappaRT, KappaRS,  #ifdef    ALLOW_NONHYDROSTATIC
304       I       myTime,myIter,myThid)  C--       Apply OBC to W if in N-H mode
305           ENDIF            IF (useOBCS.AND.nonHydrostatic) THEN
306                CALL OBCS_APPLY_W( bi, bj, k, wVel, myThid )
307  C--      Recompute density after mixing            ENDIF
308  #ifdef  INCLUDE_FIND_RHO_CALL  #endif    /* ALLOW_NONHYDROSTATIC */
309           CALL FIND_RHO(  #endif    /* ALLOW_OBCS */
310       I      bi, bj, iMin, iMax, jMin, jMax, k, k, eosType,  
311       O      rhoKm1,  C--       Calculate gradients of potential density for isoneutral
312       I      myThid )  C         slope terms (e.g. GM/Redi tensor or IVDC diffusivity)
313  #endif  c         IF ( k.GT.1 .AND. (useGMRedi.OR.ivdc_kappa.NE.0.) ) THEN
314          ENDIF            IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.) ) THEN
315    #ifdef ALLOW_AUTODIFF_TAMC
316  C--     Calculate buoyancy  CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
317          CALL CALC_BUOYANCY(  CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
318       I      bi,bj,iMin,iMax,jMin,jMax,k,rhoKm1,  #endif /* ALLOW_AUTODIFF_TAMC */
319       O      buoyKm1,              CALL FIND_RHO(
320       I      myThid )       I        bi, bj, iMin, iMax, jMin, jMax, k, k, eosType,
321         I        theta, salt,
322  C--     Integrate hydrostatic balance for phiHyd with BC of       O        rhoK,
 C--     phiHyd(z=0)=0  
         CALL CALC_PHI_HYD(  
      I      bi,bj,iMin,iMax,jMin,jMax,k,buoyKm1,buoyKm1,  
      U      phiHyd,  
      I      myThid )  
   
 #ifdef ALLOW_GMREDI  
         IF ( useGMRedi ) THEN  
         CALL GRAD_SIGMA(  
      I            bi, bj, iMin, iMax, jMin, jMax, k,  
      I            rhoKm1, rhoKm1, rhoKm1,  
      O            sigmaX, sigmaY, sigmaR,  
      I            myThid )  
         ELSE  
          DO j=1-OLy,sNy+OLy  
           DO i=1-OLx,sNx+OLx  
            sigmaX(i,j,k) = 0. _d 0  
            sigmaY(i,j,k) = 0. _d 0  
            sigmaR(i,j,k) = 0. _d 0  
           ENDDO  
          ENDDO  
         ENDIF  
 #endif  
   
 C--     Start of downward loop  
         DO k=2,Nr  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
          kkey = (ikey-1)*(Nr-2+1) + (k-2) + 1  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
          BOTTOM_LAYER = k .EQ. Nr  
   
 #ifdef DO_PIPELINED_CORRECTION_STEP  
          IF ( .NOT. BOTTOM_LAYER ) THEN  
 C--       Update fields in layer below according to tendency terms  
           CALL CORRECTION_STEP(  
      I         bi,bj,iMin,iMax,jMin,jMax,k+1,  
      I         etaSurfX,etaSurfY,myTime,myThid)  
 #ifdef ALLOW_OBCS  
           IF (openBoundaries) THEN  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE uvel (:,:,k,bi,bj) = comlev1_bibj_k  
 CADJ &   , key = kkey, byte = isbyte  
 CADJ STORE vvel (:,:,k,bi,bj) = comlev1_bibj_k  
 CADJ &   , key = kkey, byte = isbyte  
 CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k  
 CADJ &   , key = kkey, byte = isbyte  
 CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj_k  
 CADJ &   , key = kkey, byte = isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
              CALL APPLY_OBCS1( bi, bj, k+1, myThid )  
           END IF  
 #endif  
          ENDIF  
 #endif /* DO_PIPELINED_CORRECTION_STEP */  
   
 C--      Density of k level (below W(k)) reference to k level  
 #ifdef  INCLUDE_FIND_RHO_CALL  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k  
 CADJ &   , key = kkey, byte = isbyte  
 CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj_k  
 CADJ &   , key = kkey, byte = isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
          CALL FIND_RHO(  
      I      bi, bj, iMin, iMax, jMin, jMax,  k, k, eosType,  
      O      rhoK,  
      I      myThid )  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 cph(   storing not necessary  
 cphCADJ STORE rhoK(:,:)  = comlev1_bibj_k, key = kkey, byte = isbyte  
 cph)  
 #endif /* ALLOW_AUTODIFF_TAMC */  
 #endif  
   
          IF (.NOT. BOTTOM_LAYER) THEN  
   
 C--       Check static stability with layer below and mix as needed.  
 C--       Density of k+1 level (below W(k+1)) reference to k level.  
 #ifdef  INCLUDE_FIND_RHO_CALL  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE theta(:,:,k+1,bi,bj) = comlev1_bibj_k  
 CADJ &   , key = kkey, byte = isbyte  
 CADJ STORE salt (:,:,k+1,bi,bj) = comlev1_bibj_k  
 CADJ &   , key = kkey, byte = isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
           CALL FIND_RHO(  
      I       bi, bj, iMin, iMax, jMin, jMax,  k+1, k, eosType,  
      O       rhoKp1,  
      I       myThid )  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE rhoKp1(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
 #endif  
   
 #ifdef  INCLUDE_CONVECT_CALL  
           CALL CONVECT(  
      I        bi,bj,iMin,iMax,jMin,jMax,k+1,rhoK,rhoKp1,  
      U        ConvectCount,  
      I        myTime,myIter,myThid)  
   
 #endif  
   
 C--      Implicit Vertical Diffusion for Convection  
          IF (ivdc_kappa.NE.0.) THEN  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE rhoKm1(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
             CALL CALC_IVDC(  
      I       bi,bj,iMin,iMax,jMin,jMax,k+1,rhoKm1,rhoKp1,  
      U       ConvectCount, KappaRT, KappaRS,  
      I       myTime,myIter,myThid)  
          END IF  
   
 C--       Recompute density after mixing  
 #ifdef  INCLUDE_FIND_RHO_CALL  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k  
 CADJ &   , key = kkey, byte = isbyte  
 CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj_k  
 CADJ &   , key = kkey, byte = isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
           CALL FIND_RHO(  
      I       bi, bj, iMin, iMax, jMin, jMax, k, k, eosType,  
      O       rhoK,  
      I       myThid )  
 #endif  
   
 C--            IF (.NOT. BOTTOM_LAYER) ends here  
          ENDIF  
   
 C--      Calculate buoyancy  
          CALL CALC_BUOYANCY(  
      I       bi,bj,iMin,iMax,jMin,jMax,k,rhoK,  
      O       buoyK,  
      I       myThid )  
   
 C--      Integrate hydrostatic balance for phiHyd with BC of  
 C--      phiHyd(z=0)=0  
          CALL CALC_PHI_HYD(  
      I        bi,bj,iMin,iMax,jMin,jMax,k,buoyKm1,buoyK,  
      U        phiHyd,  
323       I        myThid )       I        myThid )
324                IF (k.GT.1) THEN
 #ifdef  INCLUDE_FIND_RHO_CALL  
 C--      Calculate iso-neutral slopes for the GM/Redi parameterisation  
   
325  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
326  CADJ STORE theta(:,:,k-1,bi,bj) = comlev1_bibj_k  CADJ STORE theta(:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
327  CADJ &   , key = kkey, byte = isbyte  CADJ STORE salt (:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
 CADJ STORE salt (:,:,k-1,bi,bj) = comlev1_bibj_k  
 CADJ &   , key = kkey, byte = isbyte  
328  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
329                 CALL FIND_RHO(
          CALL FIND_RHO(  
330       I        bi, bj, iMin, iMax, jMin, jMax, k-1, k, eosType,       I        bi, bj, iMin, iMax, jMin, jMax, k-1, k, eosType,
331       O        rhoTmp,       I        theta, salt,
332         O        rhoKm1,
333       I        myThid )       I        myThid )
334  #endif              ENDIF
335                CALL GRAD_SIGMA(
   
 #ifdef ALLOW_GMREDI  
          IF ( useGMRedi ) THEN  
          CALL GRAD_SIGMA(  
336       I             bi, bj, iMin, iMax, jMin, jMax, k,       I             bi, bj, iMin, iMax, jMin, jMax, k,
337       I             rhoK, rhotmp, rhoK,       I             rhoK, rhoKm1, rhoK,
338       O             sigmaX, sigmaY, sigmaR,       O             sigmaX, sigmaY, sigmaR,
339       I             myThid )       I             myThid )
340           ELSE            ENDIF
           DO j=1-OLy,sNy+OLy  
            DO i=1-OLx,sNx+OLx  
             sigmaX(i,j,k) = 0. _d 0  
             sigmaY(i,j,k) = 0. _d 0  
             sigmaR(i,j,k) = 0. _d 0  
            ENDDO  
           ENDDO  
          ENDIF  
 #endif  
341    
342           DO J=jMin,jMax  C--       Implicit Vertical Diffusion for Convection
343            DO I=iMin,iMax  c ==> should use sigmaR !!!
344  #ifdef  INCLUDE_FIND_RHO_CALL            IF (k.GT.1 .AND. ivdc_kappa.NE.0.) THEN
345             rhoKm1 (I,J) = rhoK(I,J)              CALL CALC_IVDC(
346  #endif       I        bi, bj, iMin, iMax, jMin, jMax, k,
347             buoyKm1(I,J) = buoyK(I,J)       I        rhoKm1, rhoK,
348            ENDDO       U        ConvectCount, KappaRT, KappaRS,
349           ENDDO       I        myTime, myIter, myThid)
350              ENDIF
351    
352  C--     end of k loop  C--     end of diagnostic k loop (Nr:1)
353          ENDDO          ENDDO
354    
355  C     Determines forcing terms based on external fields  #ifdef ALLOW_AUTODIFF_TAMC
356  C     relaxation terms, etc.  cph avoids recomputation of integrate_for_w
357        CALL EXTERNAL_FORCING_SURF(  CADJ STORE wvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
358    #endif /* ALLOW_AUTODIFF_TAMC */
359    
360    #ifdef  ALLOW_OBCS
361    C--     Calculate future values on open boundaries
362            IF (useOBCS) THEN
363              CALL OBCS_CALC( bi, bj, myTime+deltaT,
364         I            uVel, vVel, wVel, theta, salt,
365         I            myThid )
366            ENDIF
367    #endif  /* ALLOW_OBCS */
368    
369    C--     Determines forcing terms based on external fields
370    C       relaxation terms, etc.
371            CALL EXTERNAL_FORCING_SURF(
372       I             bi, bj, iMin, iMax, jMin, jMax,       I             bi, bj, iMin, iMax, jMin, jMax,
373       I             myThid )       I             myThid )
   
374  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
375    cph needed for KPP
376    CADJ STORE surfacetendencyU(:,:,bi,bj)
377    CADJ &     = comlev1_bibj, key=ikey, byte=isbyte
378    CADJ STORE surfacetendencyV(:,:,bi,bj)
379    CADJ &     = comlev1_bibj, key=ikey, byte=isbyte
380    CADJ STORE surfacetendencyS(:,:,bi,bj)
381    CADJ &     = comlev1_bibj, key=ikey, byte=isbyte
382    CADJ STORE surfacetendencyT(:,:,bi,bj)
383    CADJ &     = comlev1_bibj, key=ikey, byte=isbyte
384    #endif /* ALLOW_AUTODIFF_TAMC */
385    
386  CADJ STORE surfacetendencyu(:,:,bi,bj)  #ifdef  ALLOW_GMREDI
 CADJ &   , surfacetendencyv(:,:,bi,bj)  
 CADJ &   , surfacetendencys(:,:,bi,bj)  
 CADJ &   , surfacetendencyt(:,:,bi,bj)  
 CADJ &                        = comlev1_bibj, key=ikey, byte=isbyte  
387    
388  # ifdef ALLOW_GMREDI  #ifdef ALLOW_AUTODIFF_TAMC
389  CADJ STORE sigmaX(:,:,:) = comlev1, key=ikey, byte=isbyte  CADJ STORE sigmaX(:,:,:) = comlev1, key=ikey, byte=isbyte
390  CADJ STORE sigmaY(:,:,:) = comlev1, key=ikey, byte=isbyte  CADJ STORE sigmaY(:,:,:) = comlev1, key=ikey, byte=isbyte
391  CADJ STORE sigmaR(:,:,:) = comlev1, key=ikey, byte=isbyte  CADJ STORE sigmaR(:,:,:) = comlev1, key=ikey, byte=isbyte
 # endif /* ALLOW_GMREDI */  
   
392  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
393    C--     Calculate iso-neutral slopes for the GM/Redi parameterisation
 #ifdef ALLOW_GMREDI  
394          IF (useGMRedi) THEN          IF (useGMRedi) THEN
395            DO k=1, Nr            DO k=1,Nr
396              CALL GMREDI_CALC_TENSOR(              CALL GMREDI_CALC_TENSOR(
397       I             bi, bj, iMin, iMax, jMin, jMax, k,       I             bi, bj, iMin, iMax, jMin, jMax, k,
398       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
# Line 661  CADJ STORE sigmaR(:,:,:) = comlev1, key= Line 408  CADJ STORE sigmaR(:,:,:) = comlev1, key=
408            ENDDO            ENDDO
409  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
410          ENDIF          ENDIF
 #endif  
411    
412  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
413  CADJ STORE KappaRT(:,:,:)     = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE Kwx(:,:,:,bi,bj)   = comlev1_bibj, key=ikey, byte=isbyte
414  CADJ STORE KappaRS(:,:,:)     = comlev1_bibj, key=ikey, byte=isbyte  CADJ STORE Kwy(:,:,:,bi,bj)   = comlev1_bibj, key=ikey, byte=isbyte
415    CADJ STORE Kwz(:,:,:,bi,bj)   = comlev1_bibj, key=ikey, byte=isbyte
 #ifdef ALLOW_GMREDI  
 C-- R.G. We need to define a new tape since Kw use mythid instead of bi,bj  
 CADJ STORE Kwx(:,:,:,myThid)  = comlev1_bibj, key=ikey, byte=isbyte  
 CADJ STORE Kwy(:,:,:,myThid)  = comlev1_bibj, key=ikey, byte=isbyte  
 CADJ STORE Kwz(:,:,:,myThid)  = comlev1_bibj, key=ikey, byte=isbyte  
 #endif  
   
 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  
   
 C--     dummy initialization to break data flow because  
 C--     calc_div_ghat has a condition for initialization  
         DO J=jMin,jMax  
            DO I=iMin,iMax  
               cg2d_b(i,j,bi,bj) = 0.0  
            ENDDO  
         ENDDO  
416  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
417    
418  #ifdef ALLOW_KPP  #endif  /* ALLOW_GMREDI */
 C--   Compute KPP mixing coefficients  
         IF (useKPP) THEN  
419    
420            CALL TIMER_START('KPP_CALC               [DYNAMICS]', myThid)  #ifdef  ALLOW_KPP
421    C--     Compute KPP mixing coefficients
422            IF (useKPP) THEN
423            CALL KPP_CALC(            CALL KPP_CALC(
424       I                  bi, bj, myTime, myThid )       I                  bi, bj, myTime, myThid )
           CALL TIMER_STOP ('KPP_CALC               [DYNAMICS]', myThid)  
   
425  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
426          ELSE          ELSE
427            DO j=1-OLy,sNy+OLy            CALL KPP_CALC_DUMMY(
428              DO i=1-OLx,sNx+OLx       I                  bi, bj, myTime, myThid )
               KPPhbl (i,j,bi,bj) = 1.0  
               KPPfrac(i,j,bi,bj) = 0.0  
               DO k = 1,Nr  
                  KPPghat   (i,j,k,bi,bj) = 0.0  
                  KPPviscAz (i,j,k,bi,bj) = viscAz  
                  KPPdiffKzT(i,j,k,bi,bj) = diffKzT  
                  KPPdiffKzS(i,j,k,bi,bj) = diffKzS  
               ENDDO  
             ENDDO  
           ENDDO  
429  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
430          ENDIF          ENDIF
431    
# Line 723  CADJ &   , KPPfrac   (:,:  ,bi,bj) Line 438  CADJ &   , KPPfrac   (:,:  ,bi,bj)
438  CADJ &                 = comlev1_bibj, key=ikey, byte=isbyte  CADJ &                 = comlev1_bibj, key=ikey, byte=isbyte
439  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
440    
441  #endif /* ALLOW_KPP */  #endif  /* ALLOW_KPP */
442    
443  C--     Start of upward loop  #ifdef ALLOW_AUTODIFF_TAMC
444          DO k = Nr, 1, -1  CADJ STORE KappaRT(:,:,:)     = comlev1_bibj, key = ikey, byte = isbyte
445    CADJ STORE KappaRS(:,:,:)     = comlev1_bibj, key = ikey, byte = isbyte
446    CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
447    CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
448    CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
449    CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
450    #ifdef ALLOW_PASSIVE_TRACER
451    CADJ STORE tr1  (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
452    #endif
453    #endif /* ALLOW_AUTODIFF_TAMC */
454    
455  C--      km1    Points to level above k (=k-1)  #ifdef ALLOW_AIM
456  C--      kup    Cycles through 1,2 to point to layer above  C       AIM - atmospheric intermediate model, physics package code.
457  C--      kDown  Cycles through 2,1 to point to current layer  C note(jmc) : phiHyd=0 at this point but is not really used in Molteni Physics
458            IF ( useAIM ) THEN
459           km1  =max(1,k-1)           CALL TIMER_START('AIM_DO_ATMOS_PHYS      [DYNAMICS]', myThid)
460           kup  =1+MOD(k+1,2)           CALL AIM_DO_ATMOS_PHYSICS( phiHyd, bi, bj, myTime, myThid )
461           kDown=1+MOD(k,2)           CALL TIMER_STOP ('AIM_DO_ATMOS_PHYS      [DYNAMICS]', myThid)
462            ENDIF
463           iMin = 1-OLx+2  #endif /* ALLOW_AIM */
          iMax = sNx+OLx-1  
          jMin = 1-OLy+2  
          jMax = sNy+OLy-1  
464    
 #ifdef ALLOW_AUTODIFF_TAMC  
          kkey = (ikey-1)*(Nr-1+1) + (k-1) + 1  
 CADJ STORE rvel  (:,:,kdown)  = comlev1_bibj_k, key=kkey, byte=isbyte  
 CADJ STORE rTrans(:,:)        = comlev1_bibj_k, key=kkey, byte=isbyte  
 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 */  
465    
466  C--      Get temporary terms used by tendency routines  C--     Start of thermodynamics loop
467           CALL CALC_COMMON_FACTORS (          DO k=Nr,1,-1
468       I        bi,bj,iMin,iMax,jMin,jMax,k,km1,kup,kDown,  #ifdef ALLOW_AUTODIFF_TAMC
469       O        xA,yA,uTrans,vTrans,rTrans,rVel,maskC,maskUp,  C? Patrick Is this formula correct?
470       I        myThid)  cph Yes, but I rewrote it.
471    cph Also, the KappaR? need the index and subscript k!
472             kkey = (ikey-1)*Nr + k
473    #endif /* ALLOW_AUTODIFF_TAMC */
474    
475    C--       km1    Points to level above k (=k-1)
476    C--       kup    Cycles through 1,2 to point to layer above
477    C--       kDown  Cycles through 2,1 to point to current layer
478    
479              km1  = MAX(1,k-1)
480              kup  = 1+MOD(k+1,2)
481              kDown= 1+MOD(k,2)
482    
483              iMin = 1-OLx
484              iMax = sNx+OLx
485              jMin = 1-OLy
486              jMax = sNy+OLy
487    
488    C--       Get temporary terms used by tendency routines
489              CALL CALC_COMMON_FACTORS (
490         I         bi,bj,iMin,iMax,jMin,jMax,k,
491         O         xA,yA,uTrans,vTrans,rTrans,maskUp,
492         I         myThid)
493    
494  #ifdef ALLOW_OBCS  #ifdef ALLOW_AUTODIFF_TAMC
495          IF (openBoundaries) THEN  CADJ STORE KappaRT(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte
496           CALL APPLY_OBCS3( bi, bj, k, kup, rTrans, rVel, myThid )  CADJ STORE KappaRS(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte
497          ENDIF  #endif /* ALLOW_AUTODIFF_TAMC */
 #endif  
498    
499  #ifdef  INCLUDE_CALC_DIFFUSIVITY_CALL  #ifdef  INCLUDE_CALC_DIFFUSIVITY_CALL
500  C--      Calculate the total vertical diffusivity  C--      Calculate the total vertical diffusivity
501           CALL CALC_DIFFUSIVITY(           CALL CALC_DIFFUSIVITY(
502       I        bi,bj,iMin,iMax,jMin,jMax,k,       I        bi,bj,iMin,iMax,jMin,jMax,k,
503       I        maskC,maskUp,       I        maskUp,
504       O        KappaRT,KappaRS,KappaRU,KappaRV,       O        KappaRT,KappaRS,KappaRU,KappaRV,
505       I        myThid)       I        myThid)
506  #endif  #endif
507  C--      Calculate accelerations in the momentum equations  
508           IF ( momStepping ) THEN            iMin = 1-OLx+2
509            CALL CALC_MOM_RHS(            iMax = sNx+OLx-1
510       I         bi,bj,iMin,iMax,jMin,jMax,k,km1,kup,kDown,            jMin = 1-OLy+2
511       I         xA,yA,uTrans,vTrans,rTrans,rVel,maskC,            jMax = sNy+OLy-1
512       I         phiHyd,KappaRU,KappaRV,  
513       U         aTerm,xTerm,cTerm,mTerm,pTerm,  C--      Calculate active tracer tendencies (gT,gS,...)
514       U         fZon, fMer, fVerU, fVerV,  C        and step forward storing result in gTnm1, gSnm1, etc.
      I         myTime, myThid)  
 #ifdef ALLOW_AUTODIFF_TAMC  
 #ifdef INCLUDE_CD_CODE  
          ELSE  
             DO j=1-OLy,sNy+OLy  
                DO i=1-OLx,sNx+OLx  
                   guCD(i,j,k,bi,bj) = 0.0  
                   gvCD(i,j,k,bi,bj) = 0.0  
                END DO  
             END DO  
 #endif  
 #endif /* ALLOW_AUTODIFF_TAMC */  
          ENDIF  
 C--      Calculate active tracer tendencies  
515           IF ( tempStepping ) THEN           IF ( tempStepping ) THEN
516            CALL CALC_GT(             CALL CALC_GT(
517       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
518       I         xA,yA,uTrans,vTrans,rTrans,maskUp,maskC,       I         xA,yA,uTrans,vTrans,rTrans,maskUp,
519       I         KappaRT,       I         KappaRT,
520       U         aTerm,xTerm,fZon,fMer,fVerT,       U         fVerT,
521       I         myTime, myThid)       I         myTime, myThid)
522               tauAB = 0.5d0 + abEps
523               CALL TIMESTEP_TRACER(
524         I         bi,bj,iMin,iMax,jMin,jMax,k,tauAB,
525         I         theta, gT,
526         U         gTnm1,
527         I         myIter, myThid)
528           ENDIF           ENDIF
529           IF ( saltStepping ) THEN           IF ( saltStepping ) THEN
530            CALL CALC_GS(             CALL CALC_GS(
531       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
532       I         xA,yA,uTrans,vTrans,rTrans,maskUp,maskC,       I         xA,yA,uTrans,vTrans,rTrans,maskUp,
533       I         KappaRS,       I         KappaRS,
534       U         aTerm,xTerm,fZon,fMer,fVerS,       U         fVerS,
535       I         myTime, myThid)       I         myTime, myThid)
536               tauAB = 0.5d0 + abEps
537               CALL TIMESTEP_TRACER(
538         I         bi,bj,iMin,iMax,jMin,jMax,k,tauAB,
539         I         salt, gS,
540         U         gSnm1,
541         I         myIter, myThid)
542           ENDIF           ENDIF
543  #ifdef ALLOW_OBCS  #ifdef ALLOW_PASSIVE_TRACER
544  C--      Calculate future values on open boundaries           IF ( tr1Stepping ) THEN
545           IF (openBoundaries) THEN             CALL CALC_GTR1(
546  Caja      CALL CYCLE_OBCS( k, bi, bj, myThid )       I         bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
547            CALL SET_OBCS( k, bi, bj, myTime+deltaTclock, myThid )       I         xA,yA,uTrans,vTrans,rTrans,maskUp,
548         I         KappaRT,
549         U         fVerTr1,
550         I         myTime, myThid)
551               tauAB = 0.5d0 + abEps
552               CALL TIMESTEP_TRACER(
553         I         bi,bj,iMin,iMax,jMin,jMax,k,tauAB,
554         I         Tr1, gTr1,
555         U         gTr1NM1,
556         I         myIter, myThid)
557           ENDIF           ENDIF
558  #endif  #endif
 C--      Prediction step (step forward all model variables)  
          CALL TIMESTEP(  
      I       bi,bj,iMin,iMax,jMin,jMax,k,  
      I       myIter, myThid)  
 #ifdef ALLOW_OBCS  
 C--      Apply open boundary conditions  
          IF (openBoundaries) THEN  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE gunm1(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  
 CADJ STORE gvnm1(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  
 CADJ STORE gwnm1(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
559    
560              CALL APPLY_OBCS2( bi, bj, k, myThid )  #ifdef   ALLOW_OBCS
561    C--      Apply open boundary conditions
562             IF (useOBCS) THEN
563               CALL OBCS_APPLY_TS( bi, bj, k, gTnm1, gSnm1, myThid )
564           END IF           END IF
565  #endif  #endif   /* ALLOW_OBCS */
566    
567  C--      Freeze water  C--      Freeze water
568           IF (allowFreezing) THEN           IF (allowFreezing) THEN
569  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
570  CADJ STORE gTNm1(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE gTNm1(:,:,k,bi,bj) = comlev1_bibj_k
571    CADJ &   , key = kkey, byte = isbyte
572  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
573              CALL FREEZE( bi, bj, iMin, iMax, jMin, jMax, k, myThid )              CALL FREEZE( bi, bj, iMin, iMax, jMin, jMax, k, myThid )
574           END IF           END IF
575    
576  #ifdef DIVG_IN_DYNAMICS  C--     end of thermodynamic k loop (Nr:1)
 C--      Diagnose barotropic divergence of predicted fields  
          CALL CALC_DIV_GHAT(  
      I       bi,bj,iMin,iMax,jMin,jMax,k,  
      I       xA,yA,  
      I       myThid)  
 #endif /* DIVG_IN_DYNAMICS */  
   
 C--      Cumulative diagnostic calculations (ie. time-averaging)  
 #ifdef INCLUDE_DIAGNOSTICS_INTERFACE_CODE  
          IF (taveFreq.GT.0.) THEN  
           CALL DO_TIME_AVERAGES(  
      I                           myTime, myIter, bi, bj, k, kup, kDown,  
      I                           rVel, ConvectCount,  
      I                           myThid )  
          ENDIF  
 #endif  
   
   
 C--     k loop  
577          ENDDO          ENDDO
578    
579    
580  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
581             maximpl = 6  C? Patrick? What about this one?
582             iikey = (ikey-1)*maximpl  cph Keys iikey and idkey don't seem to be needed
583    cph since storing occurs on different tape for each
584    cph impldiff call anyways.
585    cph Thus, common block comlev1_impl isn't needed either.
586    cph Storing below needed in the case useGMREDI.
587            iikey = (ikey-1)*maximpl
588  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
589    
590  C--     Implicit diffusion  C--     Implicit diffusion
# Line 875  CADJ STORE gTNm1(:,:,:,bi,bj) = comlev1_ Line 597  CADJ STORE gTNm1(:,:,:,bi,bj) = comlev1_
597  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
598              CALL IMPLDIFF(              CALL IMPLDIFF(
599       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
600       I         deltaTtracer, KappaRT,recip_HFacC,       I         deltaTtracer, KappaRT, recip_HFacC,
601       U         gTNm1,       U         gTNm1,
602       I         myThid )       I         myThid )
603           END IF           ENDIF
604    
605           IF (saltStepping) THEN           IF (saltStepping) THEN
606  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 887  CADJ STORE gSNm1(:,:,:,bi,bj) = comlev1_ Line 609  CADJ STORE gSNm1(:,:,:,bi,bj) = comlev1_
609  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
610              CALL IMPLDIFF(              CALL IMPLDIFF(
611       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
612       I         deltaTtracer, KappaRS,recip_HFacC,       I         deltaTtracer, KappaRS, recip_HFacC,
613       U         gSNm1,       U         gSNm1,
614       I         myThid )       I         myThid )
615             ENDIF
616    
617    #ifdef ALLOW_PASSIVE_TRACER
618             IF (tr1Stepping) THEN
619    #ifdef ALLOW_AUTODIFF_TAMC
620    CADJ STORE gTr1Nm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte
621    #endif /* ALLOW_AUTODIFF_TAMC */
622              CALL IMPLDIFF(
623         I      bi, bj, iMin, iMax, jMin, jMax,
624         I      deltaTtracer, KappaRT, recip_HFacC,
625         U      gTr1Nm1,
626         I      myThid )
627             ENDIF
628    #endif
629    
630    #ifdef   ALLOW_OBCS
631    C--      Apply open boundary conditions
632             IF (useOBCS) THEN
633               DO K=1,Nr
634                 CALL OBCS_APPLY_TS( bi, bj, k, gTnm1, gSnm1, myThid )
635               ENDDO
636           END IF           END IF
637    #endif   /* ALLOW_OBCS */
638    
639  C--     implicitDiffusion  C--     End If implicitDiffusion
640          ENDIF          ENDIF
641    
642  C--     Implicit viscosity  C--     Start computation of dynamics
643          IF (implicitViscosity) THEN          iMin = 1-OLx+2
644            iMax = sNx+OLx-1
645            jMin = 1-OLy+2
646            jMax = sNy+OLy-1
647    
648    C--     Explicit part of the Surface Potentiel Gradient (add in TIMESTEP)
649    C       (note: this loop will be replaced by CALL CALC_GRAD_ETA)
650            IF (implicSurfPress.NE.1.) THEN
651              CALL CALC_GRAD_PHI_SURF(
652         I         bi,bj,iMin,iMax,jMin,jMax,
653         I         etaN,
654         O         phiSurfX,phiSurfY,
655         I         myThid )                        
656            ENDIF
657    
658           IF (momStepping) THEN  C--     Start of dynamics loop
659  #ifdef ALLOW_AUTODIFF_TAMC          DO k=1,Nr
660           idkey = iikey + 3  
661    C--       km1    Points to level above k (=k-1)
662    C--       kup    Cycles through 1,2 to point to layer above
663    C--       kDown  Cycles through 2,1 to point to current layer
664    
665              km1  = MAX(1,k-1)
666              kup  = 1+MOD(k+1,2)
667              kDown= 1+MOD(k,2)
668    
669    C--      Integrate hydrostatic balance for phiHyd with BC of
670    C        phiHyd(z=0)=0
671    C        distinguishe between Stagger and Non Stagger time stepping
672             IF (staggerTimeStep) THEN
673               CALL CALC_PHI_HYD(
674         I        bi,bj,iMin,iMax,jMin,jMax,k,
675         I        gTnm1, gSnm1,
676         U        phiHyd,
677         I        myThid )
678             ELSE
679               CALL CALC_PHI_HYD(
680         I        bi,bj,iMin,iMax,jMin,jMax,k,
681         I        theta, salt,
682         U        phiHyd,
683         I        myThid )
684             ENDIF
685    
686    C--      Calculate accelerations in the momentum equations (gU, gV, ...)
687    C        and step forward storing the result in gUnm1, gVnm1, etc...
688             IF ( momStepping ) THEN
689               CALL CALC_MOM_RHS(
690         I         bi,bj,iMin,iMax,jMin,jMax,k,kup,kDown,
691         I         phiHyd,KappaRU,KappaRV,
692         U         fVerU, fVerV,
693         I         myTime, myThid)
694               CALL TIMESTEP(
695         I         bi,bj,iMin,iMax,jMin,jMax,k,
696         I         phiHyd, phiSurfX, phiSurfY,
697         I         myIter, myThid)
698    
699    #ifdef   ALLOW_OBCS
700    C--      Apply open boundary conditions
701             IF (useOBCS) THEN
702               CALL OBCS_APPLY_UV( bi, bj, k, gUnm1, gVnm1, myThid )
703             END IF
704    #endif   /* ALLOW_OBCS */
705    
706    #ifdef   ALLOW_AUTODIFF_TAMC
707    #ifdef   INCLUDE_CD_CODE
708             ELSE
709               DO j=1-OLy,sNy+OLy
710                 DO i=1-OLx,sNx+OLx
711                   guCD(i,j,k,bi,bj) = 0.0
712                   gvCD(i,j,k,bi,bj) = 0.0
713                 END DO
714               END DO
715    #endif   /* INCLUDE_CD_CODE */
716    #endif   /* ALLOW_AUTODIFF_TAMC */
717             ENDIF
718    
719    
720    C--     end of dynamics k loop (1:Nr)
721            ENDDO
722    
723    
724    
725    C--     Implicit viscosity
726            IF (implicitViscosity.AND.momStepping) THEN
727    #ifdef    ALLOW_AUTODIFF_TAMC
728              idkey = iikey + 3
729  CADJ STORE gUNm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte  CADJ STORE gUNm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte
730  #endif /* ALLOW_AUTODIFF_TAMC */  #endif    /* ALLOW_AUTODIFF_TAMC */
731            CALL IMPLDIFF(            CALL IMPLDIFF(
732       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
733       I         deltaTmom, KappaRU,recip_HFacW,       I         deltaTmom, KappaRU,recip_HFacW,
734       U         gUNm1,       U         gUNm1,
735       I         myThid )       I         myThid )
736  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef    ALLOW_AUTODIFF_TAMC
737           idkey = iikey + 4            idkey = iikey + 4
738  CADJ STORE gVNm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte  CADJ STORE gVNm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte
739  #endif /* ALLOW_AUTODIFF_TAMC */  #endif    /* ALLOW_AUTODIFF_TAMC */
740            CALL IMPLDIFF(            CALL IMPLDIFF(
741       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
742       I         deltaTmom, KappaRV,recip_HFacS,       I         deltaTmom, KappaRV,recip_HFacS,
743       U         gVNm1,       U         gVNm1,
744       I         myThid )       I         myThid )
745    
746  #ifdef INCLUDE_CD_CODE  #ifdef   ALLOW_OBCS
747    C--      Apply open boundary conditions
748             IF (useOBCS) THEN
749               DO K=1,Nr
750                 CALL OBCS_APPLY_UV( bi, bj, k, gUnm1, gVnm1, myThid )
751               ENDDO
752             END IF
753    #endif   /* ALLOW_OBCS */
754    
755  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef    INCLUDE_CD_CODE
756           idkey = iikey + 5  #ifdef    ALLOW_AUTODIFF_TAMC
757              idkey = iikey + 5
758  CADJ STORE vVelD(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte  CADJ STORE vVelD(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte
759  #endif /* ALLOW_AUTODIFF_TAMC */  #endif    /* ALLOW_AUTODIFF_TAMC */
760            CALL IMPLDIFF(            CALL IMPLDIFF(
761       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
762       I         deltaTmom, KappaRU,recip_HFacW,       I         deltaTmom, KappaRU,recip_HFacW,
763       U         vVelD,       U         vVelD,
764       I         myThid )       I         myThid )
765  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef    ALLOW_AUTODIFF_TAMC
766          idkey = iikey + 6            idkey = iikey + 6
767  CADJ STORE uVelD(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte  CADJ STORE uVelD(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte
768  #endif /* ALLOW_AUTODIFF_TAMC */  #endif    /* ALLOW_AUTODIFF_TAMC */
769            CALL IMPLDIFF(            CALL IMPLDIFF(
770       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
771       I         deltaTmom, KappaRV,recip_HFacS,       I         deltaTmom, KappaRV,recip_HFacS,
772       U         uVelD,       U         uVelD,
773       I         myThid )       I         myThid )
774    #endif    /* INCLUDE_CD_CODE */
775  #endif  C--     End If implicitViscosity.AND.momStepping
   
 C--      momStepping  
          ENDIF  
   
 C--     implicitViscosity  
776          ENDIF          ENDIF
777    
778    Cjmc : add for phiHyd output <- but not working if multi tile per CPU
779    c       IF ( DIFFERENT_MULTIPLE(dumpFreq,myTime+deltaTClock,myTime)
780    c    &  .AND. buoyancyRelation .eq. 'ATMOSPHERIC' ) THEN
781    c         WRITE(suff,'(I10.10)') myIter+1
782    c         CALL WRITE_FLD_XYZ_RL('PH.',suff,phiHyd,myIter+1,myThid)
783    c       ENDIF
784    Cjmc(end)
785    
786    #ifdef ALLOW_TIMEAVE
787            IF (taveFreq.GT.0.) THEN
788              CALL TIMEAVE_CUMUL_1T(phiHydtave, phiHyd, Nr,
789         I                              deltaTclock, bi, bj, myThid)
790              IF (ivdc_kappa.NE.0.) THEN
791                CALL TIMEAVE_CUMULATE(ConvectCountTave, ConvectCount, Nr,
792         I                              deltaTclock, bi, bj, myThid)
793              ENDIF
794            ENDIF
795    #endif /* ALLOW_TIMEAVE */
796    
797         ENDDO         ENDDO
798        ENDDO        ENDDO
799    
800    #ifndef EXCLUDE_DEBUGMODE
801          If (debugMode) THEN
802           CALL DEBUG_STATS_RL(1,EtaN,'EtaN (DYNAMICS)',myThid)
803           CALL DEBUG_STATS_RL(Nr,uVel,'Uvel (DYNAMICS)',myThid)
804           CALL DEBUG_STATS_RL(Nr,vVel,'Vvel (DYNAMICS)',myThid)
805           CALL DEBUG_STATS_RL(Nr,wVel,'Wvel (DYNAMICS)',myThid)
806           CALL DEBUG_STATS_RL(Nr,theta,'Theta (DYNAMICS)',myThid)
807           CALL DEBUG_STATS_RL(Nr,salt,'Salt (DYNAMICS)',myThid)
808           CALL DEBUG_STATS_RL(Nr,Gu,'Gu (DYNAMICS)',myThid)
809           CALL DEBUG_STATS_RL(Nr,Gv,'Gv (DYNAMICS)',myThid)
810           CALL DEBUG_STATS_RL(Nr,Gt,'Gt (DYNAMICS)',myThid)
811           CALL DEBUG_STATS_RL(Nr,Gs,'Gs (DYNAMICS)',myThid)
812           CALL DEBUG_STATS_RL(Nr,GuNm1,'GuNm1 (DYNAMICS)',myThid)
813           CALL DEBUG_STATS_RL(Nr,GvNm1,'GvNm1 (DYNAMICS)',myThid)
814           CALL DEBUG_STATS_RL(Nr,GtNm1,'GtNm1 (DYNAMICS)',myThid)
815           CALL DEBUG_STATS_RL(Nr,GsNm1,'GsNm1 (DYNAMICS)',myThid)
816          ENDIF
817    #endif
818    
819        RETURN        RETURN
820        END        END
Legend:
Removed from v.1.57  
changed lines
  Added in v.1.74
  ViewVC Help
Powered by ViewVC 1.1.22