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