/[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.21 - (show annotations) (download)
Tue Jul 20 17:46:38 2004 UTC (19 years, 10 months ago) by adcroft
Branch: MAIN
CVS Tags: checkpoint54d_post, checkpoint54e_post, checkpoint55, checkpoint54f_post
Changes since 1.20: +2 -3 lines
Replaced CPP_OPTIONS.h with MOM_VECINV_OPTIONS
 - this was causing a call to diagnostics to not happen

1 C $Header: /u/gcmpack/MITgcm/pkg/mom_vecinv/mom_vecinv.F,v 1.20 2004/06/02 13:23:55 adcroft 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 I myTime, myIter, myThid)
11 C /==========================================================\
12 C | S/R MOM_VECINV |
13 C | o Form the right hand-side of the momentum equation. |
14 C |==========================================================|
15 C | Terms are evaluated one layer at a time working from |
16 C | the bottom to the top. The vertically integrated |
17 C | barotropic flow tendency term is evluated by summing the |
18 C | tendencies. |
19 C | Notes: |
20 C | We have not sorted out an entirely satisfactory formula |
21 C | for the diffusion equation bc with lopping. The present |
22 C | form produces a diffusive flux that does not scale with |
23 C | open-area. Need to do something to solidfy this and to |
24 C | deal "properly" with thin walls. |
25 C \==========================================================/
26 IMPLICIT NONE
27
28 C == Global variables ==
29 #include "SIZE.h"
30 #include "DYNVARS.h"
31 #include "EEPARAMS.h"
32 #include "PARAMS.h"
33 #include "GRID.h"
34 #ifdef ALLOW_TIMEAVE
35 #include "TIMEAVE_STATV.h"
36 #endif
37
38 C == Routine arguments ==
39 C fVerU - Flux of momentum in the vertical
40 C fVerV direction out of the upper face of a cell K
41 C ( flux into the cell above ).
42 C dPhiHydX,Y :: Gradient (X & Y dir.) of Hydrostatic Potential
43 C bi, bj, iMin, iMax, jMin, jMax - Range of points for which calculation
44 C results will be set.
45 C kUp, kDown - Index for upper and lower layers.
46 C myThid - Instance number for this innvocation of CALC_MOM_RHS
47 _RL dPhiHydX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
48 _RL dPhiHydY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
49 _RL KappaRU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
50 _RL KappaRV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
51 _RL fVerU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
52 _RL fVerV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
53 INTEGER kUp,kDown
54 _RL myTime
55 INTEGER myIter
56 INTEGER myThid
57 INTEGER bi,bj,iMin,iMax,jMin,jMax
58
59 #ifdef ALLOW_MOM_VECINV
60
61 C == Functions ==
62 LOGICAL DIFFERENT_MULTIPLE
63 EXTERNAL DIFFERENT_MULTIPLE
64
65 C == Local variables ==
66 _RL aF (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
67 _RL vF (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
68 _RL vrF (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
69 _RL uCf (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
70 _RL vCf (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
71 _RL mT (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
72 _RL pF (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
73 _RL del2u(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
74 _RL del2v(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
75 _RL tension(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
76 _RL strain(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
77 _RS hFacZ(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
78 _RS r_hFacZ(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
79 _RS xA(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
80 _RS yA(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 _RL uDiss(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
86 _RL vDiss(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
87 C I,J,K - Loop counters
88 INTEGER i,j,k
89 C rVelMaskOverride - Factor for imposing special surface boundary conditions
90 C ( set according to free-surface condition ).
91 C hFacROpen - Lopped cell factos used tohold fraction of open
92 C hFacRClosed and closed cell wall.
93 _RL rVelMaskOverride
94 C xxxFac - On-off tracer parameters used for switching terms off.
95 _RL uDudxFac
96 _RL AhDudxFac
97 _RL A4DuxxdxFac
98 _RL vDudyFac
99 _RL AhDudyFac
100 _RL A4DuyydyFac
101 _RL rVelDudrFac
102 _RL ArDudrFac
103 _RL fuFac
104 _RL phxFac
105 _RL mtFacU
106 _RL uDvdxFac
107 _RL AhDvdxFac
108 _RL A4DvxxdxFac
109 _RL vDvdyFac
110 _RL AhDvdyFac
111 _RL A4DvyydyFac
112 _RL rVelDvdrFac
113 _RL ArDvdrFac
114 _RL fvFac
115 _RL phyFac
116 _RL vForcFac
117 _RL mtFacV
118 _RL wVelBottomOverride
119 LOGICAL bottomDragTerms
120 LOGICAL writeDiag
121 _RL KE(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
122 _RL omega3(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
123 _RL vort3(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
124 _RL hDiv(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
125
126 #ifdef ALLOW_AUTODIFF_TAMC
127 C-- only the kDown part of fverU/V is set in this subroutine
128 C-- the kUp is still required
129 C-- In the case of mom_fluxform Kup is set as well
130 C-- (at least in part)
131 fVerU(1,1,kUp) = fVerU(1,1,kUp)
132 fVerV(1,1,kUp) = fVerV(1,1,kUp)
133 #endif
134
135 rVelMaskOverride=1.
136 IF ( k .EQ. 1 ) rVelMaskOverride=freeSurfFac
137 wVelBottomOverride=1.
138 IF (k.EQ.Nr) wVelBottomOverride=0.
139 writeDiag = DIFFERENT_MULTIPLE(diagFreq, myTime,
140 & myTime-deltaTClock)
141
142 C Initialise intermediate terms
143 DO J=1-OLy,sNy+OLy
144 DO I=1-OLx,sNx+OLx
145 aF(i,j) = 0.
146 vF(i,j) = 0.
147 vrF(i,j) = 0.
148 uCf(i,j) = 0.
149 vCf(i,j) = 0.
150 mT(i,j) = 0.
151 pF(i,j) = 0.
152 del2u(i,j) = 0.
153 del2v(i,j) = 0.
154 dStar(i,j) = 0.
155 zStar(i,j) = 0.
156 uDiss(i,j) = 0.
157 vDiss(i,j) = 0.
158 vort3(i,j) = 0.
159 omega3(i,j) = 0.
160 ke(i,j) = 0.
161 #ifdef ALLOW_AUTODIFF_TAMC
162 strain(i,j) = 0. _d 0
163 tension(i,j) = 0. _d 0
164 #endif
165 ENDDO
166 ENDDO
167
168 C-- Term by term tracer parmeters
169 C o U momentum equation
170 uDudxFac = afFacMom*1.
171 AhDudxFac = vfFacMom*1.
172 A4DuxxdxFac = vfFacMom*1.
173 vDudyFac = afFacMom*1.
174 AhDudyFac = vfFacMom*1.
175 A4DuyydyFac = vfFacMom*1.
176 rVelDudrFac = afFacMom*1.
177 ArDudrFac = vfFacMom*1.
178 mTFacU = mtFacMom*1.
179 fuFac = cfFacMom*1.
180 phxFac = pfFacMom*1.
181 C o V momentum equation
182 uDvdxFac = afFacMom*1.
183 AhDvdxFac = vfFacMom*1.
184 A4DvxxdxFac = vfFacMom*1.
185 vDvdyFac = afFacMom*1.
186 AhDvdyFac = vfFacMom*1.
187 A4DvyydyFac = vfFacMom*1.
188 rVelDvdrFac = afFacMom*1.
189 ArDvdrFac = vfFacMom*1.
190 mTFacV = mtFacMom*1.
191 fvFac = cfFacMom*1.
192 phyFac = pfFacMom*1.
193 vForcFac = foFacMom*1.
194
195 IF ( no_slip_bottom
196 & .OR. bottomDragQuadratic.NE.0.
197 & .OR. bottomDragLinear.NE.0.) THEN
198 bottomDragTerms=.TRUE.
199 ELSE
200 bottomDragTerms=.FALSE.
201 ENDIF
202
203 C-- with stagger time stepping, grad Phi_Hyp is directly incoporated in TIMESTEP
204 IF (staggerTimeStep) THEN
205 phxFac = 0.
206 phyFac = 0.
207 ENDIF
208
209 C-- Calculate open water fraction at vorticity points
210 CALL MOM_CALC_HFACZ(bi,bj,k,hFacZ,r_hFacZ,myThid)
211
212 C---- Calculate common quantities used in both U and V equations
213 C Calculate tracer cell face open areas
214 DO j=1-OLy,sNy+OLy
215 DO i=1-OLx,sNx+OLx
216 xA(i,j) = _dyG(i,j,bi,bj)
217 & *drF(k)*_hFacW(i,j,k,bi,bj)
218 yA(i,j) = _dxG(i,j,bi,bj)
219 & *drF(k)*_hFacS(i,j,k,bi,bj)
220 ENDDO
221 ENDDO
222
223 C Make local copies of horizontal flow field
224 DO j=1-OLy,sNy+OLy
225 DO i=1-OLx,sNx+OLx
226 uFld(i,j) = uVel(i,j,k,bi,bj)
227 vFld(i,j) = vVel(i,j,k,bi,bj)
228 ENDDO
229 ENDDO
230
231 C note (jmc) : Dissipation and Vort3 advection do not necesary
232 C use the same maskZ (and hFacZ) => needs 2 call(s)
233 c CALL MOM_VI_HFACZ_DISS(bi,bj,k,hFacZ,r_hFacZ,myThid)
234
235 CALL MOM_CALC_KE(bi,bj,k,2,uFld,vFld,KE,myThid)
236
237 CALL MOM_CALC_HDIV(bi,bj,k,2,uFld,vFld,hDiv,myThid)
238
239 CALL MOM_CALC_RELVORT3(bi,bj,k,uFld,vFld,hFacZ,vort3,myThid)
240
241 IF (useAbsVorticity)
242 & CALL MOM_CALC_ABSVORT3(bi,bj,k,vort3,omega3,myThid)
243
244 IF (momViscosity) THEN
245 C Calculate del^2 u and del^2 v for bi-harmonic term
246 IF (viscA4.NE.0.
247 & .OR. viscA4Grid.NE.0.
248 & .OR. viscC4leith.NE.0.
249 & ) THEN
250 CALL MOM_VI_DEL2UV(bi,bj,k,hDiv,vort3,hFacZ,
251 O del2u,del2v,
252 & myThid)
253 CALL MOM_CALC_HDIV(bi,bj,k,2,del2u,del2v,dStar,myThid)
254 CALL MOM_CALC_RELVORT3(
255 & bi,bj,k,del2u,del2v,hFacZ,zStar,myThid)
256 ENDIF
257 C Calculate dissipation terms for U and V equations
258 C in terms of vorticity and divergence
259 IF (viscAh.NE.0. .OR. viscA4.NE.0.
260 & .OR. viscAhGrid.NE.0. .OR. viscA4Grid.NE.0.
261 & .OR. viscC2leith.NE.0. .OR. viscC4leith.NE.0.
262 & ) THEN
263 CALL MOM_VI_HDISSIP(bi,bj,k,hDiv,vort3,hFacZ,dStar,zStar,
264 O uDiss,vDiss,
265 & myThid)
266 ENDIF
267 C or in terms of tension and strain
268 IF (viscAstrain.NE.0. .OR. viscAtension.NE.0.) THEN
269 CALL MOM_CALC_TENSION(bi,bj,k,uFld,vFld,
270 O tension,
271 I myThid)
272 CALL MOM_CALC_STRAIN(bi,bj,k,uFld,vFld,hFacZ,
273 O strain,
274 I myThid)
275 CALL MOM_HDISSIP(bi,bj,k,
276 I tension,strain,hFacZ,viscAtension,viscAstrain,
277 O uDiss,vDiss,
278 I myThid)
279 ENDIF
280 ENDIF
281
282 C- Return to standard hfacZ (min-4) and mask vort3 accordingly:
283 c CALL MOM_VI_MASK_VORT3(bi,bj,k,hFacZ,r_hFacZ,vort3,myThid)
284
285 C---- Zonal momentum equation starts here
286
287 C-- Vertical flux (fVer is at upper face of "u" cell)
288
289 C Eddy component of vertical flux (interior component only) -> vrF
290 IF (momViscosity.AND..NOT.implicitViscosity)
291 & CALL MOM_U_RVISCFLUX(bi,bj,k,uVel,KappaRU,vrF,myThid)
292
293 C Combine fluxes
294 DO j=jMin,jMax
295 DO i=iMin,iMax
296 fVerU(i,j,kDown) = ArDudrFac*vrF(i,j)
297 ENDDO
298 ENDDO
299
300 C-- Tendency is minus divergence of the fluxes + coriolis + pressure term
301 DO j=2-Oly,sNy+Oly-1
302 DO i=2-Olx,sNx+Olx-1
303 gU(i,j,k,bi,bj) = uDiss(i,j)
304 & -_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)
305 & *recip_rAw(i,j,bi,bj)
306 & *(
307 & +fVerU(i,j,kUp)*rkFac - fVerU(i,j,kDown)*rkFac
308 & )
309 & - phxFac*dPhiHydX(i,j)
310 ENDDO
311 ENDDO
312
313 C-- No-slip and drag BCs appear as body forces in cell abutting topography
314 IF (momViscosity.AND.no_slip_sides) THEN
315 C- No-slip BCs impose a drag at walls...
316 CALL MOM_U_SIDEDRAG(bi,bj,k,uFld,del2u,hFacZ,vF,myThid)
317 DO j=jMin,jMax
318 DO i=iMin,iMax
319 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+vF(i,j)
320 ENDDO
321 ENDDO
322 ENDIF
323
324 C- No-slip BCs impose a drag at bottom
325 IF (momViscosity.AND.bottomDragTerms) THEN
326 CALL MOM_U_BOTTOMDRAG(bi,bj,k,uFld,KE,KappaRU,vF,myThid)
327 DO j=jMin,jMax
328 DO i=iMin,iMax
329 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+vF(i,j)
330 ENDDO
331 ENDDO
332 ENDIF
333
334 C-- Metric terms for curvilinear grid systems
335 c IF (usingSphericalPolarMTerms) THEN
336 C o Spherical polar grid metric terms
337 c CALL MOM_U_METRIC_NH(bi,bj,k,uFld,wVel,mT,myThid)
338 c DO j=jMin,jMax
339 c DO i=iMin,iMax
340 c gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+mTFacU*mT(i,j)
341 c ENDDO
342 c ENDDO
343 c ENDIF
344
345 C---- Meridional momentum equation starts here
346
347 C-- Vertical flux (fVer is at upper face of "v" cell)
348
349 C Eddy component of vertical flux (interior component only) -> vrF
350 IF (momViscosity.AND..NOT.implicitViscosity)
351 & CALL MOM_V_RVISCFLUX(bi,bj,k,vVel,KappaRV,vrf,myThid)
352
353 C Combine fluxes -> fVerV
354 DO j=jMin,jMax
355 DO i=iMin,iMax
356 fVerV(i,j,kDown) = ArDvdrFac*vrF(i,j)
357 ENDDO
358 ENDDO
359
360 C-- Tendency is minus divergence of the fluxes + coriolis + pressure term
361 DO j=jMin,jMax
362 DO i=iMin,iMax
363 gV(i,j,k,bi,bj) = vDiss(i,j)
364 & -_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
365 & *recip_rAs(i,j,bi,bj)
366 & *(
367 & +fVerV(i,j,kUp)*rkFac - fVerV(i,j,kDown)*rkFac
368 & )
369 & - phyFac*dPhiHydY(i,j)
370 ENDDO
371 ENDDO
372
373 C-- No-slip and drag BCs appear as body forces in cell abutting topography
374 IF (momViscosity.AND.no_slip_sides) THEN
375 C- No-slip BCs impose a drag at walls...
376 CALL MOM_V_SIDEDRAG(bi,bj,k,vFld,del2v,hFacZ,vF,myThid)
377 DO j=jMin,jMax
378 DO i=iMin,iMax
379 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vF(i,j)
380 ENDDO
381 ENDDO
382 ENDIF
383 C- No-slip BCs impose a drag at bottom
384 IF (momViscosity.AND.bottomDragTerms) THEN
385 CALL MOM_V_BOTTOMDRAG(bi,bj,k,vFld,KE,KappaRV,vF,myThid)
386 DO j=jMin,jMax
387 DO i=iMin,iMax
388 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vF(i,j)
389 ENDDO
390 ENDDO
391 ENDIF
392
393 C-- Metric terms for curvilinear grid systems
394 c IF (usingSphericalPolarMTerms) THEN
395 C o Spherical polar grid metric terms
396 c CALL MOM_V_METRIC_NH(bi,bj,k,vFld,wVel,mT,myThid)
397 c DO j=jMin,jMax
398 c DO i=iMin,iMax
399 c gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+mTFacV*mT(i,j)
400 c ENDDO
401 c ENDDO
402 c ENDIF
403
404 C-- Horizontal Coriolis terms
405 IF (useCoriolis .AND. .NOT.useCDscheme
406 & .AND. .NOT. useAbsVorticity) THEN
407 CALL MOM_VI_CORIOLIS(bi,bj,k,uFld,vFld,hFacZ,r_hFacZ,
408 & uCf,vCf,myThid)
409 DO j=jMin,jMax
410 DO i=iMin,iMax
411 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
412 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
413 ENDDO
414 ENDDO
415 IF ( writeDiag ) THEN
416 CALL WRITE_LOCAL_RL('fV','I10',1,uCf,bi,bj,k,myIter,myThid)
417 CALL WRITE_LOCAL_RL('fU','I10',1,vCf,bi,bj,k,myIter,myThid)
418 ENDIF
419 ENDIF
420
421 IF (momAdvection) THEN
422 C-- Horizontal advection of relative vorticity
423 IF (useAbsVorticity) THEN
424 CALL MOM_VI_U_CORIOLIS(bi,bj,K,vFld,omega3,hFacZ,r_hFacZ,
425 & uCf,myThid)
426 ELSE
427 CALL MOM_VI_U_CORIOLIS(bi,bj,k,vFld,vort3,hFacZ,r_hFacZ,
428 & uCf,myThid)
429 ENDIF
430 c CALL MOM_VI_U_CORIOLIS_C4(bi,bj,K,vFld,vort3,r_hFacZ,uCf,myThid)
431 DO j=jMin,jMax
432 DO i=iMin,iMax
433 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
434 ENDDO
435 ENDDO
436 IF (useAbsVorticity) THEN
437 CALL MOM_VI_V_CORIOLIS(bi,bj,K,uFld,omega3,hFacZ,r_hFacZ,
438 & vCf,myThid)
439 ELSE
440 CALL MOM_VI_V_CORIOLIS(bi,bj,k,uFld,vort3,hFacZ,r_hFacZ,
441 & vCf,myThid)
442 ENDIF
443 c CALL MOM_VI_V_CORIOLIS_C4(bi,bj,K,uFld,vort3,r_hFacZ,vCf,myThid)
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 CALL WRITE_LOCAL_RL('zV','I10',1,uCf,bi,bj,k,myIter,myThid)
452 CALL WRITE_LOCAL_RL('zU','I10',1,vCf,bi,bj,k,myIter,myThid)
453 ENDIF
454 #ifdef ALLOW_TIMEAVE
455 #ifndef HRCUBE
456 IF (taveFreq.GT.0.) THEN
457 CALL TIMEAVE_CUMUL_1K1T(uZetatave,vCf,deltaTClock,
458 & Nr, k, bi, bj, myThid)
459 CALL TIMEAVE_CUMUL_1K1T(vZetatave,uCf,deltaTClock,
460 & Nr, k, bi, bj, myThid)
461 ENDIF
462 #endif /* ALLOW_TIMEAVE */
463 #endif /* ndef HRCUBE */
464
465 C-- Vertical shear terms (-w*du/dr & -w*dv/dr)
466 IF ( .NOT. momImplVertAdv ) THEN
467 CALL MOM_VI_U_VERTSHEAR(bi,bj,K,uVel,wVel,uCf,myThid)
468 DO j=jMin,jMax
469 DO i=iMin,iMax
470 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
471 ENDDO
472 ENDDO
473 CALL MOM_VI_V_VERTSHEAR(bi,bj,K,vVel,wVel,vCf,myThid)
474 DO j=jMin,jMax
475 DO i=iMin,iMax
476 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
477 ENDDO
478 ENDDO
479 ENDIF
480
481 C-- Bernoulli term
482 CALL MOM_VI_U_GRAD_KE(bi,bj,K,KE,uCf,myThid)
483 DO j=jMin,jMax
484 DO i=iMin,iMax
485 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
486 ENDDO
487 ENDDO
488 CALL MOM_VI_V_GRAD_KE(bi,bj,K,KE,vCf,myThid)
489 DO j=jMin,jMax
490 DO i=iMin,iMax
491 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
492 ENDDO
493 ENDDO
494 IF ( writeDiag ) THEN
495 CALL WRITE_LOCAL_RL('KEx','I10',1,uCf,bi,bj,k,myIter,myThid)
496 CALL WRITE_LOCAL_RL('KEy','I10',1,vCf,bi,bj,k,myIter,myThid)
497 ENDIF
498
499 C-- end if momAdvection
500 ENDIF
501
502 C-- Set du/dt & dv/dt on boundaries to zero
503 DO j=jMin,jMax
504 DO i=iMin,iMax
505 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)*_maskW(i,j,k,bi,bj)
506 gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)*_maskS(i,j,k,bi,bj)
507 ENDDO
508 ENDDO
509
510
511 IF ( writeDiag ) THEN
512 CALL WRITE_LOCAL_RL('Ds','I10',1,strain,bi,bj,k,myIter,myThid)
513 CALL WRITE_LOCAL_RL('Dt','I10',1,tension,bi,bj,k,myIter,myThid)
514 CALL WRITE_LOCAL_RL('Du','I10',1,uDiss,bi,bj,k,myIter,myThid)
515 CALL WRITE_LOCAL_RL('Dv','I10',1,vDiss,bi,bj,k,myIter,myThid)
516 CALL WRITE_LOCAL_RL('Z3','I10',1,vort3,bi,bj,k,myIter,myThid)
517 CALL WRITE_LOCAL_RL('W3','I10',1,omega3,bi,bj,k,myIter,myThid)
518 CALL WRITE_LOCAL_RL('KE','I10',1,KE,bi,bj,k,myIter,myThid)
519 CALL WRITE_LOCAL_RL('D','I10',1,hdiv,bi,bj,k,myIter,myThid)
520 ENDIF
521
522 #endif /* ALLOW_MOM_VECINV */
523
524 RETURN
525 END

  ViewVC Help
Powered by ViewVC 1.1.22