/[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.58 by adcroft, Fri Feb 2 21:04:48 2001 UTC revision 1.67 by heimbach, Mon May 14 21:46:17 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"
# Line 33  C     == Global variables === Line 33  C     == Global variables ===
33  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
34  # include "tamc.h"  # include "tamc.h"
35  # include "tamc_keys.h"  # include "tamc_keys.h"
36    # include "FFIELDS.h"
37    # ifdef ALLOW_KPP
38    #  include "KPP.h"
39    # endif
40    # ifdef ALLOW_GMREDI
41    #  include "GMREDI.h"
42    # endif
43  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
44    
45  #ifdef ALLOW_KPP  #ifdef ALLOW_TIMEAVE
46  # include "KPP.h"  #include "TIMEAVE_STATV.h"
47  #endif  #endif
48    
49  C     == Routine arguments ==  C     == Routine arguments ==
# Line 51  C     == Local variables Line 58  C     == Local variables
58  C     xA, yA                 - Per block temporaries holding face areas  C     xA, yA                 - Per block temporaries holding face areas
59  C     uTrans, vTrans, rTrans - Per block temporaries holding flow  C     uTrans, vTrans, rTrans - Per block temporaries holding flow
60  C                              transport  C                              transport
61  C     rVel                     o uTrans: Zonal transport  C                              o uTrans: Zonal transport
62  C                              o vTrans: Meridional transport  C                              o vTrans: Meridional transport
63  C                              o rTrans: Vertical transport  C                              o rTrans: Vertical transport
 C                              o rVel:   Vertical velocity at upper and  
 C                                        lower cell faces.  
64  C     maskC,maskUp             o maskC: land/water mask for tracer cells  C     maskC,maskUp             o maskC: land/water mask for tracer cells
65  C                              o maskUp: land/water mask for W points  C                              o maskUp: land/water mask for W points
66  C     fVer[STUV]               o fVer: Vertical flux term - note fVer  C     fVer[STUV]               o fVer: Vertical flux term - note fVer
# Line 65  C                                      v Line 70  C                                      v
70  C     rhoK, rhoKM1   - Density at current level, and level above  C     rhoK, rhoKM1   - Density 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 84  C                      index into fVerTe Line 88  C                      index into fVerTe
88        _RL uTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL uTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
89        _RL vTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL vTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
90        _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)  
91        _RS maskC   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _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)
93        _RL fVerT   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerT   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
# Line 94  C                      index into fVerTe Line 97  C                      index into fVerTe
97        _RL phiHyd  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL phiHyd  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
98        _RL rhokm1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rhokm1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
99        _RL rhok    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rhok    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
100          _RL phiSurfX(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
101          _RL phiSurfY(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
102        _RL KappaRT (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL KappaRT (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
103        _RL KappaRS (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL KappaRS (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
104        _RL KappaRU (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL KappaRU (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
# Line 102  C                      index into fVerTe Line 107  C                      index into fVerTe
107        _RL sigmaY  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaY  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
108        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
109    
110  C This is currently also used by IVDC and Diagnostics  C This is currently used by IVDC and Diagnostics
 C #ifdef INCLUDE_CONVECT_CALL  
111        _RL ConvectCount (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL ConvectCount (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
 C #endif  
112    
113        INTEGER iMin, iMax        INTEGER iMin, iMax
114        INTEGER jMin, jMax        INTEGER jMin, jMax
# Line 113  C #endif Line 116  C #endif
116        INTEGER i, j        INTEGER i, j
117        INTEGER k, km1, kup, kDown        INTEGER k, km1, kup, kDown
118    
119  #ifdef ALLOW_AUTODIFF_TAMC  Cjmc : add for phiHyd output <- but not working if multi tile per CPU
120        INTEGER    isbyte  c     CHARACTER*(MAX_LEN_MBUF) suff
121        PARAMETER( isbyte = 4 )  c     LOGICAL  DIFFERENT_MULTIPLE
122    c     EXTERNAL DIFFERENT_MULTIPLE
123        INTEGER act1, act2, act3, act4  Cjmc(end)
124        INTEGER max1, max2, max3  
       INTEGER iikey, kkey  
       INTEGER maximpl  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
125  C---    The algorithm...  C---    The algorithm...
126  C  C
127  C       "Correction Step"  C       "Correction Step"
# Line 137  C       "Calculation of Gs" Line 136  C       "Calculation of Gs"
136  C       ===================  C       ===================
137  C       This is where all the accelerations and tendencies (ie.  C       This is where all the accelerations and tendencies (ie.
138  C       physics, parameterizations etc...) are calculated  C       physics, parameterizations etc...) are calculated
 C         rVel = sum_r ( div. u[n] )  
139  C         rho = rho ( theta[n], salt[n] )  C         rho = rho ( theta[n], salt[n] )
140  C         b   = b(rho, theta)  C         b   = b(rho, theta)
141  C         K31 = K31 ( rho )  C         K31 = K31 ( rho )
142  C         Gu[n] = Gu( u[n], v[n], rVel, b, ... )  C         Gu[n] = Gu( u[n], v[n], wVel, b, ... )
143  C         Gv[n] = Gv( u[n], v[n], rVel, b, ... )  C         Gv[n] = Gv( u[n], v[n], wVel, b, ... )
144  C         Gt[n] = Gt( theta[n], u[n], v[n], rVel, K31, ... )  C         Gt[n] = Gt( theta[n], u[n], v[n], wVel, K31, ... )
145  C         Gs[n] = Gs( salt[n], u[n], v[n], rVel, K31, ... )  C         Gs[n] = Gs( salt[n], u[n], v[n], wVel, K31, ... )
146  C  C
147  C       "Time-stepping" or "Prediction"  C       "Time-stepping" or "Prediction"
148  C       ================================  C       ================================
# Line 195  C     uninitialised but inert locations. Line 193  C     uninitialised but inert locations.
193          rhoKM1 (i,j) = 0. _d 0          rhoKM1 (i,j) = 0. _d 0
194          rhok   (i,j) = 0. _d 0          rhok   (i,j) = 0. _d 0
195          maskC  (i,j) = 0. _d 0          maskC  (i,j) = 0. _d 0
196            phiSurfX(i,j) = 0. _d 0
197            phiSurfY(i,j) = 0. _d 0
198         ENDDO         ENDDO
199        ENDDO        ENDDO
200    
# Line 208  CHPF$ INDEPENDENT Line 208  CHPF$ INDEPENDENT
208    
209  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
210  C--    HPF directive to help TAMC  C--    HPF directive to help TAMC
211  CHPF$  INDEPENDENT, NEW (rTrans,rVel,fVerT,fVerS,fVerU,fVerV  CHPF$  INDEPENDENT, NEW (rTrans,fVerT,fVerS,fVerU,fVerV
212  CHPF$&                  ,phiHyd,utrans,vtrans,maskc,xA,yA  CHPF$&                  ,phiHyd,utrans,vtrans,maskc,xA,yA
213  CHPF$&                  ,KappaRT,KappaRS,KappaRU,KappaRV  CHPF$&                  ,KappaRT,KappaRS,KappaRU,KappaRV
214  CHPF$&                  )  CHPF$&                  )
# Line 237  C--     Set up work arrays that need val Line 237  C--     Set up work arrays that need val
237          DO j=1-OLy,sNy+OLy          DO j=1-OLy,sNy+OLy
238           DO i=1-OLx,sNx+OLx           DO i=1-OLx,sNx+OLx
239            rTrans(i,j)   = 0. _d 0            rTrans(i,j)   = 0. _d 0
           rVel  (i,j,1) = 0. _d 0  
           rVel  (i,j,2) = 0. _d 0  
240            fVerT (i,j,1) = 0. _d 0            fVerT (i,j,1) = 0. _d 0
241            fVerT (i,j,2) = 0. _d 0            fVerT (i,j,2) = 0. _d 0
242            fVerS (i,j,1) = 0. _d 0            fVerS (i,j,1) = 0. _d 0
# Line 253  C--     Set up work arrays that need val Line 251  C--     Set up work arrays that need val
251          DO k=1,Nr          DO k=1,Nr
252           DO j=1-OLy,sNy+OLy           DO j=1-OLy,sNy+OLy
253            DO i=1-OLx,sNx+OLx            DO i=1-OLx,sNx+OLx
254  #ifdef INCLUDE_CONVECT_CALL  C This is currently also used by IVDC and Diagnostics
255             ConvectCount(i,j,k) = 0.             ConvectCount(i,j,k) = 0.
 #endif  
256             KappaRT(i,j,k) = 0. _d 0             KappaRT(i,j,k) = 0. _d 0
257             KappaRS(i,j,k) = 0. _d 0             KappaRS(i,j,k) = 0. _d 0
258            ENDDO            ENDDO
# Line 268  C--     Set up work arrays that need val Line 265  C--     Set up work arrays that need val
265          jMax = sNy+OLy          jMax = sNy+OLy
266    
267    
268    #ifdef ALLOW_AUTODIFF_TAMC
269    CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
270    CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
271    CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
272    CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
273    #endif /* ALLOW_AUTODIFF_TAMC */
274    
275  C--     Start of diagnostic loop  C--     Start of diagnostic loop
276          DO k=Nr,1,-1          DO k=Nr,1,-1
277    
278  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
279  C? Patrick, is this formula correct now that we change the loop range?  C? Patrick, is this formula correct now that we change the loop range?
280  C? Do we still need this?  C? Do we still need this?
281           kkey = (ikey-1)*(Nr-2+1) + (k-2) + 1  cph kkey formula corrected.
282    cph Needed for rhok, rhokm1, in the case useGMREDI.
283             kkey = (ikey-1)*Nr + k
284    CADJ STORE rhokm1(:,:) = comlev1_bibj_k ,       key=kkey, byte=isbyte
285    CADJ STORE rhok  (:,:) = comlev1_bibj_k ,       key=kkey, byte=isbyte
286  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
287    
288  C--       Integrate continuity vertically for vertical velocity  C--       Integrate continuity vertically for vertical velocity
# Line 285  C--       Integrate continuity verticall Line 293  C--       Integrate continuity verticall
293    
294  #ifdef    ALLOW_OBCS  #ifdef    ALLOW_OBCS
295  #ifdef    ALLOW_NONHYDROSTATIC  #ifdef    ALLOW_NONHYDROSTATIC
296  C--       Calculate future values on open boundaries  C--       Apply OBC to W if in N-H mode
297            IF (useOBCS.AND.nonHydrostatic) THEN            IF (useOBCS.AND.nonHydrostatic) THEN
298              CALL OBCS_APPLY_W( bi, bj, k, wVel, myThid )              CALL OBCS_APPLY_W( bi, bj, k, wVel, myThid )
299            ENDIF            ENDIF
# Line 296  C--       Calculate gradients of potenti Line 304  C--       Calculate gradients of potenti
304  C         slope terms (e.g. GM/Redi tensor or IVDC diffusivity)  C         slope terms (e.g. GM/Redi tensor or IVDC diffusivity)
305  c         IF ( k.GT.1 .AND. (useGMRedi.OR.ivdc_kappa.NE.0.) ) THEN  c         IF ( k.GT.1 .AND. (useGMRedi.OR.ivdc_kappa.NE.0.) ) THEN
306            IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.) ) THEN            IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.) ) THEN
307    #ifdef ALLOW_AUTODIFF_TAMC
308    CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
309    CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
310    #endif /* ALLOW_AUTODIFF_TAMC */
311              CALL FIND_RHO(              CALL FIND_RHO(
312       I        bi, bj, iMin, iMax, jMin, jMax, k, k, eosType,       I        bi, bj, iMin, iMax, jMin, jMax, k, k, eosType,
313       I        theta, salt,       I        theta, salt,
314       O        rhoK,       O        rhoK,
315       I        myThid )       I        myThid )
316              IF (k.GT.1) CALL FIND_RHO(              IF (k.GT.1) THEN
317    #ifdef ALLOW_AUTODIFF_TAMC
318    CADJ STORE theta(:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
319    CADJ STORE salt (:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
320    #endif /* ALLOW_AUTODIFF_TAMC */
321                 CALL FIND_RHO(
322       I        bi, bj, iMin, iMax, jMin, jMax, k-1, k, eosType,       I        bi, bj, iMin, iMax, jMin, jMax, k-1, k, eosType,
323       I        theta, salt,       I        theta, salt,
324       O        rhoKm1,       O        rhoKm1,
325       I        myThid )       I        myThid )
326                ENDIF
327              CALL GRAD_SIGMA(              CALL GRAD_SIGMA(
328       I             bi, bj, iMin, iMax, jMin, jMax, k,       I             bi, bj, iMin, iMax, jMin, jMax, k,
329       I             rhoK, rhoKm1, rhoK,       I             rhoK, rhoKm1, rhoK,
# Line 321  c ==> should use sigmaR !!! Line 339  c ==> should use sigmaR !!!
339       I        rhoKm1, rhoK,       I        rhoKm1, rhoK,
340       U        ConvectCount, KappaRT, KappaRS,       U        ConvectCount, KappaRT, KappaRS,
341       I        myTime, myIter, myThid)       I        myTime, myIter, myThid)
342            END IF            ENDIF
343    
344  C--     end of diagnostic k loop (Nr:1)  C--     end of diagnostic k loop (Nr:1)
345          ENDDO          ENDDO
346    
347    #ifdef ALLOW_AUTODIFF_TAMC
348    cph avoids recomputation of integrate_for_w
349    CADJ STORE wvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte
350    #endif /* ALLOW_AUTODIFF_TAMC */
351    
352  #ifdef  ALLOW_OBCS  #ifdef  ALLOW_OBCS
353  C--     Calculate future values on open boundaries  C--     Calculate future values on open boundaries
354          IF (useOBCS) THEN          IF (useOBCS) THEN
# Line 340  C       relaxation terms, etc. Line 363  C       relaxation terms, etc.
363          CALL EXTERNAL_FORCING_SURF(          CALL EXTERNAL_FORCING_SURF(
364       I             bi, bj, iMin, iMax, jMin, jMax,       I             bi, bj, iMin, iMax, jMin, jMax,
365       I             myThid )       I             myThid )
366    #ifdef ALLOW_AUTODIFF_TAMC
367    cph needed for KPP
368    CADJ STORE surfacetendencyU(:,:,bi,bj)
369    CADJ &     = comlev1_bibj, key=ikey, byte=isbyte
370    CADJ STORE surfacetendencyV(:,:,bi,bj)
371    CADJ &     = comlev1_bibj, key=ikey, byte=isbyte
372    CADJ STORE surfacetendencyS(:,:,bi,bj)
373    CADJ &     = comlev1_bibj, key=ikey, byte=isbyte
374    CADJ STORE surfacetendencyT(:,:,bi,bj)
375    CADJ &     = comlev1_bibj, key=ikey, byte=isbyte
376    #endif /* ALLOW_AUTODIFF_TAMC */
377    
378  #ifdef  ALLOW_GMREDI  #ifdef  ALLOW_GMREDI
379    
380    #ifdef ALLOW_AUTODIFF_TAMC
381    CADJ STORE sigmaX(:,:,:) = comlev1, key=ikey, byte=isbyte
382    CADJ STORE sigmaY(:,:,:) = comlev1, key=ikey, byte=isbyte
383    CADJ STORE sigmaR(:,:,:) = comlev1, key=ikey, byte=isbyte
384    #endif /* ALLOW_AUTODIFF_TAMC */
385  C--     Calculate iso-neutral slopes for the GM/Redi parameterisation  C--     Calculate iso-neutral slopes for the GM/Redi parameterisation
386          IF (useGMRedi) THEN          IF (useGMRedi) THEN
387            DO k=1,Nr            DO k=1,Nr
# Line 360  C--     Calculate iso-neutral slopes for Line 400  C--     Calculate iso-neutral slopes for
400            ENDDO            ENDDO
401  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
402          ENDIF          ENDIF
403    
404    #ifdef ALLOW_AUTODIFF_TAMC
405    CADJ STORE Kwx(:,:,:,bi,bj)   = comlev1_bibj, key=ikey, byte=isbyte
406    CADJ STORE Kwy(:,:,:,bi,bj)   = comlev1_bibj, key=ikey, byte=isbyte
407    CADJ STORE Kwz(:,:,:,bi,bj)   = comlev1_bibj, key=ikey, byte=isbyte
408    #endif /* ALLOW_AUTODIFF_TAMC */
409    
410  #endif  /* ALLOW_GMREDI */  #endif  /* ALLOW_GMREDI */
411    
412  #ifdef  ALLOW_KPP  #ifdef  ALLOW_KPP
# Line 367  C--     Compute KPP mixing coefficients Line 414  C--     Compute KPP mixing coefficients
414          IF (useKPP) THEN          IF (useKPP) THEN
415            CALL KPP_CALC(            CALL KPP_CALC(
416       I                  bi, bj, myTime, myThid )       I                  bi, bj, myTime, myThid )
417    #ifdef ALLOW_AUTODIFF_TAMC
418            ELSE
419              CALL KPP_CALC_DUMMY(
420         I                  bi, bj, myTime, myThid )
421    #endif /* ALLOW_AUTODIFF_TAMC */
422          ENDIF          ENDIF
423    
424    #ifdef ALLOW_AUTODIFF_TAMC
425    CADJ STORE KPPghat   (:,:,:,bi,bj)
426    CADJ &   , KPPviscAz (:,:,:,bi,bj)
427    CADJ &   , KPPdiffKzT(:,:,:,bi,bj)
428    CADJ &   , KPPdiffKzS(:,:,:,bi,bj)
429    CADJ &   , KPPfrac   (:,:  ,bi,bj)
430    CADJ &                 = comlev1_bibj, key=ikey, byte=isbyte
431    #endif /* ALLOW_AUTODIFF_TAMC */
432    
433  #endif  /* ALLOW_KPP */  #endif  /* ALLOW_KPP */
434    
435  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 392  C note(jmc) : phiHyd=0 at this point but Line 454  C note(jmc) : phiHyd=0 at this point but
454    
455  C--     Start of thermodynamics loop  C--     Start of thermodynamics loop
456          DO k=Nr,1,-1          DO k=Nr,1,-1
457    #ifdef ALLOW_AUTODIFF_TAMC
458    C? Patrick Is this formula correct?
459    cph Yes, but I rewrote it.
460    cph Also, the KappaR? need the index and subscript k!
461             kkey = (ikey-1)*Nr + k
462    #endif /* ALLOW_AUTODIFF_TAMC */
463    
464  C--       km1    Points to level above k (=k-1)  C--       km1    Points to level above k (=k-1)
465  C--       kup    Cycles through 1,2 to point to layer above  C--       kup    Cycles through 1,2 to point to layer above
# Line 406  C--       kDown  Cycles through 2,1 to p Line 474  C--       kDown  Cycles through 2,1 to p
474            jMin = 1-OLy+2            jMin = 1-OLy+2
475            jMax = sNy+OLy-1            jMax = sNy+OLy-1
476    
 #ifdef ALLOW_AUTODIFF_TAMC  
 CPatrick Is this formula correct?  
          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(:,:,:)    = comlev1_bibj_k, key = kkey, byte = isbyte  
 CADJ STORE KappaRS(:,:,:)    = comlev1_bibj_k, key = kkey, byte = isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
477  C--      Get temporary terms used by tendency routines  C--      Get temporary terms used by tendency routines
478           CALL CALC_COMMON_FACTORS (           CALL CALC_COMMON_FACTORS (
479       I        bi,bj,iMin,iMax,jMin,jMax,k,km1,kup,kDown,       I        bi,bj,iMin,iMax,jMin,jMax,k,km1,kup,kDown,
480       O        xA,yA,uTrans,vTrans,rTrans,rVel,maskC,maskUp,       O        xA,yA,uTrans,vTrans,rTrans,maskC,maskUp,
481       I        myThid)       I        myThid)
482    
483    #ifdef ALLOW_AUTODIFF_TAMC
484    CADJ STORE KappaRT(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte
485    CADJ STORE KappaRS(:,:,k)    = comlev1_bibj_k, key=kkey, byte=isbyte
486    #endif /* ALLOW_AUTODIFF_TAMC */
487    
488  #ifdef  INCLUDE_CALC_DIFFUSIVITY_CALL  #ifdef  INCLUDE_CALC_DIFFUSIVITY_CALL
489  C--      Calculate the total vertical diffusivity  C--      Calculate the total vertical diffusivity
490           CALL CALC_DIFFUSIVITY(           CALL CALC_DIFFUSIVITY(
# Line 480  C--     end of thermodynamic k loop (Nr: Line 544  C--     end of thermodynamic k loop (Nr:
544    
545    
546  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
547  CPatrick? What about this one?  C? Patrick? What about this one?
548             maximpl = 6  cph Keys iikey and idkey don't seem to be needed
549             iikey = (ikey-1)*maximpl  cph since storing occurs on different tape for each
550    cph impldiff call anyways.
551    cph Thus, common block comlev1_impl isn't needed either.
552    cph Storing below needed in the case useGMREDI.
553            iikey = (ikey-1)*maximpl
554  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
555    
556  C--     Implicit diffusion  C--     Implicit diffusion
557          IF (implicitDiffusion) THEN          IF (implicitDiffusion) THEN
558    
559            IF (tempStepping) THEN           IF (tempStepping) THEN
560  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
561              idkey = iikey + 1              idkey = iikey + 1
562    CADJ STORE gTNm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte
563  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
564              CALL IMPLDIFF(              CALL IMPLDIFF(
565       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
# Line 502  C--     Implicit diffusion Line 571  C--     Implicit diffusion
571           IF (saltStepping) THEN           IF (saltStepping) THEN
572  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
573           idkey = iikey + 2           idkey = iikey + 2
574    CADJ STORE gSNm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte
575  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
576              CALL IMPLDIFF(              CALL IMPLDIFF(
577       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
# Line 522  C--      Apply open boundary conditions Line 592  C--      Apply open boundary conditions
592  C--     End If implicitDiffusion  C--     End If implicitDiffusion
593          ENDIF          ENDIF
594    
595    C--     Start computation of dynamics
596            iMin = 1-OLx+2
597            iMax = sNx+OLx-1
598            jMin = 1-OLy+2
599            jMax = sNy+OLy-1
600    
601    C--     Explicit part of the Surface Potentiel Gradient (add in TIMESTEP)
602    C       (note: this loop will be replaced by CALL CALC_GRAD_ETA)
603            IF (implicSurfPress.NE.1.) THEN
604              CALL CALC_GRAD_PHI_SURF(
605         I         bi,bj,iMin,iMax,jMin,jMax,
606         I         etaN,
607         O         phiSurfX,phiSurfY,
608         I         myThid )                        
609            ENDIF
610    
611  C--     Start of dynamics loop  C--     Start of dynamics loop
612          DO k=1,Nr          DO k=1,Nr
# Line 535  C--       kDown  Cycles through 2,1 to p Line 619  C--       kDown  Cycles through 2,1 to p
619            kup  = 1+MOD(k+1,2)            kup  = 1+MOD(k+1,2)
620            kDown= 1+MOD(k,2)            kDown= 1+MOD(k,2)
621    
           iMin = 1-OLx+2  
           iMax = sNx+OLx-1  
           jMin = 1-OLy+2  
           jMax = sNy+OLy-1  
   
622  C--      Integrate hydrostatic balance for phiHyd with BC of  C--      Integrate hydrostatic balance for phiHyd with BC of
623  C        phiHyd(z=0)=0  C        phiHyd(z=0)=0
624  C        distinguishe between Stagger and Non Stagger time stepping  C        distinguishe between Stagger and Non Stagger time stepping
# Line 566  C        and step forward storing the re Line 645  C        and step forward storing the re
645       U         fVerU, fVerV,       U         fVerU, fVerV,
646       I         myTime, myThid)       I         myTime, myThid)
647             CALL TIMESTEP(             CALL TIMESTEP(
648       I         bi,bj,iMin,iMax,jMin,jMax,k,phiHyd,       I         bi,bj,iMin,iMax,jMin,jMax,k,
649         I         phiHyd, phiSurfX, phiSurfY,
650       I         myIter, myThid)       I         myIter, myThid)
651    
652  #ifdef   ALLOW_OBCS  #ifdef   ALLOW_OBCS
# Line 599  C--     Implicit viscosity Line 679  C--     Implicit viscosity
679          IF (implicitViscosity.AND.momStepping) THEN          IF (implicitViscosity.AND.momStepping) THEN
680  #ifdef    ALLOW_AUTODIFF_TAMC  #ifdef    ALLOW_AUTODIFF_TAMC
681            idkey = iikey + 3            idkey = iikey + 3
682    CADJ STORE gUNm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte
683  #endif    /* ALLOW_AUTODIFF_TAMC */  #endif    /* ALLOW_AUTODIFF_TAMC */
684            CALL IMPLDIFF(            CALL IMPLDIFF(
685       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
# Line 607  C--     Implicit viscosity Line 688  C--     Implicit viscosity
688       I         myThid )       I         myThid )
689  #ifdef    ALLOW_AUTODIFF_TAMC  #ifdef    ALLOW_AUTODIFF_TAMC
690            idkey = iikey + 4            idkey = iikey + 4
691    CADJ STORE gVNm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte
692  #endif    /* ALLOW_AUTODIFF_TAMC */  #endif    /* ALLOW_AUTODIFF_TAMC */
693            CALL IMPLDIFF(            CALL IMPLDIFF(
694       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
# Line 626  C--      Apply open boundary conditions Line 708  C--      Apply open boundary conditions
708  #ifdef    INCLUDE_CD_CODE  #ifdef    INCLUDE_CD_CODE
709  #ifdef    ALLOW_AUTODIFF_TAMC  #ifdef    ALLOW_AUTODIFF_TAMC
710            idkey = iikey + 5            idkey = iikey + 5
711    CADJ STORE vVelD(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte
712  #endif    /* ALLOW_AUTODIFF_TAMC */  #endif    /* ALLOW_AUTODIFF_TAMC */
713            CALL IMPLDIFF(            CALL IMPLDIFF(
714       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
# Line 634  C--      Apply open boundary conditions Line 717  C--      Apply open boundary conditions
717       I         myThid )       I         myThid )
718  #ifdef    ALLOW_AUTODIFF_TAMC  #ifdef    ALLOW_AUTODIFF_TAMC
719            idkey = iikey + 6            idkey = iikey + 6
720    CADJ STORE uVelD(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte
721  #endif    /* ALLOW_AUTODIFF_TAMC */  #endif    /* ALLOW_AUTODIFF_TAMC */
722            CALL IMPLDIFF(            CALL IMPLDIFF(
723       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
# Line 644  C--      Apply open boundary conditions Line 728  C--      Apply open boundary conditions
728  C--     End If implicitViscosity.AND.momStepping  C--     End If implicitViscosity.AND.momStepping
729          ENDIF          ENDIF
730    
731    Cjmc : add for phiHyd output <- but not working if multi tile per CPU
732    c       IF ( DIFFERENT_MULTIPLE(dumpFreq,myTime+deltaTClock,myTime)
733    c    &  .AND. buoyancyRelation .eq. 'ATMOSPHERIC' ) THEN
734    c         WRITE(suff,'(I10.10)') myIter+1
735    c         CALL WRITE_FLD_XYZ_RL('PH.',suff,phiHyd,myIter+1,myThid)
736    c       ENDIF
737    Cjmc(end)
738    
739    #ifdef ALLOW_TIMEAVE
740            IF (taveFreq.GT.0.) THEN
741              CALL TIMEAVE_CUMULATE(phiHydtave, phiHyd, Nr,
742         I                              deltaTclock, bi, bj, myThid)
743              IF (ivdc_kappa.NE.0.) THEN
744                CALL TIMEAVE_CUMULATE(ConvectCountTave, ConvectCount, Nr,
745         I                              deltaTclock, bi, bj, myThid)
746              ENDIF
747            ENDIF
748    #endif /* ALLOW_TIMEAVE */
749    
750         ENDDO         ENDDO
751        ENDDO        ENDDO
752    
Legend:
Removed from v.1.58  
changed lines
  Added in v.1.67
  ViewVC Help
Powered by ViewVC 1.1.22