/[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.127 - (hide annotations) (download)
Mon Apr 22 02:38:07 2013 UTC (11 years, 1 month ago) by jmc
Branch: MAIN
CVS Tags: checkpoint64g
Changes since 1.126: +5 -15 lines 
move forcing adjustment (balancing surface forcing) in specific S/R
where might also be applied exchanges (if needed) to forcing arrays.
1 jmc 1.127 C $Header: /u/gcmpack/MITgcm/model/src/do_oceanic_phys.F,v 1.126 2013/04/13 20:47:18 heimbach 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     #include "TIMEAVE_STATV.h"
46     #endif
47 jmc 1.1
48     #ifdef ALLOW_AUTODIFF_TAMC
49 heimbach 1.80 # include "AUTODIFF_MYFIELDS.h"
50 jmc 1.1 # include "tamc.h"
51     # include "tamc_keys.h"
52     # include "FFIELDS.h"
53 heimbach 1.54 # include "SURFACE.h"
54 jmc 1.1 # include "EOS.h"
55     # ifdef ALLOW_KPP
56     # include "KPP.h"
57     # endif
58 gforget 1.91 # ifdef ALLOW_GGL90
59     # include "GGL90.h"
60     # endif
61 jmc 1.1 # ifdef ALLOW_GMREDI
62     # include "GMREDI.h"
63     # endif
64     # ifdef ALLOW_EBM
65     # include "EBM.h"
66     # endif
67 jmc 1.29 # ifdef ALLOW_EXF
68     # include "ctrl.h"
69 jmc 1.40 # include "EXF_FIELDS.h"
70 jmc 1.29 # ifdef ALLOW_BULKFORMULAE
71 jmc 1.40 # include "EXF_CONSTANTS.h"
72 jmc 1.29 # endif
73     # endif
74     # ifdef ALLOW_SEAICE
75 jmc 1.114 # include "SEAICE_SIZE.h"
76 jmc 1.29 # include "SEAICE.h"
77 gforget 1.124 # include "SEAICE_PARAMS.h"
78 jmc 1.29 # endif
79 heimbach 1.105 # ifdef ALLOW_THSICE
80     # include "THSICE_VARS.h"
81     # endif
82 heimbach 1.75 # ifdef ALLOW_SALT_PLUME
83     # include "SALT_PLUME.h"
84     # endif
85 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
86    
87     C !INPUT/OUTPUT PARAMETERS:
88     C == Routine arguments ==
89 jmc 1.18 C myTime :: Current time in simulation
90     C myIter :: Current iteration number in simulation
91     C myThid :: Thread number for this instance of the routine.
92 jmc 1.1 _RL myTime
93     INTEGER myIter
94     INTEGER myThid
95    
96     C !LOCAL VARIABLES:
97     C == Local variables
98 jmc 1.47 C rhoK, rhoKm1 :: Density at current level, and level above
99 jmc 1.18 C iMin, iMax :: Ranges and sub-block indices on which calculations
100 jmc 1.1 C jMin, jMax are applied.
101 jmc 1.18 C bi, bj :: tile indices
102     C i,j,k :: loop indices
103 jmc 1.47 _RL rhoKp1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
104     _RL rhoKm1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
105 jmc 1.1 _RL sigmaX (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
106     _RL sigmaY (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
107     _RL sigmaR (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
108     INTEGER iMin, iMax
109     INTEGER jMin, jMax
110     INTEGER bi, bj
111 jmc 1.18 INTEGER i, j, k
112 jmc 1.17 INTEGER doDiagsRho
113     #ifdef ALLOW_DIAGNOSTICS
114     LOGICAL DIAGNOSTICS_IS_ON
115     EXTERNAL DIAGNOSTICS_IS_ON
116     #endif /* ALLOW_DIAGNOSTICS */
117 jmc 1.1
118     CEOP
119    
120 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
121     C-- dummy statement to end declaration part
122     itdkey = 1
123     #endif /* ALLOW_AUTODIFF_TAMC */
124    
125 jmc 1.1 #ifdef ALLOW_DEBUG
126 jmc 1.96 IF (debugMode) CALL DEBUG_ENTER('DO_OCEANIC_PHYS',myThid)
127 jmc 1.1 #endif
128 jmc 1.36
129 jmc 1.17 doDiagsRho = 0
130     #ifdef ALLOW_DIAGNOSTICS
131     IF ( useDiagnostics .AND. fluidIsWater ) THEN
132 jmc 1.110 IF ( DIAGNOSTICS_IS_ON('MXLDEPTH',myThid) )
133 jmc 1.71 & doDiagsRho = doDiagsRho + 1
134     IF ( DIAGNOSTICS_IS_ON('DRHODR ',myThid) )
135     & doDiagsRho = doDiagsRho + 2
136 jmc 1.110 IF ( DIAGNOSTICS_IS_ON('WdRHO_P ',myThid) )
137 jmc 1.71 & doDiagsRho = doDiagsRho + 4
138 jmc 1.110 IF ( DIAGNOSTICS_IS_ON('WdRHOdP ',myThid) )
139     & doDiagsRho = doDiagsRho + 8
140 jmc 1.17 ENDIF
141     #endif /* ALLOW_DIAGNOSTICS */
142    
143 jmc 1.82 #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 heimbach 1.100 # 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 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('OBCS_CALC',myThid)
153 heimbach 1.100 # endif
154 jmc 1.120 CALL OBCS_CALC( myTime+deltaTClock, myIter+1,
155 jahn 1.83 I uVel, vVel, wVel, theta, salt, myThid )
156 jmc 1.82 ENDIF
157     #endif /* ALLOW_OBCS */
158 jmc 1.69
159 gforget 1.87 #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 dimitri 1.113 #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 jmc 1.121 #ifndef OLD_THSICE_CALL_SEQUENCE
182     #if (defined ALLOW_THSICE) && !(defined ALLOW_ATM2D)
183     IF ( useThSIce .AND. fluidIsWater ) THEN
184 jmc 1.122 # ifdef ALLOW_AUTODIFF_TAMC
185     CADJ STORE uice,vice = comlev1, key = ikey_dynamics,
186     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 heimbach 1.125 CADJ STORE hocemxl = comlev1, key = ikey_dynamics,
196     CADJ & kind = isbyte
197 heimbach 1.126 CADJ STORE icflxsw, snowprc = comlev1, key = ikey_dynamics,
198     CADJ & kind = isbyte
199 jmc 1.122 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 heimbach 1.126 # endif /* ALLOW_AUTODIFF_TAMC */
214 jmc 1.121 # 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 jmc 1.29 #ifdef ALLOW_SEAICE
227 gforget 1.124 # 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 jmc 1.29 IF ( useSEAICE ) THEN
256 heimbach 1.62 # ifdef ALLOW_AUTODIFF_TAMC
257 heimbach 1.65 cph-adj-test(
258 heimbach 1.81 CADJ STORE area = comlev1, key=ikey_dynamics, kind=isbyte
259     CADJ STORE hsnow = comlev1, key=ikey_dynamics, kind=isbyte
260 heimbach 1.88 CADJ STORE heff = comlev1, key=ikey_dynamics, kind=isbyte
261 heimbach 1.81 CADJ STORE empmr,qsw,theta = comlev1, key = ikey_dynamics,
262 heimbach 1.77 CADJ & kind = isbyte
263 heimbach 1.65 cph-adj-test)
264 heimbach 1.77 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 heimbach 1.34 cph# ifdef EXF_READ_EVAP
269 heimbach 1.77 CADJ STORE evap = comlev1, key = ikey_dynamics,
270     CADJ & kind = isbyte
271 heimbach 1.34 cph# endif
272 heimbach 1.77 CADJ STORE uvel,vvel = comlev1, key = ikey_dynamics,
273     CADJ & kind = isbyte
274 heimbach 1.95 # ifdef SEAICE_CGRID
275 heimbach 1.93 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 heimbach 1.62 # ifdef SEAICE_ALLOW_DYNAMICS
281 heimbach 1.77 CADJ STORE uice = comlev1, key = ikey_dynamics,
282     CADJ & kind = isbyte
283     CADJ STORE vice = comlev1, key = ikey_dynamics,
284     CADJ & kind = isbyte
285 heimbach 1.62 # ifdef SEAICE_ALLOW_EVP
286 heimbach 1.77 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 heimbach 1.62 # endif
293     # endif
294 heimbach 1.104 cph# ifdef SEAICE_SALINITY
295 heimbach 1.77 CADJ STORE salt = comlev1, key = ikey_dynamics,
296     CADJ & kind = isbyte
297 heimbach 1.104 cph# endif
298 heimbach 1.62 # ifdef ATMOSPHERIC_LOADING
299 heimbach 1.77 CADJ STORE pload = comlev1, key = ikey_dynamics,
300     CADJ & kind = isbyte
301     CADJ STORE siceload = comlev1, key = ikey_dynamics,
302     CADJ & kind = isbyte
303 heimbach 1.62 # endif
304     # ifdef NONLIN_FRSURF
305 heimbach 1.77 CADJ STORE recip_hfacc = comlev1, key = ikey_dynamics,
306     CADJ & kind = isbyte
307 heimbach 1.62 # endif
308 heimbach 1.78 # ifdef ANNUAL_BALANCE
309     CADJ STORE balance_itcount = comlev1, key = ikey_dynamics,
310     CADJ & kind = isbyte
311     # endif /* ANNUAL_BALANCE */
312 heimbach 1.126 # ifdef ALLOW_THSICE
313     CADJ STORE fu, fv = comlev1, key = ikey_dynamics,
314     CADJ & kind = isbyte
315     # endif
316     # endif /* ALLOW_AUTODIFF_TAMC */
317 heimbach 1.62 # ifdef ALLOW_DEBUG
318 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('SEAICE_MODEL',myThid)
319 heimbach 1.62 # endif
320 jmc 1.29 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 heimbach 1.62 # ifdef ALLOW_COST
324 heimbach 1.57 CALL SEAICE_COST_SENSI ( myTime, myIter, myThid )
325 heimbach 1.62 # endif
326 heimbach 1.35 ENDIF
327 jmc 1.29 #endif /* ALLOW_SEAICE */
328    
329 heimbach 1.64 #ifdef ALLOW_AUTODIFF_TAMC
330 heimbach 1.77 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 heimbach 1.64 # ifdef ALLOW_SEAICE
335 heimbach 1.77 CADJ STORE area = comlev1, key = ikey_dynamics,
336     CADJ & kind = isbyte
337 heimbach 1.64 # endif
338     #endif
339    
340 jmc 1.121 #ifdef OLD_THSICE_CALL_SEQUENCE
341 jscott 1.30 #if (defined ALLOW_THSICE) && !(defined ALLOW_ATM2D)
342 jmc 1.14 IF ( useThSIce .AND. fluidIsWater ) THEN
343 heimbach 1.101 # 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 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('THSICE_MAIN',myThid)
358 heimbach 1.101 # endif
359 jmc 1.5 C-- Step forward Therm.Sea-Ice variables
360     C and modify forcing terms including effects from ice
361     CALL TIMER_START('THSICE_MAIN [DO_OCEANIC_PHYS]', myThid)
362     CALL THSICE_MAIN( myTime, myIter, myThid )
363     CALL TIMER_STOP( 'THSICE_MAIN [DO_OCEANIC_PHYS]', myThid)
364     ENDIF
365     #endif /* ALLOW_THSICE */
366 jmc 1.121 #endif /* OLD_THSICE_CALL_SEQUENCE */
367 jmc 1.5
368 mlosch 1.21 #ifdef ALLOW_SHELFICE
369 heimbach 1.92 # ifdef ALLOW_AUTODIFF_TAMC
370     CADJ STORE salt, theta = comlev1, key = ikey_dynamics,
371     CADJ & kind = isbyte
372     # endif
373 mlosch 1.21 IF ( useShelfIce .AND. fluidIsWater ) THEN
374     #ifdef ALLOW_DEBUG
375 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('SHELFICE_THERMODYNAMICS',myThid)
376 mlosch 1.21 #endif
377 jmc 1.47 C compute temperature and (virtual) salt flux at the
378 mlosch 1.21 C shelf-ice ocean interface
379     CALL TIMER_START('SHELFICE_THERMODYNAMICS [DO_OCEANIC_PHYS]',
380     & myThid)
381     CALL SHELFICE_THERMODYNAMICS( myTime, myIter, myThid )
382     CALL TIMER_STOP( 'SHELFICE_THERMODYNAMICS [DO_OCEANIC_PHYS]',
383     & myThid)
384     ENDIF
385     #endif /* ALLOW_SHELFICE */
386    
387 dimitri 1.85 #ifdef ALLOW_ICEFRONT
388     IF ( useICEFRONT .AND. fluidIsWater ) THEN
389     #ifdef ALLOW_DEBUG
390 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('ICEFRONT_THERMODYNAMICS',myThid)
391 dimitri 1.85 #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 jmc 1.112 #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 jmc 1.5 C-- Freeze water at the surface
409 heimbach 1.104 IF ( allowFreezing ) THEN
410 jmc 1.5 #ifdef ALLOW_AUTODIFF_TAMC
411 heimbach 1.77 CADJ STORE theta = comlev1, key = ikey_dynamics,
412     CADJ & kind = isbyte
413 jmc 1.5 #endif
414     CALL FREEZE_SURFACE( myTime, myIter, myThid )
415     ENDIF
416    
417 jmc 1.28 #ifdef ALLOW_OCN_COMPON_INTERF
418 jmc 1.5 C-- Apply imported data (from coupled interface) to forcing fields
419 jmc 1.28 C jmc: do not know precisely where to put this call (bf or af thSIce ?)
420 jmc 1.36 IF ( useCoupler ) THEN
421 jmc 1.5 CALL OCN_APPLY_IMPORT( .TRUE., myTime, myIter, myThid )
422 jmc 1.36 ENDIF
423 jmc 1.28 #endif /* ALLOW_OCN_COMPON_INTERF */
424 jmc 1.5
425 jmc 1.127 C-- Balance and Apply exchanges to surface forcing
426     IF ( fluidIsWater ) THEN
427     CALL EXTERNAL_FORCING_ADJUST( myTime, myIter, myThid )
428     ENDIF
429 jmc 1.25
430 jmc 1.1 #ifdef ALLOW_AUTODIFF_TAMC
431     C-- HPF directive to help TAMC
432     CHPF$ INDEPENDENT
433 jmc 1.106 #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 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
438     DO bj=myByLo(myThid),myByHi(myThid)
439     #ifdef ALLOW_AUTODIFF_TAMC
440 heimbach 1.15 C-- HPF directive to help TAMC
441     CHPF$ INDEPENDENT
442 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
443     DO bi=myBxLo(myThid),myBxHi(myThid)
444    
445     #ifdef ALLOW_AUTODIFF_TAMC
446     act1 = bi - myBxLo(myThid)
447     max1 = myBxHi(myThid) - myBxLo(myThid) + 1
448     act2 = bj - myByLo(myThid)
449     max2 = myByHi(myThid) - myByLo(myThid) + 1
450     act3 = myThid - 1
451     max3 = nTx*nTy
452     act4 = ikey_dynamics - 1
453     itdkey = (act1 + 1) + act2*max1
454     & + act3*max1*max2
455     & + act4*max1*max2*max3
456 jmc 1.74 #endif /* ALLOW_AUTODIFF_TAMC */
457 jmc 1.1
458     C-- Set up work arrays with valid (i.e. not NaN) values
459     C These inital values do not alter the numerical results. They
460     C just ensure that all memory references are to valid floating
461     C point numbers. This prevents spurious hardware signals due to
462     C uninitialised but inert locations.
463    
464 jmc 1.69 #ifdef ALLOW_AUTODIFF_TAMC
465 jmc 1.1 DO j=1-OLy,sNy+OLy
466     DO i=1-OLx,sNx+OLx
467 jmc 1.69 rhoKm1 (i,j) = 0. _d 0
468 jmc 1.47 rhoKp1 (i,j) = 0. _d 0
469 jmc 1.1 ENDDO
470     ENDDO
471 jmc 1.69 #endif /* ALLOW_AUTODIFF_TAMC */
472 jmc 1.1
473     DO k=1,Nr
474     DO j=1-OLy,sNy+OLy
475     DO i=1-OLx,sNx+OLx
476 jmc 1.69 C This is currently used by GMRedi, IVDC, MXL-depth and Diagnostics
477 jmc 1.1 sigmaX(i,j,k) = 0. _d 0
478     sigmaY(i,j,k) = 0. _d 0
479     sigmaR(i,j,k) = 0. _d 0
480     #ifdef ALLOW_AUTODIFF_TAMC
481     cph all the following init. are necessary for TAF
482     cph although some of these are re-initialised later.
483 heimbach 1.109 rhoInSitu(i,j,k,bi,bj) = 0.
484 jmc 1.1 IVDConvCount(i,j,k,bi,bj) = 0.
485     # ifdef ALLOW_GMREDI
486     Kwx(i,j,k,bi,bj) = 0. _d 0
487     Kwy(i,j,k,bi,bj) = 0. _d 0
488     Kwz(i,j,k,bi,bj) = 0. _d 0
489     # ifdef GM_NON_UNITY_DIAGONAL
490     Kux(i,j,k,bi,bj) = 0. _d 0
491     Kvy(i,j,k,bi,bj) = 0. _d 0
492     # endif
493     # ifdef GM_EXTRA_DIAGONAL
494     Kuz(i,j,k,bi,bj) = 0. _d 0
495     Kvz(i,j,k,bi,bj) = 0. _d 0
496     # endif
497     # ifdef GM_BOLUS_ADVEC
498     GM_PsiX(i,j,k,bi,bj) = 0. _d 0
499     GM_PsiY(i,j,k,bi,bj) = 0. _d 0
500     # endif
501     # ifdef GM_VISBECK_VARIABLE_K
502     VisbeckK(i,j,bi,bj) = 0. _d 0
503     # endif
504     # endif /* ALLOW_GMREDI */
505 heimbach 1.42 # 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 gforget 1.91 # 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 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
515     ENDDO
516     ENDDO
517     ENDDO
518    
519     iMin = 1-OLx
520     iMax = sNx+OLx
521     jMin = 1-OLy
522     jMax = sNy+OLy
523    
524     #ifdef ALLOW_AUTODIFF_TAMC
525 jmc 1.96 CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
526 heimbach 1.77 CADJ & kind = isbyte
527 jmc 1.96 CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
528 heimbach 1.77 CADJ & kind = isbyte
529 heimbach 1.12 CADJ STORE totphihyd(:,:,:,bi,bj)
530 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
531 heimbach 1.77 CADJ & kind = isbyte
532 heimbach 1.10 # ifdef ALLOW_KPP
533 jmc 1.96 CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
534 heimbach 1.77 CADJ & kind = isbyte
535 jmc 1.96 CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
536 heimbach 1.77 CADJ & kind = isbyte
537 heimbach 1.10 # endif
538 heimbach 1.115 # ifdef ALLOW_SALT_PLUME
539     CADJ STORE saltplumedepth(:,:,bi,bj) = comlev1_bibj, key=itdkey,
540     CADJ & kind = isbyte
541     # endif
542 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
543    
544 jmc 1.71 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
547 jmc 1.110 IF (debugMode) CALL DEBUG_CALL('FIND_RHO_2D (xNr)',myThid)
548 jmc 1.71 #endif
549     #ifdef ALLOW_AUTODIFF_TAMC
550 jmc 1.110 IF ( fluidIsWater ) THEN
551 jmc 1.71 #endif /* ALLOW_AUTODIFF_TAMC */
552 jmc 1.69 #ifdef ALLOW_DOWN_SLOPE
553 jmc 1.110 IF ( useDOWN_SLOPE ) THEN
554     DO k=1,Nr
555 jmc 1.69 CALL DWNSLP_CALC_RHO(
556     I theta, salt,
557 jmc 1.71 O rhoInSitu(1-OLx,1-OLy,k,bi,bj),
558 jmc 1.69 I k, bi, bj, myTime, myIter, myThid )
559 jmc 1.110 ENDDO
560     ENDIF
561 jmc 1.71 #endif /* ALLOW_DOWN_SLOPE */
562 dimitri 1.107 #ifdef ALLOW_BBL
563 jmc 1.110 IF ( useBBL ) THEN
564 dimitri 1.108 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 jmc 1.110 DO k=Nr,1,-1
568 dimitri 1.107 CALL BBL_CALC_RHO(
569     I theta, salt,
570     O rhoInSitu,
571     I k, bi, bj, myTime, myIter, myThid )
572    
573 jmc 1.110 ENDDO
574     ENDIF
575 dimitri 1.107 #endif /* ALLOW_BBL */
576 jmc 1.110 IF ( .NOT. ( useDOWN_SLOPE .OR. useBBL ) ) THEN
577     DO k=1,Nr
578 jmc 1.71 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 jmc 1.110 ENDDO
585     ENDIF
586 jmc 1.74 #ifdef ALLOW_AUTODIFF_TAMC
587 jmc 1.110 ELSE
588 jmc 1.74 C- fluid is not water:
589 jmc 1.110 DO k=1,Nr
590 jmc 1.74 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 jmc 1.110 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
602    
603     C-- Start of diagnostic loop
604     DO k=Nr,1,-1
605    
606     #ifdef ALLOW_AUTODIFF_TAMC
607     C? Patrick, is this formula correct now that we change the loop range?
608     C? Do we still need this?
609     cph kkey formula corrected.
610     cph Needed for rhoK, rhoKm1, in the case useGMREDI.
611     kkey = (itdkey-1)*Nr + k
612 jmc 1.74 #endif /* ALLOW_AUTODIFF_TAMC */
613 jmc 1.69
614 jmc 1.110 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 jmc 1.1 C-- Calculate gradients of potential density for isoneutral
622     C slope terms (e.g. GM/Redi tensor or IVDC diffusivity)
623 jmc 1.17 IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.)
624 jmc 1.116 & .OR. usePP81 .OR. useMY82 .OR. useGGL90
625 dimitri 1.61 & .OR. useSALT_PLUME .OR. doDiagsRho.GE.1 ) THEN
626 jmc 1.1 IF (k.GT.1) THEN
627     #ifdef ALLOW_AUTODIFF_TAMC
628 jmc 1.96 CADJ STORE theta(:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey,
629 heimbach 1.77 CADJ & kind = isbyte
630 jmc 1.96 CADJ STORE salt (:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey,
631 heimbach 1.77 CADJ & kind = isbyte
632 jmc 1.96 CADJ STORE rhokm1 (bi,bj) = comlev1_bibj_k, key=kkey,
633 heimbach 1.77 CADJ & kind = isbyte
634 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
635 jmc 1.68 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 jmc 1.1 ENDIF
642     #ifdef ALLOW_DEBUG
643 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('GRAD_SIGMA',myThid)
644 jmc 1.1 #endif
645 heimbach 1.31 cph Avoid variable aliasing for adjoint !!!
646     DO j=jMin,jMax
647     DO i=iMin,iMax
648 jmc 1.71 rhoKp1(i,j) = rhoInSitu(i,j,k,bi,bj)
649 heimbach 1.31 ENDDO
650     ENDDO
651 jmc 1.1 CALL GRAD_SIGMA(
652     I bi, bj, iMin, iMax, jMin, jMax, k,
653 jmc 1.71 I rhoInSitu(1-OLx,1-OLy,k,bi,bj), rhoKm1, rhoKp1,
654 jmc 1.1 O sigmaX, sigmaY, sigmaR,
655     I myThid )
656 gforget 1.66 #ifdef ALLOW_AUTODIFF_TAMC
657 jmc 1.69 #ifdef GMREDI_WITH_STABLE_ADJOINT
658 gforget 1.66 cgf zero out adjoint fields to stabilize pkg/gmredi adjoint
659     cgf -> cuts adjoint dependency from slope to state
660 jmc 1.69 CALL ZERO_ADJ_LOC( Nr, sigmaX, myThid)
661     CALL ZERO_ADJ_LOC( Nr, sigmaY, myThid)
662     CALL ZERO_ADJ_LOC( Nr, sigmaR, myThid)
663 gforget 1.66 #endif
664     #endif /* ALLOW_AUTODIFF_TAMC */
665 jmc 1.1 ENDIF
666    
667     C-- Implicit Vertical Diffusion for Convection
668     IF (k.GT.1 .AND. ivdc_kappa.NE.0.) THEN
669     #ifdef ALLOW_DEBUG
670 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('CALC_IVDC',myThid)
671 jmc 1.1 #endif
672     CALL CALC_IVDC(
673     I bi, bj, iMin, iMax, jMin, jMax, k,
674 mlosch 1.111 I sigmaR,
675 jmc 1.1 I myTime, myIter, myThid)
676     ENDIF
677    
678 jmc 1.17 #ifdef ALLOW_DIAGNOSTICS
679 jmc 1.110 IF ( doDiagsRho.GE.4 ) THEN
680     CALL DIAGS_RHO_L( doDiagsRho, k, bi, bj,
681     I rhoInSitu(1-OLx,1-OLy,1,bi,bj),
682 jmc 1.74 I rhoKm1, wVel,
683 jmc 1.71 I myTime, myIter, myThid )
684 jmc 1.17 ENDIF
685     #endif
686    
687 jmc 1.1 C-- end of diagnostic k loop (Nr:1)
688     ENDDO
689    
690 heimbach 1.57 #ifdef ALLOW_AUTODIFF_TAMC
691 jmc 1.69 CADJ STORE IVDConvCount(:,:,:,bi,bj)
692 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
693 heimbach 1.77 CADJ & kind = isbyte
694 heimbach 1.57 #endif
695    
696 jmc 1.47 C-- Diagnose Mixed Layer Depth:
697 jmc 1.110 IF ( useGMRedi .OR. MOD(doDiagsRho,2).EQ.1 ) THEN
698 jmc 1.71 CALL CALC_OCE_MXLAYER(
699     I rhoInSitu(1-OLx,1-OLy,1,bi,bj), sigmaR,
700     I bi, bj, myTime, myIter, myThid )
701 jmc 1.47 ENDIF
702 heimbach 1.53
703 dimitri 1.52 #ifdef ALLOW_SALT_PLUME
704 dimitri 1.61 IF ( useSALT_PLUME ) THEN
705 jmc 1.71 CALL SALT_PLUME_CALC_DEPTH(
706     I rhoInSitu(1-OLx,1-OLy,1,bi,bj), sigmaR,
707     I bi, bj, myTime, myIter, myThid )
708 dimitri 1.60 ENDIF
709 dimitri 1.61 #endif /* ALLOW_SALT_PLUME */
710    
711 jmc 1.8 #ifdef ALLOW_DIAGNOSTICS
712 jmc 1.71 IF ( MOD(doDiagsRho,4).GE.2 ) THEN
713 jmc 1.16 CALL DIAGNOSTICS_FILL (sigmaR, 'DRHODR ', 0, Nr,
714     & 2, bi, bj, myThid)
715 jmc 1.8 ENDIF
716 dimitri 1.61 #endif /* ALLOW_DIAGNOSTICS */
717 jmc 1.8
718 jmc 1.1 C-- Determines forcing terms based on external fields
719     C relaxation terms, etc.
720     #ifdef ALLOW_DEBUG
721 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('EXTERNAL_FORCING_SURF',myThid)
722 jmc 1.1 #endif
723 heimbach 1.23 #ifdef ALLOW_AUTODIFF_TAMC
724     CADJ STORE EmPmR(:,:,bi,bj)
725 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
726 heimbach 1.77 CADJ & kind = isbyte
727 heimbach 1.26 # ifdef EXACT_CONSERV
728 heimbach 1.23 CADJ STORE PmEpR(:,:,bi,bj)
729 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
730 heimbach 1.77 CADJ & kind = isbyte
731 heimbach 1.26 # endif
732 heimbach 1.27 # ifdef NONLIN_FRSURF
733     CADJ STORE hFac_surfC(:,:,bi,bj)
734 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
735 heimbach 1.77 CADJ & kind = isbyte
736 heimbach 1.27 CADJ STORE recip_hFacC(:,:,:,bi,bj)
737 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
738 heimbach 1.77 CADJ & kind = isbyte
739 heimbach 1.97 # 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 heimbach 1.23 #endif
747 jmc 1.36 CALL EXTERNAL_FORCING_SURF(
748 jmc 1.1 I bi, bj, iMin, iMax, jMin, jMax,
749     I myTime, myIter, myThid )
750 heimbach 1.27 #ifdef ALLOW_AUTODIFF_TAMC
751     # ifdef EXACT_CONSERV
752     cph-test
753     cphCADJ STORE PmEpR(:,:,bi,bj)
754 jmc 1.96 cphCADJ & = comlev1_bibj, key=itdkey,
755 heimbach 1.77 cphCADJ & kind = isbyte
756 heimbach 1.27 # endif
757     #endif
758 jmc 1.1
759     #ifdef ALLOW_AUTODIFF_TAMC
760     cph needed for KPP
761 jmc 1.4 CADJ STORE surfaceForcingU(:,:,bi,bj)
762 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
763 heimbach 1.77 CADJ & kind = isbyte
764 jmc 1.4 CADJ STORE surfaceForcingV(:,:,bi,bj)
765 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
766 heimbach 1.77 CADJ & kind = isbyte
767 jmc 1.4 CADJ STORE surfaceForcingS(:,:,bi,bj)
768 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
769 heimbach 1.77 CADJ & kind = isbyte
770 jmc 1.4 CADJ STORE surfaceForcingT(:,:,bi,bj)
771 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
772 heimbach 1.77 CADJ & kind = isbyte
773 jmc 1.4 CADJ STORE surfaceForcingTice(:,:,bi,bj)
774 jmc 1.96 CADJ & = comlev1_bibj, key=itdkey,
775 heimbach 1.77 CADJ & kind = isbyte
776 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
777    
778     #ifdef ALLOW_KPP
779     C-- Compute KPP mixing coefficients
780     IF (useKPP) THEN
781     #ifdef ALLOW_DEBUG
782 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('KPP_CALC',myThid)
783 jmc 1.1 #endif
784 dfer 1.76 CALL TIMER_START('KPP_CALC [DO_OCEANIC_PHYS]', myThid)
785 jmc 1.1 CALL KPP_CALC(
786 jmc 1.44 I bi, bj, myTime, myIter, myThid )
787 dfer 1.76 CALL TIMER_STOP ('KPP_CALC [DO_OCEANIC_PHYS]', myThid)
788 jmc 1.1 #ifdef ALLOW_AUTODIFF_TAMC
789     ELSE
790     CALL KPP_CALC_DUMMY(
791 jmc 1.44 I bi, bj, myTime, myIter, myThid )
792 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
793     ENDIF
794    
795     #endif /* ALLOW_KPP */
796    
797 mlosch 1.6 #ifdef ALLOW_PP81
798     C-- Compute PP81 mixing coefficients
799     IF (usePP81) THEN
800     #ifdef ALLOW_DEBUG
801 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('PP81_CALC',myThid)
802 mlosch 1.6 #endif
803     CALL PP81_CALC(
804 jmc 1.116 I bi, bj, sigmaR, myTime, myIter, myThid )
805 mlosch 1.6 ENDIF
806     #endif /* ALLOW_PP81 */
807    
808     #ifdef ALLOW_MY82
809     C-- Compute MY82 mixing coefficients
810     IF (useMY82) THEN
811     #ifdef ALLOW_DEBUG
812 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('MY82_CALC',myThid)
813 mlosch 1.6 #endif
814     CALL MY82_CALC(
815 jmc 1.116 I bi, bj, sigmaR, myTime, myIter, myThid )
816 mlosch 1.6 ENDIF
817     #endif /* ALLOW_MY82 */
818    
819 mlosch 1.9 #ifdef ALLOW_GGL90
820 gforget 1.91 #ifdef ALLOW_AUTODIFF_TAMC
821     CADJ STORE GGL90TKE (:,:,:,bi,bj) = comlev1_bibj, key=itdkey,
822     CADJ & kind = isbyte
823     #endif /* ALLOW_AUTODIFF_TAMC */
824 mlosch 1.9 C-- Compute GGL90 mixing coefficients
825     IF (useGGL90) THEN
826     #ifdef ALLOW_DEBUG
827 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('GGL90_CALC',myThid)
828 mlosch 1.9 #endif
829 dfer 1.76 CALL TIMER_START('GGL90_CALC [DO_OCEANIC_PHYS]', myThid)
830 mlosch 1.9 CALL GGL90_CALC(
831 jmc 1.116 I bi, bj, sigmaR, myTime, myIter, myThid )
832 dfer 1.76 CALL TIMER_STOP ('GGL90_CALC [DO_OCEANIC_PHYS]', myThid)
833 mlosch 1.9 ENDIF
834     #endif /* ALLOW_GGL90 */
835    
836 jmc 1.20 #ifdef ALLOW_TIMEAVE
837 jmc 1.36 IF ( taveFreq.GT. 0. _d 0 ) THEN
838 jmc 1.20 CALL TIMEAVE_SURF_FLUX( bi, bj, myTime, myIter, myThid)
839     ENDIF
840     IF (taveFreq.GT.0. .AND. ivdc_kappa.NE.0.) THEN
841     CALL TIMEAVE_CUMULATE(ConvectCountTave, IVDConvCount,
842 jmc 1.120 I Nr, deltaTClock, bi, bj, myThid)
843 jmc 1.20 ENDIF
844     #endif /* ALLOW_TIMEAVE */
845    
846 jmc 1.69 #ifdef ALLOW_GMREDI
847 jmc 1.47 #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 jmc 1.86 cph but I have avoided the #ifdef for now, in case more things change
853 jmc 1.96 CADJ STORE sigmaX(:,:,:) = comlev1_bibj, key=itdkey,
854 heimbach 1.77 CADJ & kind = isbyte
855 jmc 1.96 CADJ STORE sigmaY(:,:,:) = comlev1_bibj, key=itdkey,
856 heimbach 1.77 CADJ & kind = isbyte
857 jmc 1.96 CADJ STORE sigmaR(:,:,:) = comlev1_bibj, key=itdkey,
858 heimbach 1.77 CADJ & kind = isbyte
859 jmc 1.47 # 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 jmc 1.96 IF (debugMode) CALL DEBUG_CALL('GMREDI_CALC_TENSOR',myThid)
866 jmc 1.47 #endif
867     CALL GMREDI_CALC_TENSOR(
868 jmc 1.51 I iMin, iMax, jMin, jMax,
869 jmc 1.47 I sigmaX, sigmaY, sigmaR,
870 jmc 1.51 I bi, bj, myTime, myIter, myThid )
871 jmc 1.47 #ifdef ALLOW_AUTODIFF_TAMC
872     ELSE
873     CALL GMREDI_CALC_TENSOR_DUMMY(
874 jmc 1.51 I iMin, iMax, jMin, jMax,
875 jmc 1.47 I sigmaX, sigmaY, sigmaR,
876 jmc 1.51 I bi, bj, myTime, myIter, myThid )
877 jmc 1.47 #endif /* ALLOW_AUTODIFF_TAMC */
878     ENDIF
879 jmc 1.69 #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 jmc 1.71 I bi, bj, kSurfC, rhoInSitu,
887 jmc 1.69 I myTime, myIter, myThid )
888     ELSE
889     CALL DWNSLP_CALC_FLOW(
890 jmc 1.71 I bi, bj, kLowC, rhoInSitu,
891 jmc 1.69 I myTime, myIter, myThid )
892     ENDIF
893     ENDIF
894     #endif /* ALLOW_DOWN_SLOPE */
895 jmc 1.47
896 jmc 1.106 C-- end bi,bj loops.
897     ENDDO
898     ENDDO
899    
900 gforget 1.117 #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 gforget 1.119 CALL BALANCE_RELAX( myTime, myIter, myThid )
910 gforget 1.117 #endif /* ALLOW_BALANCE_RELAX */
911    
912 jmc 1.98 #ifndef ALLOW_AUTODIFF_TAMC
913     C--- if fluid Is Water: end
914 jmc 1.106 ENDIF
915 jmc 1.98 #endif
916    
917 dimitri 1.107 #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 dimitri 1.94 #ifdef ALLOW_MYPACKAGE
925 jmc 1.106 IF ( useMYPACKAGE ) THEN
926     CALL MYPACKAGE_CALC_RHS(
927     I myTime, myIter, myThid )
928     ENDIF
929 dimitri 1.94 #endif /* ALLOW_MYPACKAGE */
930    
931 jmc 1.99 #ifdef ALLOW_GMREDI
932     IF ( useGMRedi ) THEN
933     CALL GMREDI_DO_EXCH( myTime, myIter, myThid )
934     ENDIF
935     #endif /* ALLOW_GMREDI */
936    
937 jmc 1.98 #ifdef ALLOW_KPP
938 jmc 1.45 IF (useKPP) THEN
939     CALL KPP_DO_EXCH( myThid )
940     ENDIF
941 jmc 1.98 #endif /* ALLOW_KPP */
942 jmc 1.45
943 jmc 1.18 #ifdef ALLOW_DIAGNOSTICS
944     IF ( fluidIsWater .AND. useDiagnostics ) THEN
945 jmc 1.74 CALL DIAGS_RHO_G(
946 jmc 1.110 I rhoInSitu, uVel, vVel, wVel,
947 jmc 1.71 I myTime, myIter, myThid )
948 jmc 1.18 CALL DIAGS_OCEANIC_SURF_FLUX( myTime, myIter, myThid )
949     ENDIF
950 jmc 1.19 IF ( ivdc_kappa.NE.0 .AND. useDiagnostics ) THEN
951 jmc 1.71 CALL DIAGNOSTICS_FILL( IVDConvCount, 'CONVADJ ',
952     & 0, Nr, 0, 1, 1, myThid )
953 jmc 1.19 ENDIF
954 jmc 1.18 #endif
955    
956 gforget 1.123 #ifdef ALLOW_ECCO
957     CALL ECCO_PHYS(mythid)
958     #endif
959    
960 jmc 1.1 #ifdef ALLOW_DEBUG
961 jmc 1.96 IF (debugMode) CALL DEBUG_LEAVE('DO_OCEANIC_PHYS',myThid)
962 jmc 1.1 #endif
963    
964     RETURN
965     END
  ViewVC Help
Powered by ViewVC 1.1.22