/[MITgcm]/MITgcm/pkg/mom_vecinv/mom_vecinv.F
ViewVC logotype

Contents of /MITgcm/pkg/mom_vecinv/mom_vecinv.F

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


Revision 1.47 - (show annotations) (download)
Fri Sep 16 19:32:20 2005 UTC (18 years, 8 months ago) by baylor
Branch: MAIN
Changes since 1.46: +18 -26 lines
Move calculation of viscosities to separate file mom_common/mom_calc_visc.F.  This allows Leith, LeithD,
and Smagorinsky to be used simultaneously, and soon will allow them to be used in strain-tension,
fluxform, and nonhydrostatic modes.  Also, introduce viscC4Smag, a biharmonic Smagorinsky
viscosity as in Griffies and Hallberg.

1 C $Header: /u/gcmpack/MITgcm/pkg/mom_vecinv/mom_vecinv.F,v 1.46 2005/09/04 19:29:03 jmc Exp $
2 C $Name: $
3
4 #include "MOM_VECINV_OPTIONS.h"
5
6 SUBROUTINE MOM_VECINV(
7 I bi,bj,iMin,iMax,jMin,jMax,k,kUp,kDown,
8 I KappaRU, KappaRV,
9 U fVerU, fVerV,
10 O guDiss, gvDiss,
11 I myTime, myIter, myThid)
12 C /==========================================================\
13 C | S/R MOM_VECINV |
14 C | o Form the right hand-side of the momentum equation. |
15 C |==========================================================|
16 C | Terms are evaluated one layer at a time working from |
17 C | the bottom to the top. The vertically integrated |
18 C | barotropic flow tendency term is evluated by summing the |
19 C | tendencies. |
20 C | Notes: |
21 C | We have not sorted out an entirely satisfactory formula |
22 C | for the diffusion equation bc with lopping. The present |
23 C | form produces a diffusive flux that does not scale with |
24 C | open-area. Need to do something to solidfy this and to |
25 C | deal "properly" with thin walls. |
26 C \==========================================================/
27 IMPLICIT NONE
28
29 C == Global variables ==
30 #include "SIZE.h"
31 #include "DYNVARS.h"
32 #include "EEPARAMS.h"
33 #include "PARAMS.h"
34 #ifdef ALLOW_MNC
35 #include "MNC_PARAMS.h"
36 #endif
37 #include "GRID.h"
38 #ifdef ALLOW_TIMEAVE
39 #include "TIMEAVE_STATV.h"
40 #endif
41
42 C == Routine arguments ==
43 C fVerU :: Flux of momentum in the vertical direction, out of the upper
44 C fVerV :: face of a cell K ( flux into the cell above ).
45 C guDiss :: dissipation tendency (all explicit terms), u component
46 C gvDiss :: dissipation tendency (all explicit terms), v component
47 C bi, bj, iMin, iMax, jMin, jMax - Range of points for which calculation
48 C results will be set.
49 C kUp, kDown - Index for upper and lower layers.
50 C myThid - Instance number for this innvocation of CALC_MOM_RHS
51 _RL KappaRU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
52 _RL KappaRV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
53 _RL fVerU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
54 _RL fVerV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
55 _RL guDiss(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
56 _RL gvDiss(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
57 INTEGER kUp,kDown
58 _RL myTime
59 INTEGER myIter
60 INTEGER myThid
61 INTEGER bi,bj,iMin,iMax,jMin,jMax
62
63 #ifdef ALLOW_MOM_VECINV
64
65 C == Functions ==
66 LOGICAL DIFFERENT_MULTIPLE
67 EXTERNAL DIFFERENT_MULTIPLE
68
69 C == Local variables ==
70 _RL vF (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
71 _RL vrF (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
72 _RL uCf (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
73 _RL vCf (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
74 c _RL mT (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
75 _RL del2u(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
76 _RL del2v(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
77 _RL tension(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
78 _RL strain(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
79 _RS hFacZ(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
80 _RS r_hFacZ(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
81 _RL uFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
82 _RL vFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
83 _RL dStar(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
84 _RL zStar(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
85 C I,J,K - Loop counters
86 INTEGER i,j,k
87 C xxxFac - On-off tracer parameters used for switching terms off.
88 _RL ArDudrFac
89 c _RL mtFacU
90 _RL ArDvdrFac
91 c _RL mtFacV
92 LOGICAL bottomDragTerms
93 LOGICAL writeDiag
94 _RL KE(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
95 _RL omega3(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
96 _RL vort3(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
97 _RL hDiv(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
98
99 #ifdef ALLOW_MNC
100 INTEGER offsets(9)
101 #endif
102
103 #ifdef ALLOW_AUTODIFF_TAMC
104 C-- only the kDown part of fverU/V is set in this subroutine
105 C-- the kUp is still required
106 C-- In the case of mom_fluxform Kup is set as well
107 C-- (at least in part)
108 fVerU(1,1,kUp) = fVerU(1,1,kUp)
109 fVerV(1,1,kUp) = fVerV(1,1,kUp)
110 #endif
111
112 writeDiag = DIFFERENT_MULTIPLE(diagFreq, myTime, deltaTClock)
113
114 #ifdef ALLOW_MNC
115 IF (useMNC .AND. snapshot_mnc .AND. writeDiag) THEN
116 IF ((bi .EQ. 1).AND.(bj .EQ. 1).AND.(k .EQ. 1)) THEN
117 CALL MNC_CW_SET_UDIM('mom_vi', -1, myThid)
118 CALL MNC_CW_RL_W_S('D','mom_vi',0,0,'T',myTime,myThid)
119 CALL MNC_CW_SET_UDIM('mom_vi', 0, myThid)
120 CALL MNC_CW_I_W_S('I','mom_vi',0,0,'iter',myIter,myThid)
121 ENDIF
122 DO i = 1,9
123 offsets(i) = 0
124 ENDDO
125 offsets(3) = k
126 C write(*,*) 'offsets = ',(offsets(i),i=1,9)
127 ENDIF
128 #endif /* ALLOW_MNC */
129
130 C Initialise intermediate terms
131 DO J=1-OLy,sNy+OLy
132 DO I=1-OLx,sNx+OLx
133 vF(i,j) = 0.
134 vrF(i,j) = 0.
135 uCf(i,j) = 0.
136 vCf(i,j) = 0.
137 c mT(i,j) = 0.
138 del2u(i,j) = 0.
139 del2v(i,j) = 0.
140 dStar(i,j) = 0.
141 zStar(i,j) = 0.
142 guDiss(i,j)= 0.
143 gvDiss(i,j)= 0.
144 vort3(i,j) = 0.
145 omega3(i,j)= 0.
146 ke(i,j) = 0.
147 #ifdef ALLOW_AUTODIFF_TAMC
148 strain(i,j) = 0. _d 0
149 tension(i,j) = 0. _d 0
150 #endif
151 ENDDO
152 ENDDO
153
154 C-- Term by term tracer parmeters
155 C o U momentum equation
156 ArDudrFac = vfFacMom*1.
157 c mTFacU = mtFacMom*1.
158 C o V momentum equation
159 ArDvdrFac = vfFacMom*1.
160 c mTFacV = mtFacMom*1.
161
162 IF ( no_slip_bottom
163 & .OR. bottomDragQuadratic.NE.0.
164 & .OR. bottomDragLinear.NE.0.) THEN
165 bottomDragTerms=.TRUE.
166 ELSE
167 bottomDragTerms=.FALSE.
168 ENDIF
169
170 C-- Calculate open water fraction at vorticity points
171 CALL MOM_CALC_HFACZ(bi,bj,k,hFacZ,r_hFacZ,myThid)
172
173 C Make local copies of horizontal flow field
174 DO j=1-OLy,sNy+OLy
175 DO i=1-OLx,sNx+OLx
176 uFld(i,j) = uVel(i,j,k,bi,bj)
177 vFld(i,j) = vVel(i,j,k,bi,bj)
178 ENDDO
179 ENDDO
180
181 C note (jmc) : Dissipation and Vort3 advection do not necesary
182 C use the same maskZ (and hFacZ) => needs 2 call(s)
183 c CALL MOM_VI_HFACZ_DISS(bi,bj,k,hFacZ,r_hFacZ,myThid)
184
185 CALL MOM_CALC_KE(bi,bj,k,2,uFld,vFld,KE,myThid)
186
187 CALL MOM_CALC_HDIV(bi,bj,k,2,uFld,vFld,hDiv,myThid)
188
189 CALL MOM_CALC_RELVORT3(bi,bj,k,uFld,vFld,hFacZ,vort3,myThid)
190
191 CALL MOM_CALC_TENSION(bi,bj,k,uFld,vFld,tension,myThid)
192
193 CALL MOM_CALC_STRAIN(bi,bj,k,uFld,vFld,hFacZ,strain,myThid)
194
195 IF (useAbsVorticity)
196 & CALL MOM_CALC_ABSVORT3(bi,bj,k,vort3,omega3,myThid)
197
198 IF (momViscosity) THEN
199 C Calculate del^2 u and del^2 v for bi-harmonic term
200 IF ( (viscA4.NE.0. .AND. no_slip_sides)
201 & .OR. viscA4D.NE.0. .OR. viscA4Z.NE.0.
202 & .OR. viscA4Grid.NE.0.
203 & .OR. viscC4leith.NE.0.
204 & .OR. viscC4leithD.NE.0.
205 & .OR. viscC4smag.NE.0.
206 & ) THEN
207 CALL MOM_VI_DEL2UV(bi,bj,k,hDiv,vort3,hFacZ,
208 O del2u,del2v,
209 & myThid)
210 ENDIF
211
212 C Calculate dissipation terms for U and V equations
213
214 C in terms of tension and strain
215 IF (useStrainTensionVisc) THEN
216 CALL MOM_HDISSIP(bi,bj,k,hDiv,vort3,tension,strain,KE,
217 I hFacZ,
218 O guDiss,gvDiss,
219 I myThid)
220 ELSE
221 C in terms of vorticity and divergence
222 CALL MOM_VI_HDISSIP(bi,bj,k,hDiv,vort3,tension,strain,KE,
223 I hFacZ,dStar,zStar,
224 O guDiss,gvDiss,
225 & myThid)
226 ENDIF
227 ENDIF
228
229 C- Return to standard hfacZ (min-4) and mask vort3 accordingly:
230 c CALL MOM_VI_MASK_VORT3(bi,bj,k,hFacZ,r_hFacZ,vort3,myThid)
231
232 C---- Zonal momentum equation starts here
233
234 C-- Vertical flux (fVer is at upper face of "u" cell)
235
236 C Eddy component of vertical flux (interior component only) -> vrF
237 IF (momViscosity.AND..NOT.implicitViscosity) THEN
238 CALL MOM_U_RVISCFLUX(bi,bj,k+1,uVel,KappaRU,vrF,myThid)
239
240 C Combine fluxes
241 DO j=jMin,jMax
242 DO i=iMin,iMax
243 fVerU(i,j,kDown) = ArDudrFac*vrF(i,j)
244 ENDDO
245 ENDDO
246
247 C-- Tendency is minus divergence of the fluxes
248 DO j=2-Oly,sNy+Oly-1
249 DO i=2-Olx,sNx+Olx-1
250 guDiss(i,j) = guDiss(i,j)
251 & -_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)
252 & *recip_rAw(i,j,bi,bj)
253 & *(
254 & fVerU(i,j,kDown) - fVerU(i,j,kUp)
255 & )*rkSign
256 ENDDO
257 ENDDO
258 ENDIF
259
260 C-- No-slip and drag BCs appear as body forces in cell abutting topography
261 IF (momViscosity.AND.no_slip_sides) THEN
262 C- No-slip BCs impose a drag at walls...
263 CALL MOM_U_SIDEDRAG(bi,bj,k,uFld,del2u,hFacZ,vF,myThid)
264 DO j=jMin,jMax
265 DO i=iMin,iMax
266 guDiss(i,j) = guDiss(i,j)+vF(i,j)
267 ENDDO
268 ENDDO
269 ENDIF
270
271 C- No-slip BCs impose a drag at bottom
272 IF (momViscosity.AND.bottomDragTerms) THEN
273 CALL MOM_U_BOTTOMDRAG(bi,bj,k,uFld,KE,KappaRU,vF,myThid)
274 DO j=jMin,jMax
275 DO i=iMin,iMax
276 guDiss(i,j) = guDiss(i,j)+vF(i,j)
277 ENDDO
278 ENDDO
279 ENDIF
280
281 C-- Metric terms for curvilinear grid systems
282 c IF (usingSphericalPolarMTerms) THEN
283 C o Spherical polar grid metric terms
284 c CALL MOM_U_METRIC_NH(bi,bj,k,uFld,wVel,mT,myThid)
285 c DO j=jMin,jMax
286 c DO i=iMin,iMax
287 c gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+mTFacU*mT(i,j)
288 c ENDDO
289 c ENDDO
290 c ENDIF
291
292 C---- Meridional momentum equation starts here
293
294 C-- Vertical flux (fVer is at upper face of "v" cell)
295
296 C Eddy component of vertical flux (interior component only) -> vrF
297 IF (momViscosity.AND..NOT.implicitViscosity) THEN
298 CALL MOM_V_RVISCFLUX(bi,bj,k+1,vVel,KappaRV,vrF,myThid)
299
300 C Combine fluxes -> fVerV
301 DO j=jMin,jMax
302 DO i=iMin,iMax
303 fVerV(i,j,kDown) = ArDvdrFac*vrF(i,j)
304 ENDDO
305 ENDDO
306
307 C-- Tendency is minus divergence of the fluxes
308 DO j=jMin,jMax
309 DO i=iMin,iMax
310 gvDiss(i,j) = gvDiss(i,j)
311 & -_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
312 & *recip_rAs(i,j,bi,bj)
313 & *(
314 & fVerV(i,j,kDown) - fVerV(i,j,kUp)
315 & )*rkSign
316 ENDDO
317 ENDDO
318 ENDIF
319
320 C-- No-slip and drag BCs appear as body forces in cell abutting topography
321 IF (momViscosity.AND.no_slip_sides) THEN
322 C- No-slip BCs impose a drag at walls...
323 CALL MOM_V_SIDEDRAG(bi,bj,k,vFld,del2v,hFacZ,vF,myThid)
324 DO j=jMin,jMax
325 DO i=iMin,iMax
326 gvDiss(i,j) = gvDiss(i,j)+vF(i,j)
327 ENDDO
328 ENDDO
329 ENDIF
330 C- No-slip BCs impose a drag at bottom
331 IF (momViscosity.AND.bottomDragTerms) THEN
332 CALL MOM_V_BOTTOMDRAG(bi,bj,k,vFld,KE,KappaRV,vF,myThid)
333 DO j=jMin,jMax
334 DO i=iMin,iMax
335 gvDiss(i,j) = gvDiss(i,j)+vF(i,j)
336 ENDDO
337 ENDDO
338 ENDIF
339
340 C-- Metric terms for curvilinear grid systems
341 c IF (usingSphericalPolarMTerms) THEN
342 C o Spherical polar grid metric terms
343 c CALL MOM_V_METRIC_NH(bi,bj,k,vFld,wVel,mT,myThid)
344 c DO j=jMin,jMax
345 c DO i=iMin,iMax
346 c gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+mTFacV*mT(i,j)
347 c ENDDO
348 c ENDDO
349 c ENDIF
350
351 C-- Horizontal Coriolis terms
352 c IF (useCoriolis .AND. .NOT.useCDscheme
353 c & .AND. .NOT. useAbsVorticity) THEN
354 C- jmc: change it to keep the Coriolis terms when useAbsVorticity=T & momAdvection=F
355 IF ( useCoriolis .AND.
356 & .NOT.( useCDscheme .OR. useAbsVorticity.AND.momAdvection )
357 & ) THEN
358 IF (useAbsVorticity) THEN
359 CALL MOM_VI_U_CORIOLIS(bi,bj,K,vFld,omega3,hFacZ,r_hFacZ,
360 & uCf,myThid)
361 CALL MOM_VI_V_CORIOLIS(bi,bj,K,uFld,omega3,hFacZ,r_hFacZ,
362 & vCf,myThid)
363 ELSE
364 CALL MOM_VI_CORIOLIS(bi,bj,k,uFld,vFld,hFacZ,r_hFacZ,
365 & uCf,vCf,myThid)
366 ENDIF
367 DO j=jMin,jMax
368 DO i=iMin,iMax
369 gU(i,j,k,bi,bj) = uCf(i,j)
370 gV(i,j,k,bi,bj) = vCf(i,j)
371 ENDDO
372 ENDDO
373
374 IF ( writeDiag ) THEN
375 IF (snapshot_mdsio) THEN
376 CALL WRITE_LOCAL_RL('fV','I10',1,uCf,bi,bj,k,myIter,myThid)
377 CALL WRITE_LOCAL_RL('fU','I10',1,vCf,bi,bj,k,myIter,myThid)
378 ENDIF
379 #ifdef ALLOW_MNC
380 IF (useMNC .AND. snapshot_mnc) THEN
381 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'fV', uCf,
382 & offsets, myThid)
383 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'fU', vCf,
384 & offsets, myThid)
385 ENDIF
386 #endif /* ALLOW_MNC */
387 ENDIF
388 #ifdef ALLOW_DIAGNOSTICS
389 IF ( useDiagnostics ) THEN
390 CALL DIAGNOSTICS_FILL(uCf,'Um_Cori ',k,1,2,bi,bj,myThid)
391 CALL DIAGNOSTICS_FILL(vCf,'Vm_Cori ',k,1,2,bi,bj,myThid)
392 ENDIF
393 #endif /* ALLOW_DIAGNOSTICS */
394
395 ELSE
396 DO j=jMin,jMax
397 DO i=iMin,iMax
398 gU(i,j,k,bi,bj) = 0. _d 0
399 gV(i,j,k,bi,bj) = 0. _d 0
400 ENDDO
401 ENDDO
402 ENDIF
403
404 IF (momAdvection) THEN
405 C-- Horizontal advection of relative (or absolute) vorticity
406 IF (highOrderVorticity.AND.useAbsVorticity) THEN
407 CALL MOM_VI_U_CORIOLIS_C4(bi,bj,k,vFld,omega3,r_hFacZ,
408 & uCf,myThid)
409 ELSEIF (highOrderVorticity) THEN
410 CALL MOM_VI_U_CORIOLIS_C4(bi,bj,k,vFld,vort3, r_hFacZ,
411 & uCf,myThid)
412 ELSEIF (useAbsVorticity) THEN
413 CALL MOM_VI_U_CORIOLIS(bi,bj,K,vFld,omega3,hFacZ,r_hFacZ,
414 & uCf,myThid)
415 ELSE
416 CALL MOM_VI_U_CORIOLIS(bi,bj,k,vFld,vort3, hFacZ,r_hFacZ,
417 & uCf,myThid)
418 ENDIF
419 DO j=jMin,jMax
420 DO i=iMin,iMax
421 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
422 ENDDO
423 ENDDO
424 IF (highOrderVorticity.AND.useAbsVorticity) THEN
425 CALL MOM_VI_V_CORIOLIS_C4(bi,bj,K,uFld,omega3,r_hFacZ,
426 & vCf,myThid)
427 ELSEIF (highOrderVorticity) THEN
428 CALL MOM_VI_V_CORIOLIS_C4(bi,bj,K,uFld,vort3, r_hFacZ,
429 & vCf,myThid)
430 ELSEIF (useAbsVorticity) THEN
431 CALL MOM_VI_V_CORIOLIS(bi,bj,K,uFld,omega3,hFacZ,r_hFacZ,
432 & vCf,myThid)
433 ELSE
434 CALL MOM_VI_V_CORIOLIS(bi,bj,k,uFld,vort3, hFacZ,r_hFacZ,
435 & vCf,myThid)
436 ENDIF
437 DO j=jMin,jMax
438 DO i=iMin,iMax
439 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
440 ENDDO
441 ENDDO
442
443 IF ( writeDiag ) THEN
444 IF (snapshot_mdsio) THEN
445 CALL WRITE_LOCAL_RL('zV','I10',1,uCf,bi,bj,k,myIter,myThid)
446 CALL WRITE_LOCAL_RL('zU','I10',1,vCf,bi,bj,k,myIter,myThid)
447 ENDIF
448 #ifdef ALLOW_MNC
449 IF (useMNC .AND. snapshot_mnc) THEN
450 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'zV', uCf,
451 & offsets, myThid)
452 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'zU', vCf,
453 & offsets, myThid)
454 ENDIF
455 #endif /* ALLOW_MNC */
456 ENDIF
457
458 #ifdef ALLOW_TIMEAVE
459 IF (taveFreq.GT.0.) THEN
460 CALL TIMEAVE_CUMUL_1K1T(uZetatave,vCf,deltaTClock,
461 & Nr, k, bi, bj, myThid)
462 CALL TIMEAVE_CUMUL_1K1T(vZetatave,uCf,deltaTClock,
463 & Nr, k, bi, bj, myThid)
464 ENDIF
465 #endif /* ALLOW_TIMEAVE */
466 #ifdef ALLOW_DIAGNOSTICS
467 IF ( useDiagnostics ) THEN
468 CALL DIAGNOSTICS_FILL(uCf,'Um_AdvZ3',k,1,2,bi,bj,myThid)
469 CALL DIAGNOSTICS_FILL(vCf,'Vm_AdvZ3',k,1,2,bi,bj,myThid)
470 ENDIF
471 #endif /* ALLOW_DIAGNOSTICS */
472
473 C-- Vertical shear terms (-w*du/dr & -w*dv/dr)
474 IF ( .NOT. momImplVertAdv ) THEN
475 CALL MOM_VI_U_VERTSHEAR(bi,bj,K,uVel,wVel,uCf,myThid)
476 DO j=jMin,jMax
477 DO i=iMin,iMax
478 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
479 ENDDO
480 ENDDO
481 CALL MOM_VI_V_VERTSHEAR(bi,bj,K,vVel,wVel,vCf,myThid)
482 DO j=jMin,jMax
483 DO i=iMin,iMax
484 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
485 ENDDO
486 ENDDO
487 #ifdef ALLOW_DIAGNOSTICS
488 IF ( useDiagnostics ) THEN
489 CALL DIAGNOSTICS_FILL(uCf,'Um_AdvRe',k,1,2,bi,bj,myThid)
490 CALL DIAGNOSTICS_FILL(vCf,'Vm_AdvRe',k,1,2,bi,bj,myThid)
491 ENDIF
492 #endif /* ALLOW_DIAGNOSTICS */
493 ENDIF
494
495 C-- Bernoulli term
496 CALL MOM_VI_U_GRAD_KE(bi,bj,K,KE,uCf,myThid)
497 DO j=jMin,jMax
498 DO i=iMin,iMax
499 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
500 ENDDO
501 ENDDO
502 CALL MOM_VI_V_GRAD_KE(bi,bj,K,KE,vCf,myThid)
503 DO j=jMin,jMax
504 DO i=iMin,iMax
505 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
506 ENDDO
507 ENDDO
508 IF ( writeDiag ) THEN
509 IF (snapshot_mdsio) THEN
510 CALL WRITE_LOCAL_RL('KEx','I10',1,uCf,bi,bj,k,myIter,myThid)
511 CALL WRITE_LOCAL_RL('KEy','I10',1,vCf,bi,bj,k,myIter,myThid)
512 ENDIF
513 #ifdef ALLOW_MNC
514 IF (useMNC .AND. snapshot_mnc) THEN
515 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'KEx', uCf,
516 & offsets, myThid)
517 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'KEy', vCf,
518 & offsets, myThid)
519 ENDIF
520 #endif /* ALLOW_MNC */
521 ENDIF
522
523 C-- end if momAdvection
524 ENDIF
525
526 C-- Set du/dt & dv/dt on boundaries to zero
527 DO j=jMin,jMax
528 DO i=iMin,iMax
529 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)*_maskW(i,j,k,bi,bj)
530 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)*_maskS(i,j,k,bi,bj)
531 ENDDO
532 ENDDO
533
534 #ifdef ALLOW_DEBUG
535 IF ( debugLevel .GE. debLevB
536 & .AND. k.EQ.4 .AND. myIter.EQ.nIter0
537 & .AND. nPx.EQ.1 .AND. nPy.EQ.1
538 & .AND. useCubedSphereExchange ) THEN
539 CALL DEBUG_CS_CORNER_UV( ' uDiss,vDiss from MOM_VECINV',
540 & guDiss,gvDiss, k, standardMessageUnit,bi,bj,myThid )
541 ENDIF
542 #endif /* ALLOW_DEBUG */
543
544 IF ( writeDiag ) THEN
545 IF (snapshot_mdsio) THEN
546 CALL WRITE_LOCAL_RL('Ds','I10',1,strain,bi,bj,k,myIter,myThid)
547 CALL WRITE_LOCAL_RL('Dt','I10',1,tension,bi,bj,k,myIter,
548 & myThid)
549 CALL WRITE_LOCAL_RL('Du','I10',1,guDiss,bi,bj,k,myIter,myThid)
550 CALL WRITE_LOCAL_RL('Dv','I10',1,gvDiss,bi,bj,k,myIter,myThid)
551 CALL WRITE_LOCAL_RL('Z3','I10',1,vort3,bi,bj,k,myIter,myThid)
552 CALL WRITE_LOCAL_RL('W3','I10',1,omega3,bi,bj,k,myIter,myThid)
553 CALL WRITE_LOCAL_RL('KE','I10',1,KE,bi,bj,k,myIter,myThid)
554 CALL WRITE_LOCAL_RL('D','I10',1,hDiv,bi,bj,k,myIter,myThid)
555 ENDIF
556 #ifdef ALLOW_MNC
557 IF (useMNC .AND. snapshot_mnc) THEN
558 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Ds',strain,
559 & offsets, myThid)
560 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Dt',tension,
561 & offsets, myThid)
562 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Du',guDiss,
563 & offsets, myThid)
564 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Dv',gvDiss,
565 & offsets, myThid)
566 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Z3',vort3,
567 & offsets, myThid)
568 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'W3',omega3,
569 & offsets, myThid)
570 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'KE',KE,
571 & offsets, myThid)
572 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'D', hDiv,
573 & offsets, myThid)
574 ENDIF
575 #endif /* ALLOW_MNC */
576 ENDIF
577
578 #ifdef ALLOW_DIAGNOSTICS
579 IF ( useDiagnostics ) THEN
580 CALL DIAGNOSTICS_FILL(KE, 'momKE ',k,1,2,bi,bj,myThid)
581 CALL DIAGNOSTICS_FILL(hDiv, 'momHDiv ',k,1,2,bi,bj,myThid)
582 CALL DIAGNOSTICS_FILL(vort3, 'momVort3',k,1,2,bi,bj,myThid)
583 CALL DIAGNOSTICS_FILL(gU(1-Olx,1-Oly,k,bi,bj),
584 & 'Um_Advec',k,1,2,bi,bj,myThid)
585 CALL DIAGNOSTICS_FILL(gV(1-Olx,1-Oly,k,bi,bj),
586 & 'Vm_Advec',k,1,2,bi,bj,myThid)
587 IF (momViscosity) THEN
588 CALL DIAGNOSTICS_FILL(guDiss,'Um_Diss ',k,1,2,bi,bj,myThid)
589 CALL DIAGNOSTICS_FILL(gvDiss,'Vm_Diss ',k,1,2,bi,bj,myThid)
590 ENDIF
591 ENDIF
592 #endif /* ALLOW_DIAGNOSTICS */
593
594 #endif /* ALLOW_MOM_VECINV */
595
596 RETURN
597 END

  ViewVC Help
Powered by ViewVC 1.1.22