/[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.95 by heimbach, Thu Oct 28 14:36:47 2010 UTC revision 1.132 by jmc, Tue Jul 9 16:01:01 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 39  C     == Global variables === Line 42  C     == Global variables ===
42  #include "GRID.h"  #include "GRID.h"
43  #include "DYNVARS.h"  #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_OFFLINE
48  #include "FFIELDS.h"  # include "OFFLINE_SWITCH.h"
49  #endif  #endif
50    
51  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 72  C     == Global variables === Line 75  C     == Global variables ===
75  #  endif  #  endif
76  # endif  # endif
77  # ifdef ALLOW_SEAICE  # ifdef ALLOW_SEAICE
78    #  include "SEAICE_SIZE.h"
79  #  include "SEAICE.h"  #  include "SEAICE.h"
80    #  include "SEAICE_PARAMS.h"
81    # endif
82    # ifdef ALLOW_THSICE
83    #  include "THSICE_VARS.h"
84  # endif  # endif
85  # ifdef ALLOW_SALT_PLUME  # ifdef ALLOW_SALT_PLUME
86  #  include "SALT_PLUME.h"  #  include "SALT_PLUME.h"
# Line 105  C     i,j,k         :: loop indices Line 113  C     i,j,k         :: loop indices
113        INTEGER bi, bj        INTEGER bi, bj
114        INTEGER i, j, k        INTEGER i, j, k
115        INTEGER doDiagsRho        INTEGER doDiagsRho
116          LOGICAL calcGMRedi, calcKPP, calcConvect
117  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
118        LOGICAL  DIAGNOSTICS_IS_ON        LOGICAL  DIAGNOSTICS_IS_ON
119        EXTERNAL DIAGNOSTICS_IS_ON        EXTERNAL DIAGNOSTICS_IS_ON
# Line 118  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('WRHOMASS',myThid) )          IF ( DIAGNOSTICS_IS_ON('MXLDEPTH',myThid) )
137       &       doDiagsRho = doDiagsRho + 1       &       doDiagsRho = doDiagsRho + 1
138          IF ( DIAGNOSTICS_IS_ON('DRHODR  ',myThid) )          IF ( DIAGNOSTICS_IS_ON('DRHODR  ',myThid) )
139       &       doDiagsRho = doDiagsRho + 2       &       doDiagsRho = doDiagsRho + 2
140          IF ( DIAGNOSTICS_IS_ON('MXLDEPTH',myThid) )          IF ( DIAGNOSTICS_IS_ON('WdRHO_P ',myThid) )
141       &       doDiagsRho = doDiagsRho + 4       &       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_ADDFLUID        calcGMRedi  = useGMRedi
148        DO bj=myByLo(myThid),myByHi(myThid)        calcKPP     = useKPP
149         DO bi=myBxLo(myThid),myBxHi(myThid)        calcConvect = ivdc_kappa.NE.0.
150          DO K=1,Nr  #ifdef ALLOW_OFFLINE
151           DO j=1-OLy,sNy+OLy        IF ( useOffLine ) THEN
152            DO i=1-OLx,sNx+OLx          calcGMRedi = useGMRedi .AND. .NOT.offlineLoadGMRedi
153             addmass(I,J,K,bi,bj) =  0. _d 0          calcKPP    = useKPP    .AND. .NOT.offlineLoadKPP
154            ENDDO          calcConvect=calcConvect.AND. .NOT.offlineLoadConvec
155           ENDDO        ENDIF
156          ENDDO  #endif /* ALLOW_OFFLINE */
        ENDDO  
       ENDDO  
 #endif /* ALLOW_ADDFLUID */  
157    
158  #ifdef  ALLOW_OBCS  #ifdef  ALLOW_OBCS
159        IF (useOBCS) THEN        IF (useOBCS) THEN
160  C--   Calculate future values on open boundaries  C--   Calculate future values on open boundaries
161  C--   moved before SEAICE_MODEL call since SEAICE_MODEL needs seaice-obcs fields  C--   moved before SEAICE_MODEL call since SEAICE_MODEL needs seaice-obcs fields
162  #ifdef ALLOW_DEBUG  # ifdef ALLOW_AUTODIFF_TAMC
163         IF ( debugLevel .GE. debLevB )  CADJ STORE theta = comlev1, key=ikey_dynamics, kind=isbyte
164       &     CALL DEBUG_CALL('OBCS_CALC',myThid)  CADJ STORE salt  = comlev1, key=ikey_dynamics, kind=isbyte
165  #endif  # endif
166         CALL OBCS_CALC( myTime+deltaTclock, myIter+1,  # ifdef ALLOW_DEBUG
167           IF (debugMode) CALL DEBUG_CALL('OBCS_CALC',myThid)
168    # endif
169           CALL OBCS_CALC( myTime+deltaTClock, myIter+1,
170       I                 uVel, vVel, wVel, theta, salt, myThid )       I                 uVel, vVel, wVel, theta, salt, myThid )
171        ENDIF        ENDIF
172  #endif  /* ALLOW_OBCS */  #endif  /* ALLOW_OBCS */
# Line 176  C--   moved before SEAICE_MODEL call sin Line 186  C--   moved before SEAICE_MODEL call sin
186  # endif  # endif
187  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
188    
189    #ifdef ALLOW_FRAZIL
190          IF ( useFRAZIL ) THEN
191    C--   Freeze water in the ocean interior and let it rise to the surface
192           CALL FRAZIL_CALC_RHS( myTime, myIter, myThid )
193          ENDIF
194    #endif /* ALLOW_FRAZIL */
195    
196    #ifndef OLD_THSICE_CALL_SEQUENCE
197    #if (defined ALLOW_THSICE) && !(defined ALLOW_ATM2D)
198          IF ( useThSIce .AND. fluidIsWater ) THEN
199    # ifdef ALLOW_AUTODIFF_TAMC
200    CADJ STORE uice,vice         = comlev1, key = ikey_dynamics,
201    CADJ &     kind = isbyte
202    CADJ STORE iceMask,iceHeight = comlev1, key = ikey_dynamics,
203    CADJ &     kind = isbyte
204    CADJ STORE snowHeight, Tsrf  = comlev1, key = ikey_dynamics,
205    CADJ &     kind = isbyte
206    CADJ STORE Qice1, Qice2      = comlev1, key = ikey_dynamics,
207    CADJ &     kind = isbyte
208    CADJ STORE sHeating, snowAge = comlev1, key = ikey_dynamics,
209    CADJ &     kind = isbyte
210    CADJ STORE hocemxl = comlev1, key = ikey_dynamics,
211    CADJ &     kind = isbyte
212    CADJ STORE icflxsw = comlev1, key = ikey_dynamics,
213    CADJ &     kind = isbyte
214    CADJ STORE salt,theta        = comlev1, key = ikey_dynamics,
215    CADJ &     kind = isbyte
216    CADJ STORE uvel,vvel         = comlev1, key = ikey_dynamics,
217    CADJ &     kind = isbyte
218    CADJ STORE qnet,qsw, empmr   = comlev1, key = ikey_dynamics,
219    CADJ &     kind = isbyte
220    CADJ STORE atemp,aqh,precip  = comlev1, key = ikey_dynamics,
221    CADJ &     kind = isbyte
222    CADJ STORE swdown,lwdown     = comlev1, key = ikey_dynamics,
223    CADJ &     kind = isbyte
224    #  ifdef NONLIN_FRSURF
225    CADJ STORE hFac_surfC       = comlev1, key = ikey_dynamics,
226    CADJ &     kind = isbyte
227    #  endif
228    # endif /* ALLOW_AUTODIFF_TAMC */
229    # ifdef ALLOW_DEBUG
230            IF (debugMode) CALL DEBUG_CALL('THSICE_MAIN',myThid)
231    # endif
232    C--     Step forward Therm.Sea-Ice variables
233    C       and modify forcing terms including effects from ice
234            CALL TIMER_START('THSICE_MAIN     [DO_OCEANIC_PHYS]', myThid)
235            CALL THSICE_MAIN( myTime, myIter, myThid )
236            CALL TIMER_STOP( 'THSICE_MAIN     [DO_OCEANIC_PHYS]', myThid)
237          ENDIF
238    #endif /* ALLOW_THSICE */
239    #endif /* ndef OLD_THSICE_CALL_SEQUENCE */
240    
241    #ifdef ALLOW_SEAICE
242    # ifdef ALLOW_AUTODIFF_TAMC
243    CADJ STORE area   = comlev1, key=ikey_dynamics, kind=isbyte
244    CADJ STORE fu,fv  = comlev1, key=ikey_dynamics, kind=isbyte
245    CADJ STORE qnet   = comlev1, key=ikey_dynamics, kind=isbyte
246    CADJ STORE qsw    = comlev1, key=ikey_dynamics, kind=isbyte
247    CADJ STORE theta  = comlev1, key=ikey_dynamics, kind=isbyte
248    CADJ STORE salt   = comlev1, key=ikey_dynamics, kind=isbyte
249    #if (defined ALLOW_EXF) && (defined ALLOW_ATM_TEMP)
250    CADJ STORE evap   = comlev1, key=ikey_dynamics, kind=isbyte
251    #endif
252          IF ( .NOT.useSEAICE ) THEN
253          IF ( SEAICEadjMODE .EQ. -1 ) THEN
254            CALL SEAICE_FAKE( myTime, myIter, myThid )
255          ENDIF
256          ENDIF
257    CADJ STORE area   = comlev1, key=ikey_dynamics, kind=isbyte
258    CADJ STORE fu,fv  = comlev1, key=ikey_dynamics, kind=isbyte
259    CADJ STORE qnet   = comlev1, key=ikey_dynamics, kind=isbyte
260    CADJ STORE qsw    = comlev1, key=ikey_dynamics, kind=isbyte
261    CADJ STORE theta  = comlev1, key=ikey_dynamics, kind=isbyte
262    CADJ STORE salt   = comlev1, key=ikey_dynamics, kind=isbyte
263    #if (defined ALLOW_EXF) && (defined ALLOW_ATM_TEMP)
264    CADJ STORE evap   = comlev1, key=ikey_dynamics, kind=isbyte
265    #endif
266    # endif /* ALLOW_AUTODIFF_TAMC */
267    #endif /* ALLOW_SEAICE */
268    
269  #ifdef ALLOW_SEAICE  #ifdef ALLOW_SEAICE
270        IF ( useSEAICE ) THEN        IF ( useSEAICE ) THEN
271  # ifdef ALLOW_AUTODIFF_TAMC  # ifdef ALLOW_AUTODIFF_TAMC
# Line 183  cph-adj-test( Line 273  cph-adj-test(
273  CADJ STORE area   = comlev1, key=ikey_dynamics, kind=isbyte  CADJ STORE area   = comlev1, key=ikey_dynamics, kind=isbyte
274  CADJ STORE hsnow  = comlev1, key=ikey_dynamics, kind=isbyte  CADJ STORE hsnow  = comlev1, key=ikey_dynamics, kind=isbyte
275  CADJ STORE heff   = comlev1, key=ikey_dynamics, kind=isbyte  CADJ STORE heff   = comlev1, key=ikey_dynamics, kind=isbyte
276  CADJ STORE empmr,qsw,theta   = comlev1, key = ikey_dynamics,  CADJ STORE tices  = comlev1, key=ikey_dynamics, kind=isbyte
277  CADJ &     kind = isbyte  CADJ STORE empmr, qnet  = comlev1, key=ikey_dynamics, kind=isbyte
278    CADJ STORE qsw,saltflux = comlev1, key=ikey_dynamics, kind=isbyte
279    CADJ STORE fu, fv = comlev1, key=ikey_dynamics, kind=isbyte
280    cCADJ STORE theta  = comlev1, key=ikey_dynamics, kind=isbyte
281    cCADJ STORE salt   = comlev1, key=ikey_dynamics, kind=isbyte
282  cph-adj-test)  cph-adj-test)
283    c#ifdef ALLOW_EXF
284  CADJ STORE atemp,aqh,precip    = comlev1, key = ikey_dynamics,  CADJ STORE atemp,aqh,precip    = comlev1, key = ikey_dynamics,
285  CADJ &     kind = isbyte  CADJ &     kind = isbyte
286  CADJ STORE swdown,lwdown       = comlev1, key = ikey_dynamics,  CADJ STORE swdown,lwdown       = comlev1, key = ikey_dynamics,
287  CADJ &     kind = isbyte  CADJ &     kind = isbyte
 cph# ifdef EXF_READ_EVAP  
288  CADJ STORE evap                = comlev1, key = ikey_dynamics,  CADJ STORE evap                = comlev1, key = ikey_dynamics,
289  CADJ &     kind = isbyte  CADJ &     kind = isbyte
290  cph# endif  CADJ STORE uwind,vwind         = comlev1, key = ikey_dynamics,
291    CADJ &     kind = isbyte
292    c#endif
293  CADJ STORE uvel,vvel           = comlev1, key = ikey_dynamics,  CADJ STORE uvel,vvel           = comlev1, key = ikey_dynamics,
294  CADJ &     kind = isbyte  CADJ &     kind = isbyte
295  #  ifdef SEAICE_CGRID  #  ifdef SEAICE_CGRID
# Line 207  CADJ STORE uice                = comlev1 Line 303  CADJ STORE uice                = comlev1
303  CADJ &     kind = isbyte  CADJ &     kind = isbyte
304  CADJ STORE vice                = comlev1, key = ikey_dynamics,  CADJ STORE vice                = comlev1, key = ikey_dynamics,
305  CADJ &     kind = isbyte  CADJ &     kind = isbyte
306    CADJ STORE dwatn               = comlev1, key = ikey_dynamics,
307    CADJ &     kind = isbyte
308  #   ifdef SEAICE_ALLOW_EVP  #   ifdef SEAICE_ALLOW_EVP
309  CADJ STORE seaice_sigma1       = comlev1, key = ikey_dynamics,  CADJ STORE seaice_sigma1       = comlev1, key = ikey_dynamics,
310  CADJ &     kind = isbyte  CADJ &     kind = isbyte
# Line 216  CADJ STORE seaice_sigma12      = comlev1 Line 314  CADJ STORE seaice_sigma12      = comlev1
314  CADJ &     kind = isbyte  CADJ &     kind = isbyte
315  #   endif  #   endif
316  #  endif  #  endif
317  #  ifdef SEAICE_SALINITY  #  ifdef SEAICE_VARIABLE_SALINITY
318  CADJ STORE salt                = comlev1, key = ikey_dynamics,  CADJ STORE hsalt               = comlev1, key = ikey_dynamics,
319  CADJ &     kind = isbyte  CADJ &     kind = isbyte
320  #  endif  #  endif
321  #  ifdef ATMOSPHERIC_LOADING  #  ifdef ATMOSPHERIC_LOADING
# Line 234  CADJ &     kind = isbyte Line 332  CADJ &     kind = isbyte
332  CADJ STORE balance_itcount     = comlev1, key = ikey_dynamics,  CADJ STORE balance_itcount     = comlev1, key = ikey_dynamics,
333  CADJ &     kind = isbyte  CADJ &     kind = isbyte
334  #  endif /* ANNUAL_BALANCE */  #  endif /* ANNUAL_BALANCE */
335  # endif  #  ifdef ALLOW_THSICE
336    C-- store thSIce vars before advection (called from SEAICE_MODEL) update them:
337    CADJ STORE iceMask,iceHeight = comlev1, key = ikey_dynamics,
338    CADJ &     kind = isbyte
339    CADJ STORE snowHeight, hOceMxL = comlev1, key = ikey_dynamics,
340    CADJ &     kind = isbyte
341    CADJ STORE Qice1, Qice2  = comlev1, key = ikey_dynamics,
342    CADJ &     kind = isbyte
343    #  endif /* ALLOW_THSICE */
344    # endif /* ALLOW_AUTODIFF_TAMC */
345  # ifdef ALLOW_DEBUG  # ifdef ALLOW_DEBUG
346          IF ( debugLevel .GE. debLevB )          IF (debugMode) CALL DEBUG_CALL('SEAICE_MODEL',myThid)
      &    CALL DEBUG_CALL('SEAICE_MODEL',myThid)  
347  # endif  # endif
348          CALL TIMER_START('SEAICE_MODEL    [DO_OCEANIC_PHYS]', myThid)          CALL TIMER_START('SEAICE_MODEL    [DO_OCEANIC_PHYS]', myThid)
349          CALL SEAICE_MODEL( myTime, myIter, myThid )          CALL SEAICE_MODEL( myTime, myIter, myThid )
# Line 259  CADJ &     kind = isbyte Line 365  CADJ &     kind = isbyte
365  # endif  # endif
366  #endif  #endif
367    
368    #ifdef OLD_THSICE_CALL_SEQUENCE
369  #if (defined ALLOW_THSICE) && !(defined ALLOW_ATM2D)  #if (defined ALLOW_THSICE) && !(defined ALLOW_ATM2D)
370        IF ( useThSIce .AND. fluidIsWater ) THEN        IF ( useThSIce .AND. fluidIsWater ) THEN
371  #ifdef ALLOW_DEBUG  # ifdef ALLOW_AUTODIFF_TAMC
372          IF ( debugLevel .GE. debLevB )  cph(
373       &    CALL DEBUG_CALL('THSICE_MAIN',myThid)  #  ifdef NONLIN_FRSURF
374  #endif  CADJ STORE uice,vice        = comlev1, key = ikey_dynamics,
375    CADJ &     kind = isbyte
376    CADJ STORE salt,theta       = comlev1, key = ikey_dynamics,
377    CADJ &     kind = isbyte
378    CADJ STORE qnet,qsw, empmr  = comlev1, key = ikey_dynamics,
379    CADJ &     kind = isbyte
380    CADJ STORE hFac_surfC       = comlev1, key = ikey_dynamics,
381    CADJ &     kind = isbyte
382    #  endif
383    # endif
384    # ifdef ALLOW_DEBUG
385            IF (debugMode) CALL DEBUG_CALL('THSICE_MAIN',myThid)
386    # endif
387  C--     Step forward Therm.Sea-Ice variables  C--     Step forward Therm.Sea-Ice variables
388  C       and modify forcing terms including effects from ice  C       and modify forcing terms including effects from ice
389          CALL TIMER_START('THSICE_MAIN     [DO_OCEANIC_PHYS]', myThid)          CALL TIMER_START('THSICE_MAIN     [DO_OCEANIC_PHYS]', myThid)
# Line 272  C       and modify forcing terms includi Line 391  C       and modify forcing terms includi
391          CALL TIMER_STOP( 'THSICE_MAIN     [DO_OCEANIC_PHYS]', myThid)          CALL TIMER_STOP( 'THSICE_MAIN     [DO_OCEANIC_PHYS]', myThid)
392        ENDIF        ENDIF
393  #endif /* ALLOW_THSICE */  #endif /* ALLOW_THSICE */
394    #endif /* OLD_THSICE_CALL_SEQUENCE */
395    
396  #ifdef ALLOW_SHELFICE  #ifdef ALLOW_SHELFICE
397  # ifdef ALLOW_AUTODIFF_TAMC  # ifdef ALLOW_AUTODIFF_TAMC
# Line 280  CADJ &     kind = isbyte Line 400  CADJ &     kind = isbyte
400  # endif  # endif
401        IF ( useShelfIce .AND. fluidIsWater ) THEN        IF ( useShelfIce .AND. fluidIsWater ) THEN
402  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
403          IF ( debugLevel .GE. debLevB )         IF (debugMode) CALL DEBUG_CALL('SHELFICE_THERMODYNAMICS',myThid)
      &    CALL DEBUG_CALL('SHELFICE_THERMODYNAMICS',myThid)  
404  #endif  #endif
405  C     compute temperature and (virtual) salt flux at the  C     compute temperature and (virtual) salt flux at the
406  C     shelf-ice ocean interface  C     shelf-ice ocean interface
# Line 296  C     shelf-ice ocean interface Line 415  C     shelf-ice ocean interface
415  #ifdef ALLOW_ICEFRONT  #ifdef ALLOW_ICEFRONT
416        IF ( useICEFRONT .AND. fluidIsWater ) THEN        IF ( useICEFRONT .AND. fluidIsWater ) THEN
417  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
418          IF ( debugLevel .GE. debLevB )         IF (debugMode) CALL DEBUG_CALL('ICEFRONT_THERMODYNAMICS',myThid)
      &    CALL DEBUG_CALL('ICEFRONT_THERMODYNAMICS',myThid)  
419  #endif  #endif
420  C     compute temperature and (virtual) salt flux at the  C     compute temperature and (virtual) salt flux at the
421  C     ice-front ocean interface  C     ice-front ocean interface
# Line 309  C     ice-front ocean interface Line 427  C     ice-front ocean interface
427        ENDIF        ENDIF
428  #endif /* ALLOW_ICEFRONT */  #endif /* ALLOW_ICEFRONT */
429    
430    #ifdef ALLOW_SALT_PLUME
431          IF ( useSALT_PLUME ) THEN
432              CALL SALT_PLUME_DO_EXCH( myTime, myIter, myThid )
433          ENDIF
434    #endif /* ALLOW_SALT_PLUME */
435    
436  C--   Freeze water at the surface  C--   Freeze water at the surface
437          IF ( allowFreezing ) THEN
438  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
439  CADJ STORE theta = comlev1, key = ikey_dynamics,  CADJ STORE theta = comlev1, key = ikey_dynamics,
440  CADJ &     kind = isbyte  CADJ &     kind = isbyte
441  #endif  #endif
442        IF ( allowFreezing ) THEN          CALL FREEZE_SURFACE( myTime, myIter, myThid )
         CALL FREEZE_SURFACE(  myTime, myIter, myThid )  
443        ENDIF        ENDIF
444    
445  #ifdef ALLOW_OCN_COMPON_INTERF  #ifdef ALLOW_OCN_COMPON_INTERF
# Line 326  C jmc: do not know precisely where to pu Line 450  C jmc: do not know precisely where to pu
450        ENDIF        ENDIF
451  #endif /* ALLOW_OCN_COMPON_INTERF */  #endif /* ALLOW_OCN_COMPON_INTERF */
452    
453  #ifdef ALLOW_BALANCE_FLUXES        iMin = 1-OLx
454  C     balance fluxes        iMax = sNx+OLx
455        IF ( balanceEmPmR )        jMin = 1-OLy
456       &      CALL REMOVE_MEAN_RS( 1, EmPmR, maskInC, maskInC, rA, drF,        jMax = sNy+OLy
457       &        'EmPmR', myTime, myThid )  
458        IF ( balanceQnet )  C---  Determines forcing terms based on external fields
459       &      CALL REMOVE_MEAN_RS( 1, Qnet,  maskInC, maskInC, rA, drF,  C     relaxation terms, etc.
460       &        'Qnet ', myTime, myThid )  #ifdef ALLOW_DEBUG
461  #endif /* ALLOW_BALANCE_FLUXES */        IF (debugMode) CALL DEBUG_CALL('EXTERNAL_FORCING_SURF',myThid)
462    #endif
463    #ifdef ALLOW_AUTODIFF_TAMC
464    CADJ STORE salt, theta = comlev1, key = ikey_dynamics,
465    CADJ &     kind = isbyte
466    #else  /* ALLOW_AUTODIFF_TAMC */
467    C--   if fluid is not water, by-pass surfaceForcing, find_rho, gmredi
468    C     and all vertical mixing schemes, but keep OBCS_CALC
469          IF ( fluidIsWater ) THEN
470    #endif /* ALLOW_AUTODIFF_TAMC */
471            CALL EXTERNAL_FORCING_SURF(
472         I             iMin, iMax, jMin, jMax,
473         I             myTime, myIter, myThid )
474    
475  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
476  C--   HPF directive to help TAMC  C--   HPF directive to help TAMC
# Line 358  CHPF$ INDEPENDENT Line 494  CHPF$ INDEPENDENT
494            itdkey = (act1 + 1) + act2*max1            itdkey = (act1 + 1) + act2*max1
495       &                      + act3*max1*max2       &                      + act3*max1*max2
496       &                      + act4*max1*max2*max3       &                      + act4*max1*max2*max3
 #else  /* ALLOW_AUTODIFF_TAMC */  
 C     if fluid is not water, by-pass find_rho, gmredi, surfaceForcing  
 C     and all vertical mixing schemes, but keep OBCS_CALC  
         IF ( fluidIsWater ) THEN  
497  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
498    
499  C--   Set up work arrays with valid (i.e. not NaN) values  C--   Set up work arrays with valid (i.e. not NaN) values
# Line 369  C     These inital values do not alter t Line 501  C     These inital values do not alter t
501  C     just ensure that all memory references are to valid floating  C     just ensure that all memory references are to valid floating
502  C     point numbers. This prevents spurious hardware signals due to  C     point numbers. This prevents spurious hardware signals due to
503  C     uninitialised but inert locations.  C     uninitialised but inert locations.
504            DO k=1,Nr
505             DO j=1-OLy,sNy+OLy
506              DO i=1-OLx,sNx+OLx
507    C This is currently used by GMRedi, IVDC, MXL-depth  and Diagnostics
508               sigmaX(i,j,k) = 0. _d 0
509               sigmaY(i,j,k) = 0. _d 0
510               sigmaR(i,j,k) = 0. _d 0
511              ENDDO
512             ENDDO
513            ENDDO
514    
515  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
516          DO j=1-OLy,sNy+OLy          DO j=1-OLy,sNy+OLy
# Line 377  C     uninitialised but inert locations. Line 519  C     uninitialised but inert locations.
519            rhoKp1 (i,j)   = 0. _d 0            rhoKp1 (i,j)   = 0. _d 0
520           ENDDO           ENDDO
521          ENDDO          ENDDO
522  #endif /* ALLOW_AUTODIFF_TAMC */  cph all the following init. are necessary for TAF
523    cph although some of these are re-initialised later.
524            DO k=1,Nr
525             DO j=1-OLy,sNy+OLy
526              DO i=1-OLx,sNx+OLx
527               rhoInSitu(i,j,k,bi,bj) = 0.
528    # ifdef ALLOW_GGL90
529               GGL90viscArU(i,j,k,bi,bj)  = 0. _d 0
530               GGL90viscArV(i,j,k,bi,bj)  = 0. _d 0
531               GGL90diffKr(i,j,k,bi,bj)  = 0. _d 0
532    # endif /* ALLOW_GGL90 */
533              ENDDO
534             ENDDO
535            ENDDO
536    #ifdef ALLOW_OFFLINE
537           IF ( calcConvect ) THEN
538    #endif
539          DO k=1,Nr          DO k=1,Nr
540           DO j=1-OLy,sNy+OLy           DO j=1-OLy,sNy+OLy
541            DO i=1-OLx,sNx+OLx            DO i=1-OLx,sNx+OLx
 C This is currently used by GMRedi, IVDC, MXL-depth  and Diagnostics  
            sigmaX(i,j,k) = 0. _d 0  
            sigmaY(i,j,k) = 0. _d 0  
            sigmaR(i,j,k) = 0. _d 0  
 #ifdef ALLOW_AUTODIFF_TAMC  
 cph all the following init. are necessary for TAF  
 cph although some of these are re-initialised later.  
 c          rhoInSitu(i,j,k,bi,bj) = 0.  
542             IVDConvCount(i,j,k,bi,bj) = 0.             IVDConvCount(i,j,k,bi,bj) = 0.
543              ENDDO
544             ENDDO
545            ENDDO
546    #ifdef ALLOW_OFFLINE
547           ENDIF
548           IF ( calcGMRedi ) THEN
549    #endif
550  # ifdef ALLOW_GMREDI  # ifdef ALLOW_GMREDI
551            DO k=1,Nr
552             DO j=1-OLy,sNy+OLy
553              DO i=1-OLx,sNx+OLx
554             Kwx(i,j,k,bi,bj)  = 0. _d 0             Kwx(i,j,k,bi,bj)  = 0. _d 0
555             Kwy(i,j,k,bi,bj)  = 0. _d 0             Kwy(i,j,k,bi,bj)  = 0. _d 0
556             Kwz(i,j,k,bi,bj)  = 0. _d 0             Kwz(i,j,k,bi,bj)  = 0. _d 0
# Line 410  c          rhoInSitu(i,j,k,bi,bj) = 0. Line 569  c          rhoInSitu(i,j,k,bi,bj) = 0.
569  #  ifdef GM_VISBECK_VARIABLE_K  #  ifdef GM_VISBECK_VARIABLE_K
570             VisbeckK(i,j,bi,bj)   = 0. _d 0             VisbeckK(i,j,bi,bj)   = 0. _d 0
571  #  endif  #  endif
572              ENDDO
573             ENDDO
574            ENDDO
575  # endif /* ALLOW_GMREDI */  # endif /* ALLOW_GMREDI */
576    #ifdef ALLOW_OFFLINE
577           ENDIF
578           IF ( calcKPP ) THEN
579    #endif
580  # ifdef ALLOW_KPP  # ifdef ALLOW_KPP
581            DO k=1,Nr
582             DO j=1-OLy,sNy+OLy
583              DO i=1-OLx,sNx+OLx
584             KPPdiffKzS(i,j,k,bi,bj)  = 0. _d 0             KPPdiffKzS(i,j,k,bi,bj)  = 0. _d 0
585             KPPdiffKzT(i,j,k,bi,bj)  = 0. _d 0             KPPdiffKzT(i,j,k,bi,bj)  = 0. _d 0
 # endif /* ALLOW_KPP */  
 # ifdef ALLOW_GGL90  
            GGL90viscArU(i,j,k,bi,bj)  = 0. _d 0  
            GGL90viscArV(i,j,k,bi,bj)  = 0. _d 0  
            GGL90diffKr(i,j,k,bi,bj)  = 0. _d 0  
 # endif /* ALLOW_GGL90 */  
 #endif /* ALLOW_AUTODIFF_TAMC */  
586            ENDDO            ENDDO
587           ENDDO           ENDDO
588          ENDDO          ENDDO
589    # endif /* ALLOW_KPP */
590          iMin = 1-OLx  #ifdef ALLOW_OFFLINE
591          iMax = sNx+OLx         ENDIF
592          jMin = 1-OLy  #endif
593          jMax = sNy+OLy  #endif /* ALLOW_AUTODIFF_TAMC */
594    
595  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
596  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey,  CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
597  CADJ &     kind = isbyte  CADJ &     kind = isbyte
598  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,  CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
599  CADJ &     kind = isbyte  CADJ &     kind = isbyte
600  CADJ STORE totphihyd(:,:,:,bi,bj)  CADJ STORE totphihyd(:,:,:,bi,bj)
601  CADJ &     = comlev1_bibj, key=itdkey,  CADJ &     = comlev1_bibj, key=itdkey,
602  CADJ &     kind = isbyte  CADJ &     kind = isbyte
603  # ifdef ALLOW_KPP  # ifdef ALLOW_KPP
604  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,  CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
605  CADJ &     kind = isbyte  CADJ &     kind = isbyte
606  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,  CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
607    CADJ &     kind = isbyte
608    # endif
609    # ifdef ALLOW_SALT_PLUME
610    CADJ STORE saltplumedepth(:,:,bi,bj) = comlev1_bibj, key=itdkey,
611  CADJ &     kind = isbyte  CADJ &     kind = isbyte
612  # endif  # endif
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
 #ifdef ALLOW_DEBUG  
         IF ( debugLevel .GE. debLevB )  
      &    CALL DEBUG_MSG('ENTERING UPWARD K LOOP',myThid)  
 #endif  
   
 C--     Start of diagnostic loop  
         DO k=Nr,1,-1  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 C? Patrick, is this formula correct now that we change the loop range?  
 C? Do we still need this?  
 cph kkey formula corrected.  
 cph Needed for rhoK, rhoKm1, in the case useGMREDI.  
           kkey = (itdkey-1)*Nr + k  
613  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
614    
615  C--   Always compute density (stored in common block) here; even when it is not  C--   Always compute density (stored in common block) here; even when it is not
616  C     needed here, will be used anyway in calc_phi_hyd (data flow easier this way)  C     needed here, will be used anyway in calc_phi_hyd (data flow easier this way)
617  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
618            IF ( debugLevel .GE. debLevB )          IF (debugMode) CALL DEBUG_CALL('FIND_RHO_2D (xNr)',myThid)
      &       CALL DEBUG_CALL('FIND_RHO_2D',myThid)  
619  #endif  #endif
620  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
621            IF ( fluidIsWater ) THEN          IF ( fluidIsWater ) THEN
 CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey,  
 CADJ &     kind = isbyte  
 CADJ STORE salt(:,:,k,bi,bj)  = comlev1_bibj_k, key=kkey,  
 CADJ &     kind = isbyte  
622  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
623  #ifdef ALLOW_DOWN_SLOPE  #ifdef ALLOW_DOWN_SLOPE
624            IF ( useDOWN_SLOPE ) THEN           IF ( useDOWN_SLOPE ) THEN
625               DO k=1,Nr
626              CALL DWNSLP_CALC_RHO(              CALL DWNSLP_CALC_RHO(
627       I                  theta, salt,       I                  theta, salt,
628       O                  rhoInSitu(1-OLx,1-OLy,k,bi,bj),       O                  rhoInSitu(1-OLx,1-OLy,k,bi,bj),
629       I                  k, bi, bj, myTime, myIter, myThid )       I                  k, bi, bj, myTime, myIter, myThid )
630            ELSE             ENDDO
631             ENDIF
632  #endif /* ALLOW_DOWN_SLOPE */  #endif /* ALLOW_DOWN_SLOPE */
633    #ifdef ALLOW_BBL
634             IF ( useBBL ) THEN
635    C     pkg/bbl requires in-situ bbl density for depths equal to and deeper than the bbl.
636    C     To reduce computation and storage requirement, these densities are stored in the
637    C     dry grid boxes of rhoInSitu.  See BBL_CALC_RHO for details.
638               DO k=Nr,1,-1
639                CALL BBL_CALC_RHO(
640         I                  theta, salt,
641         O                  rhoInSitu,
642         I                  k, bi, bj, myTime, myIter, myThid )
643    
644               ENDDO
645             ENDIF
646    #endif /* ALLOW_BBL */
647             IF ( .NOT. ( useDOWN_SLOPE .OR. useBBL ) ) THEN
648               DO k=1,Nr
649              CALL FIND_RHO_2D(              CALL FIND_RHO_2D(
650       I                iMin, iMax, jMin, jMax, k,       I                iMin, iMax, jMin, jMax, k,
651       I                theta(1-OLx,1-OLy,k,bi,bj),       I                theta(1-OLx,1-OLy,k,bi,bj),
652       I                salt (1-OLx,1-OLy,k,bi,bj),       I                salt (1-OLx,1-OLy,k,bi,bj),
653       O                rhoInSitu(1-OLx,1-OLy,k,bi,bj),       O                rhoInSitu(1-OLx,1-OLy,k,bi,bj),
654       I                k, bi, bj, myThid )       I                k, bi, bj, myThid )
655  #ifdef ALLOW_DOWN_SLOPE             ENDDO
656            ENDIF           ENDIF
 #endif /* ALLOW_DOWN_SLOPE */  
657  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
658            ELSE          ELSE
659  C-        fluid is not water:  C-        fluid is not water:
660              DO k=1,Nr
661             DO j=1-OLy,sNy+OLy             DO j=1-OLy,sNy+OLy
662              DO i=1-OLx,sNx+OLx              DO i=1-OLx,sNx+OLx
663                rhoInSitu(i,j,k,bi,bj) = 0.                rhoInSitu(i,j,k,bi,bj) = 0.
664              ENDDO              ENDDO
665             ENDDO             ENDDO
666            ENDIF            ENDDO
667            ENDIF
668    #endif /* ALLOW_AUTODIFF_TAMC */
669    
670    #ifdef ALLOW_DEBUG
671            IF (debugMode) CALL DEBUG_MSG('ENTERING UPWARD K LOOP',myThid)
672    #endif
673    
674    C--     Start of diagnostic loop
675            DO k=Nr,1,-1
676    
677    #ifdef ALLOW_AUTODIFF_TAMC
678    C? Patrick, is this formula correct now that we change the loop range?
679    C? Do we still need this?
680    cph kkey formula corrected.
681    cph Needed for rhoK, rhoKm1, in the case useGMREDI.
682              kkey = (itdkey-1)*Nr + k
683  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
684    
685    c#ifdef ALLOW_AUTODIFF_TAMC
686    cCADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey,
687    cCADJ &     kind = isbyte
688    cCADJ STORE salt(:,:,k,bi,bj)  = comlev1_bibj_k, key=kkey,
689    cCADJ &     kind = isbyte
690    c#endif /* ALLOW_AUTODIFF_TAMC */
691    
692  C--       Calculate gradients of potential density for isoneutral  C--       Calculate gradients of potential density for isoneutral
693  C         slope terms (e.g. GM/Redi tensor or IVDC diffusivity)  C         slope terms (e.g. GM/Redi tensor or IVDC diffusivity)
694            IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.)            IF ( calcGMRedi .OR. (k.GT.1 .AND. calcConvect)
695         &         .OR. usePP81 .OR. useMY82 .OR. useGGL90
696       &         .OR. useSALT_PLUME .OR. doDiagsRho.GE.1 ) THEN       &         .OR. useSALT_PLUME .OR. doDiagsRho.GE.1 ) THEN
697              IF (k.GT.1) THEN              IF (k.GT.1) THEN
698  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
699  CADJ STORE theta(:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey,  CADJ STORE theta(:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey,
700  CADJ &     kind = isbyte  CADJ &     kind = isbyte
701  CADJ STORE salt (:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey,  CADJ STORE salt (:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey,
702  CADJ &     kind = isbyte  CADJ &     kind = isbyte
703  CADJ STORE rhokm1 (bi,bj)       = comlev1_bibj_k, key=kkey,  CADJ STORE rhokm1 (bi,bj)       = comlev1_bibj_k, key=kkey,
704  CADJ &     kind = isbyte  CADJ &     kind = isbyte
705  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
706               CALL FIND_RHO_2D(               CALL FIND_RHO_2D(
# Line 524  CADJ &     kind = isbyte Line 711  CADJ &     kind = isbyte
711       I                 k-1, bi, bj, myThid )       I                 k-1, bi, bj, myThid )
712              ENDIF              ENDIF
713  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
714              IF ( debugLevel .GE. debLevB )              IF (debugMode) CALL DEBUG_CALL('GRAD_SIGMA',myThid)
      &       CALL DEBUG_CALL('GRAD_SIGMA',myThid)  
715  #endif  #endif
716  cph Avoid variable aliasing for adjoint !!!  cph Avoid variable aliasing for adjoint !!!
717              DO j=jMin,jMax              DO j=jMin,jMax
# Line 550  cgf -> cuts adjoint dependency from slop Line 736  cgf -> cuts adjoint dependency from slop
736            ENDIF            ENDIF
737    
738  C--       Implicit Vertical Diffusion for Convection  C--       Implicit Vertical Diffusion for Convection
739  c ==> should use sigmaR !!!            IF (k.GT.1 .AND. calcConvect) THEN
           IF (k.GT.1 .AND. ivdc_kappa.NE.0.) THEN  
740  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
741              IF ( debugLevel .GE. debLevB )              IF (debugMode) CALL DEBUG_CALL('CALC_IVDC',myThid)
      &       CALL DEBUG_CALL('CALC_IVDC',myThid)  
742  #endif  #endif
743              CALL CALC_IVDC(              CALL CALC_IVDC(
744       I        bi, bj, iMin, iMax, jMin, jMax, k,       I        bi, bj, iMin, iMax, jMin, jMax, k,
745       I        rhoKm1, rhoInSitu(1-OLx,1-OLy,k,bi,bj),       I        sigmaR,
746       I        myTime, myIter, myThid)       I        myTime, myIter, myThid)
747            ENDIF            ENDIF
748    
749  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
750            IF ( MOD(doDiagsRho,2).EQ.1 ) THEN            IF ( doDiagsRho.GE.4 ) THEN
751              CALL DIAGS_RHO_L( k, bi, bj,              CALL DIAGS_RHO_L( doDiagsRho, k, bi, bj,
752       I                        rhoInSitu(1-OLx,1-OLy,k,bi,bj),       I                        rhoInSitu(1-OLx,1-OLy,1,bi,bj),
753       I                        rhoKm1, wVel,       I                        rhoKm1, wVel,
754       I                        myTime, myIter, myThid )       I                        myTime, myIter, myThid )
755            ENDIF            ENDIF
# Line 576  C--     end of diagnostic k loop (Nr:1) Line 760  C--     end of diagnostic k loop (Nr:1)
760    
761  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
762  CADJ STORE IVDConvCount(:,:,:,bi,bj)  CADJ STORE IVDConvCount(:,:,:,bi,bj)
763  CADJ &     = comlev1_bibj, key=itdkey,  CADJ &     = comlev1_bibj, key=itdkey,
764  CADJ &     kind = isbyte  CADJ &     kind = isbyte
765  #endif  #endif
766    
767  C--     Diagnose Mixed Layer Depth:  C--     Diagnose Mixed Layer Depth:
768          IF ( useGMRedi .OR. doDiagsRho.GE.4 ) THEN          IF ( calcGMRedi .OR. MOD(doDiagsRho,2).EQ.1 ) THEN
769            CALL CALC_OCE_MXLAYER(            CALL CALC_OCE_MXLAYER(
770       I              rhoInSitu(1-OLx,1-OLy,1,bi,bj), sigmaR,       I              rhoInSitu(1-OLx,1-OLy,1,bi,bj), sigmaR,
771       I              bi, bj, myTime, myIter, myThid )       I              bi, bj, myTime, myIter, myThid )
# Line 602  C--     Diagnose Mixed Layer Depth: Line 786  C--     Diagnose Mixed Layer Depth:
786          ENDIF          ENDIF
787  #endif /* ALLOW_DIAGNOSTICS */  #endif /* ALLOW_DIAGNOSTICS */
788    
789  C--     Determines forcing terms based on external fields  C--    This is where EXTERNAL_FORCING_SURF(bi,bj) used to be called;
790  C       relaxation terms, etc.  C      now called earlier, before bi,bj loop.
 #ifdef ALLOW_DEBUG  
         IF ( debugLevel .GE. debLevB )  
      &    CALL DEBUG_CALL('EXTERNAL_FORCING_SURF',myThid)  
 #endif  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE EmPmR(:,:,bi,bj)  
 CADJ &     = comlev1_bibj, key=itdkey,  
 CADJ &     kind = isbyte  
 # ifdef EXACT_CONSERV  
 CADJ STORE PmEpR(:,:,bi,bj)  
 CADJ &     = comlev1_bibj, key=itdkey,  
 CADJ &     kind = isbyte  
 # endif  
 # ifdef NONLIN_FRSURF  
 CADJ STORE hFac_surfC(:,:,bi,bj)  
 CADJ &     = comlev1_bibj, key=itdkey,  
 CADJ &     kind = isbyte  
 CADJ STORE recip_hFacC(:,:,:,bi,bj)  
 CADJ &     = comlev1_bibj, key=itdkey,  
 CADJ &     kind = isbyte  
 # endif  
 #endif  
         CALL EXTERNAL_FORCING_SURF(  
      I             bi, bj, iMin, iMax, jMin, jMax,  
      I             myTime, myIter, myThid )  
 #ifdef ALLOW_AUTODIFF_TAMC  
 # ifdef EXACT_CONSERV  
 cph-test  
 cphCADJ STORE PmEpR(:,:,bi,bj)  
 cphCADJ &     = comlev1_bibj, key=itdkey,  
 cphCADJ &     kind = isbyte  
 # endif  
 #endif  
791    
792  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
793  cph needed for KPP  cph needed for KPP
794  CADJ STORE surfaceForcingU(:,:,bi,bj)  CADJ STORE surfaceForcingU(:,:,bi,bj)
795  CADJ &     = comlev1_bibj, key=itdkey,  CADJ &     = comlev1_bibj, key=itdkey,
796  CADJ &     kind = isbyte  CADJ &     kind = isbyte
797  CADJ STORE surfaceForcingV(:,:,bi,bj)  CADJ STORE surfaceForcingV(:,:,bi,bj)
798  CADJ &     = comlev1_bibj, key=itdkey,  CADJ &     = comlev1_bibj, key=itdkey,
799  CADJ &     kind = isbyte  CADJ &     kind = isbyte
800  CADJ STORE surfaceForcingS(:,:,bi,bj)  CADJ STORE surfaceForcingS(:,:,bi,bj)
801  CADJ &     = comlev1_bibj, key=itdkey,  CADJ &     = comlev1_bibj, key=itdkey,
802  CADJ &     kind = isbyte  CADJ &     kind = isbyte
803  CADJ STORE surfaceForcingT(:,:,bi,bj)  CADJ STORE surfaceForcingT(:,:,bi,bj)
804  CADJ &     = comlev1_bibj, key=itdkey,  CADJ &     = comlev1_bibj, key=itdkey,
805  CADJ &     kind = isbyte  CADJ &     kind = isbyte
806  CADJ STORE surfaceForcingTice(:,:,bi,bj)  CADJ STORE surfaceForcingTice(:,:,bi,bj)
807  CADJ &     = comlev1_bibj, key=itdkey,  CADJ &     = comlev1_bibj, key=itdkey,
808  CADJ &     kind = isbyte  CADJ &     kind = isbyte
809  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
810    
811  #ifdef  ALLOW_KPP  #ifdef  ALLOW_KPP
812  C--     Compute KPP mixing coefficients  C--     Compute KPP mixing coefficients
813          IF (useKPP) THEN          IF ( calcKPP ) THEN
814  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
815            IF ( debugLevel .GE. debLevB )            IF (debugMode) CALL DEBUG_CALL('KPP_CALC',myThid)
      &     CALL DEBUG_CALL('KPP_CALC',myThid)  
816  #endif  #endif
817            CALL TIMER_START('KPP_CALC [DO_OCEANIC_PHYS]', myThid)            CALL TIMER_START('KPP_CALC [DO_OCEANIC_PHYS]', myThid)
818            CALL KPP_CALC(            CALL KPP_CALC(
819       I                  bi, bj, myTime, myIter, myThid )       I                  bi, bj, myTime, myIter, myThid )
820            CALL TIMER_STOP ('KPP_CALC [DO_OCEANIC_PHYS]', myThid)            CALL TIMER_STOP ('KPP_CALC [DO_OCEANIC_PHYS]', myThid)
821  #ifdef ALLOW_AUTODIFF_TAMC  #if (defined ALLOW_AUTODIFF_TAMC) && !(defined ALLOW_OFFLINE)
822          ELSE          ELSE
823            CALL KPP_CALC_DUMMY(            CALL KPP_CALC_DUMMY(
824       I                  bi, bj, myTime, myIter, myThid )       I                  bi, bj, myTime, myIter, myThid )
825  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC and not ALLOW_OFFLINE */
826          ENDIF          ENDIF
   
827  #endif  /* ALLOW_KPP */  #endif  /* ALLOW_KPP */
828    
829  #ifdef  ALLOW_PP81  #ifdef  ALLOW_PP81
830  C--     Compute PP81 mixing coefficients  C--     Compute PP81 mixing coefficients
831          IF (usePP81) THEN          IF (usePP81) THEN
832  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
833            IF ( debugLevel .GE. debLevB )            IF (debugMode) CALL DEBUG_CALL('PP81_CALC',myThid)
      &     CALL DEBUG_CALL('PP81_CALC',myThid)  
834  #endif  #endif
835            CALL PP81_CALC(            CALL PP81_CALC(
836       I                  bi, bj, myTime, myThid )       I                     bi, bj, sigmaR, myTime, myIter, myThid )
837          ENDIF          ENDIF
838  #endif /* ALLOW_PP81 */  #endif /* ALLOW_PP81 */
839    
# Line 693  C--     Compute PP81 mixing coefficients Line 841  C--     Compute PP81 mixing coefficients
841  C--     Compute MY82 mixing coefficients  C--     Compute MY82 mixing coefficients
842          IF (useMY82) THEN          IF (useMY82) THEN
843  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
844            IF ( debugLevel .GE. debLevB )            IF (debugMode) CALL DEBUG_CALL('MY82_CALC',myThid)
      &     CALL DEBUG_CALL('MY82_CALC',myThid)  
845  #endif  #endif
846            CALL MY82_CALC(            CALL MY82_CALC(
847       I                  bi, bj, myTime, myThid )       I                     bi, bj, sigmaR, myTime, myIter, myThid )
848          ENDIF          ENDIF
849  #endif /* ALLOW_MY82 */  #endif /* ALLOW_MY82 */
850    
# Line 709  CADJ &     kind = isbyte Line 856  CADJ &     kind = isbyte
856  C--     Compute GGL90 mixing coefficients  C--     Compute GGL90 mixing coefficients
857          IF (useGGL90) THEN          IF (useGGL90) THEN
858  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
859            IF ( debugLevel .GE. debLevB )            IF (debugMode) CALL DEBUG_CALL('GGL90_CALC',myThid)
      &     CALL DEBUG_CALL('GGL90_CALC',myThid)  
860  #endif  #endif
861            CALL TIMER_START('GGL90_CALC [DO_OCEANIC_PHYS]', myThid)            CALL TIMER_START('GGL90_CALC [DO_OCEANIC_PHYS]', myThid)
862            CALL GGL90_CALC(            CALL GGL90_CALC(
863       I                  bi, bj, myTime, myIter, myThid )       I                     bi, bj, sigmaR, myTime, myIter, myThid )
864            CALL TIMER_STOP ('GGL90_CALC [DO_OCEANIC_PHYS]', myThid)            CALL TIMER_STOP ('GGL90_CALC [DO_OCEANIC_PHYS]', myThid)
865          ENDIF          ENDIF
866  #endif /* ALLOW_GGL90 */  #endif /* ALLOW_GGL90 */
# Line 723  C--     Compute GGL90 mixing coefficient Line 869  C--     Compute GGL90 mixing coefficient
869          IF ( taveFreq.GT. 0. _d 0 ) THEN          IF ( taveFreq.GT. 0. _d 0 ) THEN
870            CALL TIMEAVE_SURF_FLUX( bi, bj, myTime, myIter, myThid)            CALL TIMEAVE_SURF_FLUX( bi, bj, myTime, myIter, myThid)
871          ENDIF          ENDIF
872          IF (taveFreq.GT.0. .AND. ivdc_kappa.NE.0.) THEN          IF ( taveFreq.GT.0. .AND. calcConvect ) THEN
873            CALL TIMEAVE_CUMULATE(ConvectCountTave, IVDConvCount,            CALL TIMEAVE_CUMULATE(ConvectCountTave, IVDConvCount,
874       I                           Nr, deltaTclock, bi, bj, myThid)       I                           Nr, deltaTClock, bi, bj, myThid)
875          ENDIF          ENDIF
876  #endif /* ALLOW_TIMEAVE */  #endif /* ALLOW_TIMEAVE */
877    
# Line 736  cph storing here is needed only for one Line 882  cph storing here is needed only for one
882  cph define GM_BOLUS_ADVEC  cph define GM_BOLUS_ADVEC
883  cph keep it although TAF says you dont need to.  cph keep it although TAF says you dont need to.
884  cph but I have avoided the #ifdef for now, in case more things change  cph but I have avoided the #ifdef for now, in case more things change
885  CADJ STORE sigmaX(:,:,:)        = comlev1_bibj, key=itdkey,  CADJ STORE sigmaX(:,:,:)        = comlev1_bibj, key=itdkey,
886  CADJ &     kind = isbyte  CADJ &     kind = isbyte
887  CADJ STORE sigmaY(:,:,:)        = comlev1_bibj, key=itdkey,  CADJ STORE sigmaY(:,:,:)        = comlev1_bibj, key=itdkey,
888  CADJ &     kind = isbyte  CADJ &     kind = isbyte
889  CADJ STORE sigmaR(:,:,:)        = comlev1_bibj, key=itdkey,  CADJ STORE sigmaR(:,:,:)        = comlev1_bibj, key=itdkey,
890  CADJ &     kind = isbyte  CADJ &     kind = isbyte
891  # endif  # endif
892  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
893    
894  C--     Calculate iso-neutral slopes for the GM/Redi parameterisation  C--     Calculate iso-neutral slopes for the GM/Redi parameterisation
895          IF (useGMRedi) THEN          IF ( calcGMRedi ) THEN
896  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
897            IF ( debugLevel .GE. debLevB )            IF (debugMode) CALL DEBUG_CALL('GMREDI_CALC_TENSOR',myThid)
      &     CALL DEBUG_CALL('GMREDI_CALC_TENSOR',myThid)  
898  #endif  #endif
899            CALL GMREDI_CALC_TENSOR(            CALL GMREDI_CALC_TENSOR(
900       I             iMin, iMax, jMin, jMax,       I             iMin, iMax, jMin, jMax,
901       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
902       I             bi, bj, myTime, myIter, myThid )       I             bi, bj, myTime, myIter, myThid )
903  #ifdef ALLOW_AUTODIFF_TAMC  #if (defined ALLOW_AUTODIFF_TAMC) && !(defined ALLOW_OFFLINE)
904          ELSE          ELSE
905            CALL GMREDI_CALC_TENSOR_DUMMY(            CALL GMREDI_CALC_TENSOR_DUMMY(
906       I             iMin, iMax, jMin, jMax,       I             iMin, iMax, jMin, jMax,
907       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
908       I             bi, bj, myTime, myIter, myThid )       I             bi, bj, myTime, myIter, myThid )
909  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC and not ALLOW_OFFLINE */
910          ENDIF          ENDIF
911  #endif /* ALLOW_GMREDI */  #endif /* ALLOW_GMREDI */
912    
# Line 780  C--     Calculate Downsloping Flow for D Line 925  C--     Calculate Downsloping Flow for D
925          ENDIF          ENDIF
926  #endif /* ALLOW_DOWN_SLOPE */  #endif /* ALLOW_DOWN_SLOPE */
927    
928  #ifdef ALLOW_MYPACKAGE  C--   end bi,bj loops.
929          IF ( useMYPACKAGE ) THEN         ENDDO
930           CALL MYPACKAGE_CALC_RHS(        ENDDO
      I        bi, bj, myTime, myIter, myThid )  
         ENDIF  
 #endif /* ALLOW_MYPACKAGE */  
931    
932  #ifndef ALLOW_AUTODIFF_TAMC  #ifndef ALLOW_AUTODIFF_TAMC
933  C---  if fluid Is Water: end  C---  if fluid Is Water: end
934          ENDIF        ENDIF
935  #endif  #endif
936    
937  C--   end bi,bj loops.  #ifdef ALLOW_BBL
938         ENDDO        IF ( useBBL ) THEN
939        ENDDO         CALL BBL_CALC_RHS(
940         I                          myTime, myIter, myThid )
941          ENDIF
942    #endif /* ALLOW_BBL */
943    
944  #ifdef  ALLOW_KPP  #ifdef ALLOW_MYPACKAGE
945        IF (useKPP) THEN        IF ( useMYPACKAGE ) THEN
946           CALL MYPACKAGE_CALC_RHS(
947         I                          myTime, myIter, myThid )
948          ENDIF
949    #endif /* ALLOW_MYPACKAGE */
950    
951    #ifdef ALLOW_GMREDI
952          IF ( calcGMRedi ) THEN
953            CALL GMREDI_DO_EXCH( myTime, myIter, myThid )
954          ENDIF
955    #endif /* ALLOW_GMREDI */
956    
957    #ifdef ALLOW_KPP
958          IF ( calcKPP ) THEN
959          CALL KPP_DO_EXCH( myThid )          CALL KPP_DO_EXCH( myThid )
960        ENDIF        ENDIF
961  #endif  /* ALLOW_KPP */  #endif /* ALLOW_KPP */
962    
963  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
964        IF ( fluidIsWater .AND. useDiagnostics ) THEN        IF ( fluidIsWater .AND. useDiagnostics ) THEN
965          CALL DIAGS_RHO_G(          CALL DIAGS_RHO_G(
966       I                    rhoInSitu, uVel, vVel,       I                    rhoInSitu, uVel, vVel, wVel,
967       I                    myTime, myIter, myThid )       I                    myTime, myIter, myThid )
968          CALL DIAGS_OCEANIC_SURF_FLUX( myTime, myIter, myThid )          CALL DIAGS_OCEANIC_SURF_FLUX( myTime, myIter, myThid )
969        ENDIF        ENDIF
970        IF ( ivdc_kappa.NE.0 .AND. useDiagnostics ) THEN        IF ( calcConvect .AND. useDiagnostics ) THEN
971          CALL DIAGNOSTICS_FILL( IVDConvCount, 'CONVADJ ',          CALL DIAGNOSTICS_FILL( IVDConvCount, 'CONVADJ ',
972       &                               0, Nr, 0, 1, 1, myThid )       &                               0, Nr, 0, 1, 1, myThid )
973        ENDIF        ENDIF
974  #endif  #endif
975    
976    #ifdef ALLOW_ECCO
977          CALL ECCO_PHYS( myThid )
978    #endif
979    
980  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
981        IF ( debugLevel .GE. debLevB )        IF (debugMode) CALL DEBUG_LEAVE('DO_OCEANIC_PHYS',myThid)
      &     CALL DEBUG_LEAVE('DO_OCEANIC_PHYS',myThid)  
982  #endif  #endif
983    
984        RETURN        RETURN
Legend:
Removed from v.1.95  
changed lines
  Added in v.1.132
  ViewVC Help
Powered by ViewVC 1.1.22