/[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.52 - (show annotations) (download)
Wed Sep 28 15:53:20 2005 UTC (18 years, 7 months ago) by jmc
Branch: MAIN
Changes since 1.51: +9 -9 lines
use selectKEscheme in vecinv from and make vertshear consistent.

1 C $Header: /u/gcmpack/MITgcm/pkg/mom_vecinv/mom_vecinv.F,v 1.51 2005/09/27 13:38:42 baylor 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 _RL viscAh_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
99 _RL viscAh_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
100 _RL viscA4_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
101 _RL viscA4_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
102 LOGICAL harmonic,biharmonic,useVariableViscosity
103
104 #ifdef ALLOW_MNC
105 INTEGER offsets(9)
106 #endif
107
108 #ifdef ALLOW_AUTODIFF_TAMC
109 C-- only the kDown part of fverU/V is set in this subroutine
110 C-- the kUp is still required
111 C-- In the case of mom_fluxform Kup is set as well
112 C-- (at least in part)
113 fVerU(1,1,kUp) = fVerU(1,1,kUp)
114 fVerV(1,1,kUp) = fVerV(1,1,kUp)
115 #endif
116
117 writeDiag = DIFFERENT_MULTIPLE(diagFreq, myTime, deltaTClock)
118
119 #ifdef ALLOW_MNC
120 IF (useMNC .AND. snapshot_mnc .AND. writeDiag) THEN
121 IF ((bi .EQ. 1).AND.(bj .EQ. 1).AND.(k .EQ. 1)) THEN
122 CALL MNC_CW_SET_UDIM('mom_vi', -1, myThid)
123 CALL MNC_CW_RL_W_S('D','mom_vi',0,0,'T',myTime,myThid)
124 CALL MNC_CW_SET_UDIM('mom_vi', 0, myThid)
125 CALL MNC_CW_I_W_S('I','mom_vi',0,0,'iter',myIter,myThid)
126 ENDIF
127 DO i = 1,9
128 offsets(i) = 0
129 ENDDO
130 offsets(3) = k
131 C write(*,*) 'offsets = ',(offsets(i),i=1,9)
132 ENDIF
133 #endif /* ALLOW_MNC */
134
135 C Initialise intermediate terms
136 DO J=1-OLy,sNy+OLy
137 DO I=1-OLx,sNx+OLx
138 vF(i,j) = 0.
139 vrF(i,j) = 0.
140 uCf(i,j) = 0.
141 vCf(i,j) = 0.
142 c mT(i,j) = 0.
143 del2u(i,j) = 0.
144 del2v(i,j) = 0.
145 dStar(i,j) = 0.
146 zStar(i,j) = 0.
147 guDiss(i,j)= 0.
148 gvDiss(i,j)= 0.
149 vort3(i,j) = 0.
150 omega3(i,j)= 0.
151 ke(i,j) = 0.
152 viscAh_Z(i,j) = 0.
153 viscAh_D(i,j) = 0.
154 viscA4_Z(i,j) = 0.
155 viscA4_D(i,j) = 0.
156
157 #ifdef ALLOW_AUTODIFF_TAMC
158 strain(i,j) = 0. _d 0
159 tension(i,j) = 0. _d 0
160 #endif
161 ENDDO
162 ENDDO
163
164 C-- Term by term tracer parmeters
165 C o U momentum equation
166 ArDudrFac = vfFacMom*1.
167 c mTFacU = mtFacMom*1.
168 C o V momentum equation
169 ArDvdrFac = vfFacMom*1.
170 c mTFacV = mtFacMom*1.
171
172 IF ( no_slip_bottom
173 & .OR. bottomDragQuadratic.NE.0.
174 & .OR. bottomDragLinear.NE.0.) THEN
175 bottomDragTerms=.TRUE.
176 ELSE
177 bottomDragTerms=.FALSE.
178 ENDIF
179
180 C-- Calculate open water fraction at vorticity points
181 CALL MOM_CALC_HFACZ(bi,bj,k,hFacZ,r_hFacZ,myThid)
182
183 C Make local copies of horizontal flow field
184 DO j=1-OLy,sNy+OLy
185 DO i=1-OLx,sNx+OLx
186 uFld(i,j) = uVel(i,j,k,bi,bj)
187 vFld(i,j) = vVel(i,j,k,bi,bj)
188 ENDDO
189 ENDDO
190
191 C note (jmc) : Dissipation and Vort3 advection do not necesary
192 C use the same maskZ (and hFacZ) => needs 2 call(s)
193 c CALL MOM_VI_HFACZ_DISS(bi,bj,k,hFacZ,r_hFacZ,myThid)
194
195 CALL MOM_CALC_KE(bi,bj,k,selectKEscheme,uFld,vFld,KE,myThid)
196
197 CALL MOM_CALC_HDIV(bi,bj,k,2,uFld,vFld,hDiv,myThid)
198
199 CALL MOM_CALC_RELVORT3(bi,bj,k,uFld,vFld,hFacZ,vort3,myThid)
200
201 IF (useAbsVorticity)
202 & CALL MOM_CALC_ABSVORT3(bi,bj,k,vort3,omega3,myThid)
203
204 IF (momViscosity) THEN
205 C Calculate Viscosities
206 CALL MOM_CALC_TENSION(bi,bj,k,uFld,vFld,tension,myThid)
207
208 CALL MOM_CALC_STRAIN(bi,bj,k,uFld,vFld,hFacZ,strain,myThid)
209
210 CALL MOM_CALC_VISC(
211 I bi,bj,k,
212 O viscAh_Z,viscAh_D,viscA4_Z,viscA4_D,
213 O harmonic,biharmonic,useVariableViscosity,
214 I hDiv,vort3,tension,strain,KE,hfacZ,
215 I myThid)
216
217 C Calculate del^2 u and del^2 v for bi-harmonic term
218 IF (biharmonic) THEN
219 CALL MOM_VI_DEL2UV(bi,bj,k,hDiv,vort3,hFacZ,
220 O del2u,del2v,
221 & myThid)
222 CALL MOM_CALC_HDIV(bi,bj,k,2,del2u,del2v,dStar,myThid)
223 CALL MOM_CALC_RELVORT3(bi,bj,k,
224 & del2u,del2v,hFacZ,zStar,myThid)
225 ENDIF
226
227 C Calculate dissipation terms for U and V equations
228
229 C in terms of tension and strain
230 IF (useStrainTensionVisc) THEN
231 CALL MOM_HDISSIP(bi,bj,k,hDiv,vort3,tension,strain,KE,
232 I hFacZ,
233 I viscAh_Z,viscAh_D,viscA4_Z,viscA4_D,
234 I harmonic,biharmonic,useVariableViscosity,
235 O guDiss,gvDiss,
236 I myThid)
237 ELSE
238 C in terms of vorticity and divergence
239 CALL MOM_VI_HDISSIP(bi,bj,k,hDiv,vort3,tension,strain,KE,
240 I hFacZ,dStar,zStar,
241 I viscAh_Z,viscAh_D,viscA4_Z,viscA4_D,
242 I harmonic,biharmonic,useVariableViscosity,
243 O guDiss,gvDiss,
244 & myThid)
245 ENDIF
246 ENDIF
247
248 C- Return to standard hfacZ (min-4) and mask vort3 accordingly:
249 c CALL MOM_VI_MASK_VORT3(bi,bj,k,hFacZ,r_hFacZ,vort3,myThid)
250
251 C---- Zonal momentum equation starts here
252
253 C-- Vertical flux (fVer is at upper face of "u" cell)
254
255 C Eddy component of vertical flux (interior component only) -> vrF
256 IF (momViscosity.AND..NOT.implicitViscosity) THEN
257 CALL MOM_U_RVISCFLUX(bi,bj,k+1,uVel,KappaRU,vrF,myThid)
258
259 C Combine fluxes
260 DO j=jMin,jMax
261 DO i=iMin,iMax
262 fVerU(i,j,kDown) = ArDudrFac*vrF(i,j)
263 ENDDO
264 ENDDO
265
266 C-- Tendency is minus divergence of the fluxes
267 DO j=2-Oly,sNy+Oly-1
268 DO i=2-Olx,sNx+Olx-1
269 guDiss(i,j) = guDiss(i,j)
270 & -_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)
271 & *recip_rAw(i,j,bi,bj)
272 & *(
273 & fVerU(i,j,kDown) - fVerU(i,j,kUp)
274 & )*rkSign
275 ENDDO
276 ENDDO
277 ENDIF
278
279 C-- No-slip and drag BCs appear as body forces in cell abutting topography
280 IF (momViscosity.AND.no_slip_sides) THEN
281 C- No-slip BCs impose a drag at walls...
282 CALL MOM_U_SIDEDRAG(
283 I bi,bj,k,
284 I uFld, del2u, hFacZ,
285 I viscAh_Z,viscA4_Z,
286 I harmonic,biharmonic,useVariableViscosity,
287 O vF,
288 I myThid)
289 DO j=jMin,jMax
290 DO i=iMin,iMax
291 guDiss(i,j) = guDiss(i,j)+vF(i,j)
292 ENDDO
293 ENDDO
294 ENDIF
295
296 C- No-slip BCs impose a drag at bottom
297 IF (momViscosity.AND.bottomDragTerms) THEN
298 CALL MOM_U_BOTTOMDRAG(bi,bj,k,uFld,KE,KappaRU,vF,myThid)
299 DO j=jMin,jMax
300 DO i=iMin,iMax
301 guDiss(i,j) = guDiss(i,j)+vF(i,j)
302 ENDDO
303 ENDDO
304 ENDIF
305
306 C-- Metric terms for curvilinear grid systems
307 c IF (usingSphericalPolarMTerms) THEN
308 C o Spherical polar grid metric terms
309 c CALL MOM_U_METRIC_NH(bi,bj,k,uFld,wVel,mT,myThid)
310 c DO j=jMin,jMax
311 c DO i=iMin,iMax
312 c gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+mTFacU*mT(i,j)
313 c ENDDO
314 c ENDDO
315 c ENDIF
316
317 C---- Meridional momentum equation starts here
318
319 C-- Vertical flux (fVer is at upper face of "v" cell)
320
321 C Eddy component of vertical flux (interior component only) -> vrF
322 IF (momViscosity.AND..NOT.implicitViscosity) THEN
323 CALL MOM_V_RVISCFLUX(bi,bj,k+1,vVel,KappaRV,vrF,myThid)
324
325 C Combine fluxes -> fVerV
326 DO j=jMin,jMax
327 DO i=iMin,iMax
328 fVerV(i,j,kDown) = ArDvdrFac*vrF(i,j)
329 ENDDO
330 ENDDO
331
332 C-- Tendency is minus divergence of the fluxes
333 DO j=jMin,jMax
334 DO i=iMin,iMax
335 gvDiss(i,j) = gvDiss(i,j)
336 & -_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
337 & *recip_rAs(i,j,bi,bj)
338 & *(
339 & fVerV(i,j,kDown) - fVerV(i,j,kUp)
340 & )*rkSign
341 ENDDO
342 ENDDO
343 ENDIF
344
345 C-- No-slip and drag BCs appear as body forces in cell abutting topography
346 IF (momViscosity.AND.no_slip_sides) THEN
347 C- No-slip BCs impose a drag at walls...
348 CALL MOM_V_SIDEDRAG(
349 I bi,bj,k,
350 I vFld, del2v, hFacZ,
351 I viscAh_Z,viscA4_Z,
352 I harmonic,biharmonic,useVariableViscosity,
353 O vF,
354 I myThid)
355 DO j=jMin,jMax
356 DO i=iMin,iMax
357 gvDiss(i,j) = gvDiss(i,j)+vF(i,j)
358 ENDDO
359 ENDDO
360 ENDIF
361 C- No-slip BCs impose a drag at bottom
362 IF (momViscosity.AND.bottomDragTerms) THEN
363 CALL MOM_V_BOTTOMDRAG(bi,bj,k,vFld,KE,KappaRV,vF,myThid)
364 DO j=jMin,jMax
365 DO i=iMin,iMax
366 gvDiss(i,j) = gvDiss(i,j)+vF(i,j)
367 ENDDO
368 ENDDO
369 ENDIF
370
371 C-- Metric terms for curvilinear grid systems
372 c IF (usingSphericalPolarMTerms) THEN
373 C o Spherical polar grid metric terms
374 c CALL MOM_V_METRIC_NH(bi,bj,k,vFld,wVel,mT,myThid)
375 c DO j=jMin,jMax
376 c DO i=iMin,iMax
377 c gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+mTFacV*mT(i,j)
378 c ENDDO
379 c ENDDO
380 c ENDIF
381
382 C-- Horizontal Coriolis terms
383 c IF (useCoriolis .AND. .NOT.useCDscheme
384 c & .AND. .NOT. useAbsVorticity) THEN
385 C- jmc: change it to keep the Coriolis terms when useAbsVorticity=T & momAdvection=F
386 IF ( useCoriolis .AND.
387 & .NOT.( useCDscheme .OR. useAbsVorticity.AND.momAdvection )
388 & ) THEN
389 IF (useAbsVorticity) THEN
390 CALL MOM_VI_U_CORIOLIS(bi,bj,K,vFld,omega3,hFacZ,r_hFacZ,
391 & uCf,myThid)
392 CALL MOM_VI_V_CORIOLIS(bi,bj,K,uFld,omega3,hFacZ,r_hFacZ,
393 & vCf,myThid)
394 ELSE
395 CALL MOM_VI_CORIOLIS(bi,bj,k,uFld,vFld,hFacZ,r_hFacZ,
396 & uCf,vCf,myThid)
397 ENDIF
398 DO j=jMin,jMax
399 DO i=iMin,iMax
400 gU(i,j,k,bi,bj) = uCf(i,j)
401 gV(i,j,k,bi,bj) = vCf(i,j)
402 ENDDO
403 ENDDO
404
405 IF ( writeDiag ) THEN
406 IF (snapshot_mdsio) THEN
407 CALL WRITE_LOCAL_RL('fV','I10',1,uCf,bi,bj,k,myIter,myThid)
408 CALL WRITE_LOCAL_RL('fU','I10',1,vCf,bi,bj,k,myIter,myThid)
409 ENDIF
410 #ifdef ALLOW_MNC
411 IF (useMNC .AND. snapshot_mnc) THEN
412 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'fV', uCf,
413 & offsets, myThid)
414 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'fU', vCf,
415 & offsets, myThid)
416 ENDIF
417 #endif /* ALLOW_MNC */
418 ENDIF
419 #ifdef ALLOW_DIAGNOSTICS
420 IF ( useDiagnostics ) THEN
421 CALL DIAGNOSTICS_FILL(uCf,'Um_Cori ',k,1,2,bi,bj,myThid)
422 CALL DIAGNOSTICS_FILL(vCf,'Vm_Cori ',k,1,2,bi,bj,myThid)
423 ENDIF
424 #endif /* ALLOW_DIAGNOSTICS */
425
426 ELSE
427 DO j=jMin,jMax
428 DO i=iMin,iMax
429 gU(i,j,k,bi,bj) = 0. _d 0
430 gV(i,j,k,bi,bj) = 0. _d 0
431 ENDDO
432 ENDDO
433 ENDIF
434
435 IF (momAdvection) THEN
436 C-- Horizontal advection of relative (or absolute) vorticity
437 IF (highOrderVorticity.AND.useAbsVorticity) THEN
438 CALL MOM_VI_U_CORIOLIS_C4(bi,bj,k,vFld,omega3,r_hFacZ,
439 & uCf,myThid)
440 ELSEIF (highOrderVorticity) THEN
441 CALL MOM_VI_U_CORIOLIS_C4(bi,bj,k,vFld,vort3, r_hFacZ,
442 & uCf,myThid)
443 ELSEIF (useAbsVorticity) THEN
444 CALL MOM_VI_U_CORIOLIS(bi,bj,K,vFld,omega3,hFacZ,r_hFacZ,
445 & uCf,myThid)
446 ELSE
447 CALL MOM_VI_U_CORIOLIS(bi,bj,k,vFld,vort3, hFacZ,r_hFacZ,
448 & uCf,myThid)
449 ENDIF
450 DO j=jMin,jMax
451 DO i=iMin,iMax
452 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
453 ENDDO
454 ENDDO
455 IF (highOrderVorticity.AND.useAbsVorticity) THEN
456 CALL MOM_VI_V_CORIOLIS_C4(bi,bj,K,uFld,omega3,r_hFacZ,
457 & vCf,myThid)
458 ELSEIF (highOrderVorticity) THEN
459 CALL MOM_VI_V_CORIOLIS_C4(bi,bj,K,uFld,vort3, r_hFacZ,
460 & vCf,myThid)
461 ELSEIF (useAbsVorticity) THEN
462 CALL MOM_VI_V_CORIOLIS(bi,bj,K,uFld,omega3,hFacZ,r_hFacZ,
463 & vCf,myThid)
464 ELSE
465 CALL MOM_VI_V_CORIOLIS(bi,bj,k,uFld,vort3, hFacZ,r_hFacZ,
466 & vCf,myThid)
467 ENDIF
468 DO j=jMin,jMax
469 DO i=iMin,iMax
470 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
471 ENDDO
472 ENDDO
473
474 IF ( writeDiag ) THEN
475 IF (snapshot_mdsio) THEN
476 CALL WRITE_LOCAL_RL('zV','I10',1,uCf,bi,bj,k,myIter,myThid)
477 CALL WRITE_LOCAL_RL('zU','I10',1,vCf,bi,bj,k,myIter,myThid)
478 ENDIF
479 #ifdef ALLOW_MNC
480 IF (useMNC .AND. snapshot_mnc) THEN
481 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'zV', uCf,
482 & offsets, myThid)
483 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'zU', vCf,
484 & offsets, myThid)
485 ENDIF
486 #endif /* ALLOW_MNC */
487 ENDIF
488
489 #ifdef ALLOW_TIMEAVE
490 IF (taveFreq.GT.0.) THEN
491 CALL TIMEAVE_CUMUL_1K1T(uZetatave,vCf,deltaTClock,
492 & Nr, k, bi, bj, myThid)
493 CALL TIMEAVE_CUMUL_1K1T(vZetatave,uCf,deltaTClock,
494 & Nr, k, bi, bj, myThid)
495 ENDIF
496 #endif /* ALLOW_TIMEAVE */
497 #ifdef ALLOW_DIAGNOSTICS
498 IF ( useDiagnostics ) THEN
499 CALL DIAGNOSTICS_FILL(uCf,'Um_AdvZ3',k,1,2,bi,bj,myThid)
500 CALL DIAGNOSTICS_FILL(vCf,'Vm_AdvZ3',k,1,2,bi,bj,myThid)
501 ENDIF
502 #endif /* ALLOW_DIAGNOSTICS */
503
504 C-- Vertical shear terms (-w*du/dr & -w*dv/dr)
505 IF ( .NOT. momImplVertAdv ) THEN
506 CALL MOM_VI_U_VERTSHEAR(bi,bj,K,uVel,wVel,uCf,myThid)
507 DO j=jMin,jMax
508 DO i=iMin,iMax
509 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
510 ENDDO
511 ENDDO
512 CALL MOM_VI_V_VERTSHEAR(bi,bj,K,vVel,wVel,vCf,myThid)
513 DO j=jMin,jMax
514 DO i=iMin,iMax
515 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
516 ENDDO
517 ENDDO
518 #ifdef ALLOW_DIAGNOSTICS
519 IF ( useDiagnostics ) THEN
520 CALL DIAGNOSTICS_FILL(uCf,'Um_AdvRe',k,1,2,bi,bj,myThid)
521 CALL DIAGNOSTICS_FILL(vCf,'Vm_AdvRe',k,1,2,bi,bj,myThid)
522 ENDIF
523 #endif /* ALLOW_DIAGNOSTICS */
524 ENDIF
525
526 C-- Bernoulli term
527 CALL MOM_VI_U_GRAD_KE(bi,bj,K,KE,uCf,myThid)
528 DO j=jMin,jMax
529 DO i=iMin,iMax
530 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
531 ENDDO
532 ENDDO
533 CALL MOM_VI_V_GRAD_KE(bi,bj,K,KE,vCf,myThid)
534 DO j=jMin,jMax
535 DO i=iMin,iMax
536 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
537 ENDDO
538 ENDDO
539 IF ( writeDiag ) THEN
540 IF (snapshot_mdsio) THEN
541 CALL WRITE_LOCAL_RL('KEx','I10',1,uCf,bi,bj,k,myIter,myThid)
542 CALL WRITE_LOCAL_RL('KEy','I10',1,vCf,bi,bj,k,myIter,myThid)
543 ENDIF
544 #ifdef ALLOW_MNC
545 IF (useMNC .AND. snapshot_mnc) THEN
546 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'KEx', uCf,
547 & offsets, myThid)
548 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'KEy', vCf,
549 & offsets, myThid)
550 ENDIF
551 #endif /* ALLOW_MNC */
552 ENDIF
553
554 C-- end if momAdvection
555 ENDIF
556
557 C-- Set du/dt & dv/dt on boundaries to zero
558 DO j=jMin,jMax
559 DO i=iMin,iMax
560 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)*_maskW(i,j,k,bi,bj)
561 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)*_maskS(i,j,k,bi,bj)
562 ENDDO
563 ENDDO
564
565 #ifdef ALLOW_DEBUG
566 IF ( debugLevel .GE. debLevB
567 & .AND. k.EQ.4 .AND. myIter.EQ.nIter0
568 & .AND. nPx.EQ.1 .AND. nPy.EQ.1
569 & .AND. useCubedSphereExchange ) THEN
570 CALL DEBUG_CS_CORNER_UV( ' uDiss,vDiss from MOM_VECINV',
571 & guDiss,gvDiss, k, standardMessageUnit,bi,bj,myThid )
572 ENDIF
573 #endif /* ALLOW_DEBUG */
574
575 IF ( writeDiag ) THEN
576 IF (snapshot_mdsio) THEN
577 CALL WRITE_LOCAL_RL('Ds','I10',1,strain,bi,bj,k,myIter,myThid)
578 CALL WRITE_LOCAL_RL('Dt','I10',1,tension,bi,bj,k,myIter,
579 & myThid)
580 CALL WRITE_LOCAL_RL('Du','I10',1,guDiss,bi,bj,k,myIter,myThid)
581 CALL WRITE_LOCAL_RL('Dv','I10',1,gvDiss,bi,bj,k,myIter,myThid)
582 CALL WRITE_LOCAL_RL('Z3','I10',1,vort3,bi,bj,k,myIter,myThid)
583 CALL WRITE_LOCAL_RL('W3','I10',1,omega3,bi,bj,k,myIter,myThid)
584 CALL WRITE_LOCAL_RL('KE','I10',1,KE,bi,bj,k,myIter,myThid)
585 CALL WRITE_LOCAL_RL('D','I10',1,hDiv,bi,bj,k,myIter,myThid)
586 ENDIF
587 #ifdef ALLOW_MNC
588 IF (useMNC .AND. snapshot_mnc) THEN
589 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Ds',strain,
590 & offsets, myThid)
591 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Dt',tension,
592 & offsets, myThid)
593 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Du',guDiss,
594 & offsets, myThid)
595 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Dv',gvDiss,
596 & offsets, myThid)
597 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Z3',vort3,
598 & offsets, myThid)
599 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'W3',omega3,
600 & offsets, myThid)
601 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'KE',KE,
602 & offsets, myThid)
603 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'D', hDiv,
604 & offsets, myThid)
605 ENDIF
606 #endif /* ALLOW_MNC */
607 ENDIF
608
609 #ifdef ALLOW_DIAGNOSTICS
610 IF ( useDiagnostics ) THEN
611 CALL DIAGNOSTICS_FILL(KE, 'momKE ',k,1,2,bi,bj,myThid)
612 CALL DIAGNOSTICS_FILL(hDiv, 'momHDiv ',k,1,2,bi,bj,myThid)
613 CALL DIAGNOSTICS_FILL(vort3, 'momVort3',k,1,2,bi,bj,myThid)
614 CALL DIAGNOSTICS_FILL(gU(1-Olx,1-Oly,k,bi,bj),
615 & 'Um_Advec',k,1,2,bi,bj,myThid)
616 CALL DIAGNOSTICS_FILL(gV(1-Olx,1-Oly,k,bi,bj),
617 & 'Vm_Advec',k,1,2,bi,bj,myThid)
618 IF (momViscosity) THEN
619 CALL DIAGNOSTICS_FILL(strain,'Strain ',k,1,2,bi,bj,myThid)
620 CALL DIAGNOSTICS_FILL(tension,'Tension ',k,1,2,bi,bj,myThid)
621 CALL DIAGNOSTICS_FILL(guDiss,'Um_Diss ',k,1,2,bi,bj,myThid)
622 CALL DIAGNOSTICS_FILL(gvDiss,'Vm_Diss ',k,1,2,bi,bj,myThid)
623 ENDIF
624 ENDIF
625 #endif /* ALLOW_DIAGNOSTICS */
626
627 #endif /* ALLOW_MOM_VECINV */
628
629 RETURN
630 END

  ViewVC Help
Powered by ViewVC 1.1.22