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

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

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

Revision 1.71 - (hide annotations) (download)
Fri Jul 2 15:50:24 2004 UTC (19 years, 11 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint54, checkpoint54a_pre, checkpoint53g_post
Changes since 1.70: +1 -3 lines 
Keeping up with JMC's latest modifs.
1 heimbach 1.71 C $Header: /u/gcmpack/MITgcm/model/src/thermodynamics.F,v 1.70 2004/07/01 21:45:11 jmc Exp $
2 adcroft 1.1 C $Name: $
3    
4 edhill 1.51 #include "PACKAGES_CONFIG.h"
5 adcroft 1.1 #include "CPP_OPTIONS.h"
6 jmc 1.60 #ifdef ALLOW_PTRACERS
7     # include "PTRACERS_OPTIONS.h"
8     #endif
9 edhill 1.51
10 jmc 1.21 #ifdef ALLOW_AUTODIFF_TAMC
11     # ifdef ALLOW_GMREDI
12     # include "GMREDI_OPTIONS.h"
13     # endif
14     # ifdef ALLOW_KPP
15     # include "KPP_OPTIONS.h"
16 heimbach 1.42 # endif
17 jmc 1.21 #endif /* ALLOW_AUTODIFF_TAMC */
18 adcroft 1.1
19 cnh 1.9 CBOP
20     C !ROUTINE: THERMODYNAMICS
21     C !INTERFACE:
22 adcroft 1.1 SUBROUTINE THERMODYNAMICS(myTime, myIter, myThid)
23 cnh 1.9 C !DESCRIPTION: \bv
24     C *==========================================================*
25     C | SUBROUTINE THERMODYNAMICS
26     C | o Controlling routine for the prognostic part of the
27     C | thermo-dynamics.
28     C *===========================================================
29     C | The algorithm...
30     C |
31     C | "Correction Step"
32     C | =================
33     C | Here we update the horizontal velocities with the surface
34     C | pressure such that the resulting flow is either consistent
35     C | with the free-surface evolution or the rigid-lid:
36     C | U[n] = U* + dt x d/dx P
37     C | V[n] = V* + dt x d/dy P
38     C |
39     C | "Calculation of Gs"
40     C | ===================
41     C | This is where all the accelerations and tendencies (ie.
42     C | physics, parameterizations etc...) are calculated
43     C | rho = rho ( theta[n], salt[n] )
44     C | b = b(rho, theta)
45     C | K31 = K31 ( rho )
46     C | Gu[n] = Gu( u[n], v[n], wVel, b, ... )
47     C | Gv[n] = Gv( u[n], v[n], wVel, b, ... )
48     C | Gt[n] = Gt( theta[n], u[n], v[n], wVel, K31, ... )
49     C | Gs[n] = Gs( salt[n], u[n], v[n], wVel, K31, ... )
50     C |
51     C | "Time-stepping" or "Prediction"
52     C | ================================
53     C | The models variables are stepped forward with the appropriate
54     C | time-stepping scheme (currently we use Adams-Bashforth II)
55     C | - For momentum, the result is always *only* a "prediction"
56     C | in that the flow may be divergent and will be "corrected"
57     C | later with a surface pressure gradient.
58     C | - Normally for tracers the result is the new field at time
59     C | level [n+1} *BUT* in the case of implicit diffusion the result
60     C | is also *only* a prediction.
61     C | - We denote "predictors" with an asterisk (*).
62     C | U* = U[n] + dt x ( 3/2 Gu[n] - 1/2 Gu[n-1] )
63     C | V* = V[n] + dt x ( 3/2 Gv[n] - 1/2 Gv[n-1] )
64     C | theta[n+1] = theta[n] + dt x ( 3/2 Gt[n] - 1/2 atG[n-1] )
65     C | salt[n+1] = salt[n] + dt x ( 3/2 Gt[n] - 1/2 atG[n-1] )
66     C | With implicit diffusion:
67     C | theta* = theta[n] + dt x ( 3/2 Gt[n] - 1/2 atG[n-1] )
68     C | salt* = salt[n] + dt x ( 3/2 Gt[n] - 1/2 atG[n-1] )
69     C | (1 + dt * K * d_zz) theta[n] = theta*
70     C | (1 + dt * K * d_zz) salt[n] = salt*
71     C |
72     C *==========================================================*
73     C \ev
74    
75     C !USES:
76 adcroft 1.1 IMPLICIT NONE
77     C == Global variables ===
78     #include "SIZE.h"
79     #include "EEPARAMS.h"
80     #include "PARAMS.h"
81     #include "DYNVARS.h"
82     #include "GRID.h"
83 adcroft 1.4 #include "GAD.h"
84 adcroft 1.1 #ifdef ALLOW_PASSIVE_TRACER
85     #include "TR1.h"
86     #endif
87 jmc 1.45 #ifdef ALLOW_PTRACERS
88     #include "PTRACERS.h"
89     #endif
90 heimbach 1.42 #ifdef ALLOW_TIMEAVE
91     #include "TIMEAVE_STATV.h"
92     #endif
93    
94 adcroft 1.1 #ifdef ALLOW_AUTODIFF_TAMC
95     # include "tamc.h"
96     # include "tamc_keys.h"
97     # include "FFIELDS.h"
98 heimbach 1.30 # include "EOS.h"
99 adcroft 1.1 # ifdef ALLOW_KPP
100     # include "KPP.h"
101 heimbach 1.67 # endif
102     # ifdef ALLOW_GMREDI
103     # include "GMREDI.h"
104 adcroft 1.1 # endif
105 heimbach 1.68 # ifdef ALLOW_EBM
106     # include "EBM.h"
107 adcroft 1.1 # endif
108     #endif /* ALLOW_AUTODIFF_TAMC */
109    
110 cnh 1.9 C !INPUT/OUTPUT PARAMETERS:
111 adcroft 1.1 C == Routine arguments ==
112     C myTime - Current time in simulation
113     C myIter - Current iteration number in simulation
114     C myThid - Thread number for this instance of the routine.
115     _RL myTime
116     INTEGER myIter
117     INTEGER myThid
118    
119 cnh 1.9 C !LOCAL VARIABLES:
120 adcroft 1.1 C == Local variables
121     C xA, yA - Per block temporaries holding face areas
122     C uTrans, vTrans, rTrans - Per block temporaries holding flow
123     C transport
124     C o uTrans: Zonal transport
125     C o vTrans: Meridional transport
126     C o rTrans: Vertical transport
127 jmc 1.64 C rTransKp1 o vertical volume transp. at interface k+1
128 adcroft 1.1 C maskUp o maskUp: land/water mask for W points
129     C fVer[STUV] o fVer: Vertical flux term - note fVer
130     C is "pipelined" in the vertical
131     C so we need an fVer for each
132     C variable.
133     C rhoK, rhoKM1 - Density at current level, and level above
134     C KappaRT, - Total diffusion in vertical for T and S.
135     C KappaRS (background + spatially varying, isopycnal term).
136 jmc 1.39 C useVariableK = T when vertical diffusion is not constant
137 adcroft 1.1 C iMin, iMax - Ranges and sub-block indices on which calculations
138     C jMin, jMax are applied.
139     C bi, bj
140     C k, kup, - Index for layer above and below. kup and kDown
141     C kDown, km1 are switched with layer to be the appropriate
142     C index into fVerTerm.
143     _RS xA (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
144     _RS yA (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
145     _RL uTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
146     _RL vTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
147     _RL rTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
148 jmc 1.64 _RL rTransKp1(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
149 adcroft 1.1 _RS maskUp (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
150     _RL fVerT (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
151     _RL fVerS (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
152 heimbach 1.55 #ifdef ALLOW_PASSIVE_TRACER
153 adcroft 1.1 _RL fVerTr1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
154 heimbach 1.55 #endif
155     #ifdef ALLOW_PTRACERS
156     _RL fVerP (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2,PTRACERS_num)
157     #endif
158 adcroft 1.1 _RL rhokm1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
159     _RL rhok (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
160     _RL KappaRT (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
161     _RL KappaRS (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
162     _RL sigmaX (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
163     _RL sigmaY (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
164     _RL sigmaR (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
165 jmc 1.64 _RL kp1Msk
166 jmc 1.39 LOGICAL useVariableK
167 adcroft 1.1 INTEGER iMin, iMax
168     INTEGER jMin, jMax
169     INTEGER bi, bj
170     INTEGER i, j
171     INTEGER k, km1, kup, kDown
172 heimbach 1.55 INTEGER iTracer, ip
173 adcroft 1.1
174 cnh 1.9 CEOP
175 adcroft 1.40
176 edhill 1.56 #ifdef ALLOW_DEBUG
177 heimbach 1.43 IF ( debugLevel .GE. debLevB )
178 jmc 1.63 & CALL DEBUG_ENTER('THERMODYNAMICS',myThid)
179 adcroft 1.40 #endif
180 adcroft 1.1
181     #ifdef ALLOW_AUTODIFF_TAMC
182     C-- dummy statement to end declaration part
183     ikey = 1
184 heimbach 1.30 itdkey = 1
185 adcroft 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
186    
187     #ifdef ALLOW_AUTODIFF_TAMC
188     C-- HPF directive to help TAMC
189     CHPF$ INDEPENDENT
190     #endif /* ALLOW_AUTODIFF_TAMC */
191    
192     DO bj=myByLo(myThid),myByHi(myThid)
193    
194     #ifdef ALLOW_AUTODIFF_TAMC
195     C-- HPF directive to help TAMC
196 heimbach 1.2 CHPF$ INDEPENDENT, NEW (rTrans,fVerT,fVerS
197 jmc 1.37 CHPF$& ,utrans,vtrans,xA,yA
198 heimbach 1.2 CHPF$& ,KappaRT,KappaRS
199 adcroft 1.1 CHPF$& )
200     #endif /* ALLOW_AUTODIFF_TAMC */
201    
202     DO bi=myBxLo(myThid),myBxHi(myThid)
203    
204     #ifdef ALLOW_AUTODIFF_TAMC
205     act1 = bi - myBxLo(myThid)
206     max1 = myBxHi(myThid) - myBxLo(myThid) + 1
207     act2 = bj - myByLo(myThid)
208     max2 = myByHi(myThid) - myByLo(myThid) + 1
209     act3 = myThid - 1
210     max3 = nTx*nTy
211     act4 = ikey_dynamics - 1
212 heimbach 1.30 itdkey = (act1 + 1) + act2*max1
213 adcroft 1.1 & + act3*max1*max2
214     & + act4*max1*max2*max3
215     #endif /* ALLOW_AUTODIFF_TAMC */
216    
217 heimbach 1.41 C-- Set up work arrays with valid (i.e. not NaN) values
218     C These inital values do not alter the numerical results. They
219     C just ensure that all memory references are to valid floating
220     C point numbers. This prevents spurious hardware signals due to
221     C uninitialised but inert locations.
222    
223 adcroft 1.1 DO j=1-OLy,sNy+OLy
224     DO i=1-OLx,sNx+OLx
225 heimbach 1.41 xA(i,j) = 0. _d 0
226     yA(i,j) = 0. _d 0
227     uTrans(i,j) = 0. _d 0
228     vTrans(i,j) = 0. _d 0
229     rhok (i,j) = 0. _d 0
230     rhoKM1 (i,j) = 0. _d 0
231 adcroft 1.1 rTrans (i,j) = 0. _d 0
232 jmc 1.64 rTransKp1(i,j) = 0. _d 0
233 adcroft 1.1 fVerT (i,j,1) = 0. _d 0
234     fVerT (i,j,2) = 0. _d 0
235     fVerS (i,j,1) = 0. _d 0
236     fVerS (i,j,2) = 0. _d 0
237 heimbach 1.55 #ifdef ALLOW_PASSIVE_TRACER
238 adcroft 1.1 fVerTr1(i,j,1) = 0. _d 0
239     fVerTr1(i,j,2) = 0. _d 0
240 heimbach 1.55 #endif
241     #ifdef ALLOW_PTRACERS
242     DO ip=1,PTRACERS_num
243     fVerP (i,j,1,ip) = 0. _d 0
244     fVerP (i,j,2,ip) = 0. _d 0
245     ENDDO
246     #endif
247 adcroft 1.1 ENDDO
248     ENDDO
249    
250     DO k=1,Nr
251     DO j=1-OLy,sNy+OLy
252     DO i=1-OLx,sNx+OLx
253     C This is currently also used by IVDC and Diagnostics
254 heimbach 1.20 sigmaX(i,j,k) = 0. _d 0
255     sigmaY(i,j,k) = 0. _d 0
256     sigmaR(i,j,k) = 0. _d 0
257 heimbach 1.30 KappaRT(i,j,k) = 0. _d 0
258     KappaRS(i,j,k) = 0. _d 0
259 jmc 1.45 C- tracer tendency needs to be set to zero (moved here from gad_calc_rhs):
260 heimbach 1.30 gT(i,j,k,bi,bj) = 0. _d 0
261     gS(i,j,k,bi,bj) = 0. _d 0
262     # ifdef ALLOW_PASSIVE_TRACER
263 edhill 1.51 ceh3 needs an IF ( use PASSIVE_TRACER) THEN
264 heimbach 1.5 gTr1(i,j,k,bi,bj) = 0. _d 0
265 heimbach 1.30 # endif
266 heimbach 1.42 # ifdef ALLOW_PTRACERS
267 edhill 1.51 ceh3 this should have an IF ( usePTRACERS ) THEN
268 heimbach 1.42 DO iTracer=1,PTRACERS_numInUse
269     gPTr(i,j,k,bi,bj,itracer) = 0. _d 0
270     ENDDO
271     # endif
272 jmc 1.45 #ifdef ALLOW_AUTODIFF_TAMC
273     cph all the following init. are necessary for TAF
274     cph although some of these are re-initialised later.
275 heimbach 1.30 # ifdef ALLOW_GMREDI
276     Kwx(i,j,k,bi,bj) = 0. _d 0
277     Kwy(i,j,k,bi,bj) = 0. _d 0
278     Kwz(i,j,k,bi,bj) = 0. _d 0
279     # ifdef GM_NON_UNITY_DIAGONAL
280     Kux(i,j,k,bi,bj) = 0. _d 0
281     Kvy(i,j,k,bi,bj) = 0. _d 0
282     # endif
283     # ifdef GM_EXTRA_DIAGONAL
284     Kuz(i,j,k,bi,bj) = 0. _d 0
285     Kvz(i,j,k,bi,bj) = 0. _d 0
286     # endif
287     # ifdef GM_BOLUS_ADVEC
288     GM_PsiX(i,j,k,bi,bj) = 0. _d 0
289     GM_PsiY(i,j,k,bi,bj) = 0. _d 0
290 heimbach 1.36 # endif
291     # ifdef GM_VISBECK_VARIABLE_K
292     VisbeckK(i,j,bi,bj) = 0. _d 0
293 heimbach 1.30 # endif
294     # endif /* ALLOW_GMREDI */
295     #endif /* ALLOW_AUTODIFF_TAMC */
296 adcroft 1.1 ENDDO
297     ENDDO
298     ENDDO
299    
300 jmc 1.21 iMin = 1-OLx
301 adcroft 1.1 iMax = sNx+OLx
302 jmc 1.21 jMin = 1-OLy
303 adcroft 1.1 jMax = sNy+OLy
304    
305     #ifdef ALLOW_AUTODIFF_TAMC
306 heimbach 1.30 CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
307     CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
308 heimbach 1.38 CADJ STORE totphihyd
309     CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
310 heimbach 1.11 #ifdef ALLOW_KPP
311 heimbach 1.30 CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
312     CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
313 heimbach 1.11 #endif
314 adcroft 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
315    
316 edhill 1.56 #ifdef ALLOW_DEBUG
317 heimbach 1.43 IF ( debugLevel .GE. debLevB )
318     & CALL DEBUG_MSG('ENTERING UPWARD K LOOP',myThid)
319 adcroft 1.40 #endif
320    
321 adcroft 1.1 C-- Start of diagnostic loop
322     DO k=Nr,1,-1
323    
324     #ifdef ALLOW_AUTODIFF_TAMC
325     C? Patrick, is this formula correct now that we change the loop range?
326     C? Do we still need this?
327     cph kkey formula corrected.
328     cph Needed for rhok, rhokm1, in the case useGMREDI.
329 heimbach 1.30 kkey = (itdkey-1)*Nr + k
330 adcroft 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
331    
332     C-- Integrate continuity vertically for vertical velocity
333 jmc 1.27 c CALL INTEGRATE_FOR_W(
334     c I bi, bj, k, uVel, vVel,
335     c O wVel,
336     c I myThid )
337 adcroft 1.1
338     #ifdef ALLOW_OBCS
339     #ifdef ALLOW_NONHYDROSTATIC
340     C-- Apply OBC to W if in N-H mode
341 jmc 1.27 c IF (useOBCS.AND.nonHydrostatic) THEN
342     c CALL OBCS_APPLY_W( bi, bj, k, wVel, myThid )
343     c ENDIF
344 adcroft 1.1 #endif /* ALLOW_NONHYDROSTATIC */
345     #endif /* ALLOW_OBCS */
346    
347 heimbach 1.22 C-- Attention: by defining "SINGLE_LAYER_MODE" in CPP_OPTIONS.h
348     C-- MOST of THERMODYNAMICS will be disabled
349     #ifndef SINGLE_LAYER_MODE
350    
351 adcroft 1.1 C-- Calculate gradients of potential density for isoneutral
352     C slope terms (e.g. GM/Redi tensor or IVDC diffusivity)
353     c IF ( k.GT.1 .AND. (useGMRedi.OR.ivdc_kappa.NE.0.) ) THEN
354     IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.) ) THEN
355 edhill 1.56 #ifdef ALLOW_DEBUG
356 heimbach 1.43 IF ( debugLevel .GE. debLevB )
357     & CALL DEBUG_CALL('FIND_RHO',myThid)
358 adcroft 1.40 #endif
359 adcroft 1.1 #ifdef ALLOW_AUTODIFF_TAMC
360     CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
361     CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
362     #endif /* ALLOW_AUTODIFF_TAMC */
363     CALL FIND_RHO(
364 mlosch 1.26 I bi, bj, iMin, iMax, jMin, jMax, k, k,
365 adcroft 1.1 I theta, salt,
366     O rhoK,
367     I myThid )
368 heimbach 1.30
369 adcroft 1.1 IF (k.GT.1) THEN
370     #ifdef ALLOW_AUTODIFF_TAMC
371     CADJ STORE theta(:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
372     CADJ STORE salt (:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
373     #endif /* ALLOW_AUTODIFF_TAMC */
374     CALL FIND_RHO(
375 mlosch 1.26 I bi, bj, iMin, iMax, jMin, jMax, k-1, k,
376 adcroft 1.1 I theta, salt,
377     O rhoKm1,
378     I myThid )
379     ENDIF
380 edhill 1.56 #ifdef ALLOW_DEBUG
381 heimbach 1.43 IF ( debugLevel .GE. debLevB )
382     & CALL DEBUG_CALL('GRAD_SIGMA',myThid)
383 adcroft 1.40 #endif
384 adcroft 1.1 CALL GRAD_SIGMA(
385     I bi, bj, iMin, iMax, jMin, jMax, k,
386     I rhoK, rhoKm1, rhoK,
387     O sigmaX, sigmaY, sigmaR,
388     I myThid )
389     ENDIF
390    
391 heimbach 1.30 #ifdef ALLOW_AUTODIFF_TAMC
392     CADJ STORE rhok (:,:) = comlev1_bibj_k , key=kkey, byte=isbyte
393     CADJ STORE rhokm1 (:,:) = comlev1_bibj_k , key=kkey, byte=isbyte
394     #endif /* ALLOW_AUTODIFF_TAMC */
395 adcroft 1.1 C-- Implicit Vertical Diffusion for Convection
396     c ==> should use sigmaR !!!
397     IF (k.GT.1 .AND. ivdc_kappa.NE.0.) THEN
398 edhill 1.56 #ifdef ALLOW_DEBUG
399 heimbach 1.43 IF ( debugLevel .GE. debLevB )
400     & CALL DEBUG_CALL('CALC_IVDC',myThid)
401 adcroft 1.40 #endif
402 adcroft 1.1 CALL CALC_IVDC(
403     I bi, bj, iMin, iMax, jMin, jMax, k,
404     I rhoKm1, rhoK,
405     I myTime, myIter, myThid)
406     ENDIF
407    
408 heimbach 1.22 #endif /* SINGLE_LAYER_MODE */
409    
410 adcroft 1.1 C-- end of diagnostic k loop (Nr:1)
411     ENDDO
412    
413     #ifdef ALLOW_AUTODIFF_TAMC
414     cph avoids recomputation of integrate_for_w
415 heimbach 1.30 CADJ STORE wvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
416 adcroft 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
417    
418     #ifdef ALLOW_OBCS
419     C-- Calculate future values on open boundaries
420     IF (useOBCS) THEN
421 edhill 1.56 #ifdef ALLOW_DEBUG
422 heimbach 1.43 IF ( debugLevel .GE. debLevB )
423     & CALL DEBUG_CALL('OBCS_CALC',myThid)
424 adcroft 1.40 #endif
425 jmc 1.16 CALL OBCS_CALC( bi, bj, myTime+deltaT, myIter+1,
426 adcroft 1.1 I uVel, vVel, wVel, theta, salt,
427     I myThid )
428     ENDIF
429     #endif /* ALLOW_OBCS */
430    
431 heimbach 1.66 #ifndef ALLOW_AUTODIFF_TAMC
432 jmc 1.62 IF ( buoyancyRelation(1:7) .EQ. 'OCEANIC' ) THEN
433 heimbach 1.66 #endif
434 edhill 1.54 C-- Determines forcing terms based on external fields
435     C relaxation terms, etc.
436 edhill 1.56 #ifdef ALLOW_DEBUG
437 edhill 1.54 IF ( debugLevel .GE. debLevB )
438     & CALL DEBUG_CALL('EXTERNAL_FORCING_SURF',myThid)
439     #endif
440 heimbach 1.68 CALL EXTERNAL_FORCING_SURF(
441 edhill 1.54 I bi, bj, iMin, iMax, jMin, jMax,
442     I myTime, myIter, myThid )
443 heimbach 1.66 #ifndef ALLOW_AUTODIFF_TAMC
444     ENDIF
445     #endif
446 edhill 1.54
447 adcroft 1.1 #ifdef ALLOW_AUTODIFF_TAMC
448     cph needed for KPP
449     CADJ STORE surfacetendencyU(:,:,bi,bj)
450 heimbach 1.30 CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
451 adcroft 1.1 CADJ STORE surfacetendencyV(:,:,bi,bj)
452 heimbach 1.30 CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
453 adcroft 1.1 CADJ STORE surfacetendencyS(:,:,bi,bj)
454 heimbach 1.30 CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
455 adcroft 1.1 CADJ STORE surfacetendencyT(:,:,bi,bj)
456 heimbach 1.30 CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
457 heimbach 1.57 # ifdef ALLOW_SEAICE
458     CADJ STORE surfacetendencyTice(:,:,bi,bj)
459     CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
460     # endif
461 adcroft 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
462    
463 heimbach 1.22 C-- Attention: by defining "SINGLE_LAYER_MODE" in CPP_OPTIONS.h
464     C-- MOST of THERMODYNAMICS will be disabled
465     #ifndef SINGLE_LAYER_MODE
466    
467 adcroft 1.1 #ifdef ALLOW_GMREDI
468    
469 heimbach 1.22 #ifdef ALLOW_AUTODIFF_TAMC
470 heimbach 1.34 cph storing here is needed only for one GMREDI_OPTIONS:
471     cph define GM_BOLUS_ADVEC
472     cph but I've avoided the #ifdef for now, in case more things change
473     CADJ STORE sigmaX(:,:,:) = comlev1_bibj, key=itdkey, byte=isbyte
474     CADJ STORE sigmaY(:,:,:) = comlev1_bibj, key=itdkey, byte=isbyte
475     CADJ STORE sigmaR(:,:,:) = comlev1_bibj, key=itdkey, byte=isbyte
476 heimbach 1.22 #endif /* ALLOW_AUTODIFF_TAMC */
477 heimbach 1.35
478 adcroft 1.1 C-- Calculate iso-neutral slopes for the GM/Redi parameterisation
479     IF (useGMRedi) THEN
480 edhill 1.56 #ifdef ALLOW_DEBUG
481 heimbach 1.43 IF ( debugLevel .GE. debLevB )
482     & CALL DEBUG_CALL('GMREDI_CALC_TENSOR',myThid)
483 adcroft 1.40 #endif
484 jmc 1.15 CALL GMREDI_CALC_TENSOR(
485     I bi, bj, iMin, iMax, jMin, jMax,
486 adcroft 1.1 I sigmaX, sigmaY, sigmaR,
487     I myThid )
488     #ifdef ALLOW_AUTODIFF_TAMC
489     ELSE
490 jmc 1.15 CALL GMREDI_CALC_TENSOR_DUMMY(
491     I bi, bj, iMin, iMax, jMin, jMax,
492 adcroft 1.1 I sigmaX, sigmaY, sigmaR,
493     I myThid )
494     #endif /* ALLOW_AUTODIFF_TAMC */
495     ENDIF
496    
497     #ifdef ALLOW_AUTODIFF_TAMC
498 heimbach 1.30 CADJ STORE Kwx(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
499     CADJ STORE Kwy(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
500     CADJ STORE Kwz(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
501 adcroft 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
502    
503     #endif /* ALLOW_GMREDI */
504    
505     #ifdef ALLOW_KPP
506     C-- Compute KPP mixing coefficients
507     IF (useKPP) THEN
508 edhill 1.56 #ifdef ALLOW_DEBUG
509 heimbach 1.43 IF ( debugLevel .GE. debLevB )
510     & CALL DEBUG_CALL('KPP_CALC',myThid)
511 adcroft 1.40 #endif
512 adcroft 1.1 CALL KPP_CALC(
513     I bi, bj, myTime, myThid )
514     #ifdef ALLOW_AUTODIFF_TAMC
515     ELSE
516     CALL KPP_CALC_DUMMY(
517     I bi, bj, myTime, myThid )
518     #endif /* ALLOW_AUTODIFF_TAMC */
519     ENDIF
520    
521     #ifdef ALLOW_AUTODIFF_TAMC
522     CADJ STORE KPPghat (:,:,:,bi,bj)
523     CADJ & , KPPfrac (:,: ,bi,bj)
524 heimbach 1.30 CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
525 adcroft 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
526    
527     #endif /* ALLOW_KPP */
528    
529     #ifdef ALLOW_AUTODIFF_TAMC
530 heimbach 1.30 CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
531     CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
532     CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
533     CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
534 adcroft 1.1 #ifdef ALLOW_PASSIVE_TRACER
535 heimbach 1.30 CADJ STORE tr1 (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
536 adcroft 1.1 #endif
537 heimbach 1.42 #ifdef ALLOW_PTRACERS
538     cph-- moved to forward_step to avoid key computation
539     cphCADJ STORE ptracer(:,:,:,bi,bj,itracer) = comlev1_bibj,
540     cphCADJ & key=itdkey, byte=isbyte
541     #endif
542 adcroft 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
543    
544     #ifdef ALLOW_AIM
545     C AIM - atmospheric intermediate model, physics package code.
546     IF ( useAIM ) THEN
547 edhill 1.56 #ifdef ALLOW_DEBUG
548 heimbach 1.43 IF ( debugLevel .GE. debLevB )
549     & CALL DEBUG_CALL('AIM_DO_PHYSICS',myThid)
550 adcroft 1.40 #endif
551 jmc 1.33 CALL TIMER_START('AIM_DO_PHYSICS [THERMODYNAMICS]', myThid)
552     CALL AIM_DO_PHYSICS( bi, bj, myTime, myIter, myThid )
553     CALL TIMER_STOP( 'AIM_DO_PHYSICS [THERMODYNAMICS]', myThid)
554 adcroft 1.1 ENDIF
555     #endif /* ALLOW_AIM */
556 jmc 1.14
557 adcroft 1.12 #ifndef DISABLE_MULTIDIM_ADVECTION
558 adcroft 1.4 C-- Some advection schemes are better calculated using a multi-dimensional
559     C method in the absence of any other terms and, if used, is done here.
560 adcroft 1.13 C
561     C The CPP flag DISABLE_MULTIDIM_ADVECTION is currently unset in GAD_OPTIONS.h
562     C The default is to use multi-dimensinal advection for non-linear advection
563     C schemes. However, for the sake of efficiency of the adjoint it is necessary
564     C to be able to exclude this scheme to avoid excessive storage and
565     C recomputation. It *is* differentiable, if you need it.
566     C Edit GAD_OPTIONS.h and #define DISABLE_MULTIDIM_ADVECTION to
567     C disable this section of code.
568 jmc 1.24 IF (tempMultiDimAdvec) THEN
569 edhill 1.56 #ifdef ALLOW_DEBUG
570 heimbach 1.43 IF ( debugLevel .GE. debLevB )
571     & CALL DEBUG_CALL('GAD_ADVECTION',myThid)
572 adcroft 1.40 #endif
573 jmc 1.63 CALL GAD_ADVECTION(
574 jmc 1.69 I tempImplVertAdv, tempAdvScheme, tempVertAdvScheme,
575     I GAD_TEMPERATURE,
576 jmc 1.63 I uVel, vVel, wVel, theta,
577     O gT,
578     I bi,bj,myTime,myIter,myThid)
579 jmc 1.23 ENDIF
580 jmc 1.24 IF (saltMultiDimAdvec) THEN
581 edhill 1.56 #ifdef ALLOW_DEBUG
582 heimbach 1.43 IF ( debugLevel .GE. debLevB )
583     & CALL DEBUG_CALL('GAD_ADVECTION',myThid)
584 adcroft 1.40 #endif
585 jmc 1.63 CALL GAD_ADVECTION(
586 jmc 1.69 I saltImplVertAdv, saltAdvScheme, saltVertAdvScheme,
587     I GAD_SALINITY,
588 jmc 1.63 I uVel, vVel, wVel, salt,
589     O gS,
590     I bi,bj,myTime,myIter,myThid)
591 adcroft 1.6 ENDIF
592 adcroft 1.17 C Since passive tracers are configurable separately from T,S we
593     C call the multi-dimensional method for PTRACERS regardless
594     C of whether multiDimAdvection is set or not.
595     #ifdef ALLOW_PTRACERS
596     IF ( usePTRACERS ) THEN
597 edhill 1.56 #ifdef ALLOW_DEBUG
598 heimbach 1.43 IF ( debugLevel .GE. debLevB )
599     & CALL DEBUG_CALL('PTRACERS_ADVECTION',myThid)
600 adcroft 1.40 #endif
601 adcroft 1.17 CALL PTRACERS_ADVECTION( bi,bj,myIter,myTime,myThid )
602     ENDIF
603     #endif /* ALLOW_PTRACERS */
604 adcroft 1.12 #endif /* DISABLE_MULTIDIM_ADVECTION */
605 adcroft 1.1
606 edhill 1.56 #ifdef ALLOW_DEBUG
607 jmc 1.63 IF ( debugLevel .GE. debLevB )
608 heimbach 1.43 & CALL DEBUG_MSG('ENTERING DOWNWARD K LOOP',myThid)
609 adcroft 1.40 #endif
610    
611 adcroft 1.1 C-- Start of thermodynamics loop
612     DO k=Nr,1,-1
613     #ifdef ALLOW_AUTODIFF_TAMC
614     C? Patrick Is this formula correct?
615     cph Yes, but I rewrote it.
616     cph Also, the KappaR? need the index and subscript k!
617 heimbach 1.30 kkey = (itdkey-1)*Nr + k
618 adcroft 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
619    
620     C-- km1 Points to level above k (=k-1)
621     C-- kup Cycles through 1,2 to point to layer above
622     C-- kDown Cycles through 2,1 to point to current layer
623    
624     km1 = MAX(1,k-1)
625     kup = 1+MOD(k+1,2)
626     kDown= 1+MOD(k,2)
627    
628     iMin = 1-OLx
629     iMax = sNx+OLx
630     jMin = 1-OLy
631     jMax = sNy+OLy
632    
633 jmc 1.64 kp1Msk=1.
634     IF (k.EQ.Nr) kp1Msk=0.
635     DO j=1-Oly,sNy+Oly
636     DO i=1-Olx,sNx+Olx
637     rTransKp1(i,j) = kp1Msk*rTrans(i,j)
638     ENDDO
639     ENDDO
640 heimbach 1.66 #ifdef ALLOW_AUTODIFF_TAMC
641     CADJ STORE rTransKp1(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
642     #endif
643 jmc 1.64
644 adcroft 1.1 C-- Get temporary terms used by tendency routines
645     CALL CALC_COMMON_FACTORS (
646     I bi,bj,iMin,iMax,jMin,jMax,k,
647     O xA,yA,uTrans,vTrans,rTrans,maskUp,
648     I myThid)
649 jmc 1.19
650 jmc 1.64 IF (k.EQ.1) THEN
651     C- Surface interface :
652     DO j=1-Oly,sNy+Oly
653     DO i=1-Olx,sNx+Olx
654     rTrans(i,j) = 0.
655     ENDDO
656     ENDDO
657     ELSE
658     C- Interior interface :
659     DO j=1-Oly,sNy+Oly
660     DO i=1-Olx,sNx+Olx
661     rTrans(i,j) = rTrans(i,j)*maskC(i,j,k-1,bi,bj)
662     ENDDO
663     ENDDO
664     ENDIF
665    
666 jmc 1.19 #ifdef ALLOW_GMREDI
667 heimbach 1.35
668 jmc 1.19 C-- Residual transp = Bolus transp + Eulerian transp
669     IF (useGMRedi) THEN
670     CALL GMREDI_CALC_UVFLOW(
671     & uTrans, vTrans, bi, bj, k, myThid)
672     IF (K.GE.2) CALL GMREDI_CALC_WFLOW(
673     & rTrans, bi, bj, k, myThid)
674     ENDIF
675 heimbach 1.35
676 heimbach 1.66 #ifdef ALLOW_AUTODIFF_TAMC
677     CADJ STORE rTrans(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
678 heimbach 1.35 #ifdef GM_BOLUS_ADVEC
679     CADJ STORE uTrans(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
680     CADJ STORE vTrans(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
681     #endif
682     #endif /* ALLOW_AUTODIFF_TAMC */
683    
684 jmc 1.19 #endif /* ALLOW_GMREDI */
685 adcroft 1.1
686     #ifdef INCLUDE_CALC_DIFFUSIVITY_CALL
687     C-- Calculate the total vertical diffusivity
688     CALL CALC_DIFFUSIVITY(
689     I bi,bj,iMin,iMax,jMin,jMax,k,
690     I maskUp,
691 heimbach 1.2 O KappaRT,KappaRS,
692 adcroft 1.1 I myThid)
693 heimbach 1.52 # ifdef ALLOW_AUTODIFF_TAMC
694     CADJ STORE KappaRT(:,:,k) = comlev1_bibj_k, key=kkey, byte=isbyte
695     CADJ STORE KappaRS(:,:,k) = comlev1_bibj_k, key=kkey, byte=isbyte
696     # endif /* ALLOW_AUTODIFF_TAMC */
697 adcroft 1.1 #endif
698    
699     iMin = 1-OLx+2
700     iMax = sNx+OLx-1
701     jMin = 1-OLy+2
702     jMax = sNy+OLy-1
703    
704     C-- Calculate active tracer tendencies (gT,gS,...)
705     C and step forward storing result in gTnm1, gSnm1, etc.
706     IF ( tempStepping ) THEN
707     CALL CALC_GT(
708     I bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
709 jmc 1.64 I xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
710 adcroft 1.1 I KappaRT,
711     U fVerT,
712 adcroft 1.7 I myTime,myIter,myThid)
713 adcroft 1.1 CALL TIMESTEP_TRACER(
714 adcroft 1.3 I bi,bj,iMin,iMax,jMin,jMax,k,tempAdvScheme,
715 adcroft 1.1 I theta, gT,
716     I myIter, myThid)
717     ENDIF
718 jmc 1.44
719 adcroft 1.1 IF ( saltStepping ) THEN
720     CALL CALC_GS(
721     I bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
722 jmc 1.64 I xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
723 adcroft 1.1 I KappaRS,
724     U fVerS,
725 adcroft 1.7 I myTime,myIter,myThid)
726 adcroft 1.1 CALL TIMESTEP_TRACER(
727 adcroft 1.3 I bi,bj,iMin,iMax,jMin,jMax,k,saltAdvScheme,
728 adcroft 1.1 I salt, gS,
729     I myIter, myThid)
730     ENDIF
731     #ifdef ALLOW_PASSIVE_TRACER
732 edhill 1.51 ceh3 needs an IF ( usePASSIVE_TRACER ) THEN
733 adcroft 1.1 IF ( tr1Stepping ) THEN
734     CALL CALC_GTR1(
735     I bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown,
736 jmc 1.64 I xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
737 adcroft 1.1 I KappaRT,
738     U fVerTr1,
739 heimbach 1.8 I myTime,myIter,myThid)
740 adcroft 1.1 CALL TIMESTEP_TRACER(
741 adcroft 1.3 I bi,bj,iMin,iMax,jMin,jMax,k,tracerAdvScheme,
742 adcroft 1.1 I Tr1, gTr1,
743 heimbach 1.8 I myIter,myThid)
744 adcroft 1.1 ENDIF
745     #endif
746 adcroft 1.17 #ifdef ALLOW_PTRACERS
747     IF ( usePTRACERS ) THEN
748 heimbach 1.42 CALL PTRACERS_INTEGRATE(
749 adcroft 1.17 I bi,bj,k,
750 jmc 1.64 I xA,yA,uTrans,vTrans,rTrans,rTransKp1,maskUp,
751 heimbach 1.55 X fVerP, KappaRS,
752 adcroft 1.17 I myIter,myTime,myThid)
753     ENDIF
754     #endif /* ALLOW_PTRACERS */
755 adcroft 1.1
756     #ifdef ALLOW_OBCS
757     C-- Apply open boundary conditions
758     IF (useOBCS) THEN
759 adcroft 1.7 CALL OBCS_APPLY_TS( bi, bj, k, gT, gS, myThid )
760 adcroft 1.1 END IF
761     #endif /* ALLOW_OBCS */
762 edhill 1.54
763 jmc 1.59 C-- Freeze water
764     C this bit of code is left here for backward compatibility.
765     C freezing at surface level has been moved to FORWARD_STEP
766     IF ( useOldFreezing .AND. .NOT. useSEAICE
767 jmc 1.61 & .AND. .NOT.(useThSIce.AND.k.EQ.1) ) THEN
768 jmc 1.59 #ifdef ALLOW_AUTODIFF_TAMC
769     CADJ STORE gT(:,:,k,bi,bj) = comlev1_bibj_k
770     CADJ & , key = kkey, byte = isbyte
771     #endif /* ALLOW_AUTODIFF_TAMC */
772     CALL FREEZE( bi, bj, iMin, iMax, jMin, jMax, k, myThid )
773     ENDIF
774 adcroft 1.1
775     C-- end of thermodynamic k loop (Nr:1)
776     ENDDO
777 cheisey 1.31
778 adcroft 1.1
779 jmc 1.63 C-- Implicit vertical advection & diffusion
780     #ifdef INCLUDE_IMPLVERTADV_CODE
781     IF ( tempImplVertAdv ) THEN
782     CALL GAD_IMPLICIT_R(
783     I tempImplVertAdv, tempAdvScheme, GAD_TEMPERATURE,
784     I kappaRT, wVel, theta,
785     U gT,
786     I bi, bj, myTime, myIter, myThid )
787     ELSEIF ( tempStepping .AND. implicitDiffusion ) THEN
788     #else /* INCLUDE_IMPLVERTADV_CODE */
789     IF ( tempStepping .AND. implicitDiffusion ) THEN
790     #endif /* INCLUDE_IMPLVERTADV_CODE */
791 adcroft 1.1 #ifdef ALLOW_AUTODIFF_TAMC
792 heimbach 1.52 CADJ STORE KappaRT(:,:,:) = comlev1_bibj , key=itdkey, byte=isbyte
793 heimbach 1.30 CADJ STORE gT(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte
794 adcroft 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
795 jmc 1.63 CALL IMPLDIFF(
796 adcroft 1.1 I bi, bj, iMin, iMax, jMin, jMax,
797     I deltaTtracer, KappaRT, recip_HFacC,
798 adcroft 1.7 U gT,
799 adcroft 1.1 I myThid )
800 jmc 1.63 ENDIF
801 adcroft 1.1
802 jmc 1.63 #ifdef INCLUDE_IMPLVERTADV_CODE
803     IF ( saltImplVertAdv ) THEN
804     CALL GAD_IMPLICIT_R(
805     I saltImplVertAdv, saltAdvScheme, GAD_SALINITY,
806     I kappaRS, wVel, salt,
807     U gS,
808     I bi, bj, myTime, myIter, myThid )
809     ELSEIF ( saltStepping .AND. implicitDiffusion ) THEN
810     #else /* INCLUDE_IMPLVERTADV_CODE */
811     IF ( saltStepping .AND. implicitDiffusion ) THEN
812     #endif /* INCLUDE_IMPLVERTADV_CODE */
813 adcroft 1.1 #ifdef ALLOW_AUTODIFF_TAMC
814 heimbach 1.52 CADJ STORE KappaRS(:,:,:) = comlev1_bibj , key=itdkey, byte=isbyte
815 heimbach 1.30 CADJ STORE gS(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte
816 adcroft 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
817 jmc 1.63 CALL IMPLDIFF(
818 adcroft 1.1 I bi, bj, iMin, iMax, jMin, jMax,
819     I deltaTtracer, KappaRS, recip_HFacC,
820 adcroft 1.7 U gS,
821 adcroft 1.1 I myThid )
822 jmc 1.63 ENDIF
823 adcroft 1.1
824     #ifdef ALLOW_PASSIVE_TRACER
825 jmc 1.63 IF ( tr1Stepping .AND. implicitDiffusion ) THEN
826 adcroft 1.1 #ifdef ALLOW_AUTODIFF_TAMC
827 heimbach 1.30 CADJ STORE gTr1(:,:,:,bi,bj) = comlev1_bibj , key=itdkey, byte=isbyte
828 adcroft 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
829     CALL IMPLDIFF(
830     I bi, bj, iMin, iMax, jMin, jMax,
831     I deltaTtracer, KappaRT, recip_HFacC,
832 adcroft 1.7 U gTr1,
833 adcroft 1.1 I myThid )
834 jmc 1.63 ENDIF
835 adcroft 1.1 #endif
836    
837 adcroft 1.17 #ifdef ALLOW_PTRACERS
838 jmc 1.63 c #ifdef INCLUDE_IMPLVERTADV_CODE
839     c IF ( usePTRACERS .AND. ptracerImplVertAdv ) THEN
840     c ELSEIF ( usePTRACERS .AND. implicitDiffusion ) THEN
841     c #else
842     IF ( usePTRACERS .AND. implicitDiffusion ) THEN
843     C-- Vertical diffusion (implicit) for passive tracers
844 adcroft 1.17 CALL PTRACERS_IMPLDIFF( bi,bj,KappaRS,myThid )
845 jmc 1.63 ENDIF
846 adcroft 1.17 #endif /* ALLOW_PTRACERS */
847    
848 adcroft 1.1 #ifdef ALLOW_OBCS
849     C-- Apply open boundary conditions
850 jmc 1.63 IF ( ( implicitDiffusion
851     & .OR. tempImplVertAdv
852     & .OR. saltImplVertAdv
853     & ) .AND. useOBCS ) THEN
854 adcroft 1.1 DO K=1,Nr
855 adcroft 1.7 CALL OBCS_APPLY_TS( bi, bj, k, gT, gS, myThid )
856 adcroft 1.1 ENDDO
857 jmc 1.63 ENDIF
858 adcroft 1.1 #endif /* ALLOW_OBCS */
859    
860 jmc 1.39 #ifdef ALLOW_TIMEAVE
861 dimitri 1.65 #ifndef HRCUBE
862 jmc 1.39 IF (taveFreq.GT.0. .AND. ivdc_kappa.NE.0.) THEN
863 jmc 1.70 CALL TIMEAVE_CUMULATE(ConvectCountTave, IVDConvCount,
864 jmc 1.39 I Nr, deltaTclock, bi, bj, myThid)
865     ENDIF
866     useVariableK = useKPP .OR. useGMredi .OR. ivdc_kappa.NE.0.
867     IF (taveFreq.GT.0. .AND. useVariableK ) THEN
868     IF (implicitDiffusion) THEN
869     CALL TIMEAVE_CUMUL_DIF_1T(TdiffRtave, gT, kappaRT,
870     I Nr, 3, deltaTclock, bi, bj, myThid)
871     ELSE
872     CALL TIMEAVE_CUMUL_DIF_1T(TdiffRtave, theta, kappaRT,
873     I Nr, 3, deltaTclock, bi, bj, myThid)
874     ENDIF
875     ENDIF
876 dimitri 1.65 #endif /* ndef HRCUBE */
877 jmc 1.39 #endif /* ALLOW_TIMEAVE */
878    
879 heimbach 1.22 #endif /* SINGLE_LAYER_MODE */
880 adcroft 1.1
881 jmc 1.39 C-- end bi,bj loops.
882 adcroft 1.1 ENDDO
883     ENDDO
884 adcroft 1.17
885 edhill 1.56 #ifdef ALLOW_DEBUG
886 adcroft 1.17 If (debugMode) THEN
887     CALL DEBUG_STATS_RL(Nr,uVel,'Uvel (THERMODYNAMICS)',myThid)
888     CALL DEBUG_STATS_RL(Nr,vVel,'Vvel (THERMODYNAMICS)',myThid)
889     CALL DEBUG_STATS_RL(Nr,wVel,'Wvel (THERMODYNAMICS)',myThid)
890     CALL DEBUG_STATS_RL(Nr,theta,'Theta (THERMODYNAMICS)',myThid)
891     CALL DEBUG_STATS_RL(Nr,salt,'Salt (THERMODYNAMICS)',myThid)
892     CALL DEBUG_STATS_RL(Nr,Gt,'Gt (THERMODYNAMICS)',myThid)
893     CALL DEBUG_STATS_RL(Nr,Gs,'Gs (THERMODYNAMICS)',myThid)
894     CALL DEBUG_STATS_RL(Nr,GtNm1,'GtNm1 (THERMODYNAMICS)',myThid)
895     CALL DEBUG_STATS_RL(Nr,GsNm1,'GsNm1 (THERMODYNAMICS)',myThid)
896 adcroft 1.18 #ifdef ALLOW_PTRACERS
897     IF ( usePTRACERS ) THEN
898     CALL PTRACERS_DEBUG(myThid)
899     ENDIF
900     #endif /* ALLOW_PTRACERS */
901 adcroft 1.17 ENDIF
902 adcroft 1.40 #endif
903    
904 edhill 1.56 #ifdef ALLOW_DEBUG
905 heimbach 1.43 IF ( debugLevel .GE. debLevB )
906 jmc 1.63 & CALL DEBUG_LEAVE('THERMODYNAMICS',myThid)
907 adcroft 1.17 #endif
908 adcroft 1.1
909     RETURN
910     END
  ViewVC Help
Powered by ViewVC 1.1.22