/[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.40 - (show annotations) (download)
Thu Jun 9 15:57:45 2005 UTC (18 years, 11 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint57i_post
Changes since 1.39: +7 -3 lines
uncomment calls to MOM_VI_U/V_CORIOLIS_C4 (if highOrderVorticity=T)

1 C $Header: /u/gcmpack/MITgcm/pkg/mom_vecinv/mom_vecinv.F,v 1.39 2005/05/25 04:03:09 edhill 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 dPhiHydX,dPhiHydY,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 dPhiHydX,Y :: Gradient (X & Y dir.) of Hydrostatic Potential
46 C guDiss :: dissipation tendency (all explicit terms), u component
47 C gvDiss :: dissipation tendency (all explicit terms), v component
48 C bi, bj, iMin, iMax, jMin, jMax - Range of points for which calculation
49 C results will be set.
50 C kUp, kDown - Index for upper and lower layers.
51 C myThid - Instance number for this innvocation of CALC_MOM_RHS
52 _RL dPhiHydX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
53 _RL dPhiHydY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
54 _RL KappaRU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
55 _RL KappaRV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
56 _RL fVerU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
57 _RL fVerV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
58 _RL guDiss(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
59 _RL gvDiss(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
60 INTEGER kUp,kDown
61 _RL myTime
62 INTEGER myIter
63 INTEGER myThid
64 INTEGER bi,bj,iMin,iMax,jMin,jMax
65
66 #ifdef ALLOW_MOM_VECINV
67
68 C == Functions ==
69 LOGICAL DIFFERENT_MULTIPLE
70 EXTERNAL DIFFERENT_MULTIPLE
71
72 C == Local variables ==
73 _RL vF (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
74 _RL vrF (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
75 _RL uCf (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
76 _RL vCf (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
77 c _RL mT (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
78 _RL del2u(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
79 _RL del2v(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
80 _RL tension(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
81 _RL strain(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
82 _RS hFacZ(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
83 _RS r_hFacZ(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
84 _RL uFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
85 _RL vFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
86 _RL dStar(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
87 _RL zStar(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
88 C I,J,K - Loop counters
89 INTEGER i,j,k
90 C xxxFac - On-off tracer parameters used for switching terms off.
91 _RL ArDudrFac
92 _RL phxFac
93 c _RL mtFacU
94 _RL ArDvdrFac
95 _RL phyFac
96 c _RL mtFacV
97 LOGICAL bottomDragTerms
98 LOGICAL writeDiag
99 _RL KE(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
100 _RL omega3(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
101 _RL vort3(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
102 _RL hDiv(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
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 #ifdef ALLOW_AUTODIFF_TAMC
153 strain(i,j) = 0. _d 0
154 tension(i,j) = 0. _d 0
155 #endif
156 ENDDO
157 ENDDO
158
159 C-- Term by term tracer parmeters
160 C o U momentum equation
161 ArDudrFac = vfFacMom*1.
162 c mTFacU = mtFacMom*1.
163 phxFac = pfFacMom*1.
164 C o V momentum equation
165 ArDvdrFac = vfFacMom*1.
166 c mTFacV = mtFacMom*1.
167 phyFac = pfFacMom*1.
168
169 IF ( no_slip_bottom
170 & .OR. bottomDragQuadratic.NE.0.
171 & .OR. bottomDragLinear.NE.0.) THEN
172 bottomDragTerms=.TRUE.
173 ELSE
174 bottomDragTerms=.FALSE.
175 ENDIF
176
177 C-- with stagger time stepping, grad Phi_Hyp is directly incoporated in TIMESTEP
178 IF (staggerTimeStep) THEN
179 phxFac = 0.
180 phyFac = 0.
181 ENDIF
182
183 C-- Calculate open water fraction at vorticity points
184 CALL MOM_CALC_HFACZ(bi,bj,k,hFacZ,r_hFacZ,myThid)
185
186 C Make local copies of horizontal flow field
187 DO j=1-OLy,sNy+OLy
188 DO i=1-OLx,sNx+OLx
189 uFld(i,j) = uVel(i,j,k,bi,bj)
190 vFld(i,j) = vVel(i,j,k,bi,bj)
191 ENDDO
192 ENDDO
193
194 C note (jmc) : Dissipation and Vort3 advection do not necesary
195 C use the same maskZ (and hFacZ) => needs 2 call(s)
196 c CALL MOM_VI_HFACZ_DISS(bi,bj,k,hFacZ,r_hFacZ,myThid)
197
198 CALL MOM_CALC_KE(bi,bj,k,2,uFld,vFld,KE,myThid)
199
200 CALL MOM_CALC_HDIV(bi,bj,k,2,uFld,vFld,hDiv,myThid)
201
202 CALL MOM_CALC_RELVORT3(bi,bj,k,uFld,vFld,hFacZ,vort3,myThid)
203
204 IF (useAbsVorticity)
205 & CALL MOM_CALC_ABSVORT3(bi,bj,k,vort3,omega3,myThid)
206
207 IF (momViscosity) THEN
208 C Calculate del^2 u and del^2 v for bi-harmonic term
209 IF ( (viscA4.NE.0. .AND. no_slip_sides)
210 & .OR. viscA4D.NE.0. .OR. viscA4Z.NE.0.
211 & .OR. viscA4Grid.NE.0.
212 & .OR. viscC4leith.NE.0.
213 & .OR. viscC4leithD.NE.0.
214 & ) THEN
215 CALL MOM_VI_DEL2UV(bi,bj,k,hDiv,vort3,hFacZ,
216 O del2u,del2v,
217 & myThid)
218 CALL MOM_CALC_HDIV(bi,bj,k,2,del2u,del2v,dStar,myThid)
219 CALL MOM_CALC_RELVORT3(
220 & bi,bj,k,del2u,del2v,hFacZ,zStar,myThid)
221 ENDIF
222 C Calculate dissipation terms for U and V equations
223 C in terms of vorticity and divergence
224 IF ( viscAhD.NE.0. .OR. viscAhZ.NE.0.
225 & .OR. viscA4D.NE.0. .OR. viscA4Z.NE.0.
226 & .OR. viscAhGrid.NE.0. .OR. viscA4Grid.NE.0.
227 & .OR. viscC2leith.NE.0. .OR. viscC4leith.NE.0.
228 & .OR. viscC2leithD.NE.0. .OR. viscC4leithD.NE.0.
229 & ) THEN
230 CALL MOM_VI_HDISSIP(bi,bj,k,hDiv,vort3,hFacZ,dStar,zStar,
231 O guDiss,gvDiss,
232 & myThid)
233 ENDIF
234 C or in terms of tension and strain
235 IF (viscAstrain.NE.0. .OR. viscAtension.NE.0.
236 O .OR. viscC2smag.ne.0) THEN
237 CALL MOM_CALC_TENSION(bi,bj,k,uFld,vFld,
238 O tension,
239 I myThid)
240 CALL MOM_CALC_STRAIN(bi,bj,k,uFld,vFld,hFacZ,
241 O strain,
242 I myThid)
243 CALL MOM_HDISSIP(bi,bj,k,
244 I tension,strain,hFacZ,viscAtension,viscAstrain,
245 O guDiss,gvDiss,
246 I myThid)
247 ENDIF
248 ENDIF
249
250 C- Return to standard hfacZ (min-4) and mask vort3 accordingly:
251 c CALL MOM_VI_MASK_VORT3(bi,bj,k,hFacZ,r_hFacZ,vort3,myThid)
252
253 C---- Zonal momentum equation starts here
254
255 C-- Vertical flux (fVer is at upper face of "u" cell)
256
257 C Eddy component of vertical flux (interior component only) -> vrF
258 IF (momViscosity.AND..NOT.implicitViscosity) THEN
259 CALL MOM_U_RVISCFLUX(bi,bj,k,uVel,KappaRU,vrF,myThid)
260
261 C Combine fluxes
262 DO j=jMin,jMax
263 DO i=iMin,iMax
264 fVerU(i,j,kDown) = ArDudrFac*vrF(i,j)
265 ENDDO
266 ENDDO
267
268 C-- Tendency is minus divergence of the fluxes
269 DO j=2-Oly,sNy+Oly-1
270 DO i=2-Olx,sNx+Olx-1
271 guDiss(i,j) = guDiss(i,j)
272 & -_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)
273 & *recip_rAw(i,j,bi,bj)
274 & *(
275 & +fVerU(i,j,kUp)*rkFac - fVerU(i,j,kDown)*rkFac
276 & )
277 ENDDO
278 ENDDO
279 ENDIF
280
281 C-- No-slip and drag BCs appear as body forces in cell abutting topography
282 IF (momViscosity.AND.no_slip_sides) THEN
283 C- No-slip BCs impose a drag at walls...
284 CALL MOM_U_SIDEDRAG(bi,bj,k,uFld,del2u,hFacZ,vF,myThid)
285 DO j=jMin,jMax
286 DO i=iMin,iMax
287 guDiss(i,j) = guDiss(i,j)+vF(i,j)
288 ENDDO
289 ENDDO
290 ENDIF
291
292 C- No-slip BCs impose a drag at bottom
293 IF (momViscosity.AND.bottomDragTerms) THEN
294 CALL MOM_U_BOTTOMDRAG(bi,bj,k,uFld,KE,KappaRU,vF,myThid)
295 DO j=jMin,jMax
296 DO i=iMin,iMax
297 guDiss(i,j) = guDiss(i,j)+vF(i,j)
298 ENDDO
299 ENDDO
300 ENDIF
301
302 C-- Metric terms for curvilinear grid systems
303 c IF (usingSphericalPolarMTerms) THEN
304 C o Spherical polar grid metric terms
305 c CALL MOM_U_METRIC_NH(bi,bj,k,uFld,wVel,mT,myThid)
306 c DO j=jMin,jMax
307 c DO i=iMin,iMax
308 c gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+mTFacU*mT(i,j)
309 c ENDDO
310 c ENDDO
311 c ENDIF
312
313 C---- Meridional momentum equation starts here
314
315 C-- Vertical flux (fVer is at upper face of "v" cell)
316
317 C Eddy component of vertical flux (interior component only) -> vrF
318 IF (momViscosity.AND..NOT.implicitViscosity) THEN
319 CALL MOM_V_RVISCFLUX(bi,bj,k,vVel,KappaRV,vrf,myThid)
320
321 C Combine fluxes -> fVerV
322 DO j=jMin,jMax
323 DO i=iMin,iMax
324 fVerV(i,j,kDown) = ArDvdrFac*vrF(i,j)
325 ENDDO
326 ENDDO
327
328 C-- Tendency is minus divergence of the fluxes
329 DO j=jMin,jMax
330 DO i=iMin,iMax
331 gvDiss(i,j) = gvDiss(i,j)
332 & -_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
333 & *recip_rAs(i,j,bi,bj)
334 & *(
335 & +fVerV(i,j,kUp)*rkFac - fVerV(i,j,kDown)*rkFac
336 & )
337 ENDDO
338 ENDDO
339 ENDIF
340
341 C-- No-slip and drag BCs appear as body forces in cell abutting topography
342 IF (momViscosity.AND.no_slip_sides) THEN
343 C- No-slip BCs impose a drag at walls...
344 CALL MOM_V_SIDEDRAG(bi,bj,k,vFld,del2v,hFacZ,vF,myThid)
345 DO j=jMin,jMax
346 DO i=iMin,iMax
347 gvDiss(i,j) = gvDiss(i,j)+vF(i,j)
348 ENDDO
349 ENDDO
350 ENDIF
351 C- No-slip BCs impose a drag at bottom
352 IF (momViscosity.AND.bottomDragTerms) THEN
353 CALL MOM_V_BOTTOMDRAG(bi,bj,k,vFld,KE,KappaRV,vF,myThid)
354 DO j=jMin,jMax
355 DO i=iMin,iMax
356 gvDiss(i,j) = gvDiss(i,j)+vF(i,j)
357 ENDDO
358 ENDDO
359 ENDIF
360
361 C-- Metric terms for curvilinear grid systems
362 c IF (usingSphericalPolarMTerms) THEN
363 C o Spherical polar grid metric terms
364 c CALL MOM_V_METRIC_NH(bi,bj,k,vFld,wVel,mT,myThid)
365 c DO j=jMin,jMax
366 c DO i=iMin,iMax
367 c gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+mTFacV*mT(i,j)
368 c ENDDO
369 c ENDDO
370 c ENDIF
371
372 C-- Horizontal Coriolis terms
373 c IF (useCoriolis .AND. .NOT.useCDscheme
374 c & .AND. .NOT. useAbsVorticity) THEN
375 C- jmc: change it to keep the Coriolis terms when useAbsVorticity=T & momAdvection=F
376 IF ( useCoriolis .AND.
377 & .NOT.( useCDscheme .OR. useAbsVorticity.AND.momAdvection )
378 & ) THEN
379 IF (useAbsVorticity) THEN
380 CALL MOM_VI_U_CORIOLIS(bi,bj,K,vFld,omega3,hFacZ,r_hFacZ,
381 & uCf,myThid)
382 CALL MOM_VI_V_CORIOLIS(bi,bj,K,uFld,omega3,hFacZ,r_hFacZ,
383 & vCf,myThid)
384 ELSE
385 CALL MOM_VI_CORIOLIS(bi,bj,k,uFld,vFld,hFacZ,r_hFacZ,
386 & uCf,vCf,myThid)
387 ENDIF
388 DO j=jMin,jMax
389 DO i=iMin,iMax
390 gU(i,j,k,bi,bj) = uCf(i,j) - phxFac*dPhiHydX(i,j)
391 gV(i,j,k,bi,bj) = vCf(i,j) - phyFac*dPhiHydY(i,j)
392 ENDDO
393 ENDDO
394 IF ( writeDiag ) THEN
395 IF (snapshot_mdsio) THEN
396 CALL WRITE_LOCAL_RL('fV','I10',1,uCf,bi,bj,k,myIter,myThid)
397 CALL WRITE_LOCAL_RL('fU','I10',1,vCf,bi,bj,k,myIter,myThid)
398 ENDIF
399 #ifdef ALLOW_MNC
400 IF (useMNC .AND. snapshot_mnc) THEN
401 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'fV', uCf,
402 & offsets, myThid)
403 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'fU', vCf,
404 & offsets, myThid)
405 ENDIF
406 #endif /* ALLOW_MNC */
407 ENDIF
408 ELSE
409 DO j=jMin,jMax
410 DO i=iMin,iMax
411 gU(i,j,k,bi,bj) = -phxFac*dPhiHydX(i,j)
412 gV(i,j,k,bi,bj) = -phyFac*dPhiHydY(i,j)
413 ENDDO
414 ENDDO
415 ENDIF
416
417 IF (momAdvection) THEN
418 C-- Horizontal advection of relative vorticity
419 IF (useAbsVorticity) THEN
420 CALL MOM_VI_U_CORIOLIS(bi,bj,K,vFld,omega3,hFacZ,r_hFacZ,
421 & uCf,myThid)
422 ELSEIF (highOrderVorticity) THEN
423 CALL MOM_VI_U_CORIOLIS_C4(bi,bj,k,vFld,vort3,r_hFacZ,
424 & uCf,myThid)
425 ELSE
426 CALL MOM_VI_U_CORIOLIS(bi,bj,k,vFld,vort3,hFacZ,r_hFacZ,
427 & uCf,myThid)
428 ENDIF
429 DO j=jMin,jMax
430 DO i=iMin,iMax
431 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
432 ENDDO
433 ENDDO
434 IF (useAbsVorticity) THEN
435 CALL MOM_VI_V_CORIOLIS(bi,bj,K,uFld,omega3,hFacZ,r_hFacZ,
436 & vCf,myThid)
437 ELSEIF (highOrderVorticity) THEN
438 CALL MOM_VI_V_CORIOLIS_C4(bi,bj,K,uFld,vort3,r_hFacZ,
439 & vCf,myThid)
440 ELSE
441 CALL MOM_VI_V_CORIOLIS(bi,bj,k,uFld,vort3,hFacZ,r_hFacZ,
442 & vCf,myThid)
443 ENDIF
444 DO j=jMin,jMax
445 DO i=iMin,iMax
446 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
447 ENDDO
448 ENDDO
449
450 IF ( writeDiag ) THEN
451 IF (snapshot_mdsio) THEN
452 CALL WRITE_LOCAL_RL('zV','I10',1,uCf,bi,bj,k,myIter,myThid)
453 CALL WRITE_LOCAL_RL('zU','I10',1,vCf,bi,bj,k,myIter,myThid)
454 ENDIF
455 #ifdef ALLOW_MNC
456 IF (useMNC .AND. snapshot_mnc) THEN
457 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'zV', uCf,
458 & offsets, myThid)
459 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'zU', vCf,
460 & offsets, myThid)
461 ENDIF
462 #endif /* ALLOW_MNC */
463 ENDIF
464
465 #ifdef ALLOW_TIMEAVE
466 #ifndef MINIMAL_TAVE_OUTPUT
467 IF (taveFreq.GT.0.) THEN
468 CALL TIMEAVE_CUMUL_1K1T(uZetatave,vCf,deltaTClock,
469 & Nr, k, bi, bj, myThid)
470 CALL TIMEAVE_CUMUL_1K1T(vZetatave,uCf,deltaTClock,
471 & Nr, k, bi, bj, myThid)
472 ENDIF
473 #endif /* ndef MINIMAL_TAVE_OUTPUT */
474 #endif /* ALLOW_TIMEAVE */
475
476 C-- Vertical shear terms (-w*du/dr & -w*dv/dr)
477 IF ( .NOT. momImplVertAdv ) THEN
478 CALL MOM_VI_U_VERTSHEAR(bi,bj,K,uVel,wVel,uCf,myThid)
479 DO j=jMin,jMax
480 DO i=iMin,iMax
481 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
482 ENDDO
483 ENDDO
484 CALL MOM_VI_V_VERTSHEAR(bi,bj,K,vVel,wVel,vCf,myThid)
485 DO j=jMin,jMax
486 DO i=iMin,iMax
487 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
488 ENDDO
489 ENDDO
490 ENDIF
491
492 C-- Bernoulli term
493 CALL MOM_VI_U_GRAD_KE(bi,bj,K,KE,uCf,myThid)
494 DO j=jMin,jMax
495 DO i=iMin,iMax
496 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
497 ENDDO
498 ENDDO
499 CALL MOM_VI_V_GRAD_KE(bi,bj,K,KE,vCf,myThid)
500 DO j=jMin,jMax
501 DO i=iMin,iMax
502 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
503 ENDDO
504 ENDDO
505 IF ( writeDiag ) THEN
506 IF (snapshot_mdsio) THEN
507 CALL WRITE_LOCAL_RL('KEx','I10',1,uCf,bi,bj,k,myIter,myThid)
508 CALL WRITE_LOCAL_RL('KEy','I10',1,vCf,bi,bj,k,myIter,myThid)
509 ENDIF
510 #ifdef ALLOW_MNC
511 IF (useMNC .AND. snapshot_mnc) THEN
512 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'KEx', uCf,
513 & offsets, myThid)
514 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj, 'KEy', vCf,
515 & offsets, myThid)
516 ENDIF
517 #endif /* ALLOW_MNC */
518 ENDIF
519
520 C-- end if momAdvection
521 ENDIF
522
523 C-- Set du/dt & dv/dt on boundaries to zero
524 DO j=jMin,jMax
525 DO i=iMin,iMax
526 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)*_maskW(i,j,k,bi,bj)
527 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)*_maskS(i,j,k,bi,bj)
528 ENDDO
529 ENDDO
530
531 #ifdef ALLOW_DEBUG
532 IF ( debugLevel .GE. debLevB
533 & .AND. k.EQ.4 .AND. myIter.EQ.nIter0
534 & .AND. nPx.EQ.1 .AND. nPy.EQ.1
535 & .AND. useCubedSphereExchange ) THEN
536 CALL DEBUG_CS_CORNER_UV( ' uDiss,vDiss from MOM_VECINV',
537 & guDiss,gvDiss, k, standardMessageUnit,bi,bj,myThid )
538 ENDIF
539 #endif /* ALLOW_DEBUG */
540
541 IF ( writeDiag ) THEN
542 IF (snapshot_mdsio) THEN
543 CALL WRITE_LOCAL_RL('Ds','I10',1,strain,bi,bj,k,myIter,myThid)
544 CALL WRITE_LOCAL_RL('Dt','I10',1,tension,bi,bj,k,myIter,
545 & myThid)
546 CALL WRITE_LOCAL_RL('Du','I10',1,guDiss,bi,bj,k,myIter,myThid)
547 CALL WRITE_LOCAL_RL('Dv','I10',1,gvDiss,bi,bj,k,myIter,myThid)
548 CALL WRITE_LOCAL_RL('Z3','I10',1,vort3,bi,bj,k,myIter,myThid)
549 CALL WRITE_LOCAL_RL('W3','I10',1,omega3,bi,bj,k,myIter,myThid)
550 CALL WRITE_LOCAL_RL('KE','I10',1,KE,bi,bj,k,myIter,myThid)
551 CALL WRITE_LOCAL_RL('D','I10',1,hdiv,bi,bj,k,myIter,myThid)
552 ENDIF
553 #ifdef ALLOW_MNC
554 IF (useMNC .AND. snapshot_mnc) THEN
555 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Ds',strain,
556 & offsets, myThid)
557 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Dt',tension,
558 & offsets, myThid)
559 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Du',guDiss,
560 & offsets, myThid)
561 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Dv',gvDiss,
562 & offsets, myThid)
563 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'Z3',vort3,
564 & offsets, myThid)
565 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'W3',omega3,
566 & offsets, myThid)
567 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'KE',KE,
568 & offsets, myThid)
569 CALL MNC_CW_RL_W_OFFSET('D','mom_vi',bi,bj,'D', hdiv,
570 & offsets, myThid)
571 ENDIF
572 #endif /* ALLOW_MNC */
573 ENDIF
574
575 #endif /* ALLOW_MOM_VECINV */
576
577 RETURN
578 END

  ViewVC Help
Powered by ViewVC 1.1.22