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