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

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

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


Revision 1.133 - (hide annotations) (download)
Tue Sep 24 01:42:29 2013 UTC (10 years, 7 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint64q, checkpoint64p, checkpoint64s, checkpoint64r, checkpoint64u, checkpoint64t, checkpoint64o
Changes since 1.132: +3 -1 lines
also fill-up some (selected) surface diags if fluidIsAir

1 jmc 1.133 C $Header: /u/gcmpack/MITgcm/model/src/do_oceanic_phys.F,v 1.132 2013/07/09 16:01:01 jmc Exp $
2 jmc 1.1 C $Name: $
3    
4     #include "PACKAGES_CONFIG.h"
5     #include "CPP_OPTIONS.h"
6    
7     #ifdef ALLOW_AUTODIFF_TAMC
8     # ifdef ALLOW_GMREDI
9     # include "GMREDI_OPTIONS.h"
10     # endif
11     # ifdef ALLOW_KPP
12     # include "KPP_OPTIONS.h"
13     # endif
14 jmc 1.29 # ifdef ALLOW_SEAICE
15     # include "SEAICE_OPTIONS.h"
16     # endif
17 gforget 1.117 # ifdef ALLOW_EXF
18     # include "EXF_OPTIONS.h"
19     # endif
20 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
21    
22     CBOP
23     C !ROUTINE: DO_OCEANIC_PHYS
24     C !INTERFACE:
25     SUBROUTINE DO_OCEANIC_PHYS(myTime, myIter, myThid)
26     C !DESCRIPTION: \bv
27     C *==========================================================*
28 jmc 1.28 C | SUBROUTINE DO_OCEANIC_PHYS
29     C | o Controlling routine for oceanic physics and
30 jmc 1.1 C | parameterization
31     C *==========================================================*
32     C | o originally, part of S/R thermodynamics
33     C *==========================================================*
34     C \ev
35    
36     C !USES:
37     IMPLICIT NONE
38     C == Global variables ===
39     #include "SIZE.h"
40     #include "EEPARAMS.h"
41     #include "PARAMS.h"
42 jmc 1.69 #include "GRID.h"
43 jmc 1.1 #include "DYNVARS.h"
44 jmc 1.20 #ifdef ALLOW_TIMEAVE
45 jmc 1.131 # include "TIMEAVE_STATV.h"
46     #endif
47     #ifdef ALLOW_OFFLINE
48     # include "OFFLINE_SWITCH.h"
49 jmc 1.20 #endif
50 jmc 1.1
51     #ifdef ALLOW_AUTODIFF_TAMC
52 heimbach 1.80 # include "AUTODIFF_MYFIELDS.h"
53 jmc 1.1 # include "tamc.h"
54     # include "tamc_keys.h"
55     # include "FFIELDS.h"
56 heimbach 1.54 # include "SURFACE.h"
57 jmc 1.1 # include "EOS.h"
58     # ifdef ALLOW_KPP
59     # include "KPP.h"
60     # endif
61 gforget 1.91 # ifdef ALLOW_GGL90
62     # include "GGL90.h"
63     # endif
64 jmc 1.1 # ifdef ALLOW_GMREDI
65     # include "GMREDI.h"
66     # endif
67     # ifdef ALLOW_EBM
68     # include "EBM.h"
69     # endif
70 jmc 1.29 # ifdef ALLOW_EXF
71     # include "ctrl.h"
72 jmc 1.40 # include "EXF_FIELDS.h"
73 jmc 1.29 # ifdef ALLOW_BULKFORMULAE
74 jmc 1.40 # include "EXF_CONSTANTS.h"
75 jmc 1.29 # endif
76     # endif
77     # ifdef ALLOW_SEAICE
78 jmc 1.114 # include "SEAICE_SIZE.h"
79 jmc 1.29 # include "SEAICE.h"
80 gforget 1.124 # include "SEAICE_PARAMS.h"
81 jmc 1.29 # endif
82 heimbach 1.105 # ifdef ALLOW_THSICE
83     # include "THSICE_VARS.h"
84     # endif
85 heimbach 1.75 # ifdef ALLOW_SALT_PLUME
86     # include "SALT_PLUME.h"
87     # endif
88 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
89    
90     C !INPUT/OUTPUT PARAMETERS:
91     C == Routine arguments ==
92 jmc 1.18 C myTime :: Current time in simulation
93     C myIter :: Current iteration number in simulation
94     C myThid :: Thread number for this instance of the routine.
95 jmc 1.1 _RL myTime
96     INTEGER myIter
97     INTEGER myThid
98    
99     C !LOCAL VARIABLES:
100     C == Local variables
101 jmc 1.47 C rhoK, rhoKm1 :: Density at current level, and level above
102 jmc 1.18 C iMin, iMax :: Ranges and sub-block indices on which calculations
103 jmc 1.1 C jMin, jMax are applied.
104 jmc 1.18 C bi, bj :: tile indices
105     C i,j,k :: loop indices
106 jmc 1.47 _RL rhoKp1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
107     _RL rhoKm1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
108 jmc 1.1 _RL sigmaX (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
109     _RL sigmaY (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
110     _RL sigmaR (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
111     INTEGER iMin, iMax
112     INTEGER jMin, jMax
113     INTEGER bi, bj
114 jmc 1.18 INTEGER i, j, k
115 jmc 1.17 INTEGER doDiagsRho
116 jmc 1.131 LOGICAL calcGMRedi, calcKPP, calcConvect
117 jmc 1.17 #ifdef ALLOW_DIAGNOSTICS
118     LOGICAL DIAGNOSTICS_IS_ON
119     EXTERNAL DIAGNOSTICS_IS_ON
120     #endif /* ALLOW_DIAGNOSTICS */
121 jmc 1.1
122     CEOP
123    
124 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
125     C-- dummy statement to end declaration part
126     itdkey = 1
127     #endif /* ALLOW_AUTODIFF_TAMC */
128    
129 jmc 1.1 #ifdef ALLOW_DEBUG
130 jmc 1.96 IF (debugMode) CALL DEBUG_ENTER('DO_OCEANIC_PHYS',myThid)
131 jmc 1.1 #endif
132 jmc 1.36
133 jmc 1.17 doDiagsRho = 0
134     #ifdef ALLOW_DIAGNOSTICS
135     IF ( useDiagnostics .AND. fluidIsWater ) THEN
136 jmc 1.110 IF ( DIAGNOSTICS_IS_ON('MXLDEPTH',myThid) )
137 jmc 1.71 & doDiagsRho = doDiagsRho + 1
138     IF ( DIAGNOSTICS_IS_ON('DRHODR ',myThid) )
139     & doDiagsRho = doDiagsRho + 2
140 jmc 1.110 IF ( DIAGNOSTICS_IS_ON('WdRHO_P ',myThid) )
141 jmc 1.71 & doDiagsRho = doDiagsRho + 4
142 jmc 1.110 IF ( DIAGNOSTICS_IS_ON('WdRHOdP ',myThid) )
143     & doDiagsRho = doDiagsRho + 8
144 jmc 1.17 ENDIF
145     #endif /* ALLOW_DIAGNOSTICS */
146    
147 jmc 1.131 calcGMRedi = useGMRedi
148     calcKPP = useKPP
149     calcConvect = ivdc_kappa.NE.0.
150     #ifdef ALLOW_OFFLINE
151     IF ( useOffLine ) THEN
152     calcGMRedi = useGMRedi .AND. .NOT.offlineLoadGMRedi
153     calcKPP = useKPP .AND. .NOT.offlineLoadKPP
154     calcConvect=calcConvect.AND. .NOT.offlineLoadConvec
155     ENDIF
156     #endif /* ALLOW_OFFLINE */
157    
158 jmc 1.82 #ifdef ALLOW_OBCS
159     IF (useOBCS) THEN
160     C-- Calculate future values on open boundaries
161     C-- moved before SEAICE_MODEL call since SEAICE_MODEL needs seaice-obcs fields
162 heimbach 1.100 # ifdef ALLOW_AUTODIFF_TAMC
163     CADJ STORE theta = comlev1, key=ikey_dynamics, kind=isbyte
164     CADJ STORE salt = comlev1, key=ikey_dynamics, kind=isbyte
165     # endif
166     # ifdef ALLOW_DEBUG
167 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('OBCS_CALC',myThid)
168 heimbach 1.100 # endif
169 jmc 1.120 CALL OBCS_CALC( myTime+deltaTClock, myIter+1,
170 jahn 1.83 I uVel, vVel, wVel, theta, salt, myThid )
171 jmc 1.82 ENDIF
172     #endif /* ALLOW_OBCS */
173 jmc 1.69
174 gforget 1.87 #ifdef ALLOW_AUTODIFF_TAMC
175     # ifdef ALLOW_SALT_PLUME
176     DO bj=myByLo(myThid),myByHi(myThid)
177     DO bi=myBxLo(myThid),myBxHi(myThid)
178     DO j=1-OLy,sNy+OLy
179     DO i=1-OLx,sNx+OLx
180     saltPlumeDepth(i,j,bi,bj) = 0. _d 0
181     saltPlumeFlux(i,j,bi,bj) = 0. _d 0
182     ENDDO
183     ENDDO
184     ENDDO
185     ENDDO
186     # endif
187     #endif /* ALLOW_AUTODIFF_TAMC */
188    
189 dimitri 1.113 #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 jmc 1.121 #ifndef OLD_THSICE_CALL_SEQUENCE
197     #if (defined ALLOW_THSICE) && !(defined ALLOW_ATM2D)
198     IF ( useThSIce .AND. fluidIsWater ) THEN
199 jmc 1.122 # 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 heimbach 1.125 CADJ STORE hocemxl = comlev1, key = ikey_dynamics,
211     CADJ & kind = isbyte
212 jmc 1.128 CADJ STORE icflxsw = comlev1, key = ikey_dynamics,
213 heimbach 1.126 CADJ & kind = isbyte
214 jmc 1.122 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 heimbach 1.126 # endif /* ALLOW_AUTODIFF_TAMC */
229 jmc 1.121 # 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 jmc 1.29 #ifdef ALLOW_SEAICE
242 gforget 1.124 # 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
270 jmc 1.29 IF ( useSEAICE ) THEN
271 heimbach 1.62 # ifdef ALLOW_AUTODIFF_TAMC
272 heimbach 1.65 cph-adj-test(
273 heimbach 1.81 CADJ STORE area = comlev1, key=ikey_dynamics, kind=isbyte
274     CADJ STORE hsnow = comlev1, key=ikey_dynamics, kind=isbyte
275 heimbach 1.88 CADJ STORE heff = comlev1, key=ikey_dynamics, kind=isbyte
276 jmc 1.132 CADJ STORE tices = comlev1, key=ikey_dynamics, kind=isbyte
277     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 heimbach 1.65 cph-adj-test)
283 jmc 1.132 c#ifdef ALLOW_EXF
284 heimbach 1.77 CADJ STORE atemp,aqh,precip = comlev1, key = ikey_dynamics,
285     CADJ & kind = isbyte
286     CADJ STORE swdown,lwdown = comlev1, key = ikey_dynamics,
287     CADJ & kind = isbyte
288     CADJ STORE evap = comlev1, key = ikey_dynamics,
289     CADJ & kind = isbyte
290 jmc 1.132 CADJ STORE uwind,vwind = comlev1, key = ikey_dynamics,
291     CADJ & kind = isbyte
292     c#endif
293 heimbach 1.77 CADJ STORE uvel,vvel = comlev1, key = ikey_dynamics,
294     CADJ & kind = isbyte
295 heimbach 1.95 # ifdef SEAICE_CGRID
296 heimbach 1.93 CADJ STORE stressdivergencex = comlev1, key = ikey_dynamics,
297     CADJ & kind = isbyte
298     CADJ STORE stressdivergencey = comlev1, key = ikey_dynamics,
299     CADJ & kind = isbyte
300     # endif
301 heimbach 1.62 # ifdef SEAICE_ALLOW_DYNAMICS
302 heimbach 1.77 CADJ STORE uice = comlev1, key = ikey_dynamics,
303     CADJ & kind = isbyte
304     CADJ STORE vice = comlev1, key = ikey_dynamics,
305     CADJ & kind = isbyte
306 jmc 1.132 CADJ STORE dwatn = comlev1, key = ikey_dynamics,
307     CADJ & kind = isbyte
308 heimbach 1.62 # ifdef SEAICE_ALLOW_EVP
309 heimbach 1.77 CADJ STORE seaice_sigma1 = comlev1, key = ikey_dynamics,
310     CADJ & kind = isbyte
311     CADJ STORE seaice_sigma2 = comlev1, key = ikey_dynamics,
312     CADJ & kind = isbyte
313     CADJ STORE seaice_sigma12 = comlev1, key = ikey_dynamics,
314     CADJ & kind = isbyte
315 heimbach 1.62 # endif
316     # endif
317 jmc 1.132 # ifdef SEAICE_VARIABLE_SALINITY
318     CADJ STORE hsalt = comlev1, key = ikey_dynamics,
319 heimbach 1.77 CADJ & kind = isbyte
320 jmc 1.132 # endif
321 heimbach 1.62 # ifdef ATMOSPHERIC_LOADING
322 heimbach 1.77 CADJ STORE pload = comlev1, key = ikey_dynamics,
323     CADJ & kind = isbyte
324     CADJ STORE siceload = comlev1, key = ikey_dynamics,
325     CADJ & kind = isbyte
326 heimbach 1.62 # endif
327     # ifdef NONLIN_FRSURF
328 heimbach 1.77 CADJ STORE recip_hfacc = comlev1, key = ikey_dynamics,
329     CADJ & kind = isbyte
330 heimbach 1.62 # endif
331 heimbach 1.78 # ifdef ANNUAL_BALANCE
332     CADJ STORE balance_itcount = comlev1, key = ikey_dynamics,
333     CADJ & kind = isbyte
334     # endif /* ANNUAL_BALANCE */
335 heimbach 1.126 # ifdef ALLOW_THSICE
336 jmc 1.132 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 heimbach 1.126 CADJ & kind = isbyte
341 jmc 1.132 CADJ STORE Qice1, Qice2 = comlev1, key = ikey_dynamics,
342     CADJ & kind = isbyte
343     # endif /* ALLOW_THSICE */
344 heimbach 1.126 # endif /* ALLOW_AUTODIFF_TAMC */
345 heimbach 1.62 # ifdef ALLOW_DEBUG
346 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('SEAICE_MODEL',myThid)
347 heimbach 1.62 # endif
348 jmc 1.29 CALL TIMER_START('SEAICE_MODEL [DO_OCEANIC_PHYS]', myThid)
349     CALL SEAICE_MODEL( myTime, myIter, myThid )
350     CALL TIMER_STOP ('SEAICE_MODEL [DO_OCEANIC_PHYS]', myThid)
351 heimbach 1.62 # ifdef ALLOW_COST
352 heimbach 1.57 CALL SEAICE_COST_SENSI ( myTime, myIter, myThid )
353 heimbach 1.62 # endif
354 heimbach 1.35 ENDIF
355 jmc 1.29 #endif /* ALLOW_SEAICE */
356    
357 heimbach 1.64 #ifdef ALLOW_AUTODIFF_TAMC
358 heimbach 1.77 CADJ STORE sst, sss = comlev1, key = ikey_dynamics,
359     CADJ & kind = isbyte
360     CADJ STORE qsw = comlev1, key = ikey_dynamics,
361     CADJ & kind = isbyte
362 heimbach 1.64 # ifdef ALLOW_SEAICE
363 heimbach 1.77 CADJ STORE area = comlev1, key = ikey_dynamics,
364     CADJ & kind = isbyte
365 heimbach 1.64 # endif
366     #endif
367    
368 jmc 1.121 #ifdef OLD_THSICE_CALL_SEQUENCE
369 jscott 1.30 #if (defined ALLOW_THSICE) && !(defined ALLOW_ATM2D)
370 jmc 1.14 IF ( useThSIce .AND. fluidIsWater ) THEN
371 heimbach 1.101 # ifdef ALLOW_AUTODIFF_TAMC
372     cph(
373     # ifdef NONLIN_FRSURF
374     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 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('THSICE_MAIN',myThid)
386 heimbach 1.101 # endif
387 jmc 1.5 C-- Step forward Therm.Sea-Ice variables
388     C and modify forcing terms including effects from ice
389     CALL TIMER_START('THSICE_MAIN [DO_OCEANIC_PHYS]', myThid)
390     CALL THSICE_MAIN( myTime, myIter, myThid )
391     CALL TIMER_STOP( 'THSICE_MAIN [DO_OCEANIC_PHYS]', myThid)
392     ENDIF
393     #endif /* ALLOW_THSICE */
394 jmc 1.121 #endif /* OLD_THSICE_CALL_SEQUENCE */
395 jmc 1.5
396 mlosch 1.21 #ifdef ALLOW_SHELFICE
397 heimbach 1.92 # ifdef ALLOW_AUTODIFF_TAMC
398     CADJ STORE salt, theta = comlev1, key = ikey_dynamics,
399     CADJ & kind = isbyte
400     # endif
401 mlosch 1.21 IF ( useShelfIce .AND. fluidIsWater ) THEN
402     #ifdef ALLOW_DEBUG
403 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('SHELFICE_THERMODYNAMICS',myThid)
404 mlosch 1.21 #endif
405 jmc 1.47 C compute temperature and (virtual) salt flux at the
406 mlosch 1.21 C shelf-ice ocean interface
407     CALL TIMER_START('SHELFICE_THERMODYNAMICS [DO_OCEANIC_PHYS]',
408     & myThid)
409     CALL SHELFICE_THERMODYNAMICS( myTime, myIter, myThid )
410     CALL TIMER_STOP( 'SHELFICE_THERMODYNAMICS [DO_OCEANIC_PHYS]',
411     & myThid)
412     ENDIF
413     #endif /* ALLOW_SHELFICE */
414    
415 dimitri 1.85 #ifdef ALLOW_ICEFRONT
416     IF ( useICEFRONT .AND. fluidIsWater ) THEN
417     #ifdef ALLOW_DEBUG
418 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('ICEFRONT_THERMODYNAMICS',myThid)
419 dimitri 1.85 #endif
420     C compute temperature and (virtual) salt flux at the
421     C ice-front ocean interface
422     CALL TIMER_START('ICEFRONT_THERMODYNAMICS [DO_OCEANIC_PHYS]',
423     & myThid)
424     CALL ICEFRONT_THERMODYNAMICS( myTime, myIter, myThid )
425     CALL TIMER_STOP( 'ICEFRONT_THERMODYNAMICS [DO_OCEANIC_PHYS]',
426     & myThid)
427     ENDIF
428     #endif /* ALLOW_ICEFRONT */
429    
430 jmc 1.112 #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 jmc 1.5 C-- Freeze water at the surface
437 heimbach 1.104 IF ( allowFreezing ) THEN
438 jmc 1.5 #ifdef ALLOW_AUTODIFF_TAMC
439 heimbach 1.77 CADJ STORE theta = comlev1, key = ikey_dynamics,
440     CADJ & kind = isbyte
441 jmc 1.5 #endif
442 jmc 1.129 CALL FREEZE_SURFACE( myTime, myIter, myThid )
443 jmc 1.5 ENDIF
444    
445 jmc 1.28 #ifdef ALLOW_OCN_COMPON_INTERF
446 jmc 1.5 C-- Apply imported data (from coupled interface) to forcing fields
447 jmc 1.28 C jmc: do not know precisely where to put this call (bf or af thSIce ?)
448 jmc 1.36 IF ( useCoupler ) THEN
449 jmc 1.5 CALL OCN_APPLY_IMPORT( .TRUE., myTime, myIter, myThid )
450 jmc 1.36 ENDIF
451 jmc 1.28 #endif /* ALLOW_OCN_COMPON_INTERF */
452 jmc 1.5
453 jmc 1.130 iMin = 1-OLx
454     iMax = sNx+OLx
455     jMin = 1-OLy
456     jMax = sNy+OLy
457    
458     C--- Determines forcing terms based on external fields
459     C relaxation terms, etc.
460     #ifdef ALLOW_DEBUG
461     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 jmc 1.127 IF ( fluidIsWater ) THEN
470 jmc 1.130 #endif /* ALLOW_AUTODIFF_TAMC */
471     CALL EXTERNAL_FORCING_SURF(
472     I iMin, iMax, jMin, jMax,
473     I myTime, myIter, myThid )
474 jmc 1.25
475 jmc 1.1 #ifdef ALLOW_AUTODIFF_TAMC
476     C-- HPF directive to help TAMC
477     CHPF$ INDEPENDENT
478     #endif /* ALLOW_AUTODIFF_TAMC */
479     DO bj=myByLo(myThid),myByHi(myThid)
480     #ifdef ALLOW_AUTODIFF_TAMC
481 heimbach 1.15 C-- HPF directive to help TAMC
482     CHPF$ INDEPENDENT
483 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
484     DO bi=myBxLo(myThid),myBxHi(myThid)
485    
486     #ifdef ALLOW_AUTODIFF_TAMC
487     act1 = bi - myBxLo(myThid)
488     max1 = myBxHi(myThid) - myBxLo(myThid) + 1
489     act2 = bj - myByLo(myThid)
490     max2 = myByHi(myThid) - myByLo(myThid) + 1
491     act3 = myThid - 1
492     max3 = nTx*nTy
493     act4 = ikey_dynamics - 1
494     itdkey = (act1 + 1) + act2*max1
495     & + act3*max1*max2
496     & + act4*max1*max2*max3
497 jmc 1.74 #endif /* ALLOW_AUTODIFF_TAMC */
498 jmc 1.1
499     C-- Set up work arrays with valid (i.e. not NaN) values
500     C These inital values do not alter the numerical results. They
501     C just ensure that all memory references are to valid floating
502     C point numbers. This prevents spurious hardware signals due to
503     C uninitialised but inert locations.
504 jmc 1.131 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 jmc 1.1
515 jmc 1.69 #ifdef ALLOW_AUTODIFF_TAMC
516 jmc 1.1 DO j=1-OLy,sNy+OLy
517     DO i=1-OLx,sNx+OLx
518 jmc 1.69 rhoKm1 (i,j) = 0. _d 0
519 jmc 1.47 rhoKp1 (i,j) = 0. _d 0
520 jmc 1.1 ENDDO
521     ENDDO
522 jmc 1.131 cph all the following init. are necessary for TAF
523     cph although some of these are re-initialised later.
524 jmc 1.1 DO k=1,Nr
525     DO j=1-OLy,sNy+OLy
526     DO i=1-OLx,sNx+OLx
527 heimbach 1.109 rhoInSitu(i,j,k,bi,bj) = 0.
528 jmc 1.131 # 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
540     DO j=1-OLy,sNy+OLy
541     DO i=1-OLx,sNx+OLx
542 jmc 1.1 IVDConvCount(i,j,k,bi,bj) = 0.
543 jmc 1.131 ENDDO
544     ENDDO
545     ENDDO
546     #ifdef ALLOW_OFFLINE
547     ENDIF
548     IF ( calcGMRedi ) THEN
549     #endif
550 jmc 1.1 # ifdef ALLOW_GMREDI
551 jmc 1.131 DO k=1,Nr
552     DO j=1-OLy,sNy+OLy
553     DO i=1-OLx,sNx+OLx
554 jmc 1.1 Kwx(i,j,k,bi,bj) = 0. _d 0
555     Kwy(i,j,k,bi,bj) = 0. _d 0
556     Kwz(i,j,k,bi,bj) = 0. _d 0
557     # ifdef GM_NON_UNITY_DIAGONAL
558     Kux(i,j,k,bi,bj) = 0. _d 0
559     Kvy(i,j,k,bi,bj) = 0. _d 0
560     # endif
561     # ifdef GM_EXTRA_DIAGONAL
562     Kuz(i,j,k,bi,bj) = 0. _d 0
563     Kvz(i,j,k,bi,bj) = 0. _d 0
564     # endif
565     # ifdef GM_BOLUS_ADVEC
566     GM_PsiX(i,j,k,bi,bj) = 0. _d 0
567     GM_PsiY(i,j,k,bi,bj) = 0. _d 0
568     # endif
569     # ifdef GM_VISBECK_VARIABLE_K
570     VisbeckK(i,j,bi,bj) = 0. _d 0
571     # endif
572 jmc 1.131 ENDDO
573     ENDDO
574     ENDDO
575 jmc 1.1 # endif /* ALLOW_GMREDI */
576 jmc 1.131 #ifdef ALLOW_OFFLINE
577     ENDIF
578     IF ( calcKPP ) THEN
579     #endif
580 heimbach 1.42 # ifdef ALLOW_KPP
581 jmc 1.131 DO k=1,Nr
582     DO j=1-OLy,sNy+OLy
583     DO i=1-OLx,sNx+OLx
584 heimbach 1.42 KPPdiffKzS(i,j,k,bi,bj) = 0. _d 0
585     KPPdiffKzT(i,j,k,bi,bj) = 0. _d 0
586 jmc 1.1 ENDDO
587     ENDDO
588     ENDDO
589 jmc 1.131 # endif /* ALLOW_KPP */
590     #ifdef ALLOW_OFFLINE
591     ENDIF
592     #endif
593     #endif /* ALLOW_AUTODIFF_TAMC */
594 jmc 1.1
595     #ifdef ALLOW_AUTODIFF_TAMC
596 jmc 1.96 CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
597 heimbach 1.77 CADJ & kind = isbyte
598 jmc 1.96 CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
599 heimbach 1.77 CADJ & kind = isbyte
600 heimbach 1.12 CADJ STORE totphihyd(:,:,:,bi,bj)
601 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
602 heimbach 1.77 CADJ & kind = isbyte
603 heimbach 1.10 # ifdef ALLOW_KPP
604 jmc 1.96 CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
605 heimbach 1.77 CADJ & kind = isbyte
606 jmc 1.96 CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
607 heimbach 1.77 CADJ & kind = isbyte
608 heimbach 1.10 # endif
609 heimbach 1.115 # ifdef ALLOW_SALT_PLUME
610     CADJ STORE saltplumedepth(:,:,bi,bj) = comlev1_bibj, key=itdkey,
611     CADJ & kind = isbyte
612     # endif
613 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
614    
615 jmc 1.71 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)
617     #ifdef ALLOW_DEBUG
618 jmc 1.110 IF (debugMode) CALL DEBUG_CALL('FIND_RHO_2D (xNr)',myThid)
619 jmc 1.71 #endif
620     #ifdef ALLOW_AUTODIFF_TAMC
621 jmc 1.110 IF ( fluidIsWater ) THEN
622 jmc 1.71 #endif /* ALLOW_AUTODIFF_TAMC */
623 jmc 1.69 #ifdef ALLOW_DOWN_SLOPE
624 jmc 1.110 IF ( useDOWN_SLOPE ) THEN
625     DO k=1,Nr
626 jmc 1.69 CALL DWNSLP_CALC_RHO(
627     I theta, salt,
628 jmc 1.71 O rhoInSitu(1-OLx,1-OLy,k,bi,bj),
629 jmc 1.69 I k, bi, bj, myTime, myIter, myThid )
630 jmc 1.110 ENDDO
631     ENDIF
632 jmc 1.71 #endif /* ALLOW_DOWN_SLOPE */
633 dimitri 1.107 #ifdef ALLOW_BBL
634 jmc 1.110 IF ( useBBL ) THEN
635 dimitri 1.108 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 jmc 1.110 DO k=Nr,1,-1
639 dimitri 1.107 CALL BBL_CALC_RHO(
640     I theta, salt,
641     O rhoInSitu,
642     I k, bi, bj, myTime, myIter, myThid )
643    
644 jmc 1.110 ENDDO
645     ENDIF
646 dimitri 1.107 #endif /* ALLOW_BBL */
647 jmc 1.110 IF ( .NOT. ( useDOWN_SLOPE .OR. useBBL ) ) THEN
648     DO k=1,Nr
649 jmc 1.71 CALL FIND_RHO_2D(
650     I iMin, iMax, jMin, jMax, k,
651     I theta(1-OLx,1-OLy,k,bi,bj),
652     I salt (1-OLx,1-OLy,k,bi,bj),
653     O rhoInSitu(1-OLx,1-OLy,k,bi,bj),
654     I k, bi, bj, myThid )
655 jmc 1.110 ENDDO
656     ENDIF
657 jmc 1.74 #ifdef ALLOW_AUTODIFF_TAMC
658 jmc 1.110 ELSE
659 jmc 1.74 C- fluid is not water:
660 jmc 1.110 DO k=1,Nr
661 jmc 1.74 DO j=1-OLy,sNy+OLy
662     DO i=1-OLx,sNx+OLx
663     rhoInSitu(i,j,k,bi,bj) = 0.
664     ENDDO
665     ENDDO
666 jmc 1.110 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 jmc 1.74 #endif /* ALLOW_AUTODIFF_TAMC */
684 jmc 1.69
685 jmc 1.110 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 jmc 1.1 C-- Calculate gradients of potential density for isoneutral
693     C slope terms (e.g. GM/Redi tensor or IVDC diffusivity)
694 jmc 1.131 IF ( calcGMRedi .OR. (k.GT.1 .AND. calcConvect)
695 jmc 1.116 & .OR. usePP81 .OR. useMY82 .OR. useGGL90
696 dimitri 1.61 & .OR. useSALT_PLUME .OR. doDiagsRho.GE.1 ) THEN
697 jmc 1.1 IF (k.GT.1) THEN
698     #ifdef ALLOW_AUTODIFF_TAMC
699 jmc 1.96 CADJ STORE theta(:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey,
700 heimbach 1.77 CADJ & kind = isbyte
701 jmc 1.96 CADJ STORE salt (:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey,
702 heimbach 1.77 CADJ & kind = isbyte
703 jmc 1.96 CADJ STORE rhokm1 (bi,bj) = comlev1_bibj_k, key=kkey,
704 heimbach 1.77 CADJ & kind = isbyte
705 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
706 jmc 1.68 CALL FIND_RHO_2D(
707     I iMin, iMax, jMin, jMax, k,
708     I theta(1-OLx,1-OLy,k-1,bi,bj),
709     I salt (1-OLx,1-OLy,k-1,bi,bj),
710     O rhoKm1,
711     I k-1, bi, bj, myThid )
712 jmc 1.1 ENDIF
713     #ifdef ALLOW_DEBUG
714 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('GRAD_SIGMA',myThid)
715 jmc 1.1 #endif
716 heimbach 1.31 cph Avoid variable aliasing for adjoint !!!
717     DO j=jMin,jMax
718     DO i=iMin,iMax
719 jmc 1.71 rhoKp1(i,j) = rhoInSitu(i,j,k,bi,bj)
720 heimbach 1.31 ENDDO
721     ENDDO
722 jmc 1.1 CALL GRAD_SIGMA(
723     I bi, bj, iMin, iMax, jMin, jMax, k,
724 jmc 1.71 I rhoInSitu(1-OLx,1-OLy,k,bi,bj), rhoKm1, rhoKp1,
725 jmc 1.1 O sigmaX, sigmaY, sigmaR,
726     I myThid )
727 gforget 1.66 #ifdef ALLOW_AUTODIFF_TAMC
728 jmc 1.69 #ifdef GMREDI_WITH_STABLE_ADJOINT
729 gforget 1.66 cgf zero out adjoint fields to stabilize pkg/gmredi adjoint
730     cgf -> cuts adjoint dependency from slope to state
731 jmc 1.69 CALL ZERO_ADJ_LOC( Nr, sigmaX, myThid)
732     CALL ZERO_ADJ_LOC( Nr, sigmaY, myThid)
733     CALL ZERO_ADJ_LOC( Nr, sigmaR, myThid)
734 gforget 1.66 #endif
735     #endif /* ALLOW_AUTODIFF_TAMC */
736 jmc 1.1 ENDIF
737    
738     C-- Implicit Vertical Diffusion for Convection
739 jmc 1.131 IF (k.GT.1 .AND. calcConvect) THEN
740 jmc 1.1 #ifdef ALLOW_DEBUG
741 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('CALC_IVDC',myThid)
742 jmc 1.1 #endif
743     CALL CALC_IVDC(
744     I bi, bj, iMin, iMax, jMin, jMax, k,
745 mlosch 1.111 I sigmaR,
746 jmc 1.1 I myTime, myIter, myThid)
747     ENDIF
748    
749 jmc 1.17 #ifdef ALLOW_DIAGNOSTICS
750 jmc 1.110 IF ( doDiagsRho.GE.4 ) THEN
751     CALL DIAGS_RHO_L( doDiagsRho, k, bi, bj,
752     I rhoInSitu(1-OLx,1-OLy,1,bi,bj),
753 jmc 1.74 I rhoKm1, wVel,
754 jmc 1.71 I myTime, myIter, myThid )
755 jmc 1.17 ENDIF
756     #endif
757    
758 jmc 1.1 C-- end of diagnostic k loop (Nr:1)
759     ENDDO
760    
761 heimbach 1.57 #ifdef ALLOW_AUTODIFF_TAMC
762 jmc 1.69 CADJ STORE IVDConvCount(:,:,:,bi,bj)
763 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
764 heimbach 1.77 CADJ & kind = isbyte
765 heimbach 1.57 #endif
766    
767 jmc 1.47 C-- Diagnose Mixed Layer Depth:
768 jmc 1.131 IF ( calcGMRedi .OR. MOD(doDiagsRho,2).EQ.1 ) THEN
769 jmc 1.71 CALL CALC_OCE_MXLAYER(
770     I rhoInSitu(1-OLx,1-OLy,1,bi,bj), sigmaR,
771     I bi, bj, myTime, myIter, myThid )
772 jmc 1.47 ENDIF
773 heimbach 1.53
774 dimitri 1.52 #ifdef ALLOW_SALT_PLUME
775 dimitri 1.61 IF ( useSALT_PLUME ) THEN
776 jmc 1.71 CALL SALT_PLUME_CALC_DEPTH(
777     I rhoInSitu(1-OLx,1-OLy,1,bi,bj), sigmaR,
778     I bi, bj, myTime, myIter, myThid )
779 dimitri 1.60 ENDIF
780 dimitri 1.61 #endif /* ALLOW_SALT_PLUME */
781    
782 jmc 1.8 #ifdef ALLOW_DIAGNOSTICS
783 jmc 1.71 IF ( MOD(doDiagsRho,4).GE.2 ) THEN
784 jmc 1.16 CALL DIAGNOSTICS_FILL (sigmaR, 'DRHODR ', 0, Nr,
785     & 2, bi, bj, myThid)
786 jmc 1.8 ENDIF
787 dimitri 1.61 #endif /* ALLOW_DIAGNOSTICS */
788 jmc 1.8
789 jmc 1.130 C-- This is where EXTERNAL_FORCING_SURF(bi,bj) used to be called;
790     C now called earlier, before bi,bj loop.
791 jmc 1.1
792     #ifdef ALLOW_AUTODIFF_TAMC
793     cph needed for KPP
794 jmc 1.4 CADJ STORE surfaceForcingU(:,:,bi,bj)
795 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
796 heimbach 1.77 CADJ & kind = isbyte
797 jmc 1.4 CADJ STORE surfaceForcingV(:,:,bi,bj)
798 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
799 heimbach 1.77 CADJ & kind = isbyte
800 jmc 1.4 CADJ STORE surfaceForcingS(:,:,bi,bj)
801 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
802 heimbach 1.77 CADJ & kind = isbyte
803 jmc 1.4 CADJ STORE surfaceForcingT(:,:,bi,bj)
804 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
805 heimbach 1.77 CADJ & kind = isbyte
806 jmc 1.4 CADJ STORE surfaceForcingTice(:,:,bi,bj)
807 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
808 heimbach 1.77 CADJ & kind = isbyte
809 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
810    
811     #ifdef ALLOW_KPP
812     C-- Compute KPP mixing coefficients
813 jmc 1.131 IF ( calcKPP ) THEN
814 jmc 1.1 #ifdef ALLOW_DEBUG
815 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('KPP_CALC',myThid)
816 jmc 1.1 #endif
817 dfer 1.76 CALL TIMER_START('KPP_CALC [DO_OCEANIC_PHYS]', myThid)
818 jmc 1.1 CALL KPP_CALC(
819 jmc 1.44 I bi, bj, myTime, myIter, myThid )
820 dfer 1.76 CALL TIMER_STOP ('KPP_CALC [DO_OCEANIC_PHYS]', myThid)
821 jmc 1.131 #if (defined ALLOW_AUTODIFF_TAMC) && !(defined ALLOW_OFFLINE)
822 jmc 1.1 ELSE
823     CALL KPP_CALC_DUMMY(
824 jmc 1.44 I bi, bj, myTime, myIter, myThid )
825 jmc 1.131 #endif /* ALLOW_AUTODIFF_TAMC and not ALLOW_OFFLINE */
826 jmc 1.1 ENDIF
827     #endif /* ALLOW_KPP */
828    
829 mlosch 1.6 #ifdef ALLOW_PP81
830     C-- Compute PP81 mixing coefficients
831     IF (usePP81) THEN
832     #ifdef ALLOW_DEBUG
833 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('PP81_CALC',myThid)
834 mlosch 1.6 #endif
835     CALL PP81_CALC(
836 jmc 1.116 I bi, bj, sigmaR, myTime, myIter, myThid )
837 mlosch 1.6 ENDIF
838     #endif /* ALLOW_PP81 */
839    
840     #ifdef ALLOW_MY82
841     C-- Compute MY82 mixing coefficients
842     IF (useMY82) THEN
843     #ifdef ALLOW_DEBUG
844 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('MY82_CALC',myThid)
845 mlosch 1.6 #endif
846     CALL MY82_CALC(
847 jmc 1.116 I bi, bj, sigmaR, myTime, myIter, myThid )
848 mlosch 1.6 ENDIF
849     #endif /* ALLOW_MY82 */
850    
851 mlosch 1.9 #ifdef ALLOW_GGL90
852 gforget 1.91 #ifdef ALLOW_AUTODIFF_TAMC
853     CADJ STORE GGL90TKE (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
854     CADJ & kind = isbyte
855     #endif /* ALLOW_AUTODIFF_TAMC */
856 mlosch 1.9 C-- Compute GGL90 mixing coefficients
857     IF (useGGL90) THEN
858     #ifdef ALLOW_DEBUG
859 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('GGL90_CALC',myThid)
860 mlosch 1.9 #endif
861 dfer 1.76 CALL TIMER_START('GGL90_CALC [DO_OCEANIC_PHYS]', myThid)
862 mlosch 1.9 CALL GGL90_CALC(
863 jmc 1.116 I bi, bj, sigmaR, myTime, myIter, myThid )
864 dfer 1.76 CALL TIMER_STOP ('GGL90_CALC [DO_OCEANIC_PHYS]', myThid)
865 mlosch 1.9 ENDIF
866     #endif /* ALLOW_GGL90 */
867    
868 jmc 1.20 #ifdef ALLOW_TIMEAVE
869 jmc 1.36 IF ( taveFreq.GT. 0. _d 0 ) THEN
870 jmc 1.20 CALL TIMEAVE_SURF_FLUX( bi, bj, myTime, myIter, myThid)
871     ENDIF
872 jmc 1.131 IF ( taveFreq.GT.0. .AND. calcConvect ) THEN
873 jmc 1.20 CALL TIMEAVE_CUMULATE(ConvectCountTave, IVDConvCount,
874 jmc 1.120 I Nr, deltaTClock, bi, bj, myThid)
875 jmc 1.20 ENDIF
876     #endif /* ALLOW_TIMEAVE */
877    
878 jmc 1.69 #ifdef ALLOW_GMREDI
879 jmc 1.47 #ifdef ALLOW_AUTODIFF_TAMC
880     # ifndef GM_EXCLUDE_CLIPPING
881     cph storing here is needed only for one GMREDI_OPTIONS:
882     cph define GM_BOLUS_ADVEC
883     cph keep it although TAF says you dont need to.
884 jmc 1.86 cph but I have avoided the #ifdef for now, in case more things change
885 jmc 1.96 CADJ STORE sigmaX(:,:,:) = comlev1_bibj, key=itdkey,
886 heimbach 1.77 CADJ & kind = isbyte
887 jmc 1.96 CADJ STORE sigmaY(:,:,:) = comlev1_bibj, key=itdkey,
888 heimbach 1.77 CADJ & kind = isbyte
889 jmc 1.96 CADJ STORE sigmaR(:,:,:) = comlev1_bibj, key=itdkey,
890 heimbach 1.77 CADJ & kind = isbyte
891 jmc 1.47 # endif
892     #endif /* ALLOW_AUTODIFF_TAMC */
893    
894     C-- Calculate iso-neutral slopes for the GM/Redi parameterisation
895 jmc 1.131 IF ( calcGMRedi ) THEN
896 jmc 1.47 #ifdef ALLOW_DEBUG
897 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('GMREDI_CALC_TENSOR',myThid)
898 jmc 1.47 #endif
899     CALL GMREDI_CALC_TENSOR(
900 jmc 1.51 I iMin, iMax, jMin, jMax,
901 jmc 1.47 I sigmaX, sigmaY, sigmaR,
902 jmc 1.51 I bi, bj, myTime, myIter, myThid )
903 jmc 1.131 #if (defined ALLOW_AUTODIFF_TAMC) && !(defined ALLOW_OFFLINE)
904 jmc 1.47 ELSE
905     CALL GMREDI_CALC_TENSOR_DUMMY(
906 jmc 1.51 I iMin, iMax, jMin, jMax,
907 jmc 1.47 I sigmaX, sigmaY, sigmaR,
908 jmc 1.51 I bi, bj, myTime, myIter, myThid )
909 jmc 1.131 #endif /* ALLOW_AUTODIFF_TAMC and not ALLOW_OFFLINE */
910 jmc 1.47 ENDIF
911 jmc 1.69 #endif /* ALLOW_GMREDI */
912    
913     #ifdef ALLOW_DOWN_SLOPE
914     IF ( useDOWN_SLOPE ) THEN
915     C-- Calculate Downsloping Flow for Down_Slope parameterization
916     IF ( usingPCoords ) THEN
917     CALL DWNSLP_CALC_FLOW(
918 jmc 1.71 I bi, bj, kSurfC, rhoInSitu,
919 jmc 1.69 I myTime, myIter, myThid )
920     ELSE
921     CALL DWNSLP_CALC_FLOW(
922 jmc 1.71 I bi, bj, kLowC, rhoInSitu,
923 jmc 1.69 I myTime, myIter, myThid )
924     ENDIF
925     ENDIF
926     #endif /* ALLOW_DOWN_SLOPE */
927 jmc 1.47
928 jmc 1.106 C-- end bi,bj loops.
929     ENDDO
930     ENDDO
931    
932 jmc 1.98 #ifndef ALLOW_AUTODIFF_TAMC
933     C--- if fluid Is Water: end
934 jmc 1.106 ENDIF
935 jmc 1.98 #endif
936    
937 dimitri 1.107 #ifdef ALLOW_BBL
938     IF ( useBBL ) THEN
939     CALL BBL_CALC_RHS(
940     I myTime, myIter, myThid )
941     ENDIF
942     #endif /* ALLOW_BBL */
943    
944 dimitri 1.94 #ifdef ALLOW_MYPACKAGE
945 jmc 1.106 IF ( useMYPACKAGE ) THEN
946     CALL MYPACKAGE_CALC_RHS(
947     I myTime, myIter, myThid )
948     ENDIF
949 dimitri 1.94 #endif /* ALLOW_MYPACKAGE */
950    
951 jmc 1.99 #ifdef ALLOW_GMREDI
952 jmc 1.131 IF ( calcGMRedi ) THEN
953 jmc 1.99 CALL GMREDI_DO_EXCH( myTime, myIter, myThid )
954     ENDIF
955     #endif /* ALLOW_GMREDI */
956    
957 jmc 1.98 #ifdef ALLOW_KPP
958 jmc 1.131 IF ( calcKPP ) THEN
959 jmc 1.45 CALL KPP_DO_EXCH( myThid )
960     ENDIF
961 jmc 1.98 #endif /* ALLOW_KPP */
962 jmc 1.45
963 jmc 1.18 #ifdef ALLOW_DIAGNOSTICS
964     IF ( fluidIsWater .AND. useDiagnostics ) THEN
965 jmc 1.74 CALL DIAGS_RHO_G(
966 jmc 1.110 I rhoInSitu, uVel, vVel, wVel,
967 jmc 1.71 I myTime, myIter, myThid )
968 jmc 1.133 ENDIF
969     IF ( useDiagnostics ) THEN
970 jmc 1.18 CALL DIAGS_OCEANIC_SURF_FLUX( myTime, myIter, myThid )
971     ENDIF
972 jmc 1.131 IF ( calcConvect .AND. useDiagnostics ) THEN
973 jmc 1.71 CALL DIAGNOSTICS_FILL( IVDConvCount, 'CONVADJ ',
974     & 0, Nr, 0, 1, 1, myThid )
975 jmc 1.19 ENDIF
976 jmc 1.18 #endif
977    
978 gforget 1.123 #ifdef ALLOW_ECCO
979 jmc 1.129 CALL ECCO_PHYS( myThid )
980 gforget 1.123 #endif
981    
982 jmc 1.1 #ifdef ALLOW_DEBUG
983 jmc 1.96 IF (debugMode) CALL DEBUG_LEAVE('DO_OCEANIC_PHYS',myThid)
984 jmc 1.1 #endif
985    
986     RETURN
987     END

  ViewVC Help
Powered by ViewVC 1.1.22