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

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

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

revision 1.31 by heimbach, Sun Oct 22 01:11:44 2006 UTC revision 1.121 by jmc, Mon Jan 21 23:04:02 2013 UTC
# Line 14  C $Name$ Line 14  C $Name$
14  # ifdef ALLOW_SEAICE  # ifdef ALLOW_SEAICE
15  #  include "SEAICE_OPTIONS.h"  #  include "SEAICE_OPTIONS.h"
16  # endif  # endif
17    # ifdef ALLOW_EXF
18    #  include "EXF_OPTIONS.h"
19    # endif
20  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
21    
22  CBOP  CBOP
# Line 36  C     == Global variables === Line 39  C     == Global variables ===
39  #include "SIZE.h"  #include "SIZE.h"
40  #include "EEPARAMS.h"  #include "EEPARAMS.h"
41  #include "PARAMS.h"  #include "PARAMS.h"
 #include "DYNVARS.h"  
42  #include "GRID.h"  #include "GRID.h"
43    #include "DYNVARS.h"
44  #ifdef ALLOW_TIMEAVE  #ifdef ALLOW_TIMEAVE
45  #include "TIMEAVE_STATV.h"  #include "TIMEAVE_STATV.h"
46  #endif  #endif
47  #if defined (ALLOW_BALANCE_FLUXES) && !(defined ALLOW_AUTODIFF_TAMC)  #ifdef ALLOW_BALANCE_FLUXES
48  #include "FFIELDS.h"  #include "FFIELDS.h"
49  #endif  #endif
50    
51  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
52    # include "AUTODIFF_MYFIELDS.h"
53  # include "tamc.h"  # include "tamc.h"
54  # include "tamc_keys.h"  # include "tamc_keys.h"
55    #ifndef ALLOW_BALANCE_FLUXES
56  # include "FFIELDS.h"  # include "FFIELDS.h"
57    #endif
58    # include "SURFACE.h"
59  # include "EOS.h"  # include "EOS.h"
60  # ifdef ALLOW_KPP  # ifdef ALLOW_KPP
61  #  include "KPP.h"  #  include "KPP.h"
62  # endif  # endif
63    # ifdef ALLOW_GGL90
64    #  include "GGL90.h"
65    # endif
66  # ifdef ALLOW_GMREDI  # ifdef ALLOW_GMREDI
67  #  include "GMREDI.h"  #  include "GMREDI.h"
68  # endif  # endif
69  # ifdef ALLOW_EBM  # ifdef ALLOW_EBM
70  #  include "EBM.h"  #  include "EBM.h"
71  # endif  # endif
 # ifdef EXACT_CONSERV  
 #  include "SURFACE.h"  
 # endif  
72  # ifdef ALLOW_EXF  # ifdef ALLOW_EXF
73  #  include "ctrl.h"  #  include "ctrl.h"
74  #  include "exf_fields.h"  #  include "EXF_FIELDS.h"
 #  include "exf_clim_fields.h"  
75  #  ifdef ALLOW_BULKFORMULAE  #  ifdef ALLOW_BULKFORMULAE
76  #   include "exf_constants.h"  #   include "EXF_CONSTANTS.h"
77  #  endif  #  endif
78  # endif  # endif
79  # ifdef ALLOW_SEAICE  # ifdef ALLOW_SEAICE
80    #  include "SEAICE_SIZE.h"
81  #  include "SEAICE.h"  #  include "SEAICE.h"
82  # endif  # endif
83    # ifdef ALLOW_THSICE
84    #  include "THSICE_VARS.h"
85    # endif
86    # ifdef ALLOW_SALT_PLUME
87    #  include "SALT_PLUME.h"
88    # endif
89  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
90    
91  C     !INPUT/OUTPUT PARAMETERS:  C     !INPUT/OUTPUT PARAMETERS:
# Line 86  C     myThid :: Thread number for this i Line 99  C     myThid :: Thread number for this i
99    
100  C     !LOCAL VARIABLES:  C     !LOCAL VARIABLES:
101  C     == Local variables  C     == Local variables
102  C     rhoK, rhoKM1  :: Density at current level, and level above  C     rhoK, rhoKm1  :: Density at current level, and level above
103  C     iMin, iMax    :: Ranges and sub-block indices on which calculations  C     iMin, iMax    :: Ranges and sub-block indices on which calculations
104  C     jMin, jMax       are applied.  C     jMin, jMax       are applied.
105  C     bi, bj        :: tile indices  C     bi, bj        :: tile indices
106  C     i,j,k         :: loop indices  C     i,j,k         :: loop indices
107        _RL rhokp1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rhoKp1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
108        _RL rhokm1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rhoKm1  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
       _RL rhok    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
109        _RL sigmaX  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaX  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
110        _RL sigmaY  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaY  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
111        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL sigmaR  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
# Line 115  C--   dummy statement to end declaration Line 127  C--   dummy statement to end declaration
127  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
128    
129  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
130        IF ( debugLevel .GE. debLevB )        IF (debugMode) CALL DEBUG_ENTER('DO_OCEANIC_PHYS',myThid)
      &     CALL DEBUG_ENTER('DO_OCEANIC_PHYS',myThid)  
131  #endif  #endif
132    
133        doDiagsRho = 0        doDiagsRho = 0
134  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
135        IF ( useDiagnostics .AND. fluidIsWater ) THEN        IF ( useDiagnostics .AND. fluidIsWater ) THEN
136          IF ( DIAGNOSTICS_IS_ON('DRHODR  ',myThid) ) doDiagsRho = 1          IF ( DIAGNOSTICS_IS_ON('MXLDEPTH',myThid) )
137          IF ( DIAGNOSTICS_IS_ON('RHOANOSQ',myThid) .OR.       &       doDiagsRho = doDiagsRho + 1
138       &       DIAGNOSTICS_IS_ON('URHOMASS',myThid) .OR.          IF ( DIAGNOSTICS_IS_ON('DRHODR  ',myThid) )
139       &       DIAGNOSTICS_IS_ON('VRHOMASS',myThid) .OR.       &       doDiagsRho = doDiagsRho + 2
140       &       DIAGNOSTICS_IS_ON('WRHOMASS',myThid) .OR.          IF ( DIAGNOSTICS_IS_ON('WdRHO_P ',myThid) )
141       &       DIAGNOSTICS_IS_ON('WRHOMASS',myThid) ) doDiagsRho = 2       &       doDiagsRho = doDiagsRho + 4
142            IF ( DIAGNOSTICS_IS_ON('WdRHOdP ',myThid) )
143         &       doDiagsRho = doDiagsRho + 8
144        ENDIF        ENDIF
145  #endif /* ALLOW_DIAGNOSTICS */  #endif /* ALLOW_DIAGNOSTICS */
146    
147    #ifdef  ALLOW_OBCS
148          IF (useOBCS) THEN
149    C--   Calculate future values on open boundaries
150    C--   moved before SEAICE_MODEL call since SEAICE_MODEL needs seaice-obcs fields
151    # ifdef ALLOW_AUTODIFF_TAMC
152    CADJ STORE theta = comlev1, key=ikey_dynamics, kind=isbyte
153    CADJ STORE salt  = comlev1, key=ikey_dynamics, kind=isbyte
154    # endif
155    # ifdef ALLOW_DEBUG
156           IF (debugMode) CALL DEBUG_CALL('OBCS_CALC',myThid)
157    # endif
158           CALL OBCS_CALC( myTime+deltaTClock, myIter+1,
159         I                 uVel, vVel, wVel, theta, salt, myThid )
160          ENDIF
161    #endif  /* ALLOW_OBCS */
162    
163    #ifdef ALLOW_AUTODIFF_TAMC
164    # ifdef ALLOW_SALT_PLUME
165          DO bj=myByLo(myThid),myByHi(myThid)
166           DO bi=myBxLo(myThid),myBxHi(myThid)
167            DO j=1-OLy,sNy+OLy
168             DO i=1-OLx,sNx+OLx
169              saltPlumeDepth(i,j,bi,bj) = 0. _d 0
170              saltPlumeFlux(i,j,bi,bj)  = 0. _d 0
171             ENDDO
172            ENDDO
173           ENDDO
174          ENDDO
175    # endif
176    #endif /* ALLOW_AUTODIFF_TAMC */
177    
178    #ifdef ALLOW_FRAZIL
179          IF ( useFRAZIL ) THEN
180    C--   Freeze water in the ocean interior and let it rise to the surface
181           CALL FRAZIL_CALC_RHS( myTime, myIter, myThid )
182          ENDIF
183    #endif /* ALLOW_FRAZIL */
184    
185    #ifndef OLD_THSICE_CALL_SEQUENCE
186    #if (defined ALLOW_THSICE) && !(defined ALLOW_ATM2D)
187          IF ( useThSIce .AND. fluidIsWater ) THEN
188    # ifdef ALLOW_DEBUG
189            IF (debugMode) CALL DEBUG_CALL('THSICE_MAIN',myThid)
190    # endif
191    C--     Step forward Therm.Sea-Ice variables
192    C       and modify forcing terms including effects from ice
193            CALL TIMER_START('THSICE_MAIN     [DO_OCEANIC_PHYS]', myThid)
194            CALL THSICE_MAIN( myTime, myIter, myThid )
195            CALL TIMER_STOP( 'THSICE_MAIN     [DO_OCEANIC_PHYS]', myThid)
196          ENDIF
197    #endif /* ALLOW_THSICE */
198    #endif /* ndef OLD_THSICE_CALL_SEQUENCE */
199    
200  #ifdef ALLOW_SEAICE  #ifdef ALLOW_SEAICE
 C--   Call sea ice model to compute forcing/external data fields.  In  
 C     addition to computing prognostic sea-ice variables and diagnosing the  
 C     forcing/external data fields that drive the ocean model, SEAICE_MODEL  
 C     also sets theta to the freezing point under sea-ice.  The implied  
 C     surface heat flux is then stored in variable surfaceTendencyTice,  
 C     which is needed by KPP package (kpp_calc.F and kpp_transport_t.F)  
 C     to diagnose surface buoyancy fluxes and for the non-local transport  
 C     term.  Because this call precedes model thermodynamics, temperature  
 C     under sea-ice may not be "exactly" at the freezing point by the time  
 C     theta is dumped or time-averaged.  
201        IF ( useSEAICE ) THEN        IF ( useSEAICE ) THEN
202  #ifdef ALLOW_AUTODIFF_TAMC  # ifdef ALLOW_AUTODIFF_TAMC
203  CADJ STORE qnet,qsw            = comlev1, key = ikey_dynamics  cph-adj-test(
204  CADJ STORE aqh,precip,swdown   = comlev1, key = ikey_dynamics  CADJ STORE area   = comlev1, key=ikey_dynamics, kind=isbyte
205  CADJ STORE theta               = comlev1, key = ikey_dynamics  CADJ STORE hsnow  = comlev1, key=ikey_dynamics, kind=isbyte
206  # ifdef SEAICE_ALLOW_DYNAMICS  CADJ STORE heff   = comlev1, key=ikey_dynamics, kind=isbyte
207  CADJ STORE uvel,vvel           = comlev1, key = ikey_dynamics  CADJ STORE empmr,qsw,theta   = comlev1, key = ikey_dynamics,
208    CADJ &     kind = isbyte
209    cph-adj-test)
210    CADJ STORE atemp,aqh,precip    = comlev1, key = ikey_dynamics,
211    CADJ &     kind = isbyte
212    CADJ STORE swdown,lwdown       = comlev1, key = ikey_dynamics,
213    CADJ &     kind = isbyte
214    cph# ifdef EXF_READ_EVAP
215    CADJ STORE evap                = comlev1, key = ikey_dynamics,
216    CADJ &     kind = isbyte
217    cph# endif
218    CADJ STORE uvel,vvel           = comlev1, key = ikey_dynamics,
219    CADJ &     kind = isbyte
220    #  ifdef SEAICE_CGRID
221    CADJ STORE stressdivergencex   = comlev1, key = ikey_dynamics,
222    CADJ &     kind = isbyte
223    CADJ STORE stressdivergencey   = comlev1, key = ikey_dynamics,
224    CADJ &     kind = isbyte
225    #  endif
226    #  ifdef SEAICE_ALLOW_DYNAMICS
227    CADJ STORE uice                = comlev1, key = ikey_dynamics,
228    CADJ &     kind = isbyte
229    CADJ STORE vice                = comlev1, key = ikey_dynamics,
230    CADJ &     kind = isbyte
231    #   ifdef SEAICE_ALLOW_EVP
232    CADJ STORE seaice_sigma1       = comlev1, key = ikey_dynamics,
233    CADJ &     kind = isbyte
234    CADJ STORE seaice_sigma2       = comlev1, key = ikey_dynamics,
235    CADJ &     kind = isbyte
236    CADJ STORE seaice_sigma12      = comlev1, key = ikey_dynamics,
237    CADJ &     kind = isbyte
238    #   endif
239    #  endif
240    cph#  ifdef SEAICE_SALINITY
241    CADJ STORE salt                = comlev1, key = ikey_dynamics,
242    CADJ &     kind = isbyte
243    cph#  endif
244    #  ifdef ATMOSPHERIC_LOADING
245    CADJ STORE pload               = comlev1, key = ikey_dynamics,
246    CADJ &     kind = isbyte
247    CADJ STORE siceload            = comlev1, key = ikey_dynamics,
248    CADJ &     kind = isbyte
249    #  endif
250    #  ifdef NONLIN_FRSURF
251    CADJ STORE recip_hfacc         = comlev1, key = ikey_dynamics,
252    CADJ &     kind = isbyte
253    #  endif
254    #  ifdef ANNUAL_BALANCE
255    CADJ STORE balance_itcount     = comlev1, key = ikey_dynamics,
256    CADJ &     kind = isbyte
257    #  endif /* ANNUAL_BALANCE */
258    # endif
259    # ifdef ALLOW_DEBUG
260            IF (debugMode) CALL DEBUG_CALL('SEAICE_MODEL',myThid)
261  # endif  # endif
 #endif  
 #ifdef ALLOW_DEBUG  
         IF ( debugLevel .GE. debLevB )  
      &    CALL DEBUG_CALL('SEAICE_MODEL',myThid)  
 #endif  
262          CALL TIMER_START('SEAICE_MODEL    [DO_OCEANIC_PHYS]', myThid)          CALL TIMER_START('SEAICE_MODEL    [DO_OCEANIC_PHYS]', myThid)
263          CALL SEAICE_MODEL( myTime, myIter, myThid )          CALL SEAICE_MODEL( myTime, myIter, myThid )
264          CALL TIMER_STOP ('SEAICE_MODEL    [DO_OCEANIC_PHYS]', myThid)          CALL TIMER_STOP ('SEAICE_MODEL    [DO_OCEANIC_PHYS]', myThid)
265  #ifdef ALLOW_COST_ICE  # ifdef ALLOW_COST
266          CALL COST_ICE_TEST ( myTime, myIter, myThid )          CALL SEAICE_COST_SENSI ( myTime, myIter, myThid )
267  #endif  # endif
268        ENDIF        ENDIF
269  #endif /* ALLOW_SEAICE */  #endif /* ALLOW_SEAICE */
270    
271    #ifdef ALLOW_AUTODIFF_TAMC
272    CADJ STORE sst, sss           = comlev1, key = ikey_dynamics,
273    CADJ &     kind = isbyte
274    CADJ STORE qsw                = comlev1, key = ikey_dynamics,
275    CADJ &     kind = isbyte
276    # ifdef ALLOW_SEAICE
277    CADJ STORE area               = comlev1, key = ikey_dynamics,
278    CADJ &     kind = isbyte
279    # endif
280    #endif
281    
282    #ifdef OLD_THSICE_CALL_SEQUENCE
283  #if (defined ALLOW_THSICE) && !(defined ALLOW_ATM2D)  #if (defined ALLOW_THSICE) && !(defined ALLOW_ATM2D)
284        IF ( useThSIce .AND. fluidIsWater ) THEN        IF ( useThSIce .AND. fluidIsWater ) THEN
285  #ifdef ALLOW_DEBUG  # ifdef ALLOW_AUTODIFF_TAMC
286          IF ( debugLevel .GE. debLevB )  cph(
287       &    CALL DEBUG_CALL('THSICE_MAIN',myThid)  #  ifdef NONLIN_FRSURF
288  #endif  CADJ STORE uice,vice        = comlev1, key = ikey_dynamics,
289    CADJ &     kind = isbyte
290    CADJ STORE salt,theta       = comlev1, key = ikey_dynamics,
291    CADJ &     kind = isbyte
292    CADJ STORE qnet,qsw, empmr  = comlev1, key = ikey_dynamics,
293    CADJ &     kind = isbyte
294    CADJ STORE hFac_surfC       = comlev1, key = ikey_dynamics,
295    CADJ &     kind = isbyte
296    #  endif
297    # endif
298    # ifdef ALLOW_DEBUG
299            IF (debugMode) CALL DEBUG_CALL('THSICE_MAIN',myThid)
300    # endif
301  C--     Step forward Therm.Sea-Ice variables  C--     Step forward Therm.Sea-Ice variables
302  C       and modify forcing terms including effects from ice  C       and modify forcing terms including effects from ice
303          CALL TIMER_START('THSICE_MAIN     [DO_OCEANIC_PHYS]', myThid)          CALL TIMER_START('THSICE_MAIN     [DO_OCEANIC_PHYS]', myThid)
# Line 177  C       and modify forcing terms includi Line 305  C       and modify forcing terms includi
305          CALL TIMER_STOP( 'THSICE_MAIN     [DO_OCEANIC_PHYS]', myThid)          CALL TIMER_STOP( 'THSICE_MAIN     [DO_OCEANIC_PHYS]', myThid)
306        ENDIF        ENDIF
307  #endif /* ALLOW_THSICE */  #endif /* ALLOW_THSICE */
308    #endif /* OLD_THSICE_CALL_SEQUENCE */
309    
310  #ifdef ALLOW_SHELFICE  #ifdef ALLOW_SHELFICE
311    # ifdef ALLOW_AUTODIFF_TAMC
312    CADJ STORE salt, theta = comlev1, key = ikey_dynamics,
313    CADJ &     kind = isbyte
314    # endif
315        IF ( useShelfIce .AND. fluidIsWater ) THEN        IF ( useShelfIce .AND. fluidIsWater ) THEN
316  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
317          IF ( debugLevel .GE. debLevB )         IF (debugMode) CALL DEBUG_CALL('SHELFICE_THERMODYNAMICS',myThid)
      &    CALL DEBUG_CALL('SHELFICE_THERMODYNAMICS',myThid)  
318  #endif  #endif
319  C     compute temperature and (virtual) salt flux at the  C     compute temperature and (virtual) salt flux at the
320  C     shelf-ice ocean interface  C     shelf-ice ocean interface
321         CALL TIMER_START('SHELFICE_THERMODYNAMICS [DO_OCEANIC_PHYS]',         CALL TIMER_START('SHELFICE_THERMODYNAMICS [DO_OCEANIC_PHYS]',
322       &       myThid)       &       myThid)
# Line 194  C     shelf-ice ocean interface Line 326  C     shelf-ice ocean interface
326        ENDIF        ENDIF
327  #endif /* ALLOW_SHELFICE */  #endif /* ALLOW_SHELFICE */
328    
329    #ifdef ALLOW_ICEFRONT
330          IF ( useICEFRONT .AND. fluidIsWater ) THEN
331    #ifdef ALLOW_DEBUG
332           IF (debugMode) CALL DEBUG_CALL('ICEFRONT_THERMODYNAMICS',myThid)
333    #endif
334    C     compute temperature and (virtual) salt flux at the
335    C     ice-front ocean interface
336           CALL TIMER_START('ICEFRONT_THERMODYNAMICS [DO_OCEANIC_PHYS]',
337         &       myThid)
338           CALL ICEFRONT_THERMODYNAMICS( myTime, myIter, myThid )
339           CALL TIMER_STOP( 'ICEFRONT_THERMODYNAMICS [DO_OCEANIC_PHYS]',
340         &      myThid)
341          ENDIF
342    #endif /* ALLOW_ICEFRONT */
343    
344    #ifdef ALLOW_SALT_PLUME
345          IF ( useSALT_PLUME ) THEN
346              CALL SALT_PLUME_DO_EXCH( myTime, myIter, myThid )
347          ENDIF
348    #endif /* ALLOW_SALT_PLUME */
349    
350  C--   Freeze water at the surface  C--   Freeze water at the surface
351          IF ( allowFreezing ) THEN
352  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
353  CADJ STORE theta = comlev1, key = ikey_dynamics  CADJ STORE theta = comlev1, key = ikey_dynamics,
354    CADJ &     kind = isbyte
355  #endif  #endif
       IF ( allowFreezing  
      &                   .AND. .NOT. useSEAICE  
      &                   .AND. .NOT. useThSIce ) THEN  
356          CALL FREEZE_SURFACE(  myTime, myIter, myThid )          CALL FREEZE_SURFACE(  myTime, myIter, myThid )
357        ENDIF        ENDIF
358    
359  #ifdef ALLOW_OCN_COMPON_INTERF  #ifdef ALLOW_OCN_COMPON_INTERF
360  C--    Apply imported data (from coupled interface) to forcing fields  C--    Apply imported data (from coupled interface) to forcing fields
361  C jmc: do not know precisely where to put this call (bf or af thSIce ?)  C jmc: do not know precisely where to put this call (bf or af thSIce ?)
362         IF ( useCoupler ) THEN        IF ( useCoupler ) THEN
363           CALL OCN_APPLY_IMPORT( .TRUE., myTime, myIter, myThid )           CALL OCN_APPLY_IMPORT( .TRUE., myTime, myIter, myThid )
364         ENDIF        ENDIF
365  #endif /* ALLOW_OCN_COMPON_INTERF */  #endif /* ALLOW_OCN_COMPON_INTERF */
366    
367  #ifdef ALLOW_BALANCE_FLUXES  #ifdef ALLOW_BALANCE_FLUXES
368  C     balance fluxes  C     balance fluxes
369         IF ( balanceEmPmR )        IF ( balanceEmPmR .AND. (.NOT.useSeaice .OR. useThSIce) )
370       &        CALL REMOVE_MEAN_RS( 1, EmPmR, maskH, maskH, rA, drF,       &      CALL REMOVE_MEAN_RS( 1, EmPmR, maskInC, maskInC, rA, drF,
371       &        'EmPmR', myTime, myThid )       &        'EmPmR', myTime, myThid )
372         IF ( balanceQnet )        IF ( balanceQnet  .AND. (.NOT.useSeaice .OR. useThSIce) )
373       &        CALL REMOVE_MEAN_RS( 1, Qnet,  maskH, maskH, rA, drF,       &      CALL REMOVE_MEAN_RS( 1, Qnet,  maskInC, maskInC, rA, drF,
374       &        'Qnet ', myTime, myThid )       &        'Qnet ', myTime, myThid )
375  #endif /* ALLOW_BALANCE_FLUXES */  #endif /* ALLOW_BALANCE_FLUXES */
376    
377  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
378  C--   HPF directive to help TAMC  C--   HPF directive to help TAMC
379  CHPF$ INDEPENDENT  CHPF$ INDEPENDENT
380    #else  /* ALLOW_AUTODIFF_TAMC */
381    C     if fluid is not water, by-pass find_rho, gmredi, surfaceForcing
382    C     and all vertical mixing schemes, but keep OBCS_CALC
383          IF ( fluidIsWater ) THEN
384  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
385        DO bj=myByLo(myThid),myByHi(myThid)        DO bj=myByLo(myThid),myByHi(myThid)
386  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 252  C     just ensure that all memory refere Line 408  C     just ensure that all memory refere
408  C     point numbers. This prevents spurious hardware signals due to  C     point numbers. This prevents spurious hardware signals due to
409  C     uninitialised but inert locations.  C     uninitialised but inert locations.
410    
411    #ifdef ALLOW_AUTODIFF_TAMC
412          DO j=1-OLy,sNy+OLy          DO j=1-OLy,sNy+OLy
413           DO i=1-OLx,sNx+OLx           DO i=1-OLx,sNx+OLx
414            rhok   (i,j)   = 0. _d 0            rhoKm1 (i,j)   = 0. _d 0
415            rhoKM1 (i,j)   = 0. _d 0            rhoKp1 (i,j)   = 0. _d 0
           rhoKP1 (i,j)   = 0. _d 0  
416           ENDDO           ENDDO
417          ENDDO          ENDDO
418    #endif /* ALLOW_AUTODIFF_TAMC */
419    
420          DO k=1,Nr          DO k=1,Nr
421           DO j=1-OLy,sNy+OLy           DO j=1-OLy,sNy+OLy
422            DO i=1-OLx,sNx+OLx            DO i=1-OLx,sNx+OLx
423  C This is currently also used by IVDC and Diagnostics  C This is currently used by GMRedi, IVDC, MXL-depth  and Diagnostics
424             sigmaX(i,j,k) = 0. _d 0             sigmaX(i,j,k) = 0. _d 0
425             sigmaY(i,j,k) = 0. _d 0             sigmaY(i,j,k) = 0. _d 0
426             sigmaR(i,j,k) = 0. _d 0             sigmaR(i,j,k) = 0. _d 0
427  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
428  cph all the following init. are necessary for TAF  cph all the following init. are necessary for TAF
429  cph although some of these are re-initialised later.  cph although some of these are re-initialised later.
430               rhoInSitu(i,j,k,bi,bj) = 0.
431             IVDConvCount(i,j,k,bi,bj) = 0.             IVDConvCount(i,j,k,bi,bj) = 0.
432  # ifdef ALLOW_GMREDI  # ifdef ALLOW_GMREDI
433             Kwx(i,j,k,bi,bj)  = 0. _d 0             Kwx(i,j,k,bi,bj)  = 0. _d 0
# Line 291  cph although some of these are re-initia Line 449  cph although some of these are re-initia
449             VisbeckK(i,j,bi,bj)   = 0. _d 0             VisbeckK(i,j,bi,bj)   = 0. _d 0
450  #  endif  #  endif
451  # endif /* ALLOW_GMREDI */  # endif /* ALLOW_GMREDI */
452    # ifdef ALLOW_KPP
453               KPPdiffKzS(i,j,k,bi,bj)  = 0. _d 0
454               KPPdiffKzT(i,j,k,bi,bj)  = 0. _d 0
455    # endif /* ALLOW_KPP */
456    # ifdef ALLOW_GGL90
457               GGL90viscArU(i,j,k,bi,bj)  = 0. _d 0
458               GGL90viscArV(i,j,k,bi,bj)  = 0. _d 0
459               GGL90diffKr(i,j,k,bi,bj)  = 0. _d 0
460    # endif /* ALLOW_GGL90 */
461  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
462            ENDDO            ENDDO
463           ENDDO           ENDDO
# Line 302  cph although some of these are re-initia Line 469  cph although some of these are re-initia
469          jMax = sNy+OLy          jMax = sNy+OLy
470    
471  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
472  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
473  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ &     kind = isbyte
474    CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
475    CADJ &     kind = isbyte
476  CADJ STORE totphihyd(:,:,:,bi,bj)  CADJ STORE totphihyd(:,:,:,bi,bj)
477  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey,
478    CADJ &     kind = isbyte
479  # ifdef ALLOW_KPP  # ifdef ALLOW_KPP
480  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
481  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte  CADJ &     kind = isbyte
482    CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
483    CADJ &     kind = isbyte
484    # endif
485    # ifdef ALLOW_SALT_PLUME
486    CADJ STORE saltplumedepth(:,:,bi,bj) = comlev1_bibj, key=itdkey,
487    CADJ &     kind = isbyte
488  # endif  # endif
489  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
490    
491    C--   Always compute density (stored in common block) here; even when it is not
492    C     needed here, will be used anyway in calc_phi_hyd (data flow easier this way)
493  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
494          IF ( debugLevel .GE. debLevB )          IF (debugMode) CALL DEBUG_CALL('FIND_RHO_2D (xNr)',myThid)
495       &    CALL DEBUG_MSG('ENTERING UPWARD K LOOP',myThid)  #endif
496    #ifdef ALLOW_AUTODIFF_TAMC
497            IF ( fluidIsWater ) THEN
498    #endif /* ALLOW_AUTODIFF_TAMC */
499    #ifdef ALLOW_DOWN_SLOPE
500             IF ( useDOWN_SLOPE ) THEN
501               DO k=1,Nr
502                CALL DWNSLP_CALC_RHO(
503         I                  theta, salt,
504         O                  rhoInSitu(1-OLx,1-OLy,k,bi,bj),
505         I                  k, bi, bj, myTime, myIter, myThid )
506               ENDDO
507             ENDIF
508    #endif /* ALLOW_DOWN_SLOPE */
509    #ifdef ALLOW_BBL
510             IF ( useBBL ) THEN
511    C     pkg/bbl requires in-situ bbl density for depths equal to and deeper than the bbl.
512    C     To reduce computation and storage requirement, these densities are stored in the
513    C     dry grid boxes of rhoInSitu.  See BBL_CALC_RHO for details.
514               DO k=Nr,1,-1
515                CALL BBL_CALC_RHO(
516         I                  theta, salt,
517         O                  rhoInSitu,
518         I                  k, bi, bj, myTime, myIter, myThid )
519    
520               ENDDO
521             ENDIF
522    #endif /* ALLOW_BBL */
523             IF ( .NOT. ( useDOWN_SLOPE .OR. useBBL ) ) THEN
524               DO k=1,Nr
525                CALL FIND_RHO_2D(
526         I                iMin, iMax, jMin, jMax, k,
527         I                theta(1-OLx,1-OLy,k,bi,bj),
528         I                salt (1-OLx,1-OLy,k,bi,bj),
529         O                rhoInSitu(1-OLx,1-OLy,k,bi,bj),
530         I                k, bi, bj, myThid )
531               ENDDO
532             ENDIF
533    #ifdef ALLOW_AUTODIFF_TAMC
534            ELSE
535    C-        fluid is not water:
536              DO k=1,Nr
537               DO j=1-OLy,sNy+OLy
538                DO i=1-OLx,sNx+OLx
539                  rhoInSitu(i,j,k,bi,bj) = 0.
540                ENDDO
541               ENDDO
542              ENDDO
543            ENDIF
544    #endif /* ALLOW_AUTODIFF_TAMC */
545    
546    #ifdef ALLOW_DEBUG
547            IF (debugMode) CALL DEBUG_MSG('ENTERING UPWARD K LOOP',myThid)
548  #endif  #endif
549    
550  C--     Start of diagnostic loop  C--     Start of diagnostic loop
# Line 323  C--     Start of diagnostic loop Line 553  C--     Start of diagnostic loop
553  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
554  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?
555  C? Do we still need this?  C? Do we still need this?
556  cph kkey formula corrected.  cph kkey formula corrected.
557  cph Needed for rhok, rhokm1, in the case useGMREDI.  cph Needed for rhoK, rhoKm1, in the case useGMREDI.
558           kkey = (itdkey-1)*Nr + k            kkey = (itdkey-1)*Nr + k
559  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
560    
561    c#ifdef ALLOW_AUTODIFF_TAMC
562    cCADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey,
563    cCADJ &     kind = isbyte
564    cCADJ STORE salt(:,:,k,bi,bj)  = comlev1_bibj_k, key=kkey,
565    cCADJ &     kind = isbyte
566    c#endif /* ALLOW_AUTODIFF_TAMC */
567    
568  C--       Calculate gradients of potential density for isoneutral  C--       Calculate gradients of potential density for isoneutral
569  C         slope terms (e.g. GM/Redi tensor or IVDC diffusivity)  C         slope terms (e.g. GM/Redi tensor or IVDC diffusivity)
 c         IF ( k.GT.1 .AND. (useGMRedi.OR.ivdc_kappa.NE.0.) ) THEN  
570            IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.)            IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.)
571       &                   .OR. doDiagsRho.GE.1 ) THEN       &         .OR. usePP81 .OR. useMY82 .OR. useGGL90
572  #ifdef ALLOW_DEBUG       &         .OR. useSALT_PLUME .OR. doDiagsRho.GE.1 ) THEN
             IF ( debugLevel .GE. debLevB )  
      &       CALL DEBUG_CALL('FIND_RHO',myThid)  
 #endif  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  
 CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
             CALL FIND_RHO(  
      I        bi, bj, iMin, iMax, jMin, jMax, k, k,  
      I        theta, salt,  
      O        rhoK,  
      I        myThid )  
   
573              IF (k.GT.1) THEN              IF (k.GT.1) THEN
574  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
575  CADJ STORE theta(:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE theta(:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey,
576  CADJ STORE salt (:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ &     kind = isbyte
577  #endif /* ALLOW_AUTODIFF_TAMC */  CADJ STORE salt (:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey,
578               CALL FIND_RHO(  CADJ &     kind = isbyte
579       I        bi, bj, iMin, iMax, jMin, jMax, k-1, k,  CADJ STORE rhokm1 (bi,bj)       = comlev1_bibj_k, key=kkey,
580       I        theta, salt,  CADJ &     kind = isbyte
581       O        rhoKm1,  #endif /* ALLOW_AUTODIFF_TAMC */
582       I        myThid )               CALL FIND_RHO_2D(
583         I                 iMin, iMax, jMin, jMax, k,
584         I                 theta(1-OLx,1-OLy,k-1,bi,bj),
585         I                 salt (1-OLx,1-OLy,k-1,bi,bj),
586         O                 rhoKm1,
587         I                 k-1, bi, bj, myThid )
588              ENDIF              ENDIF
589  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
590              IF ( debugLevel .GE. debLevB )              IF (debugMode) CALL DEBUG_CALL('GRAD_SIGMA',myThid)
      &       CALL DEBUG_CALL('GRAD_SIGMA',myThid)  
591  #endif  #endif
592  cph Avoid variable aliasing for adjoint !!!  cph Avoid variable aliasing for adjoint !!!
593              DO j=jMin,jMax              DO j=jMin,jMax
594               DO i=iMin,iMax               DO i=iMin,iMax
595                rhoKP1(i,j) = rhoK(i,j)                rhoKp1(i,j) = rhoInSitu(i,j,k,bi,bj)
596               ENDDO               ENDDO
597              ENDDO              ENDDO
598              CALL GRAD_SIGMA(              CALL GRAD_SIGMA(
599       I             bi, bj, iMin, iMax, jMin, jMax, k,       I             bi, bj, iMin, iMax, jMin, jMax, k,
600       I             rhoK, rhoKm1, rhoKp1,       I             rhoInSitu(1-OLx,1-OLy,k,bi,bj), rhoKm1, rhoKp1,
601       O             sigmaX, sigmaY, sigmaR,       O             sigmaX, sigmaY, sigmaR,
602       I             myThid )       I             myThid )
           ENDIF  
   
603  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
604  ctest# ifndef GM_EXCLUDE_CLIPPING  #ifdef GMREDI_WITH_STABLE_ADJOINT
605  CADJ STORE rhok   (:,:) = comlev1_bibj_k ,       key=kkey, byte=isbyte  cgf zero out adjoint fields to stabilize pkg/gmredi adjoint
606  ctest# endif  cgf -> cuts adjoint dependency from slope to state
607  CADJ STORE rhokm1 (:,:) = comlev1_bibj_k ,       key=kkey, byte=isbyte              CALL ZERO_ADJ_LOC( Nr, sigmaX, myThid)
608                CALL ZERO_ADJ_LOC( Nr, sigmaY, myThid)
609                CALL ZERO_ADJ_LOC( Nr, sigmaR, myThid)
610    #endif
611  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
612              ENDIF
613    
614  C--       Implicit Vertical Diffusion for Convection  C--       Implicit Vertical Diffusion for Convection
 c ==> should use sigmaR !!!  
615            IF (k.GT.1 .AND. ivdc_kappa.NE.0.) THEN            IF (k.GT.1 .AND. ivdc_kappa.NE.0.) THEN
616  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
617              IF ( debugLevel .GE. debLevB )              IF (debugMode) CALL DEBUG_CALL('CALC_IVDC',myThid)
      &       CALL DEBUG_CALL('CALC_IVDC',myThid)  
618  #endif  #endif
619              CALL CALC_IVDC(              CALL CALC_IVDC(
620       I        bi, bj, iMin, iMax, jMin, jMax, k,       I        bi, bj, iMin, iMax, jMin, jMax, k,
621       I        rhoKm1, rhoK,       I        sigmaR,
622       I        myTime, myIter, myThid)       I        myTime, myIter, myThid)
623            ENDIF            ENDIF
624    
625  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
626            IF ( doDiagsRho.GE.2 ) THEN            IF ( doDiagsRho.GE.4 ) THEN
627              CALL DIAGS_RHO( k, bi, bj,              CALL DIAGS_RHO_L( doDiagsRho, k, bi, bj,
628       I                      rhoK, rhoKm1,       I                        rhoInSitu(1-OLx,1-OLy,1,bi,bj),
629       I                      myTime, myIter, myThid)       I                        rhoKm1, wVel,
630         I                        myTime, myIter, myThid )
631            ENDIF            ENDIF
632  #endif  #endif
633    
634  C--     end of diagnostic k loop (Nr:1)  C--     end of diagnostic k loop (Nr:1)
635          ENDDO          ENDDO
636    
637    #ifdef ALLOW_AUTODIFF_TAMC
638    CADJ STORE IVDConvCount(:,:,:,bi,bj)
639    CADJ &     = comlev1_bibj, key=itdkey,
640    CADJ &     kind = isbyte
641    #endif
642    
643    C--     Diagnose Mixed Layer Depth:
644            IF ( useGMRedi .OR. MOD(doDiagsRho,2).EQ.1 ) THEN
645              CALL CALC_OCE_MXLAYER(
646         I              rhoInSitu(1-OLx,1-OLy,1,bi,bj), sigmaR,
647         I              bi, bj, myTime, myIter, myThid )
648            ENDIF
649    
650    #ifdef ALLOW_SALT_PLUME
651            IF ( useSALT_PLUME ) THEN
652              CALL SALT_PLUME_CALC_DEPTH(
653         I              rhoInSitu(1-OLx,1-OLy,1,bi,bj), sigmaR,
654         I              bi, bj, myTime, myIter, myThid )
655            ENDIF
656    #endif /* ALLOW_SALT_PLUME */
657    
658  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
659  c       IF ( useDiagnostics .AND.          IF ( MOD(doDiagsRho,4).GE.2 ) THEN
 c    &       (useGMRedi .OR. ivdc_kappa.NE.0.) ) THEN  
         IF ( doDiagsRho.GE.1 ) THEN  
660            CALL DIAGNOSTICS_FILL (sigmaR, 'DRHODR  ', 0, Nr,            CALL DIAGNOSTICS_FILL (sigmaR, 'DRHODR  ', 0, Nr,
661       &         2, bi, bj, myThid)       &         2, bi, bj, myThid)
662          ENDIF          ENDIF
663  #endif  #endif /* ALLOW_DIAGNOSTICS */
   
 #ifdef  ALLOW_OBCS  
 C--     Calculate future values on open boundaries  
         IF (useOBCS) THEN  
 #ifdef ALLOW_DEBUG  
           IF ( debugLevel .GE. debLevB )  
      &     CALL DEBUG_CALL('OBCS_CALC',myThid)  
 #endif  
           CALL OBCS_CALC( bi, bj, myTime+deltaTclock, myIter+1,  
      I            uVel, vVel, wVel, theta, salt,  
      I            myThid )  
         ENDIF  
 #endif  /* ALLOW_OBCS */  
664    
 #ifndef ALLOW_AUTODIFF_TAMC  
         IF ( fluidIsWater ) THEN  
 #endif  
665  C--     Determines forcing terms based on external fields  C--     Determines forcing terms based on external fields
666  C       relaxation terms, etc.  C       relaxation terms, etc.
667  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
668          IF ( debugLevel .GE. debLevB )          IF (debugMode) CALL DEBUG_CALL('EXTERNAL_FORCING_SURF',myThid)
      &    CALL DEBUG_CALL('EXTERNAL_FORCING_SURF',myThid)  
669  #endif  #endif
670  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
671  CADJ STORE EmPmR(:,:,bi,bj)  CADJ STORE EmPmR(:,:,bi,bj)
672  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey,
673    CADJ &     kind = isbyte
674  # ifdef EXACT_CONSERV  # ifdef EXACT_CONSERV
675  CADJ STORE PmEpR(:,:,bi,bj)  CADJ STORE PmEpR(:,:,bi,bj)
676  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey,
677    CADJ &     kind = isbyte
678  # endif  # endif
679  # ifdef NONLIN_FRSURF  # ifdef NONLIN_FRSURF
680  CADJ STORE hFac_surfC(:,:,bi,bj)  CADJ STORE hFac_surfC(:,:,bi,bj)
681  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey,
682    CADJ &     kind = isbyte
683  CADJ STORE recip_hFacC(:,:,:,bi,bj)  CADJ STORE recip_hFacC(:,:,:,bi,bj)
684  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey,
685  # endif  CADJ &     kind = isbyte
686    #  if (defined (ALLOW_PTRACERS))
687    CADJ STORE surfaceForcingS(:,:,bi,bj)   = comlev1_bibj, key=itdkey,
688    CADJ &     kind = isbyte
689    CADJ STORE surfaceForcingT(:,:,bi,bj)   = comlev1_bibj, key=itdkey,
690    CADJ &     kind = isbyte
691    #  endif /* ALLOW_PTRACERS */
692    # endif /* NONLIN_FRSURF */
693  #endif  #endif
694           CALL EXTERNAL_FORCING_SURF(          CALL EXTERNAL_FORCING_SURF(
695       I             bi, bj, iMin, iMax, jMin, jMax,       I             bi, bj, iMin, iMax, jMin, jMax,
696       I             myTime, myIter, myThid )       I             myTime, myIter, myThid )
 #ifndef ALLOW_AUTODIFF_TAMC  
         ENDIF  
 #endif  
697  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
698  # ifdef EXACT_CONSERV  # ifdef EXACT_CONSERV
699  cph-test  cph-test
700  cphCADJ STORE PmEpR(:,:,bi,bj)  cphCADJ STORE PmEpR(:,:,bi,bj)
701  cphCADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  cphCADJ &     = comlev1_bibj, key=itdkey,
702    cphCADJ &     kind = isbyte
703  # endif  # endif
704  #endif  #endif
705    
706  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
707  cph needed for KPP  cph needed for KPP
708  CADJ STORE surfaceForcingU(:,:,bi,bj)  CADJ STORE surfaceForcingU(:,:,bi,bj)
709  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey,
710    CADJ &     kind = isbyte
711  CADJ STORE surfaceForcingV(:,:,bi,bj)  CADJ STORE surfaceForcingV(:,:,bi,bj)
712  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey,
713    CADJ &     kind = isbyte
714  CADJ STORE surfaceForcingS(:,:,bi,bj)  CADJ STORE surfaceForcingS(:,:,bi,bj)
715  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey,
716    CADJ &     kind = isbyte
717  CADJ STORE surfaceForcingT(:,:,bi,bj)  CADJ STORE surfaceForcingT(:,:,bi,bj)
718  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey,
719    CADJ &     kind = isbyte
720  CADJ STORE surfaceForcingTice(:,:,bi,bj)  CADJ STORE surfaceForcingTice(:,:,bi,bj)
721  CADJ &     = comlev1_bibj, key=itdkey, byte=isbyte  CADJ &     = comlev1_bibj, key=itdkey,
722  #endif /* ALLOW_AUTODIFF_TAMC */  CADJ &     kind = isbyte
   
 #ifdef  ALLOW_GMREDI  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 # ifndef GM_EXCLUDE_CLIPPING  
 cph storing here is needed only for one GMREDI_OPTIONS:  
 cph define GM_BOLUS_ADVEC  
 cph keep it although TAF says you dont need to.  
 cph but I've avoided the #ifdef for now, in case more things change  
 CADJ STORE sigmaX(:,:,:)        = comlev1_bibj, key=itdkey, byte=isbyte  
 CADJ STORE sigmaY(:,:,:)        = comlev1_bibj, key=itdkey, byte=isbyte  
 CADJ STORE sigmaR(:,:,:)        = comlev1_bibj, key=itdkey, byte=isbyte  
 # endif  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
 C--     Calculate iso-neutral slopes for the GM/Redi parameterisation  
         IF (useGMRedi) THEN  
 #ifdef ALLOW_DEBUG  
           IF ( debugLevel .GE. debLevB )  
      &     CALL DEBUG_CALL('GMREDI_CALC_TENSOR',myThid)  
 #endif  
           CALL GMREDI_CALC_TENSOR(  
      I             bi, bj, iMin, iMax, jMin, jMax,  
      I             sigmaX, sigmaY, sigmaR,  
      I             myThid )  
 #ifdef ALLOW_AUTODIFF_TAMC  
         ELSE  
           CALL GMREDI_CALC_TENSOR_DUMMY(  
      I             bi, bj, iMin, iMax, jMin, jMax,  
      I             sigmaX, sigmaY, sigmaR,  
      I             myThid )  
723  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
         ENDIF  
   
 #endif  /* ALLOW_GMREDI */  
724    
725  #ifdef  ALLOW_KPP  #ifdef  ALLOW_KPP
726  C--     Compute KPP mixing coefficients  C--     Compute KPP mixing coefficients
727          IF (useKPP) THEN          IF (useKPP) THEN
728  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
729            IF ( debugLevel .GE. debLevB )            IF (debugMode) CALL DEBUG_CALL('KPP_CALC',myThid)
      &     CALL DEBUG_CALL('KPP_CALC',myThid)  
730  #endif  #endif
731              CALL TIMER_START('KPP_CALC [DO_OCEANIC_PHYS]', myThid)
732            CALL KPP_CALC(            CALL KPP_CALC(
733       I                  bi, bj, myTime, myThid )       I                  bi, bj, myTime, myIter, myThid )
734              CALL TIMER_STOP ('KPP_CALC [DO_OCEANIC_PHYS]', myThid)
735  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
736          ELSE          ELSE
737            CALL KPP_CALC_DUMMY(            CALL KPP_CALC_DUMMY(
738       I                  bi, bj, myTime, myThid )       I                  bi, bj, myTime, myIter, myThid )
739  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
740          ENDIF          ENDIF
741    
# Line 535  C--     Compute KPP mixing coefficients Line 745  C--     Compute KPP mixing coefficients
745  C--     Compute PP81 mixing coefficients  C--     Compute PP81 mixing coefficients
746          IF (usePP81) THEN          IF (usePP81) THEN
747  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
748            IF ( debugLevel .GE. debLevB )            IF (debugMode) CALL DEBUG_CALL('PP81_CALC',myThid)
      &     CALL DEBUG_CALL('PP81_CALC',myThid)  
749  #endif  #endif
750            CALL PP81_CALC(            CALL PP81_CALC(
751       I                  bi, bj, myTime, myThid )       I                     bi, bj, sigmaR, myTime, myIter, myThid )
752          ENDIF          ENDIF
753  #endif /* ALLOW_PP81 */  #endif /* ALLOW_PP81 */
754    
# Line 547  C--     Compute PP81 mixing coefficients Line 756  C--     Compute PP81 mixing coefficients
756  C--     Compute MY82 mixing coefficients  C--     Compute MY82 mixing coefficients
757          IF (useMY82) THEN          IF (useMY82) THEN
758  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
759            IF ( debugLevel .GE. debLevB )            IF (debugMode) CALL DEBUG_CALL('MY82_CALC',myThid)
      &     CALL DEBUG_CALL('MY82_CALC',myThid)  
760  #endif  #endif
761            CALL MY82_CALC(            CALL MY82_CALC(
762       I                  bi, bj, myTime, myThid )       I                     bi, bj, sigmaR, myTime, myIter, myThid )
763          ENDIF          ENDIF
764  #endif /* ALLOW_MY82 */  #endif /* ALLOW_MY82 */
765    
766  #ifdef  ALLOW_GGL90  #ifdef  ALLOW_GGL90
767    #ifdef ALLOW_AUTODIFF_TAMC
768    CADJ STORE GGL90TKE (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
769    CADJ &     kind = isbyte
770    #endif /* ALLOW_AUTODIFF_TAMC */
771  C--     Compute GGL90 mixing coefficients  C--     Compute GGL90 mixing coefficients
772          IF (useGGL90) THEN          IF (useGGL90) THEN
773  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
774            IF ( debugLevel .GE. debLevB )            IF (debugMode) CALL DEBUG_CALL('GGL90_CALC',myThid)
      &     CALL DEBUG_CALL('GGL90_CALC',myThid)  
775  #endif  #endif
776              CALL TIMER_START('GGL90_CALC [DO_OCEANIC_PHYS]', myThid)
777            CALL GGL90_CALC(            CALL GGL90_CALC(
778       I                  bi, bj, myTime, myThid )       I                     bi, bj, sigmaR, myTime, myIter, myThid )
779              CALL TIMER_STOP ('GGL90_CALC [DO_OCEANIC_PHYS]', myThid)
780          ENDIF          ENDIF
781  #endif /* ALLOW_GGL90 */  #endif /* ALLOW_GGL90 */
782    
783  #ifdef ALLOW_TIMEAVE  #ifdef ALLOW_TIMEAVE
784          IF ( taveFreq.GT. 0. _d 0 .AND. fluidIsWater ) THEN          IF ( taveFreq.GT. 0. _d 0 ) THEN
785            CALL TIMEAVE_SURF_FLUX( bi, bj, myTime, myIter, myThid)            CALL TIMEAVE_SURF_FLUX( bi, bj, myTime, myIter, myThid)
786          ENDIF          ENDIF
787          IF (taveFreq.GT.0. .AND. ivdc_kappa.NE.0.) THEN          IF (taveFreq.GT.0. .AND. ivdc_kappa.NE.0.) THEN
788            CALL TIMEAVE_CUMULATE(ConvectCountTave, IVDConvCount,            CALL TIMEAVE_CUMULATE(ConvectCountTave, IVDConvCount,
789       I                           Nr, deltaTclock, bi, bj, myThid)       I                           Nr, deltaTClock, bi, bj, myThid)
790          ENDIF          ENDIF
791  #endif /* ALLOW_TIMEAVE */  #endif /* ALLOW_TIMEAVE */
792    
793    #ifdef ALLOW_GMREDI
794    #ifdef ALLOW_AUTODIFF_TAMC
795    # ifndef GM_EXCLUDE_CLIPPING
796    cph storing here is needed only for one GMREDI_OPTIONS:
797    cph define GM_BOLUS_ADVEC
798    cph keep it although TAF says you dont need to.
799    cph but I have avoided the #ifdef for now, in case more things change
800    CADJ STORE sigmaX(:,:,:)        = comlev1_bibj, key=itdkey,
801    CADJ &     kind = isbyte
802    CADJ STORE sigmaY(:,:,:)        = comlev1_bibj, key=itdkey,
803    CADJ &     kind = isbyte
804    CADJ STORE sigmaR(:,:,:)        = comlev1_bibj, key=itdkey,
805    CADJ &     kind = isbyte
806    # endif
807    #endif /* ALLOW_AUTODIFF_TAMC */
808    
809    C--     Calculate iso-neutral slopes for the GM/Redi parameterisation
810            IF (useGMRedi) THEN
811    #ifdef ALLOW_DEBUG
812              IF (debugMode) CALL DEBUG_CALL('GMREDI_CALC_TENSOR',myThid)
813    #endif
814              CALL GMREDI_CALC_TENSOR(
815         I             iMin, iMax, jMin, jMax,
816         I             sigmaX, sigmaY, sigmaR,
817         I             bi, bj, myTime, myIter, myThid )
818    #ifdef ALLOW_AUTODIFF_TAMC
819            ELSE
820              CALL GMREDI_CALC_TENSOR_DUMMY(
821         I             iMin, iMax, jMin, jMax,
822         I             sigmaX, sigmaY, sigmaR,
823         I             bi, bj, myTime, myIter, myThid )
824    #endif /* ALLOW_AUTODIFF_TAMC */
825            ENDIF
826    #endif /* ALLOW_GMREDI */
827    
828    #ifdef ALLOW_DOWN_SLOPE
829            IF ( useDOWN_SLOPE ) THEN
830    C--     Calculate Downsloping Flow for Down_Slope parameterization
831             IF ( usingPCoords ) THEN
832              CALL DWNSLP_CALC_FLOW(
833         I                bi, bj, kSurfC, rhoInSitu,
834         I                myTime, myIter, myThid )
835             ELSE
836              CALL DWNSLP_CALC_FLOW(
837         I                bi, bj, kLowC, rhoInSitu,
838         I                myTime, myIter, myThid )
839             ENDIF
840            ENDIF
841    #endif /* ALLOW_DOWN_SLOPE */
842    
843  C--   end bi,bj loops.  C--   end bi,bj loops.
844         ENDDO         ENDDO
845        ENDDO        ENDDO
846    
847    #ifdef ALLOW_BALANCE_RELAX
848    # ifdef ALLOW_AUTODIFF_TAMC
849    CADJ STORE SSSrlx = comlev1, key=ikey_dynamics, kind=isbyte
850    CADJ STORE SSSrlxTile = comlev1, key=ikey_dynamics, kind=isbyte
851    CADJ STORE SSSrlxGlob = comlev1, key=ikey_dynamics, kind=isbyte
852    CADJ STORE SSTrlx = comlev1, key=ikey_dynamics, kind=isbyte
853    CADJ STORE SSTrlxTile = comlev1, key=ikey_dynamics, kind=isbyte
854    CADJ STORE SSTrlxGlob = comlev1, key=ikey_dynamics, kind=isbyte
855    # endif /* ALLOW_AUTODIFF_TAMC */
856           CALL BALANCE_RELAX( myTime, myIter, myThid )
857    #endif /* ALLOW_BALANCE_RELAX */
858    
859    #ifndef ALLOW_AUTODIFF_TAMC
860    C---  if fluid Is Water: end
861          ENDIF
862    #endif
863    
864    #ifdef ALLOW_BBL
865          IF ( useBBL ) THEN
866           CALL BBL_CALC_RHS(
867         I                          myTime, myIter, myThid )
868          ENDIF
869    #endif /* ALLOW_BBL */
870    
871    #ifdef ALLOW_MYPACKAGE
872          IF ( useMYPACKAGE ) THEN
873           CALL MYPACKAGE_CALC_RHS(
874         I                          myTime, myIter, myThid )
875          ENDIF
876    #endif /* ALLOW_MYPACKAGE */
877    
878    #ifdef ALLOW_GMREDI
879          IF ( useGMRedi ) THEN
880            CALL GMREDI_DO_EXCH( myTime, myIter, myThid )
881          ENDIF
882    #endif /* ALLOW_GMREDI */
883    
884    #ifdef ALLOW_KPP
885          IF (useKPP) THEN
886            CALL KPP_DO_EXCH( myThid )
887          ENDIF
888    #endif /* ALLOW_KPP */
889    
890  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
891        IF ( fluidIsWater .AND. useDiagnostics ) THEN        IF ( fluidIsWater .AND. useDiagnostics ) THEN
892            CALL DIAGS_RHO_G(
893         I                    rhoInSitu, uVel, vVel, wVel,
894         I                    myTime, myIter, myThid )
895          CALL DIAGS_OCEANIC_SURF_FLUX( myTime, myIter, myThid )          CALL DIAGS_OCEANIC_SURF_FLUX( myTime, myIter, myThid )
896        ENDIF        ENDIF
897        IF ( ivdc_kappa.NE.0 .AND. useDiagnostics ) THEN        IF ( ivdc_kappa.NE.0 .AND. useDiagnostics ) THEN
898          CALL DIAGNOSTICS_FILL( IVDConvCount,'CONVADJ ',          CALL DIAGNOSTICS_FILL( IVDConvCount, 'CONVADJ ',
899       &                         0, Nr, 0, 1, 1, myThid )       &                               0, Nr, 0, 1, 1, myThid )
900        ENDIF        ENDIF
901  #endif  #endif
902    
903  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
904        IF ( debugLevel .GE. debLevB )        IF (debugMode) CALL DEBUG_LEAVE('DO_OCEANIC_PHYS',myThid)
      &     CALL DEBUG_LEAVE('DO_OCEANIC_PHYS',myThid)  
905  #endif  #endif
906    
907        RETURN        RETURN
Legend:
Removed from v.1.31  
changed lines
  Added in v.1.121
  ViewVC Help
Powered by ViewVC 1.1.22