/[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.56 - (hide annotations) (download)
Tue Feb 7 11:46:17 2006 UTC (18 years, 4 months ago) by mlosch
Branch: MAIN
CVS Tags: checkpoint58b_post, checkpoint58a_post, checkpoint58c_post
Changes since 1.55: +23 -1 lines
o add hooks for friction at water-shelfice interface

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

  ViewVC Help
Powered by ViewVC 1.1.22