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

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

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

revision 1.4 by adcroft, Mon Sep 10 01:22:48 2001 UTC revision 1.17 by adcroft, Mon Mar 4 17:26:41 2002 UTC
# Line 3  C $Name$ Line 3  C $Name$
3    
4  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
5    
6    CBOP
7    C     !ROUTINE: THERMODYNAMICS
8    C     !INTERFACE:
9        SUBROUTINE THERMODYNAMICS(myTime, myIter, myThid)        SUBROUTINE THERMODYNAMICS(myTime, myIter, myThid)
10  C     /==========================================================\  C     !DESCRIPTION: \bv
11  C     | SUBROUTINE THERMODYNAMICS                                |  C     *==========================================================*
12  C     | o Controlling routine for the prognostic part of the     |  C     | SUBROUTINE THERMODYNAMICS                                
13  C     |   thermo-dynamics.                                       |  C     | o Controlling routine for the prognostic part of the      
14  C     |==========================================================|  C     |   thermo-dynamics.                                        
15  C     \==========================================================/  C     *===========================================================
16        IMPLICIT NONE  C     | The algorithm...
17    C     |
18    C     | "Correction Step"
19    C     | =================
20    C     | Here we update the horizontal velocities with the surface
21    C     | pressure such that the resulting flow is either consistent
22    C     | with the free-surface evolution or the rigid-lid:
23    C     |   U[n] = U* + dt x d/dx P
24    C     |   V[n] = V* + dt x d/dy P
25    C     |
26    C     | "Calculation of Gs"
27    C     | ===================
28    C     | This is where all the accelerations and tendencies (ie.
29    C     | physics, parameterizations etc...) are calculated
30    C     |   rho = rho ( theta[n], salt[n] )
31    C     |   b   = b(rho, theta)
32    C     |   K31 = K31 ( rho )
33    C     |   Gu[n] = Gu( u[n], v[n], wVel, b, ... )
34    C     |   Gv[n] = Gv( u[n], v[n], wVel, b, ... )
35    C     |   Gt[n] = Gt( theta[n], u[n], v[n], wVel, K31, ... )
36    C     |   Gs[n] = Gs( salt[n], u[n], v[n], wVel, K31, ... )
37    C     |
38    C     | "Time-stepping" or "Prediction"
39    C     | ================================
40    C     | The models variables are stepped forward with the appropriate
41    C     | time-stepping scheme (currently we use Adams-Bashforth II)
42    C     | - For momentum, the result is always *only* a "prediction"
43    C     | in that the flow may be divergent and will be "corrected"
44    C     | later with a surface pressure gradient.
45    C     | - Normally for tracers the result is the new field at time
46    C     | level [n+1} *BUT* in the case of implicit diffusion the result
47    C     | is also *only* a prediction.
48    C     | - We denote "predictors" with an asterisk (*).
49    C     |   U* = U[n] + dt x ( 3/2 Gu[n] - 1/2 Gu[n-1] )
50    C     |   V* = V[n] + dt x ( 3/2 Gv[n] - 1/2 Gv[n-1] )
51    C     |   theta[n+1] = theta[n] + dt x ( 3/2 Gt[n] - 1/2 atG[n-1] )
52    C     |   salt[n+1] = salt[n] + dt x ( 3/2 Gt[n] - 1/2 atG[n-1] )
53    C     | With implicit diffusion:
54    C     |   theta* = theta[n] + dt x ( 3/2 Gt[n] - 1/2 atG[n-1] )
55    C     |   salt* = salt[n] + dt x ( 3/2 Gt[n] - 1/2 atG[n-1] )
56    C     |   (1 + dt * K * d_zz) theta[n] = theta*
57    C     |   (1 + dt * K * d_zz) salt[n] = salt*
58    C     |
59    C     *==========================================================*
60    C     \ev
61    
62    C     !USES:
63          IMPLICIT NONE
64  C     == Global variables ===  C     == Global variables ===
65  #include "SIZE.h"  #include "SIZE.h"
66  #include "EEPARAMS.h"  #include "EEPARAMS.h"
# Line 22  C     == Global variables === Line 71  C     == Global variables ===
71  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_PASSIVE_TRACER
72  #include "TR1.h"  #include "TR1.h"
73  #endif  #endif
   
74  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
75  # include "tamc.h"  # include "tamc.h"
76  # include "tamc_keys.h"  # include "tamc_keys.h"
# Line 34  C     == Global variables === Line 82  C     == Global variables ===
82  #  include "GMREDI.h"  #  include "GMREDI.h"
83  # endif  # endif
84  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
   
85  #ifdef ALLOW_TIMEAVE  #ifdef ALLOW_TIMEAVE
86  #include "TIMEAVE_STATV.h"  #include "TIMEAVE_STATV.h"
87  #endif  #endif
88    
89    C     !INPUT/OUTPUT PARAMETERS:
90  C     == Routine arguments ==  C     == Routine arguments ==
91  C     myTime - Current time in simulation  C     myTime - Current time in simulation
92  C     myIter - Current iteration number in simulation  C     myIter - Current iteration number in simulation
# Line 47  C     myThid - Thread number for this in Line 95  C     myThid - Thread number for this in
95        INTEGER myIter        INTEGER myIter
96        INTEGER myThid        INTEGER myThid
97    
98    C     !LOCAL VARIABLES:
99  C     == Local variables  C     == Local variables
100  C     xA, yA                 - Per block temporaries holding face areas  C     xA, yA                 - Per block temporaries holding face areas
101  C     uTrans, vTrans, rTrans - Per block temporaries holding flow  C     uTrans, vTrans, rTrans - Per block temporaries holding flow
# Line 94  C                      index into fVerTe Line 143  C                      index into fVerTe
143        _RL sigmaX  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaX  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
144        _RL sigmaY  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaY  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
145        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
146    C     This is currently used by IVDC and Diagnostics
 C This is currently used by IVDC and Diagnostics  
147        _RL ConvectCount (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL ConvectCount (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
   
148        INTEGER iMin, iMax        INTEGER iMin, iMax
149        INTEGER jMin, jMax        INTEGER jMin, jMax
150        INTEGER bi, bj        INTEGER bi, bj
151        INTEGER i, j        INTEGER i, j
152        INTEGER k, km1, kup, kDown        INTEGER k, km1, kup, kDown
153    
154  Cjmc : add for phiHyd output <- but not working if multi tile per CPU  CEOP
 c     CHARACTER*(MAX_LEN_MBUF) suff  
 c     LOGICAL  DIFFERENT_MULTIPLE  
 c     EXTERNAL DIFFERENT_MULTIPLE  
 Cjmc(end)  
155    
 C---    The algorithm...  
 C  
 C       "Correction Step"  
 C       =================  
 C       Here we update the horizontal velocities with the surface  
 C       pressure such that the resulting flow is either consistent  
 C       with the free-surface evolution or the rigid-lid:  
 C         U[n] = U* + dt x d/dx P  
 C         V[n] = V* + dt x d/dy P  
 C  
 C       "Calculation of Gs"  
 C       ===================  
 C       This is where all the accelerations and tendencies (ie.  
 C       physics, parameterizations etc...) are calculated  
 C         rho = rho ( theta[n], salt[n] )  
 C         b   = b(rho, theta)  
 C         K31 = K31 ( rho )  
 C         Gu[n] = Gu( u[n], v[n], wVel, b, ... )  
 C         Gv[n] = Gv( u[n], v[n], wVel, b, ... )  
 C         Gt[n] = Gt( theta[n], u[n], v[n], wVel, K31, ... )  
 C         Gs[n] = Gs( salt[n], u[n], v[n], wVel, K31, ... )  
 C  
 C       "Time-stepping" or "Prediction"  
 C       ================================  
 C       The models variables are stepped forward with the appropriate  
 C       time-stepping scheme (currently we use Adams-Bashforth II)  
 C       - For momentum, the result is always *only* a "prediction"  
 C       in that the flow may be divergent and will be "corrected"  
 C       later with a surface pressure gradient.  
 C       - Normally for tracers the result is the new field at time  
 C       level [n+1} *BUT* in the case of implicit diffusion the result  
 C       is also *only* a prediction.  
 C       - We denote "predictors" with an asterisk (*).  
 C         U* = U[n] + dt x ( 3/2 Gu[n] - 1/2 Gu[n-1] )  
 C         V* = V[n] + dt x ( 3/2 Gv[n] - 1/2 Gv[n-1] )  
 C         theta[n+1] = theta[n] + dt x ( 3/2 Gt[n] - 1/2 atG[n-1] )  
 C         salt[n+1] = salt[n] + dt x ( 3/2 Gt[n] - 1/2 atG[n-1] )  
 C       With implicit diffusion:  
 C         theta* = theta[n] + dt x ( 3/2 Gt[n] - 1/2 atG[n-1] )  
 C         salt* = salt[n] + dt x ( 3/2 Gt[n] - 1/2 atG[n-1] )  
 C         (1 + dt * K * d_zz) theta[n] = theta*  
 C         (1 + dt * K * d_zz) salt[n] = salt*  
 C---  
   
156  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
157  C--   dummy statement to end declaration part  C--   dummy statement to end declaration part
158        ikey = 1        ikey = 1
# Line 204  CHPF$&                  ) Line 203  CHPF$&                  )
203  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
204            act1 = bi - myBxLo(myThid)            act1 = bi - myBxLo(myThid)
205            max1 = myBxHi(myThid) - myBxLo(myThid) + 1            max1 = myBxHi(myThid) - myBxLo(myThid) + 1
   
206            act2 = bj - myByLo(myThid)            act2 = bj - myByLo(myThid)
207            max2 = myByHi(myThid) - myByLo(myThid) + 1            max2 = myByHi(myThid) - myByLo(myThid) + 1
   
208            act3 = myThid - 1            act3 = myThid - 1
209            max3 = nTx*nTy            max3 = nTx*nTy
   
210            act4 = ikey_dynamics - 1            act4 = ikey_dynamics - 1
   
211            ikey = (act1 + 1) + act2*max1            ikey = (act1 + 1) + act2*max1
212       &                      + act3*max1*max2       &                      + act3*max1*max2
213       &                      + act4*max1*max2*max3       &                      + act4*max1*max2*max3
# Line 238  C This is currently also used by IVDC an Line 233  C This is currently also used by IVDC an
233             ConvectCount(i,j,k) = 0.             ConvectCount(i,j,k) = 0.
234             KappaRT(i,j,k) = 0. _d 0             KappaRT(i,j,k) = 0. _d 0
235             KappaRS(i,j,k) = 0. _d 0             KappaRS(i,j,k) = 0. _d 0
236    #ifdef ALLOW_AUTODIFF_TAMC
237               gT(i,j,k,bi,bj) = 0. _d 0
238               gS(i,j,k,bi,bj) = 0. _d 0
239    #ifdef ALLOW_PASSIVE_TRACER
240               gTr1(i,j,k,bi,bj) = 0. _d 0
241    #endif
242    #endif
243            ENDDO            ENDDO
244           ENDDO           ENDDO
245          ENDDO          ENDDO
# Line 249  C This is currently also used by IVDC an Line 251  C This is currently also used by IVDC an
251    
252    
253  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
254  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
255  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
256    #ifdef ALLOW_KPP
257    CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
258    CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
259    #endif
260  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
261    
262  C--     Start of diagnostic loop  C--     Start of diagnostic loop
# Line 327  C--     end of diagnostic k loop (Nr:1) Line 333  C--     end of diagnostic k loop (Nr:1)
333    
334  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
335  cph avoids recomputation of integrate_for_w  cph avoids recomputation of integrate_for_w
336  CADJ STORE wvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  CADJ STORE wvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
337  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
338    
339  #ifdef  ALLOW_OBCS  #ifdef  ALLOW_OBCS
340  C--     Calculate future values on open boundaries  C--     Calculate future values on open boundaries
341          IF (useOBCS) THEN          IF (useOBCS) THEN
342            CALL OBCS_CALC( bi, bj, myTime+deltaT,            CALL OBCS_CALC( bi, bj, myTime+deltaT, myIter+1,
343       I            uVel, vVel, wVel, theta, salt,       I            uVel, vVel, wVel, theta, salt,
344       I            myThid )       I            myThid )
345          ENDIF          ENDIF
# Line 365  CADJ STORE sigmaR(:,:,:) = comlev1, key= Line 371  CADJ STORE sigmaR(:,:,:) = comlev1, key=
371  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
372  C--     Calculate iso-neutral slopes for the GM/Redi parameterisation  C--     Calculate iso-neutral slopes for the GM/Redi parameterisation
373          IF (useGMRedi) THEN          IF (useGMRedi) THEN
374            DO k=1,Nr            CALL GMREDI_CALC_TENSOR(
375              CALL GMREDI_CALC_TENSOR(       I             bi, bj, iMin, iMax, jMin, jMax,
      I             bi, bj, iMin, iMax, jMin, jMax, k,  
376       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
377       I             myThid )       I             myThid )
           ENDDO  
378  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
379          ELSE          ELSE
380            DO k=1, Nr            CALL GMREDI_CALC_TENSOR_DUMMY(
381              CALL GMREDI_CALC_TENSOR_DUMMY(       I             bi, bj, iMin, iMax, jMin, jMax,
      I             bi, bj, iMin, iMax, jMin, jMax, k,  
382       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
383       I             myThid )       I             myThid )
           ENDDO  
384  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
385          ENDIF          ENDIF
386    
# Line 404  C--     Compute KPP mixing coefficients Line 406  C--     Compute KPP mixing coefficients
406    
407  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
408  CADJ STORE KPPghat   (:,:,:,bi,bj)  CADJ STORE KPPghat   (:,:,:,bi,bj)
 CADJ &   , KPPviscAz (:,:,:,bi,bj)  
409  CADJ &   , KPPdiffKzT(:,:,:,bi,bj)  CADJ &   , KPPdiffKzT(:,:,:,bi,bj)
410  CADJ &   , KPPdiffKzS(:,:,:,bi,bj)  CADJ &   , KPPdiffKzS(:,:,:,bi,bj)
411  CADJ &   , KPPfrac   (:,:  ,bi,bj)  CADJ &   , KPPfrac   (:,:  ,bi,bj)
# Line 414  CADJ &                 = comlev1_bibj, k Line 415  CADJ &                 = comlev1_bibj, k
415  #endif  /* ALLOW_KPP */  #endif  /* ALLOW_KPP */
416    
417  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
418  CADJ STORE KappaRT(:,:,:)     = comlev1_bibj, key = ikey, byte = isbyte  CADJ STORE KappaRT(:,:,:)     = comlev1_bibj, key=ikey, byte=isbyte
419  CADJ STORE KappaRS(:,:,:)     = comlev1_bibj, key = ikey, byte = isbyte  CADJ STORE KappaRS(:,:,:)     = comlev1_bibj, key=ikey, byte=isbyte
420  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
421  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
422  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
423  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
424  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_PASSIVE_TRACER
425  CADJ STORE tr1  (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte  CADJ STORE tr1  (:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
426  #endif  #endif
427  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
428    
# Line 435  C note(jmc) : phiHyd=0 at this point but Line 436  C note(jmc) : phiHyd=0 at this point but
436          ENDIF          ENDIF
437  #endif /* ALLOW_AIM */  #endif /* ALLOW_AIM */
438    
439    #ifdef ALLOW_TIMEAVE
440            IF (taveFreq.GT.0. .AND. ivdc_kappa.NE.0.) THEN
441              CALL TIMEAVE_CUMULATE(ConvectCountTave, ConvectCount, Nr,
442         I                               deltaTclock, bi, bj, myThid)
443            ENDIF
444    #endif /* ALLOW_TIMEAVE */
445    
446    #ifndef DISABLE_MULTIDIM_ADVECTION
447  C--     Some advection schemes are better calculated using a multi-dimensional  C--     Some advection schemes are better calculated using a multi-dimensional
448  C       method in the absence of any other terms and, if used, is done here.  C       method in the absence of any other terms and, if used, is done here.
449    C
450    C The CPP flag DISABLE_MULTIDIM_ADVECTION is currently unset in GAD_OPTIONS.h
451    C The default is to use multi-dimensinal advection for non-linear advection
452    C schemes. However, for the sake of efficiency of the adjoint it is necessary
453    C to be able to exclude this scheme to avoid excessive storage and
454    C recomputation. It *is* differentiable, if you need it.
455    C Edit GAD_OPTIONS.h and #define DISABLE_MULTIDIM_ADVECTION to
456    C disable this section of code.
457          IF (multiDimAdvection) THEN          IF (multiDimAdvection) THEN
458           IF (tempStepping .AND.           IF (tempStepping .AND.
459       &       tempAdvScheme.NE.ENUM_CENTERED_2ND .AND.       &       tempAdvScheme.NE.ENUM_CENTERED_2ND .AND.
460       &       tempAdvScheme.NE.ENUM_UPWIND_3RD .AND.       &       tempAdvScheme.NE.ENUM_UPWIND_3RD .AND.
461       &       tempAdvScheme.NE.ENUM_CENTERED_4TH )       &       tempAdvScheme.NE.ENUM_CENTERED_4TH ) THEN
462       &   CALL GAD_ADVECTION(bi,bj,tempAdvScheme,GAD_TEMPERATURE,theta,            CALL GAD_ADVECTION(bi,bj,tempAdvScheme,GAD_TEMPERATURE,
463       U                      gT,       U                      theta,gT,
464       I                      myTime,myIter,myThid)       I                      myTime,myIter,myThid)
465             ENDIF
466           IF (saltStepping .AND.           IF (saltStepping .AND.
467       &       saltAdvScheme.NE.ENUM_CENTERED_2ND .AND.       &       saltAdvScheme.NE.ENUM_CENTERED_2ND .AND.
468       &       saltAdvScheme.NE.ENUM_UPWIND_3RD .AND.       &       saltAdvScheme.NE.ENUM_UPWIND_3RD .AND.
469       &       saltAdvScheme.NE.ENUM_CENTERED_4TH )       &       saltAdvScheme.NE.ENUM_CENTERED_4TH ) THEN
470       &   CALL GAD_ADVECTION(bi,bj,saltAdvScheme,GAD_SALINITY,salt,            CALL GAD_ADVECTION(bi,bj,saltAdvScheme,GAD_SALINITY,
471       U                      gS,       U                      salt,gS,
472       I                      myTime,myIter,myThid)       I                      myTime,myIter,myThid)
473             ENDIF
474          ENDIF          ENDIF
475    C Since passive tracers are configurable separately from T,S we
476    C call the multi-dimensional method for PTRACERS regardless
477    C of whether multiDimAdvection is set or not.
478    #ifdef ALLOW_PTRACERS
479            IF ( usePTRACERS ) THEN
480             CALL PTRACERS_ADVECTION( bi,bj,myIter,myTime,myThid )
481            ENDIF
482    #endif /* ALLOW_PTRACERS */
483    #endif /* DISABLE_MULTIDIM_ADVECTION */
484    
485  C--     Start of thermodynamics loop  C--     Start of thermodynamics loop
486          DO k=Nr,1,-1          DO k=Nr,1,-1
# Line 510  C        and step forward storing result Line 537  C        and step forward storing result
537       I         xA,yA,uTrans,vTrans,rTrans,maskUp,       I         xA,yA,uTrans,vTrans,rTrans,maskUp,
538       I         KappaRT,       I         KappaRT,
539       U         fVerT,       U         fVerT,
540       I         myTime, myThid)       I         myTime,myIter,myThid)
541             CALL TIMESTEP_TRACER(             CALL TIMESTEP_TRACER(
542       I         bi,bj,iMin,iMax,jMin,jMax,k,tempAdvScheme,       I         bi,bj,iMin,iMax,jMin,jMax,k,tempAdvScheme,
543       I         theta, gT,       I         theta, gT,
      U         gTnm1,  
544       I         myIter, myThid)       I         myIter, myThid)
545           ENDIF           ENDIF
546           IF ( saltStepping ) THEN           IF ( saltStepping ) THEN
# Line 523  C        and step forward storing result Line 549  C        and step forward storing result
549       I         xA,yA,uTrans,vTrans,rTrans,maskUp,       I         xA,yA,uTrans,vTrans,rTrans,maskUp,
550       I         KappaRS,       I         KappaRS,
551       U         fVerS,       U         fVerS,
552       I         myTime, myThid)       I         myTime,myIter,myThid)
553             CALL TIMESTEP_TRACER(             CALL TIMESTEP_TRACER(
554       I         bi,bj,iMin,iMax,jMin,jMax,k,saltAdvScheme,       I         bi,bj,iMin,iMax,jMin,jMax,k,saltAdvScheme,
555       I         salt, gS,       I         salt, gS,
      U         gSnm1,  
556       I         myIter, myThid)       I         myIter, myThid)
557           ENDIF           ENDIF
558  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_PASSIVE_TRACER
# Line 537  C        and step forward storing result Line 562  C        and step forward storing result
562       I         xA,yA,uTrans,vTrans,rTrans,maskUp,       I         xA,yA,uTrans,vTrans,rTrans,maskUp,
563       I         KappaRT,       I         KappaRT,
564       U         fVerTr1,       U         fVerTr1,
565       I         myTime, myThid)       I         myTime,myIter,myThid)
566             CALL TIMESTEP_TRACER(             CALL TIMESTEP_TRACER(
567       I         bi,bj,iMin,iMax,jMin,jMax,k,tracerAdvScheme,       I         bi,bj,iMin,iMax,jMin,jMax,k,tracerAdvScheme,
568       I         Tr1, gTr1,       I         Tr1, gTr1,
569       U         gTr1NM1,       I         myIter,myThid)
      I         myIter, myThid)  
570           ENDIF           ENDIF
571  #endif  #endif
572    #ifdef ALLOW_PTRACERS
573             IF ( usePTRACERS ) THEN
574               CALL PTRACERS_INTEGERATE(
575         I         bi,bj,k,
576         I         xA,yA,uTrans,vTrans,rTrans,maskUp,
577         X         KappaRS,
578         I         myIter,myTime,myThid)
579             ENDIF
580    #endif /* ALLOW_PTRACERS */
581    
582  #ifdef   ALLOW_OBCS  #ifdef   ALLOW_OBCS
583  C--      Apply open boundary conditions  C--      Apply open boundary conditions
584           IF (useOBCS) THEN           IF (useOBCS) THEN
585             CALL OBCS_APPLY_TS( bi, bj, k, gTnm1, gSnm1, myThid )             CALL OBCS_APPLY_TS( bi, bj, k, gT, gS, myThid )
586           END IF           END IF
587  #endif   /* ALLOW_OBCS */  #endif   /* ALLOW_OBCS */
588    
589  C--      Freeze water  C--      Freeze water
590           IF (allowFreezing) THEN           IF (allowFreezing) THEN
591  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
592  CADJ STORE gTNm1(:,:,k,bi,bj) = comlev1_bibj_k  CADJ STORE gT(:,:,k,bi,bj) = comlev1_bibj_k
593  CADJ &   , key = kkey, byte = isbyte  CADJ &   , key = kkey, byte = isbyte
594  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
595              CALL FREEZE( bi, bj, iMin, iMax, jMin, jMax, k, myThid )              CALL FREEZE( bi, bj, iMin, iMax, jMin, jMax, k, myThid )
# Line 568  C--     end of thermodynamic k loop (Nr: Line 601  C--     end of thermodynamic k loop (Nr:
601    
602  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
603  C? Patrick? What about this one?  C? Patrick? What about this one?
604  cph Keys iikey and idkey don't seem to be needed  cph Keys iikey and idkey dont seem to be needed
605  cph since storing occurs on different tape for each  cph since storing occurs on different tape for each
606  cph impldiff call anyways.  cph impldiff call anyways.
607  cph Thus, common block comlev1_impl isn't needed either.  cph Thus, common block comlev1_impl isnt needed either.
608  cph Storing below needed in the case useGMREDI.  cph Storing below needed in the case useGMREDI.
609          iikey = (ikey-1)*maximpl          iikey = (ikey-1)*maximpl
610  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
# Line 582  C--     Implicit diffusion Line 615  C--     Implicit diffusion
615           IF (tempStepping) THEN           IF (tempStepping) THEN
616  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
617              idkey = iikey + 1              idkey = iikey + 1
618  CADJ STORE gTNm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte  CADJ STORE gT(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte
619  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
620              CALL IMPLDIFF(              CALL IMPLDIFF(
621       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
622       I         deltaTtracer, KappaRT, recip_HFacC,       I         deltaTtracer, KappaRT, recip_HFacC,
623       U         gTNm1,       U         gT,
624       I         myThid )       I         myThid )
625           ENDIF           ENDIF
626    
627           IF (saltStepping) THEN           IF (saltStepping) THEN
628  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
629           idkey = iikey + 2           idkey = iikey + 2
630  CADJ STORE gSNm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte  CADJ STORE gS(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte
631  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
632              CALL IMPLDIFF(              CALL IMPLDIFF(
633       I         bi, bj, iMin, iMax, jMin, jMax,       I         bi, bj, iMin, iMax, jMin, jMax,
634       I         deltaTtracer, KappaRS, recip_HFacC,       I         deltaTtracer, KappaRS, recip_HFacC,
635       U         gSNm1,       U         gS,
636       I         myThid )       I         myThid )
637           ENDIF           ENDIF
638    
639  #ifdef ALLOW_PASSIVE_TRACER  #ifdef ALLOW_PASSIVE_TRACER
640           IF (tr1Stepping) THEN           IF (tr1Stepping) THEN
641  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
642  CADJ STORE gTr1Nm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte  CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte
643  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
644            CALL IMPLDIFF(            CALL IMPLDIFF(
645       I      bi, bj, iMin, iMax, jMin, jMax,       I      bi, bj, iMin, iMax, jMin, jMax,
646       I      deltaTtracer, KappaRT, recip_HFacC,       I      deltaTtracer, KappaRT, recip_HFacC,
647       U      gTr1Nm1,       U      gTr1,
648       I      myThid )       I      myThid )
649           ENDIF           ENDIF
650  #endif  #endif
651    
652    #ifdef ALLOW_PTRACERS
653    C Vertical diffusion (implicit) for passive tracers
654             IF ( usePTRACERS ) THEN
655               CALL PTRACERS_IMPLDIFF( bi,bj,KappaRS,myThid )
656             ENDIF
657    #endif /* ALLOW_PTRACERS */
658    
659  #ifdef   ALLOW_OBCS  #ifdef   ALLOW_OBCS
660  C--      Apply open boundary conditions  C--      Apply open boundary conditions
661           IF (useOBCS) THEN           IF (useOBCS) THEN
662             DO K=1,Nr             DO K=1,Nr
663               CALL OBCS_APPLY_TS( bi, bj, k, gTnm1, gSnm1, myThid )               CALL OBCS_APPLY_TS( bi, bj, k, gT, gS, myThid )
664             ENDDO             ENDDO
665           END IF           END IF
666  #endif   /* ALLOW_OBCS */  #endif   /* ALLOW_OBCS */
# Line 636  Ccs- Line 676  Ccs-
676        IF ( useAIM ) THEN        IF ( useAIM ) THEN
677         CALL AIM_AIM2DYN_EXCHANGES( myTime, myThid )         CALL AIM_AIM2DYN_EXCHANGES( myTime, myThid )
678        ENDIF        ENDIF
679         _EXCH_XYZ_R8(gTnm1,myThid)         _EXCH_XYZ_R8(gT,myThid)
680         _EXCH_XYZ_R8(gSnm1,myThid)         _EXCH_XYZ_R8(gS,myThid)
681  #else  #else
682        IF (staggerTimeStep.AND.useCubedSphereExchange) THEN        IF (staggerTimeStep.AND.useCubedSphereExchange) THEN
683         _EXCH_XYZ_R8(gTnm1,myThid)         _EXCH_XYZ_R8(gT,myThid)
684         _EXCH_XYZ_R8(gSnm1,myThid)         _EXCH_XYZ_R8(gS,myThid)
685        ENDIF        ENDIF
686  #endif /* ALLOW_AIM */  #endif /* ALLOW_AIM */
687    
688    #ifndef DISABLE_DEBUGMODE
689          If (debugMode) THEN
690           CALL DEBUG_STATS_RL(Nr,uVel,'Uvel (THERMODYNAMICS)',myThid)
691           CALL DEBUG_STATS_RL(Nr,vVel,'Vvel (THERMODYNAMICS)',myThid)
692           CALL DEBUG_STATS_RL(Nr,wVel,'Wvel (THERMODYNAMICS)',myThid)
693           CALL DEBUG_STATS_RL(Nr,theta,'Theta (THERMODYNAMICS)',myThid)
694           CALL DEBUG_STATS_RL(Nr,salt,'Salt (THERMODYNAMICS)',myThid)
695           CALL DEBUG_STATS_RL(Nr,Gt,'Gt (THERMODYNAMICS)',myThid)
696           CALL DEBUG_STATS_RL(Nr,Gs,'Gs (THERMODYNAMICS)',myThid)
697           CALL DEBUG_STATS_RL(Nr,GtNm1,'GtNm1 (THERMODYNAMICS)',myThid)
698           CALL DEBUG_STATS_RL(Nr,GsNm1,'GsNm1 (THERMODYNAMICS)',myThid)
699           CALL PTRACERS_DEBUG(myThid)
700          ENDIF
701    #endif
702    
703        RETURN        RETURN
704        END        END
Legend:
Removed from v.1.4  
changed lines
  Added in v.1.17
  ViewVC Help
Powered by ViewVC 1.1.22