/[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.70 - (hide annotations) (download)
Sun Aug 24 21:44:56 2008 UTC (15 years, 9 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint61c
Changes since 1.69: +2 -4 lines 
always call FREEZE_SURFACE if "allowFreezing" is set ;
(+ stop in config_check if trying to use seaice pkgs with allowFreezing=T)
1 jmc 1.70 C $Header: /u/gcmpack/MITgcm/model/src/do_oceanic_phys.F,v 1.69 2008/08/17 02:08:24 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 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
18    
19     CBOP
20     C !ROUTINE: DO_OCEANIC_PHYS
21     C !INTERFACE:
22     SUBROUTINE DO_OCEANIC_PHYS(myTime, myIter, myThid)
23     C !DESCRIPTION: \bv
24     C *==========================================================*
25 jmc 1.28 C | SUBROUTINE DO_OCEANIC_PHYS
26     C | o Controlling routine for oceanic physics and
27 jmc 1.1 C | parameterization
28     C *==========================================================*
29     C | o originally, part of S/R thermodynamics
30     C *==========================================================*
31     C \ev
32    
33     C !USES:
34     IMPLICIT NONE
35     C == Global variables ===
36     #include "SIZE.h"
37     #include "EEPARAMS.h"
38     #include "PARAMS.h"
39 jmc 1.69 #include "GRID.h"
40 jmc 1.1 #include "DYNVARS.h"
41 jmc 1.20 #ifdef ALLOW_TIMEAVE
42     #include "TIMEAVE_STATV.h"
43     #endif
44 mlosch 1.22 #if defined (ALLOW_BALANCE_FLUXES) && !(defined ALLOW_AUTODIFF_TAMC)
45     #include "FFIELDS.h"
46     #endif
47 jmc 1.1
48     #ifdef ALLOW_AUTODIFF_TAMC
49     # include "tamc.h"
50     # include "tamc_keys.h"
51     # include "FFIELDS.h"
52 heimbach 1.54 # include "SURFACE.h"
53 jmc 1.1 # include "EOS.h"
54     # ifdef ALLOW_KPP
55     # include "KPP.h"
56     # endif
57     # ifdef ALLOW_GMREDI
58     # include "GMREDI.h"
59     # endif
60     # ifdef ALLOW_EBM
61     # include "EBM.h"
62     # endif
63 jmc 1.29 # ifdef ALLOW_EXF
64     # include "ctrl.h"
65 jmc 1.40 # include "EXF_FIELDS.h"
66 jmc 1.29 # ifdef ALLOW_BULKFORMULAE
67 jmc 1.40 # include "EXF_CONSTANTS.h"
68 jmc 1.29 # endif
69     # endif
70     # ifdef ALLOW_SEAICE
71     # include "SEAICE.h"
72     # endif
73 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
74    
75     C !INPUT/OUTPUT PARAMETERS:
76     C == Routine arguments ==
77 jmc 1.18 C myTime :: Current time in simulation
78     C myIter :: Current iteration number in simulation
79     C myThid :: Thread number for this instance of the routine.
80 jmc 1.1 _RL myTime
81     INTEGER myIter
82     INTEGER myThid
83    
84     C !LOCAL VARIABLES:
85     C == Local variables
86 jmc 1.47 C rhoK, rhoKm1 :: Density at current level, and level above
87 jmc 1.18 C iMin, iMax :: Ranges and sub-block indices on which calculations
88 jmc 1.1 C jMin, jMax are applied.
89 jmc 1.18 C bi, bj :: tile indices
90     C i,j,k :: loop indices
91 jmc 1.47 _RL rhoKp1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
92     _RL rhoKm1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
93     _RL rhoK (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
94 jmc 1.1 _RL sigmaX (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
95     _RL sigmaY (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
96     _RL sigmaR (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
97     INTEGER iMin, iMax
98     INTEGER jMin, jMax
99     INTEGER bi, bj
100 jmc 1.18 INTEGER i, j, k
101 jmc 1.17 INTEGER doDiagsRho
102     #ifdef ALLOW_DIAGNOSTICS
103     LOGICAL DIAGNOSTICS_IS_ON
104     EXTERNAL DIAGNOSTICS_IS_ON
105     #endif /* ALLOW_DIAGNOSTICS */
106 jmc 1.69 #ifdef ALLOW_DOWN_SLOPE
107     _RL rhoInSitu(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
108     #endif /* ALLOW_DOWN_SLOPE */
109 jmc 1.1
110     CEOP
111    
112 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
113     C-- dummy statement to end declaration part
114     itdkey = 1
115     #endif /* ALLOW_AUTODIFF_TAMC */
116    
117 jmc 1.1 #ifdef ALLOW_DEBUG
118 jmc 1.29 IF ( debugLevel .GE. debLevB )
119     & CALL DEBUG_ENTER('DO_OCEANIC_PHYS',myThid)
120 jmc 1.1 #endif
121 jmc 1.36
122 jmc 1.17 doDiagsRho = 0
123     #ifdef ALLOW_DIAGNOSTICS
124     IF ( useDiagnostics .AND. fluidIsWater ) THEN
125     IF ( DIAGNOSTICS_IS_ON('RHOANOSQ',myThid) .OR.
126     & DIAGNOSTICS_IS_ON('URHOMASS',myThid) .OR.
127     & DIAGNOSTICS_IS_ON('VRHOMASS',myThid) .OR.
128     & DIAGNOSTICS_IS_ON('WRHOMASS',myThid) .OR.
129     & DIAGNOSTICS_IS_ON('WRHOMASS',myThid) ) doDiagsRho = 2
130 jmc 1.47 IF ( doDiagsRho.EQ.0 .AND.
131     & DIAGNOSTICS_IS_ON('MXLDEPTH',myThid) ) doDiagsRho = 1
132     IF ( doDiagsRho.EQ.0 .AND.
133     & DIAGNOSTICS_IS_ON('DRHODR ',myThid) ) doDiagsRho = 1
134 jmc 1.17 ENDIF
135     #endif /* ALLOW_DIAGNOSTICS */
136    
137 jmc 1.69
138 jmc 1.29 #ifdef ALLOW_SEAICE
139     IF ( useSEAICE ) THEN
140 heimbach 1.62 # ifdef ALLOW_AUTODIFF_TAMC
141 heimbach 1.65 cph-adj-test(
142     CADJ STORE area,empmr,qsw,theta = comlev1, key = ikey_dynamics
143     cph-adj-test)
144 heimbach 1.34 CADJ STORE atemp,aqh,precip = comlev1, key = ikey_dynamics
145     CADJ STORE swdown,lwdown = comlev1, key = ikey_dynamics
146     cph# ifdef EXF_READ_EVAP
147 heimbach 1.33 CADJ STORE evap = comlev1, key = ikey_dynamics
148 heimbach 1.34 cph# endif
149 jmc 1.29 CADJ STORE uvel,vvel = comlev1, key = ikey_dynamics
150 heimbach 1.62 # ifdef SEAICE_ALLOW_DYNAMICS
151     CADJ STORE uice = comlev1, key = ikey_dynamics
152     CADJ STORE vice = comlev1, key = ikey_dynamics
153     # ifdef SEAICE_ALLOW_EVP
154 heimbach 1.50 CADJ STORE seaice_sigma1 = comlev1, key = ikey_dynamics
155     CADJ STORE seaice_sigma2 = comlev1, key = ikey_dynamics
156     CADJ STORE seaice_sigma12 = comlev1, key = ikey_dynamics
157 heimbach 1.62 # endif
158     # endif
159     # ifdef SEAICE_SALINITY
160 heimbach 1.56 CADJ STORE salt = comlev1, key = ikey_dynamics
161 heimbach 1.62 # endif
162     # ifdef ATMOSPHERIC_LOADING
163 heimbach 1.64 CADJ STORE pload = comlev1, key = ikey_dynamics
164 heimbach 1.37 CADJ STORE siceload = comlev1, key = ikey_dynamics
165 heimbach 1.62 # endif
166     # ifdef NONLIN_FRSURF
167 heimbach 1.55 CADJ STORE recip_hfacc = comlev1, key = ikey_dynamics
168 heimbach 1.62 # endif
169 heimbach 1.55 # endif
170 heimbach 1.62 # ifdef ALLOW_DEBUG
171 jmc 1.29 IF ( debugLevel .GE. debLevB )
172     & CALL DEBUG_CALL('SEAICE_MODEL',myThid)
173 heimbach 1.62 # endif
174 jmc 1.29 CALL TIMER_START('SEAICE_MODEL [DO_OCEANIC_PHYS]', myThid)
175     CALL SEAICE_MODEL( myTime, myIter, myThid )
176     CALL TIMER_STOP ('SEAICE_MODEL [DO_OCEANIC_PHYS]', myThid)
177 heimbach 1.62 # ifdef ALLOW_COST
178 heimbach 1.57 CALL SEAICE_COST_SENSI ( myTime, myIter, myThid )
179 heimbach 1.62 # endif
180 heimbach 1.35 ENDIF
181 jmc 1.29 #endif /* ALLOW_SEAICE */
182    
183 heimbach 1.64 #ifdef ALLOW_AUTODIFF_TAMC
184     CADJ STORE sst, sss = comlev1, key = ikey_dynamics
185     CADJ STORE qsw = comlev1, key = ikey_dynamics
186     # ifdef ALLOW_SEAICE
187     CADJ STORE area = comlev1, key = ikey_dynamics
188     # endif
189     #endif
190    
191 jscott 1.30 #if (defined ALLOW_THSICE) && !(defined ALLOW_ATM2D)
192 jmc 1.14 IF ( useThSIce .AND. fluidIsWater ) THEN
193 jmc 1.5 #ifdef ALLOW_DEBUG
194     IF ( debugLevel .GE. debLevB )
195     & CALL DEBUG_CALL('THSICE_MAIN',myThid)
196     #endif
197     C-- Step forward Therm.Sea-Ice variables
198     C and modify forcing terms including effects from ice
199     CALL TIMER_START('THSICE_MAIN [DO_OCEANIC_PHYS]', myThid)
200     CALL THSICE_MAIN( myTime, myIter, myThid )
201     CALL TIMER_STOP( 'THSICE_MAIN [DO_OCEANIC_PHYS]', myThid)
202     ENDIF
203     #endif /* ALLOW_THSICE */
204    
205 mlosch 1.21 #ifdef ALLOW_SHELFICE
206     IF ( useShelfIce .AND. fluidIsWater ) THEN
207     #ifdef ALLOW_DEBUG
208     IF ( debugLevel .GE. debLevB )
209     & CALL DEBUG_CALL('SHELFICE_THERMODYNAMICS',myThid)
210     #endif
211 jmc 1.47 C compute temperature and (virtual) salt flux at the
212 mlosch 1.21 C shelf-ice ocean interface
213     CALL TIMER_START('SHELFICE_THERMODYNAMICS [DO_OCEANIC_PHYS]',
214     & myThid)
215     CALL SHELFICE_THERMODYNAMICS( myTime, myIter, myThid )
216     CALL TIMER_STOP( 'SHELFICE_THERMODYNAMICS [DO_OCEANIC_PHYS]',
217     & myThid)
218     ENDIF
219     #endif /* ALLOW_SHELFICE */
220    
221 jmc 1.5 C-- Freeze water at the surface
222     #ifdef ALLOW_AUTODIFF_TAMC
223     CADJ STORE theta = comlev1, key = ikey_dynamics
224     #endif
225 jmc 1.70 IF ( allowFreezing ) THEN
226 jmc 1.5 CALL FREEZE_SURFACE( myTime, myIter, myThid )
227     ENDIF
228    
229 jmc 1.28 #ifdef ALLOW_OCN_COMPON_INTERF
230 jmc 1.5 C-- Apply imported data (from coupled interface) to forcing fields
231 jmc 1.28 C jmc: do not know precisely where to put this call (bf or af thSIce ?)
232 jmc 1.36 IF ( useCoupler ) THEN
233 jmc 1.5 CALL OCN_APPLY_IMPORT( .TRUE., myTime, myIter, myThid )
234 jmc 1.36 ENDIF
235 jmc 1.28 #endif /* ALLOW_OCN_COMPON_INTERF */
236 jmc 1.5
237 jmc 1.25 #ifdef ALLOW_BALANCE_FLUXES
238 jmc 1.36 C balance fluxes
239     IF ( balanceEmPmR )
240 jmc 1.25 & CALL REMOVE_MEAN_RS( 1, EmPmR, maskH, maskH, rA, drF,
241     & 'EmPmR', myTime, myThid )
242 jmc 1.36 IF ( balanceQnet )
243 jmc 1.25 & CALL REMOVE_MEAN_RS( 1, Qnet, maskH, maskH, rA, drF,
244     & 'Qnet ', myTime, myThid )
245     #endif /* ALLOW_BALANCE_FLUXES */
246    
247 jmc 1.1 #ifdef ALLOW_AUTODIFF_TAMC
248     C-- HPF directive to help TAMC
249     CHPF$ INDEPENDENT
250     #endif /* ALLOW_AUTODIFF_TAMC */
251     DO bj=myByLo(myThid),myByHi(myThid)
252     #ifdef ALLOW_AUTODIFF_TAMC
253 heimbach 1.15 C-- HPF directive to help TAMC
254     CHPF$ INDEPENDENT
255 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
256     DO bi=myBxLo(myThid),myBxHi(myThid)
257    
258     #ifdef ALLOW_AUTODIFF_TAMC
259     act1 = bi - myBxLo(myThid)
260     max1 = myBxHi(myThid) - myBxLo(myThid) + 1
261     act2 = bj - myByLo(myThid)
262     max2 = myByHi(myThid) - myByLo(myThid) + 1
263     act3 = myThid - 1
264     max3 = nTx*nTy
265     act4 = ikey_dynamics - 1
266     itdkey = (act1 + 1) + act2*max1
267     & + act3*max1*max2
268     & + act4*max1*max2*max3
269 jmc 1.36 #else /* ALLOW_AUTODIFF_TAMC */
270 jmc 1.47 C if fluid is not water, by-pass find_rho, gmredi, surfaceForcing
271 jmc 1.36 C and all vertical mixing schemes, but keep OBCS_CALC
272     IF ( fluidIsWater ) THEN
273 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
274    
275     C-- Set up work arrays with valid (i.e. not NaN) values
276     C These inital values do not alter the numerical results. They
277     C just ensure that all memory references are to valid floating
278     C point numbers. This prevents spurious hardware signals due to
279     C uninitialised but inert locations.
280    
281 jmc 1.69 #ifdef ALLOW_AUTODIFF_TAMC
282 jmc 1.1 DO j=1-OLy,sNy+OLy
283     DO i=1-OLx,sNx+OLx
284 jmc 1.69 rhoKm1 (i,j) = 0. _d 0
285 jmc 1.47 rhoK (i,j) = 0. _d 0
286     rhoKp1 (i,j) = 0. _d 0
287 jmc 1.1 ENDDO
288     ENDDO
289 jmc 1.69 #endif /* ALLOW_AUTODIFF_TAMC */
290 jmc 1.1
291     DO k=1,Nr
292     DO j=1-OLy,sNy+OLy
293     DO i=1-OLx,sNx+OLx
294 jmc 1.69 C This is currently used by GMRedi, IVDC, MXL-depth and Diagnostics
295 jmc 1.1 sigmaX(i,j,k) = 0. _d 0
296     sigmaY(i,j,k) = 0. _d 0
297     sigmaR(i,j,k) = 0. _d 0
298     #ifdef ALLOW_AUTODIFF_TAMC
299     cph all the following init. are necessary for TAF
300     cph although some of these are re-initialised later.
301     IVDConvCount(i,j,k,bi,bj) = 0.
302     # ifdef ALLOW_GMREDI
303     Kwx(i,j,k,bi,bj) = 0. _d 0
304     Kwy(i,j,k,bi,bj) = 0. _d 0
305     Kwz(i,j,k,bi,bj) = 0. _d 0
306     # ifdef GM_NON_UNITY_DIAGONAL
307     Kux(i,j,k,bi,bj) = 0. _d 0
308     Kvy(i,j,k,bi,bj) = 0. _d 0
309     # endif
310     # ifdef GM_EXTRA_DIAGONAL
311     Kuz(i,j,k,bi,bj) = 0. _d 0
312     Kvz(i,j,k,bi,bj) = 0. _d 0
313     # endif
314     # ifdef GM_BOLUS_ADVEC
315     GM_PsiX(i,j,k,bi,bj) = 0. _d 0
316     GM_PsiY(i,j,k,bi,bj) = 0. _d 0
317     # endif
318     # ifdef GM_VISBECK_VARIABLE_K
319     VisbeckK(i,j,bi,bj) = 0. _d 0
320     # endif
321     # endif /* ALLOW_GMREDI */
322 heimbach 1.42 # ifdef ALLOW_KPP
323     KPPdiffKzS(i,j,k,bi,bj) = 0. _d 0
324     KPPdiffKzT(i,j,k,bi,bj) = 0. _d 0
325     # endif /* ALLOW_KPP */
326 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
327     ENDDO
328     ENDDO
329     ENDDO
330    
331     iMin = 1-OLx
332     iMax = sNx+OLx
333     jMin = 1-OLy
334     jMax = sNy+OLy
335    
336     #ifdef ALLOW_AUTODIFF_TAMC
337     CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
338     CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
339 heimbach 1.12 CADJ STORE totphihyd(:,:,:,bi,bj)
340 jmc 1.1 CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
341 heimbach 1.10 # ifdef ALLOW_KPP
342 jmc 1.1 CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
343     CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key=itdkey, byte=isbyte
344 heimbach 1.10 # endif
345 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
346    
347     #ifdef ALLOW_DEBUG
348 jmc 1.36 IF ( debugLevel .GE. debLevB )
349 jmc 1.1 & CALL DEBUG_MSG('ENTERING UPWARD K LOOP',myThid)
350     #endif
351    
352     C-- Start of diagnostic loop
353     DO k=Nr,1,-1
354    
355     #ifdef ALLOW_AUTODIFF_TAMC
356     C? Patrick, is this formula correct now that we change the loop range?
357     C? Do we still need this?
358 jmc 1.69 cph kkey formula corrected.
359 jmc 1.47 cph Needed for rhoK, rhoKm1, in the case useGMREDI.
360 jmc 1.1 kkey = (itdkey-1)*Nr + k
361     #endif /* ALLOW_AUTODIFF_TAMC */
362    
363 jmc 1.69 #ifdef ALLOW_DOWN_SLOPE
364     IF ( useDOWN_SLOPE ) THEN
365     CALL DWNSLP_CALC_RHO(
366     I theta, salt,
367     O rhoK, rhoInSitu,
368     I k, bi, bj, myTime, myIter, myThid )
369     ENDIF
370     #endif /* ALLOW_DOWN_SLOPE */
371    
372 jmc 1.1 C-- Calculate gradients of potential density for isoneutral
373     C slope terms (e.g. GM/Redi tensor or IVDC diffusivity)
374 jmc 1.17 IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.)
375 dimitri 1.61 & .OR. useSALT_PLUME .OR. doDiagsRho.GE.1 ) THEN
376 jmc 1.1 #ifdef ALLOW_DEBUG
377 jmc 1.36 IF ( debugLevel .GE. debLevB )
378 jmc 1.68 & CALL DEBUG_CALL('FIND_RHO_2D',myThid)
379 jmc 1.1 #endif
380 jmc 1.69 IF ( .NOT.useDOWN_SLOPE ) THEN
381 jmc 1.1 #ifdef ALLOW_AUTODIFF_TAMC
382     CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
383     CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
384     #endif /* ALLOW_AUTODIFF_TAMC */
385 jmc 1.69 CALL FIND_RHO_2D(
386 jmc 1.68 I iMin, iMax, jMin, jMax, k,
387     I theta(1-OLx,1-OLy,k,bi,bj),
388     I salt (1-OLx,1-OLy,k,bi,bj),
389     O rhoK,
390     I k, bi, bj, myThid )
391 jmc 1.1
392 jmc 1.69 ENDIF
393 jmc 1.1 IF (k.GT.1) THEN
394     #ifdef ALLOW_AUTODIFF_TAMC
395     CADJ STORE theta(:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
396     CADJ STORE salt (:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
397     #endif /* ALLOW_AUTODIFF_TAMC */
398 jmc 1.68 CALL FIND_RHO_2D(
399     I iMin, iMax, jMin, jMax, k,
400     I theta(1-OLx,1-OLy,k-1,bi,bj),
401     I salt (1-OLx,1-OLy,k-1,bi,bj),
402     O rhoKm1,
403     I k-1, bi, bj, myThid )
404 jmc 1.1 ENDIF
405     #ifdef ALLOW_DEBUG
406 jmc 1.36 IF ( debugLevel .GE. debLevB )
407 jmc 1.1 & CALL DEBUG_CALL('GRAD_SIGMA',myThid)
408     #endif
409 heimbach 1.31 cph Avoid variable aliasing for adjoint !!!
410     DO j=jMin,jMax
411     DO i=iMin,iMax
412 jmc 1.47 rhoKp1(i,j) = rhoK(i,j)
413 heimbach 1.31 ENDDO
414     ENDDO
415 jmc 1.1 CALL GRAD_SIGMA(
416     I bi, bj, iMin, iMax, jMin, jMax, k,
417 heimbach 1.31 I rhoK, rhoKm1, rhoKp1,
418 jmc 1.1 O sigmaX, sigmaY, sigmaR,
419     I myThid )
420 gforget 1.66 #ifdef ALLOW_AUTODIFF_TAMC
421 jmc 1.69 #ifdef GMREDI_WITH_STABLE_ADJOINT
422 gforget 1.66 cgf zero out adjoint fields to stabilize pkg/gmredi adjoint
423     cgf -> cuts adjoint dependency from slope to state
424 jmc 1.69 CALL ZERO_ADJ_LOC( Nr, sigmaX, myThid)
425     CALL ZERO_ADJ_LOC( Nr, sigmaY, myThid)
426     CALL ZERO_ADJ_LOC( Nr, sigmaR, myThid)
427 gforget 1.66 #endif
428     #endif /* ALLOW_AUTODIFF_TAMC */
429 jmc 1.1 ENDIF
430    
431     C-- Implicit Vertical Diffusion for Convection
432     c ==> should use sigmaR !!!
433     IF (k.GT.1 .AND. ivdc_kappa.NE.0.) THEN
434     #ifdef ALLOW_DEBUG
435 jmc 1.36 IF ( debugLevel .GE. debLevB )
436 jmc 1.1 & CALL DEBUG_CALL('CALC_IVDC',myThid)
437     #endif
438     CALL CALC_IVDC(
439     I bi, bj, iMin, iMax, jMin, jMax, k,
440     I rhoKm1, rhoK,
441     I myTime, myIter, myThid)
442     ENDIF
443    
444 jmc 1.17 #ifdef ALLOW_DIAGNOSTICS
445     IF ( doDiagsRho.GE.2 ) THEN
446     CALL DIAGS_RHO( k, bi, bj,
447     I rhoK, rhoKm1,
448     I myTime, myIter, myThid)
449     ENDIF
450     #endif
451    
452 jmc 1.1 C-- end of diagnostic k loop (Nr:1)
453     ENDDO
454    
455 heimbach 1.57 #ifdef ALLOW_AUTODIFF_TAMC
456 jmc 1.69 CADJ STORE IVDConvCount(:,:,:,bi,bj)
457 heimbach 1.57 CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
458     #endif
459    
460 jmc 1.47 C-- Diagnose Mixed Layer Depth:
461     IF ( useGMRedi .OR. doDiagsRho.GE.1 ) THEN
462     CALL CALC_OCE_MXLAYER( rhoK, sigmaR,
463     & bi, bj, myTime, myIter, myThid )
464     ENDIF
465 heimbach 1.53
466 dimitri 1.52 #ifdef ALLOW_SALT_PLUME
467 dimitri 1.61 IF ( useSALT_PLUME ) THEN
468 dimitri 1.60 CALL SALT_PLUME_CALC_DEPTH( rhoK, sigmaR,
469 dimitri 1.52 & bi, bj, myTime, myIter, myThid )
470 dimitri 1.60 ENDIF
471 dimitri 1.61 #endif /* ALLOW_SALT_PLUME */
472    
473 jmc 1.8 #ifdef ALLOW_DIAGNOSTICS
474 jmc 1.17 IF ( doDiagsRho.GE.1 ) THEN
475 jmc 1.16 CALL DIAGNOSTICS_FILL (sigmaR, 'DRHODR ', 0, Nr,
476     & 2, bi, bj, myThid)
477 jmc 1.8 ENDIF
478 dimitri 1.61 #endif /* ALLOW_DIAGNOSTICS */
479 jmc 1.8
480 jmc 1.1 C-- Determines forcing terms based on external fields
481     C relaxation terms, etc.
482     #ifdef ALLOW_DEBUG
483 jmc 1.36 IF ( debugLevel .GE. debLevB )
484 jmc 1.1 & CALL DEBUG_CALL('EXTERNAL_FORCING_SURF',myThid)
485     #endif
486 heimbach 1.23 #ifdef ALLOW_AUTODIFF_TAMC
487     CADJ STORE EmPmR(:,:,bi,bj)
488     CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
489 heimbach 1.26 # ifdef EXACT_CONSERV
490 heimbach 1.23 CADJ STORE PmEpR(:,:,bi,bj)
491     CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
492 heimbach 1.26 # endif
493 heimbach 1.27 # ifdef NONLIN_FRSURF
494     CADJ STORE hFac_surfC(:,:,bi,bj)
495     CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
496     CADJ STORE recip_hFacC(:,:,:,bi,bj)
497     CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
498     # endif
499 heimbach 1.23 #endif
500 jmc 1.36 CALL EXTERNAL_FORCING_SURF(
501 jmc 1.1 I bi, bj, iMin, iMax, jMin, jMax,
502     I myTime, myIter, myThid )
503 heimbach 1.27 #ifdef ALLOW_AUTODIFF_TAMC
504     # ifdef EXACT_CONSERV
505     cph-test
506     cphCADJ STORE PmEpR(:,:,bi,bj)
507     cphCADJ & = comlev1_bibj, key=itdkey, byte=isbyte
508     # endif
509     #endif
510 jmc 1.1
511     #ifdef ALLOW_AUTODIFF_TAMC
512     cph needed for KPP
513 jmc 1.4 CADJ STORE surfaceForcingU(:,:,bi,bj)
514 jmc 1.1 CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
515 jmc 1.4 CADJ STORE surfaceForcingV(:,:,bi,bj)
516 jmc 1.1 CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
517 jmc 1.4 CADJ STORE surfaceForcingS(:,:,bi,bj)
518 jmc 1.1 CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
519 jmc 1.4 CADJ STORE surfaceForcingT(:,:,bi,bj)
520 jmc 1.1 CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
521 jmc 1.4 CADJ STORE surfaceForcingTice(:,:,bi,bj)
522 jmc 1.1 CADJ & = comlev1_bibj, key=itdkey, byte=isbyte
523     #endif /* ALLOW_AUTODIFF_TAMC */
524    
525     #ifdef ALLOW_KPP
526     C-- Compute KPP mixing coefficients
527     IF (useKPP) THEN
528     #ifdef ALLOW_DEBUG
529 jmc 1.36 IF ( debugLevel .GE. debLevB )
530 jmc 1.1 & CALL DEBUG_CALL('KPP_CALC',myThid)
531     #endif
532     CALL KPP_CALC(
533 jmc 1.44 I bi, bj, myTime, myIter, myThid )
534 jmc 1.1 #ifdef ALLOW_AUTODIFF_TAMC
535     ELSE
536     CALL KPP_CALC_DUMMY(
537 jmc 1.44 I bi, bj, myTime, myIter, myThid )
538 jmc 1.1 #endif /* ALLOW_AUTODIFF_TAMC */
539     ENDIF
540    
541     #endif /* ALLOW_KPP */
542    
543 mlosch 1.6 #ifdef ALLOW_PP81
544     C-- Compute PP81 mixing coefficients
545     IF (usePP81) THEN
546     #ifdef ALLOW_DEBUG
547 jmc 1.36 IF ( debugLevel .GE. debLevB )
548 mlosch 1.6 & CALL DEBUG_CALL('PP81_CALC',myThid)
549     #endif
550     CALL PP81_CALC(
551     I bi, bj, myTime, myThid )
552     ENDIF
553     #endif /* ALLOW_PP81 */
554    
555     #ifdef ALLOW_MY82
556     C-- Compute MY82 mixing coefficients
557     IF (useMY82) THEN
558     #ifdef ALLOW_DEBUG
559 jmc 1.36 IF ( debugLevel .GE. debLevB )
560 mlosch 1.6 & CALL DEBUG_CALL('MY82_CALC',myThid)
561     #endif
562     CALL MY82_CALC(
563     I bi, bj, myTime, myThid )
564     ENDIF
565     #endif /* ALLOW_MY82 */
566    
567 mlosch 1.9 #ifdef ALLOW_GGL90
568     C-- Compute GGL90 mixing coefficients
569     IF (useGGL90) THEN
570     #ifdef ALLOW_DEBUG
571 jmc 1.36 IF ( debugLevel .GE. debLevB )
572 mlosch 1.9 & CALL DEBUG_CALL('GGL90_CALC',myThid)
573     #endif
574     CALL GGL90_CALC(
575     I bi, bj, myTime, myThid )
576     ENDIF
577     #endif /* ALLOW_GGL90 */
578    
579 jmc 1.20 #ifdef ALLOW_TIMEAVE
580 jmc 1.36 IF ( taveFreq.GT. 0. _d 0 ) THEN
581 jmc 1.20 CALL TIMEAVE_SURF_FLUX( bi, bj, myTime, myIter, myThid)
582     ENDIF
583     IF (taveFreq.GT.0. .AND. ivdc_kappa.NE.0.) THEN
584     CALL TIMEAVE_CUMULATE(ConvectCountTave, IVDConvCount,
585     I Nr, deltaTclock, bi, bj, myThid)
586     ENDIF
587     #endif /* ALLOW_TIMEAVE */
588    
589 jmc 1.69 #ifdef ALLOW_GMREDI
590 jmc 1.47 #ifdef ALLOW_AUTODIFF_TAMC
591     # ifndef GM_EXCLUDE_CLIPPING
592     cph storing here is needed only for one GMREDI_OPTIONS:
593     cph define GM_BOLUS_ADVEC
594     cph keep it although TAF says you dont need to.
595     cph but I've avoided the #ifdef for now, in case more things change
596     CADJ STORE sigmaX(:,:,:) = comlev1_bibj, key=itdkey, byte=isbyte
597     CADJ STORE sigmaY(:,:,:) = comlev1_bibj, key=itdkey, byte=isbyte
598     CADJ STORE sigmaR(:,:,:) = comlev1_bibj, key=itdkey, byte=isbyte
599     # endif
600     #endif /* ALLOW_AUTODIFF_TAMC */
601    
602     C-- Calculate iso-neutral slopes for the GM/Redi parameterisation
603     IF (useGMRedi) THEN
604     #ifdef ALLOW_DEBUG
605     IF ( debugLevel .GE. debLevB )
606     & CALL DEBUG_CALL('GMREDI_CALC_TENSOR',myThid)
607     #endif
608     CALL GMREDI_CALC_TENSOR(
609 jmc 1.51 I iMin, iMax, jMin, jMax,
610 jmc 1.47 I sigmaX, sigmaY, sigmaR,
611 jmc 1.51 I bi, bj, myTime, myIter, myThid )
612 jmc 1.47 #ifdef ALLOW_AUTODIFF_TAMC
613     ELSE
614     CALL GMREDI_CALC_TENSOR_DUMMY(
615 jmc 1.51 I iMin, iMax, jMin, jMax,
616 jmc 1.47 I sigmaX, sigmaY, sigmaR,
617 jmc 1.51 I bi, bj, myTime, myIter, myThid )
618 jmc 1.47 #endif /* ALLOW_AUTODIFF_TAMC */
619     ENDIF
620 jmc 1.69 #endif /* ALLOW_GMREDI */
621    
622     #ifdef ALLOW_DOWN_SLOPE
623     IF ( useDOWN_SLOPE ) THEN
624     C-- Calculate Downsloping Flow for Down_Slope parameterization
625     IF ( usingPCoords ) THEN
626     CALL DWNSLP_CALC_FLOW(
627     I bi, bj, kSurfC,
628     I rhoInSitu,
629     I myTime, myIter, myThid )
630     ELSE
631     CALL DWNSLP_CALC_FLOW(
632     I bi, bj, kLowC,
633     I rhoInSitu,
634     I myTime, myIter, myThid )
635     ENDIF
636     ENDIF
637     #endif /* ALLOW_DOWN_SLOPE */
638 jmc 1.47
639 jmc 1.36 #ifndef ALLOW_AUTODIFF_TAMC
640     C--- if fluid Is Water: end
641     ENDIF
642     #endif
643    
644     #ifdef ALLOW_OBCS
645     C-- Calculate future values on open boundaries
646     IF (useOBCS) THEN
647     #ifdef ALLOW_DEBUG
648     IF ( debugLevel .GE. debLevB )
649     & CALL DEBUG_CALL('OBCS_CALC',myThid)
650     #endif
651     CALL OBCS_CALC( bi, bj, myTime+deltaTclock, myIter+1,
652     I uVel, vVel, wVel, theta, salt,
653     I myThid )
654     ENDIF
655     #endif /* ALLOW_OBCS */
656    
657 jmc 1.1 C-- end bi,bj loops.
658     ENDDO
659     ENDDO
660    
661 jmc 1.45 #ifdef ALLOW_KPP
662     IF (useKPP) THEN
663     CALL KPP_DO_EXCH( myThid )
664     ENDIF
665     #endif /* ALLOW_KPP */
666    
667 jmc 1.18 #ifdef ALLOW_DIAGNOSTICS
668     IF ( fluidIsWater .AND. useDiagnostics ) THEN
669     CALL DIAGS_OCEANIC_SURF_FLUX( myTime, myIter, myThid )
670     ENDIF
671 jmc 1.19 IF ( ivdc_kappa.NE.0 .AND. useDiagnostics ) THEN
672     CALL DIAGNOSTICS_FILL( IVDConvCount,'CONVADJ ',
673     & 0, Nr, 0, 1, 1, myThid )
674     ENDIF
675 jmc 1.18 #endif
676    
677 jmc 1.1 #ifdef ALLOW_DEBUG
678 jmc 1.29 IF ( debugLevel .GE. debLevB )
679     & CALL DEBUG_LEAVE('DO_OCEANIC_PHYS',myThid)
680 jmc 1.1 #endif
681    
682     RETURN
683     END
  ViewVC Help
Powered by ViewVC 1.1.22