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

Annotation of /MITgcm/pkg/mom_vecinv/mom_vecinv.F

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


Revision 1.64 - (hide annotations) (download)
Mon Apr 25 02:03:13 2011 UTC (13 years, 2 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint62w, checkpoint62y, checkpoint62x
Changes since 1.63: +11 -1 lines
add output (function of diagFreq) of intermediate var. del2u,del2v,dStar & zStar

1 jmc 1.64 C $Header: /u/gcmpack/MITgcm/pkg/mom_vecinv/mom_vecinv.F,v 1.63 2010/03/16 00:16:50 jmc Exp $
2 adcroft 1.2 C $Name: $
3 adcroft 1.1
4 adcroft 1.21 #include "MOM_VECINV_OPTIONS.h"
5 adcroft 1.1
6 jmc 1.57 SUBROUTINE MOM_VECINV(
7 adcroft 1.1 I bi,bj,iMin,iMax,jMin,jMax,k,kUp,kDown,
8 jmc 1.43 I KappaRU, KappaRV,
9 adcroft 1.1 U fVerU, fVerV,
10 jmc 1.31 O guDiss, gvDiss,
11 jmc 1.15 I myTime, myIter, myThid)
12 adcroft 1.1 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 edhill 1.27 #ifdef ALLOW_MNC
35     #include "MNC_PARAMS.h"
36     #endif
37 adcroft 1.1 #include "GRID.h"
38 jmc 1.7 #ifdef ALLOW_TIMEAVE
39     #include "TIMEAVE_STATV.h"
40     #endif
41 heimbach 1.59 #ifdef ALLOW_AUTODIFF_TAMC
42     # include "tamc.h"
43     # include "tamc_keys.h"
44     #endif
45 adcroft 1.1
46     C == Routine arguments ==
47 jmc 1.31 C fVerU :: Flux of momentum in the vertical direction, out of the upper
48     C fVerV :: face of a cell K ( flux into the cell above ).
49     C guDiss :: dissipation tendency (all explicit terms), u component
50     C gvDiss :: dissipation tendency (all explicit terms), v component
51 adcroft 1.1 C bi, bj, iMin, iMax, jMin, jMax - Range of points for which calculation
52     C results will be set.
53     C kUp, kDown - Index for upper and lower layers.
54 jmc 1.54 C myThid :: my Thread Id number
55 adcroft 1.1 _RL KappaRU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
56     _RL KappaRV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
57     _RL fVerU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
58     _RL fVerV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
59 jmc 1.31 _RL guDiss(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
60     _RL gvDiss(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
61 adcroft 1.1 INTEGER kUp,kDown
62 jmc 1.15 _RL myTime
63 adcroft 1.2 INTEGER myIter
64 adcroft 1.1 INTEGER myThid
65     INTEGER bi,bj,iMin,iMax,jMin,jMax
66    
67 edhill 1.11 #ifdef ALLOW_MOM_VECINV
68 jmc 1.7
69 adcroft 1.2 C == Functions ==
70 jmc 1.38 LOGICAL DIFFERENT_MULTIPLE
71     EXTERNAL DIFFERENT_MULTIPLE
72 adcroft 1.2
73 adcroft 1.1 C == Local variables ==
74     _RL vF (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
75 jmc 1.54 _RL vrF(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
76     _RL uCf(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
77     _RL vCf(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
78     _RS hFacZ (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
79     _RS r_hFacZ (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
80     _RL uFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
81     _RL vFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
82     _RL del2u (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
83     _RL del2v (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
84     _RL dStar (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
85     _RL zStar (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
86     _RL tension (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
87     _RL strain (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
88     _RL KE (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
89     _RL omega3 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
90     _RL vort3 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
91     _RL hDiv (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
92     _RL viscAh_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
93     _RL viscAh_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
94     _RL viscA4_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
95     _RL viscA4_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
96     C i,j,k :: Loop counters
97 adcroft 1.1 INTEGER i,j,k
98     C xxxFac - On-off tracer parameters used for switching terms off.
99     _RL ArDudrFac
100     _RL ArDvdrFac
101 jmc 1.54 _RL sideMaskFac
102 adcroft 1.1 LOGICAL bottomDragTerms
103 jmc 1.15 LOGICAL writeDiag
104 baylor 1.50 LOGICAL harmonic,biharmonic,useVariableViscosity
105 heimbach 1.59 #ifdef ALLOW_AUTODIFF_TAMC
106     INTEGER imomkey
107     #endif
108 adcroft 1.1
109 edhill 1.25 #ifdef ALLOW_MNC
110     INTEGER offsets(9)
111 edhill 1.53 CHARACTER*(1) pf
112 edhill 1.25 #endif
113    
114 heimbach 1.9 #ifdef ALLOW_AUTODIFF_TAMC
115     C-- only the kDown part of fverU/V is set in this subroutine
116     C-- the kUp is still required
117     C-- In the case of mom_fluxform Kup is set as well
118     C-- (at least in part)
119     fVerU(1,1,kUp) = fVerU(1,1,kUp)
120     fVerV(1,1,kUp) = fVerV(1,1,kUp)
121     #endif
122    
123 heimbach 1.59 #ifdef ALLOW_AUTODIFF_TAMC
124     act0 = k - 1
125     max0 = Nr
126     act1 = bi - myBxLo(myThid)
127     max1 = myBxHi(myThid) - myBxLo(myThid) + 1
128     act2 = bj - myByLo(myThid)
129     max2 = myByHi(myThid) - myByLo(myThid) + 1
130     act3 = myThid - 1
131     max3 = nTx*nTy
132     act4 = ikey_dynamics - 1
133     imomkey = (act0 + 1)
134     & + act1*max0
135     & + act2*max0*max1
136     & + act3*max0*max1*max2
137     & + act4*max0*max1*max2*max3
138     #endif /* ALLOW_AUTODIFF_TAMC */
139    
140 jmc 1.38 writeDiag = DIFFERENT_MULTIPLE(diagFreq, myTime, deltaTClock)
141 adcroft 1.1
142 edhill 1.24 #ifdef ALLOW_MNC
143     IF (useMNC .AND. snapshot_mnc .AND. writeDiag) THEN
144 edhill 1.53 IF ( writeBinaryPrec .EQ. precFloat64 ) THEN
145     pf(1:1) = 'D'
146     ELSE
147     pf(1:1) = 'R'
148     ENDIF
149 edhill 1.25 IF ((bi .EQ. 1).AND.(bj .EQ. 1).AND.(k .EQ. 1)) THEN
150     CALL MNC_CW_SET_UDIM('mom_vi', -1, myThid)
151 edhill 1.39 CALL MNC_CW_RL_W_S('D','mom_vi',0,0,'T',myTime,myThid)
152 edhill 1.25 CALL MNC_CW_SET_UDIM('mom_vi', 0, myThid)
153 edhill 1.39 CALL MNC_CW_I_W_S('I','mom_vi',0,0,'iter',myIter,myThid)
154 edhill 1.25 ENDIF
155     DO i = 1,9
156     offsets(i) = 0
157     ENDDO
158     offsets(3) = k
159 jmc 1.61 c write(*,*) 'offsets = ',(offsets(i),i=1,9)
160 edhill 1.24 ENDIF
161     #endif /* ALLOW_MNC */
162    
163 jmc 1.61 C-- Initialise intermediate terms
164     DO j=1-OLy,sNy+OLy
165     DO i=1-OLx,sNx+OLx
166 jmc 1.31 vF(i,j) = 0.
167     vrF(i,j) = 0.
168 adcroft 1.1 uCf(i,j) = 0.
169     vCf(i,j) = 0.
170     del2u(i,j) = 0.
171     del2v(i,j) = 0.
172     dStar(i,j) = 0.
173     zStar(i,j) = 0.
174 jmc 1.31 guDiss(i,j)= 0.
175     gvDiss(i,j)= 0.
176 adcroft 1.1 vort3(i,j) = 0.
177 jmc 1.31 omega3(i,j)= 0.
178 jmc 1.54 KE(i,j) = 0.
179 jmc 1.61 C- need to initialise hDiv for MOM_VI_DEL2UV(call FILL_CS_CORNER_TR_RL)
180     hDiv(i,j) = 0.
181 baylor 1.50 viscAh_Z(i,j) = 0.
182     viscAh_D(i,j) = 0.
183     viscA4_Z(i,j) = 0.
184     viscA4_D(i,j) = 0.
185    
186 heimbach 1.8 strain(i,j) = 0. _d 0
187     tension(i,j) = 0. _d 0
188 jmc 1.60 #ifdef ALLOW_AUTODIFF_TAMC
189 heimbach 1.55 hFacZ(i,j) = 0. _d 0
190 heimbach 1.8 #endif
191 adcroft 1.1 ENDDO
192     ENDDO
193    
194     C-- Term by term tracer parmeters
195     C o U momentum equation
196     ArDudrFac = vfFacMom*1.
197     C o V momentum equation
198     ArDvdrFac = vfFacMom*1.
199    
200 jmc 1.54 C note: using standard stencil (no mask) results in under-estimating
201     C vorticity at a no-slip boundary by a factor of 2 = sideDragFactor
202     IF ( no_slip_sides ) THEN
203     sideMaskFac = sideDragFactor
204     ELSE
205     sideMaskFac = 0. _d 0
206     ENDIF
207    
208 adcroft 1.1 IF ( no_slip_bottom
209     & .OR. bottomDragQuadratic.NE.0.
210     & .OR. bottomDragLinear.NE.0.) THEN
211     bottomDragTerms=.TRUE.
212     ELSE
213     bottomDragTerms=.FALSE.
214     ENDIF
215    
216     C-- Calculate open water fraction at vorticity points
217     CALL MOM_CALC_HFACZ(bi,bj,k,hFacZ,r_hFacZ,myThid)
218    
219     C Make local copies of horizontal flow field
220     DO j=1-OLy,sNy+OLy
221     DO i=1-OLx,sNx+OLx
222     uFld(i,j) = uVel(i,j,k,bi,bj)
223     vFld(i,j) = vVel(i,j,k,bi,bj)
224     ENDDO
225     ENDDO
226    
227 jmc 1.7 C note (jmc) : Dissipation and Vort3 advection do not necesary
228     C use the same maskZ (and hFacZ) => needs 2 call(s)
229     c CALL MOM_VI_HFACZ_DISS(bi,bj,k,hFacZ,r_hFacZ,myThid)
230    
231 jmc 1.52 CALL MOM_CALC_KE(bi,bj,k,selectKEscheme,uFld,vFld,KE,myThid)
232 adcroft 1.1
233 jmc 1.54 CALL MOM_CALC_RELVORT3(bi,bj,k,uFld,vFld,hFacZ,vort3,myThid)
234 adcroft 1.1
235 jmc 1.54 IF (momViscosity) THEN
236 jmc 1.57 C-- For viscous term, compute horizontal divergence, tension & strain
237 jmc 1.54 C and mask relative vorticity (free-slip case):
238 adcroft 1.1
239 heimbach 1.59 #ifdef ALLOW_AUTODIFF_TAMC
240 jmc 1.61 CADJ STORE vort3(:,:) =
241 heimbach 1.59 CADJ & comlev1_bibj_k, key = imomkey, byte = isbyte
242     #endif
243    
244 jmc 1.54 CALL MOM_CALC_HDIV(bi,bj,k,2,uFld,vFld,hDiv,myThid)
245 adcroft 1.1
246 jmc 1.52 CALL MOM_CALC_TENSION(bi,bj,k,uFld,vFld,tension,myThid)
247    
248     CALL MOM_CALC_STRAIN(bi,bj,k,uFld,vFld,hFacZ,strain,myThid)
249    
250 jmc 1.54 C- account for no-slip / free-slip BC:
251     DO j=1-Oly,sNy+Oly
252     DO i=1-Olx,sNx+Olx
253     IF ( hFacZ(i,j).EQ.0. ) THEN
254     vort3(i,j) = sideMaskFac*vort3(i,j)
255     strain(i,j) = sideMaskFac*strain(i,j)
256     ENDIF
257     ENDDO
258     ENDDO
259    
260     C-- Calculate Viscosities
261 jmc 1.52 CALL MOM_CALC_VISC(
262 baylor 1.50 I bi,bj,k,
263     O viscAh_Z,viscAh_D,viscA4_Z,viscA4_D,
264     O harmonic,biharmonic,useVariableViscosity,
265     I hDiv,vort3,tension,strain,KE,hfacZ,
266     I myThid)
267    
268 adcroft 1.1 C Calculate del^2 u and del^2 v for bi-harmonic term
269 baylor 1.50 IF (biharmonic) THEN
270 adcroft 1.2 CALL MOM_VI_DEL2UV(bi,bj,k,hDiv,vort3,hFacZ,
271     O del2u,del2v,
272     & myThid)
273 jmc 1.48 CALL MOM_CALC_HDIV(bi,bj,k,2,del2u,del2v,dStar,myThid)
274     CALL MOM_CALC_RELVORT3(bi,bj,k,
275     & del2u,del2v,hFacZ,zStar,myThid)
276 jmc 1.64 IF ( writeDiag ) THEN
277     CALL WRITE_LOCAL_RL( 'del2u', 'I10', 1, del2u,
278     & bi,bj,k, myIter, myThid )
279     CALL WRITE_LOCAL_RL( 'del2v', 'I10', 1, del2v,
280     & bi,bj,k, myIter, myThid )
281     CALL WRITE_LOCAL_RL( 'dStar', 'I10', 1, dStar,
282     & bi,bj,k, myIter, myThid )
283     CALL WRITE_LOCAL_RL( 'zStar', 'I10', 1, zStar,
284     & bi,bj,k, myIter, myThid )
285     ENDIF
286 adcroft 1.2 ENDIF
287 baylor 1.47
288 jmc 1.54 C- Strain diagnostics:
289     IF ( writeDiag ) THEN
290     IF (snapshot_mdsio) THEN
291     CALL WRITE_LOCAL_RL('Ds','I10',1,strain,bi,bj,k,myIter,myThid)
292     ENDIF
293     #ifdef ALLOW_MNC
294     IF (useMNC .AND. snapshot_mnc) THEN
295     CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'Ds',strain,
296     & offsets, myThid)
297     ENDIF
298     #endif /* ALLOW_MNC */
299     ENDIF
300     #ifdef ALLOW_DIAGNOSTICS
301     IF ( useDiagnostics ) THEN
302     CALL DIAGNOSTICS_FILL(strain, 'Strain ',k,1,2,bi,bj,myThid)
303     ENDIF
304     #endif /* ALLOW_DIAGNOSTICS */
305    
306     C--- Calculate dissipation terms for U and V equations
307 baylor 1.47
308     C in terms of tension and strain
309     IF (useStrainTensionVisc) THEN
310 jmc 1.54 C mask strain as if free-slip since side-drag is computed separately
311     DO j=1-Oly,sNy+Oly
312     DO i=1-Olx,sNx+Olx
313     IF ( hFacZ(i,j).EQ.0. ) strain(i,j) = 0. _d 0
314     ENDDO
315     ENDDO
316 baylor 1.47 CALL MOM_HDISSIP(bi,bj,k,hDiv,vort3,tension,strain,KE,
317     I hFacZ,
318 baylor 1.50 I viscAh_Z,viscAh_D,viscA4_Z,viscA4_D,
319     I harmonic,biharmonic,useVariableViscosity,
320 baylor 1.47 O guDiss,gvDiss,
321     I myThid)
322     ELSE
323 adcroft 1.2 C in terms of vorticity and divergence
324 baylor 1.47 CALL MOM_VI_HDISSIP(bi,bj,k,hDiv,vort3,tension,strain,KE,
325     I hFacZ,dStar,zStar,
326 baylor 1.50 I viscAh_Z,viscAh_D,viscA4_Z,viscA4_D,
327     I harmonic,biharmonic,useVariableViscosity,
328 jmc 1.31 O guDiss,gvDiss,
329 jmc 1.57 & myThid)
330 adcroft 1.3 ENDIF
331 jmc 1.54 C-- if (momViscosity) end of block.
332 adcroft 1.1 ENDIF
333    
334 jmc 1.7 C- Return to standard hfacZ (min-4) and mask vort3 accordingly:
335     c CALL MOM_VI_MASK_VORT3(bi,bj,k,hFacZ,r_hFacZ,vort3,myThid)
336    
337 jmc 1.54 C--- Other dissipation terms in Zonal momentum equation
338 adcroft 1.1
339     C-- Vertical flux (fVer is at upper face of "u" cell)
340    
341     C Eddy component of vertical flux (interior component only) -> vrF
342 jmc 1.31 IF (momViscosity.AND..NOT.implicitViscosity) THEN
343 jmc 1.44 CALL MOM_U_RVISCFLUX(bi,bj,k+1,uVel,KappaRU,vrF,myThid)
344 adcroft 1.1
345     C Combine fluxes
346 jmc 1.31 DO j=jMin,jMax
347     DO i=iMin,iMax
348     fVerU(i,j,kDown) = ArDudrFac*vrF(i,j)
349     ENDDO
350 adcroft 1.1 ENDDO
351    
352 jmc 1.31 C-- Tendency is minus divergence of the fluxes
353     DO j=2-Oly,sNy+Oly-1
354     DO i=2-Olx,sNx+Olx-1
355     guDiss(i,j) = guDiss(i,j)
356 adcroft 1.1 & -_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)
357     & *recip_rAw(i,j,bi,bj)
358     & *(
359 jmc 1.42 & fVerU(i,j,kDown) - fVerU(i,j,kUp)
360     & )*rkSign
361 jmc 1.31 ENDDO
362 adcroft 1.1 ENDDO
363 jmc 1.31 ENDIF
364 adcroft 1.1
365 jmc 1.57 C-- No-slip and drag BCs appear as body forces in cell abutting topography
366 adcroft 1.1 IF (momViscosity.AND.no_slip_sides) THEN
367     C- No-slip BCs impose a drag at walls...
368 baylor 1.50 CALL MOM_U_SIDEDRAG(
369     I bi,bj,k,
370     I uFld, del2u, hFacZ,
371     I viscAh_Z,viscA4_Z,
372     I harmonic,biharmonic,useVariableViscosity,
373     O vF,
374     I myThid)
375 adcroft 1.1 DO j=jMin,jMax
376     DO i=iMin,iMax
377 jmc 1.31 guDiss(i,j) = guDiss(i,j)+vF(i,j)
378 adcroft 1.1 ENDDO
379     ENDDO
380     ENDIF
381     C- No-slip BCs impose a drag at bottom
382     IF (momViscosity.AND.bottomDragTerms) THEN
383     CALL MOM_U_BOTTOMDRAG(bi,bj,k,uFld,KE,KappaRU,vF,myThid)
384     DO j=jMin,jMax
385     DO i=iMin,iMax
386 jmc 1.31 guDiss(i,j) = guDiss(i,j)+vF(i,j)
387 adcroft 1.1 ENDDO
388     ENDDO
389     ENDIF
390 mlosch 1.56 #ifdef ALLOW_SHELFICE
391     IF (useShelfIce.AND.momViscosity.AND.bottomDragTerms) THEN
392     CALL SHELFICE_U_DRAG(bi,bj,k,uFld,KE,KappaRU,vF,myThid)
393     DO j=jMin,jMax
394     DO i=iMin,iMax
395     guDiss(i,j) = guDiss(i,j) + vF(i,j)
396     ENDDO
397     ENDDO
398     ENDIF
399     #endif /* ALLOW_SHELFICE */
400    
401 adcroft 1.1
402 jmc 1.54 C--- Other dissipation terms in Meridional momentum equation
403 adcroft 1.1
404     C-- Vertical flux (fVer is at upper face of "v" cell)
405    
406     C Eddy component of vertical flux (interior component only) -> vrF
407 jmc 1.31 IF (momViscosity.AND..NOT.implicitViscosity) THEN
408 jmc 1.44 CALL MOM_V_RVISCFLUX(bi,bj,k+1,vVel,KappaRV,vrF,myThid)
409 adcroft 1.1
410     C Combine fluxes -> fVerV
411 jmc 1.31 DO j=jMin,jMax
412     DO i=iMin,iMax
413     fVerV(i,j,kDown) = ArDvdrFac*vrF(i,j)
414     ENDDO
415 adcroft 1.1 ENDDO
416    
417 jmc 1.31 C-- Tendency is minus divergence of the fluxes
418     DO j=jMin,jMax
419     DO i=iMin,iMax
420     gvDiss(i,j) = gvDiss(i,j)
421 adcroft 1.1 & -_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
422     & *recip_rAs(i,j,bi,bj)
423     & *(
424 jmc 1.42 & fVerV(i,j,kDown) - fVerV(i,j,kUp)
425     & )*rkSign
426 jmc 1.31 ENDDO
427 adcroft 1.1 ENDDO
428 jmc 1.31 ENDIF
429 adcroft 1.1
430 jmc 1.57 C-- No-slip and drag BCs appear as body forces in cell abutting topography
431 adcroft 1.1 IF (momViscosity.AND.no_slip_sides) THEN
432     C- No-slip BCs impose a drag at walls...
433 baylor 1.50 CALL MOM_V_SIDEDRAG(
434     I bi,bj,k,
435     I vFld, del2v, hFacZ,
436     I viscAh_Z,viscA4_Z,
437     I harmonic,biharmonic,useVariableViscosity,
438     O vF,
439     I myThid)
440 adcroft 1.1 DO j=jMin,jMax
441     DO i=iMin,iMax
442 jmc 1.31 gvDiss(i,j) = gvDiss(i,j)+vF(i,j)
443 adcroft 1.1 ENDDO
444     ENDDO
445     ENDIF
446     C- No-slip BCs impose a drag at bottom
447     IF (momViscosity.AND.bottomDragTerms) THEN
448     CALL MOM_V_BOTTOMDRAG(bi,bj,k,vFld,KE,KappaRV,vF,myThid)
449     DO j=jMin,jMax
450     DO i=iMin,iMax
451 jmc 1.31 gvDiss(i,j) = gvDiss(i,j)+vF(i,j)
452 adcroft 1.1 ENDDO
453     ENDDO
454     ENDIF
455 mlosch 1.56 #ifdef ALLOW_SHELFICE
456     IF (useShelfIce.AND.momViscosity.AND.bottomDragTerms) THEN
457     CALL SHELFICE_V_DRAG(bi,bj,k,vFld,KE,KappaRU,vF,myThid)
458     DO j=jMin,jMax
459     DO i=iMin,iMax
460     gvDiss(i,j) = gvDiss(i,j) + vF(i,j)
461     ENDDO
462     ENDDO
463     ENDIF
464     #endif /* ALLOW_SHELFICE */
465    
466 adcroft 1.1
467 jmc 1.54 C- Vorticity diagnostics:
468     IF ( writeDiag ) THEN
469     IF (snapshot_mdsio) THEN
470     CALL WRITE_LOCAL_RL('Z3','I10',1,vort3, bi,bj,k,myIter,myThid)
471     ENDIF
472     #ifdef ALLOW_MNC
473     IF (useMNC .AND. snapshot_mnc) THEN
474     CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'Z3',vort3,
475     & offsets, myThid)
476     ENDIF
477     #endif /* ALLOW_MNC */
478     ENDIF
479     #ifdef ALLOW_DIAGNOSTICS
480     IF ( useDiagnostics ) THEN
481     CALL DIAGNOSTICS_FILL(vort3, 'momVort3',k,1,2,bi,bj,myThid)
482     ENDIF
483     #endif /* ALLOW_DIAGNOSTICS */
484    
485     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
486    
487     C--- Prepare for Advection & Coriolis terms:
488     C- Mask relative vorticity and calculate absolute vorticity
489     DO j=1-Oly,sNy+Oly
490     DO i=1-Olx,sNx+Olx
491     IF ( hFacZ(i,j).EQ.0. ) vort3(i,j) = 0.
492     ENDDO
493     ENDDO
494     IF (useAbsVorticity)
495     & CALL MOM_CALC_ABSVORT3(bi,bj,k,vort3,omega3,myThid)
496 adcroft 1.1
497 jmc 1.5 C-- Horizontal Coriolis terms
498 jmc 1.37 c IF (useCoriolis .AND. .NOT.useCDscheme
499     c & .AND. .NOT. useAbsVorticity) THEN
500     C- jmc: change it to keep the Coriolis terms when useAbsVorticity=T & momAdvection=F
501 jmc 1.46 IF ( useCoriolis .AND.
502 jmc 1.37 & .NOT.( useCDscheme .OR. useAbsVorticity.AND.momAdvection )
503     & ) THEN
504     IF (useAbsVorticity) THEN
505     CALL MOM_VI_U_CORIOLIS(bi,bj,K,vFld,omega3,hFacZ,r_hFacZ,
506     & uCf,myThid)
507     CALL MOM_VI_V_CORIOLIS(bi,bj,K,uFld,omega3,hFacZ,r_hFacZ,
508     & vCf,myThid)
509     ELSE
510     CALL MOM_VI_CORIOLIS(bi,bj,k,uFld,vFld,hFacZ,r_hFacZ,
511     & uCf,vCf,myThid)
512     ENDIF
513 jmc 1.5 DO j=jMin,jMax
514     DO i=iMin,iMax
515 jmc 1.43 gU(i,j,k,bi,bj) = uCf(i,j)
516     gV(i,j,k,bi,bj) = vCf(i,j)
517 jmc 1.5 ENDDO
518 adcroft 1.1 ENDDO
519 jmc 1.15 IF ( writeDiag ) THEN
520 edhill 1.24 IF (snapshot_mdsio) THEN
521     CALL WRITE_LOCAL_RL('fV','I10',1,uCf,bi,bj,k,myIter,myThid)
522     CALL WRITE_LOCAL_RL('fU','I10',1,vCf,bi,bj,k,myIter,myThid)
523     ENDIF
524     #ifdef ALLOW_MNC
525     IF (useMNC .AND. snapshot_mnc) THEN
526 edhill 1.53 CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj, 'fV', uCf,
527 edhill 1.25 & offsets, myThid)
528 edhill 1.53 CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj, 'fU', vCf,
529 edhill 1.25 & offsets, myThid)
530 edhill 1.24 ENDIF
531     #endif /* ALLOW_MNC */
532 jmc 1.15 ENDIF
533 jmc 1.46 #ifdef ALLOW_DIAGNOSTICS
534     IF ( useDiagnostics ) THEN
535     CALL DIAGNOSTICS_FILL(uCf,'Um_Cori ',k,1,2,bi,bj,myThid)
536     CALL DIAGNOSTICS_FILL(vCf,'Vm_Cori ',k,1,2,bi,bj,myThid)
537     ENDIF
538     #endif /* ALLOW_DIAGNOSTICS */
539 jmc 1.31 ELSE
540     DO j=jMin,jMax
541     DO i=iMin,iMax
542 jmc 1.43 gU(i,j,k,bi,bj) = 0. _d 0
543     gV(i,j,k,bi,bj) = 0. _d 0
544 jmc 1.31 ENDDO
545     ENDDO
546 jmc 1.5 ENDIF
547 adcroft 1.1
548 jmc 1.5 IF (momAdvection) THEN
549 jmc 1.41 C-- Horizontal advection of relative (or absolute) vorticity
550 jmc 1.62 IF ( (highOrderVorticity.OR.upwindVorticity)
551     & .AND.useAbsVorticity ) THEN
552 jmc 1.41 CALL MOM_VI_U_CORIOLIS_C4(bi,bj,k,vFld,omega3,r_hFacZ,
553 adcroft 1.20 & uCf,myThid)
554 jmc 1.62 ELSEIF ( (highOrderVorticity.OR.upwindVorticity) ) THEN
555 jmc 1.41 CALL MOM_VI_U_CORIOLIS_C4(bi,bj,k,vFld,vort3, r_hFacZ,
556     & uCf,myThid)
557 jmc 1.62 ELSEIF ( useAbsVorticity ) THEN
558 jmc 1.41 CALL MOM_VI_U_CORIOLIS(bi,bj,K,vFld,omega3,hFacZ,r_hFacZ,
559 jmc 1.40 & uCf,myThid)
560 adcroft 1.20 ELSE
561 jmc 1.41 CALL MOM_VI_U_CORIOLIS(bi,bj,k,vFld,vort3, hFacZ,r_hFacZ,
562 adcroft 1.20 & uCf,myThid)
563     ENDIF
564 jmc 1.5 DO j=jMin,jMax
565     DO i=iMin,iMax
566     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
567     ENDDO
568 adcroft 1.1 ENDDO
569 jmc 1.62 IF ( (highOrderVorticity.OR.upwindVorticity)
570     & .AND.useAbsVorticity ) THEN
571 jmc 1.41 CALL MOM_VI_V_CORIOLIS_C4(bi,bj,K,uFld,omega3,r_hFacZ,
572 adcroft 1.20 & vCf,myThid)
573 jmc 1.62 ELSEIF ( (highOrderVorticity.OR.upwindVorticity) ) THEN
574 jmc 1.41 CALL MOM_VI_V_CORIOLIS_C4(bi,bj,K,uFld,vort3, r_hFacZ,
575     & vCf,myThid)
576 jmc 1.62 ELSEIF ( useAbsVorticity ) THEN
577 jmc 1.41 CALL MOM_VI_V_CORIOLIS(bi,bj,K,uFld,omega3,hFacZ,r_hFacZ,
578 jmc 1.40 & vCf,myThid)
579 adcroft 1.20 ELSE
580 jmc 1.41 CALL MOM_VI_V_CORIOLIS(bi,bj,k,uFld,vort3, hFacZ,r_hFacZ,
581 adcroft 1.20 & vCf,myThid)
582     ENDIF
583 jmc 1.5 DO j=jMin,jMax
584     DO i=iMin,iMax
585     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
586     ENDDO
587 adcroft 1.1 ENDDO
588    
589 jmc 1.15 IF ( writeDiag ) THEN
590 edhill 1.24 IF (snapshot_mdsio) THEN
591     CALL WRITE_LOCAL_RL('zV','I10',1,uCf,bi,bj,k,myIter,myThid)
592     CALL WRITE_LOCAL_RL('zU','I10',1,vCf,bi,bj,k,myIter,myThid)
593     ENDIF
594     #ifdef ALLOW_MNC
595     IF (useMNC .AND. snapshot_mnc) THEN
596 edhill 1.53 CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj, 'zV', uCf,
597 edhill 1.25 & offsets, myThid)
598 edhill 1.53 CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj, 'zU', vCf,
599 edhill 1.25 & offsets, myThid)
600 edhill 1.24 ENDIF
601     #endif /* ALLOW_MNC */
602 jmc 1.15 ENDIF
603 edhill 1.24
604 jmc 1.7 #ifdef ALLOW_TIMEAVE
605     IF (taveFreq.GT.0.) THEN
606     CALL TIMEAVE_CUMUL_1K1T(uZetatave,vCf,deltaTClock,
607     & Nr, k, bi, bj, myThid)
608     CALL TIMEAVE_CUMUL_1K1T(vZetatave,uCf,deltaTClock,
609     & Nr, k, bi, bj, myThid)
610     ENDIF
611 dimitri 1.13 #endif /* ALLOW_TIMEAVE */
612 jmc 1.46 #ifdef ALLOW_DIAGNOSTICS
613     IF ( useDiagnostics ) THEN
614     CALL DIAGNOSTICS_FILL(uCf,'Um_AdvZ3',k,1,2,bi,bj,myThid)
615     CALL DIAGNOSTICS_FILL(vCf,'Vm_AdvZ3',k,1,2,bi,bj,myThid)
616     ENDIF
617     #endif /* ALLOW_DIAGNOSTICS */
618 jmc 1.7
619 jmc 1.5 C-- Vertical shear terms (-w*du/dr & -w*dv/dr)
620 jmc 1.12 IF ( .NOT. momImplVertAdv ) THEN
621     CALL MOM_VI_U_VERTSHEAR(bi,bj,K,uVel,wVel,uCf,myThid)
622     DO j=jMin,jMax
623     DO i=iMin,iMax
624     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
625     ENDDO
626 jmc 1.5 ENDDO
627 jmc 1.12 CALL MOM_VI_V_VERTSHEAR(bi,bj,K,vVel,wVel,vCf,myThid)
628     DO j=jMin,jMax
629     DO i=iMin,iMax
630     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
631     ENDDO
632 jmc 1.5 ENDDO
633 jmc 1.46 #ifdef ALLOW_DIAGNOSTICS
634     IF ( useDiagnostics ) THEN
635     CALL DIAGNOSTICS_FILL(uCf,'Um_AdvRe',k,1,2,bi,bj,myThid)
636     CALL DIAGNOSTICS_FILL(vCf,'Vm_AdvRe',k,1,2,bi,bj,myThid)
637     ENDIF
638     #endif /* ALLOW_DIAGNOSTICS */
639 jmc 1.12 ENDIF
640 adcroft 1.1
641     C-- Bernoulli term
642 jmc 1.5 CALL MOM_VI_U_GRAD_KE(bi,bj,K,KE,uCf,myThid)
643     DO j=jMin,jMax
644     DO i=iMin,iMax
645     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
646     ENDDO
647     ENDDO
648     CALL MOM_VI_V_GRAD_KE(bi,bj,K,KE,vCf,myThid)
649     DO j=jMin,jMax
650     DO i=iMin,iMax
651     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
652     ENDDO
653 adcroft 1.1 ENDDO
654 jmc 1.15 IF ( writeDiag ) THEN
655 edhill 1.24 IF (snapshot_mdsio) THEN
656     CALL WRITE_LOCAL_RL('KEx','I10',1,uCf,bi,bj,k,myIter,myThid)
657     CALL WRITE_LOCAL_RL('KEy','I10',1,vCf,bi,bj,k,myIter,myThid)
658     ENDIF
659     #ifdef ALLOW_MNC
660     IF (useMNC .AND. snapshot_mnc) THEN
661 edhill 1.53 CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj, 'KEx', uCf,
662 edhill 1.25 & offsets, myThid)
663 edhill 1.53 CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj, 'KEy', vCf,
664 edhill 1.25 & offsets, myThid)
665 jmc 1.54 ENDIF
666 edhill 1.24 #endif /* ALLOW_MNC */
667 jmc 1.15 ENDIF
668    
669 jmc 1.5 C-- end if momAdvection
670     ENDIF
671    
672 jmc 1.63 C-- 3.D Coriolis term (horizontal momentum, Eastward component: -fprime*w)
673 jmc 1.58 IF ( use3dCoriolis ) THEN
674 jmc 1.57 CALL MOM_U_CORIOLIS_NH(bi,bj,k,wVel,uCf,myThid)
675     DO j=jMin,jMax
676     DO i=iMin,iMax
677     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
678     ENDDO
679     ENDDO
680     IF ( usingCurvilinearGrid ) THEN
681     C- presently, non zero angleSinC array only supported with Curvilinear-Grid
682     CALL MOM_V_CORIOLIS_NH(bi,bj,k,wVel,vCf,myThid)
683     DO j=jMin,jMax
684     DO i=iMin,iMax
685     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
686     ENDDO
687     ENDDO
688     ENDIF
689     ENDIF
690    
691     C-- Non-Hydrostatic (spherical) metric terms
692     IF ( useNHMTerms ) THEN
693     CALL MOM_U_METRIC_NH(bi,bj,k,uFld,wVel,uCf,myThid)
694     DO j=jMin,jMax
695     DO i=iMin,iMax
696     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
697     ENDDO
698     ENDDO
699     CALL MOM_V_METRIC_NH(bi,bj,k,vFld,wVel,vCf,myThid)
700     DO j=jMin,jMax
701     DO i=iMin,iMax
702     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
703     ENDDO
704     ENDDO
705     ENDIF
706 jmc 1.54
707 jmc 1.5 C-- Set du/dt & dv/dt on boundaries to zero
708 adcroft 1.1 DO j=jMin,jMax
709     DO i=iMin,iMax
710 jmc 1.5 gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)*_maskW(i,j,k,bi,bj)
711     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)*_maskS(i,j,k,bi,bj)
712 adcroft 1.1 ENDDO
713     ENDDO
714 jmc 1.5
715 jmc 1.22 #ifdef ALLOW_DEBUG
716     IF ( debugLevel .GE. debLevB
717     & .AND. k.EQ.4 .AND. myIter.EQ.nIter0
718     & .AND. nPx.EQ.1 .AND. nPy.EQ.1
719     & .AND. useCubedSphereExchange ) THEN
720 jmc 1.23 CALL DEBUG_CS_CORNER_UV( ' uDiss,vDiss from MOM_VECINV',
721 jmc 1.31 & guDiss,gvDiss, k, standardMessageUnit,bi,bj,myThid )
722 jmc 1.22 ENDIF
723     #endif /* ALLOW_DEBUG */
724 adcroft 1.2
725 jmc 1.15 IF ( writeDiag ) THEN
726 edhill 1.24 IF (snapshot_mdsio) THEN
727 jmc 1.54 CALL WRITE_LOCAL_RL('W3','I10',1,omega3, bi,bj,k,myIter,myThid)
728     CALL WRITE_LOCAL_RL('KE','I10',1,KE, bi,bj,k,myIter,myThid)
729     CALL WRITE_LOCAL_RL('D', 'I10',1,hDiv, bi,bj,k,myIter,myThid)
730     CALL WRITE_LOCAL_RL('Dt','I10',1,tension,bi,bj,k,myIter,myThid)
731     CALL WRITE_LOCAL_RL('Du','I10',1,guDiss, bi,bj,k,myIter,myThid)
732     CALL WRITE_LOCAL_RL('Dv','I10',1,gvDiss, bi,bj,k,myIter,myThid)
733 edhill 1.24 ENDIF
734     #ifdef ALLOW_MNC
735     IF (useMNC .AND. snapshot_mnc) THEN
736 jmc 1.54 CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'W3',omega3,
737     & offsets, myThid)
738     CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'KE',KE,
739     & offsets, myThid)
740     CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'D', hDiv,
741 edhill 1.25 & offsets, myThid)
742 edhill 1.53 CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'Dt',tension,
743 edhill 1.25 & offsets, myThid)
744 edhill 1.53 CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'Du',guDiss,
745 edhill 1.25 & offsets, myThid)
746 edhill 1.53 CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'Dv',gvDiss,
747 edhill 1.25 & offsets, myThid)
748 edhill 1.24 ENDIF
749     #endif /* ALLOW_MNC */
750 adcroft 1.1 ENDIF
751 jmc 1.41
752 jmc 1.46 #ifdef ALLOW_DIAGNOSTICS
753     IF ( useDiagnostics ) THEN
754 jmc 1.54 CALL DIAGNOSTICS_FILL(KE, 'momKE ',k,1,2,bi,bj,myThid)
755 jmc 1.46 IF (momViscosity) THEN
756 jmc 1.54 CALL DIAGNOSTICS_FILL(hDiv, 'momHDiv ',k,1,2,bi,bj,myThid)
757 jmc 1.52 CALL DIAGNOSTICS_FILL(tension,'Tension ',k,1,2,bi,bj,myThid)
758 jmc 1.54 CALL DIAGNOSTICS_FILL(guDiss, 'Um_Diss ',k,1,2,bi,bj,myThid)
759     CALL DIAGNOSTICS_FILL(gvDiss, 'Vm_Diss ',k,1,2,bi,bj,myThid)
760 jmc 1.46 ENDIF
761 jmc 1.54 CALL DIAGNOSTICS_FILL(gU(1-Olx,1-Oly,k,bi,bj),
762     & 'Um_Advec',k,1,2,bi,bj,myThid)
763     CALL DIAGNOSTICS_FILL(gV(1-Olx,1-Oly,k,bi,bj),
764     & 'Vm_Advec',k,1,2,bi,bj,myThid)
765 jmc 1.46 ENDIF
766     #endif /* ALLOW_DIAGNOSTICS */
767    
768 edhill 1.11 #endif /* ALLOW_MOM_VECINV */
769 adcroft 1.1
770     RETURN
771     END

  ViewVC Help
Powered by ViewVC 1.1.22