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

Annotation of /MITgcm/pkg/mom_fluxform/mom_fluxform.F

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

Revision 1.24 - (hide annotations) (download)
Sun Sep 4 19:29:03 2005 UTC (18 years, 9 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint57r_post
Changes since 1.23: +62 -3 lines 
add diagnostics for (almost) each term in momentun Eq.
1 jmc 1.24 C $Header: /u/gcmpack/MITgcm/pkg/mom_fluxform/mom_fluxform.F,v 1.23 2005/07/30 22:07:00 jmc Exp $
2 adcroft 1.2 C $Name: $
3 adcroft 1.1
4 adcroft 1.3 CBOI
5     C !TITLE: pkg/mom\_advdiff
6     C !AUTHORS: adcroft@mit.edu
7 adcroft 1.4 C !INTRODUCTION: Flux-form Momentum Equations Package
8 adcroft 1.3 C
9     C Package "mom\_fluxform" provides methods for calculating explicit terms
10     C in the momentum equation cast in flux-form:
11     C \begin{eqnarray*}
12     C G^u & = & -\frac{1}{\rho} \partial_x \phi_h
13     C -\nabla \cdot {\bf v} u
14     C -fv
15     C +\frac{1}{\rho} \nabla \cdot {\bf \tau}^x
16     C + \mbox{metrics}
17     C \\
18     C G^v & = & -\frac{1}{\rho} \partial_y \phi_h
19     C -\nabla \cdot {\bf v} v
20     C +fu
21     C +\frac{1}{\rho} \nabla \cdot {\bf \tau}^y
22     C + \mbox{metrics}
23     C \end{eqnarray*}
24     C where ${\bf v}=(u,v,w)$ and $\tau$, the stress tensor, includes surface
25     C stresses as well as internal viscous stresses.
26     CEOI
27    
28 edhill 1.13 #include "MOM_FLUXFORM_OPTIONS.h"
29 adcroft 1.1
30 adcroft 1.3 CBOP
31     C !ROUTINE: MOM_FLUXFORM
32    
33     C !INTERFACE: ==========================================================
34 adcroft 1.1 SUBROUTINE MOM_FLUXFORM(
35     I bi,bj,iMin,iMax,jMin,jMax,k,kUp,kDown,
36 jmc 1.23 I KappaRU, KappaRV,
37 adcroft 1.1 U fVerU, fVerV,
38 jmc 1.23 O guDiss, gvDiss,
39     I myTime, myIter, myThid)
40 adcroft 1.3
41     C !DESCRIPTION:
42     C Calculates all the horizontal accelerations except for the implicit surface
43     C pressure gradient and implciit vertical viscosity.
44 adcroft 1.1
45 adcroft 1.3 C !USES: ===============================================================
46 adcroft 1.1 C == Global variables ==
47 adcroft 1.3 IMPLICIT NONE
48 adcroft 1.1 #include "SIZE.h"
49     #include "DYNVARS.h"
50     #include "FFIELDS.h"
51     #include "EEPARAMS.h"
52     #include "PARAMS.h"
53     #include "GRID.h"
54     #include "SURFACE.h"
55    
56 adcroft 1.3 C !INPUT PARAMETERS: ===================================================
57     C bi,bj :: tile indices
58     C iMin,iMax,jMin,jMAx :: loop ranges
59     C k :: vertical level
60     C kUp :: =1 or 2 for consecutive k
61     C kDown :: =2 or 1 for consecutive k
62     C KappaRU :: vertical viscosity
63     C KappaRV :: vertical viscosity
64     C fVerU :: vertical flux of U, 2 1/2 dim for pipe-lining
65     C fVerV :: vertical flux of V, 2 1/2 dim for pipe-lining
66 jmc 1.23 C guDiss :: dissipation tendency (all explicit terms), u component
67     C gvDiss :: dissipation tendency (all explicit terms), v component
68 jmc 1.8 C myTime :: current time
69 adcroft 1.3 C myIter :: current time-step number
70     C myThid :: thread number
71     INTEGER bi,bj,iMin,iMax,jMin,jMax
72     INTEGER k,kUp,kDown
73 adcroft 1.1 _RL KappaRU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
74     _RL KappaRV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
75     _RL fVerU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
76     _RL fVerV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
77 jmc 1.23 _RL guDiss(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
78     _RL gvDiss(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
79 jmc 1.8 _RL myTime
80 adcroft 1.2 INTEGER myIter
81 adcroft 1.1 INTEGER myThid
82    
83 adcroft 1.3 C !OUTPUT PARAMETERS: ==================================================
84     C None - updates gU() and gV() in common blocks
85    
86     C !LOCAL VARIABLES: ====================================================
87     C i,j :: loop indices
88     C vF :: viscous flux
89     C v4F :: bi-harmonic viscous flux
90     C cF :: Coriolis acceleration
91     C mT :: Metric terms
92     C fZon :: zonal fluxes
93     C fMer :: meridional fluxes
94 jmc 1.23 C fVrUp,fVrDw :: vertical viscous fluxes at interface k-1 & k
95 adcroft 1.3 INTEGER i,j
96     _RL vF(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
97     _RL v4F(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
98     _RL cF(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
99     _RL mT(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
100     _RL fZon(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
101     _RL fMer(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
102 jmc 1.23 _RL fVrUp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
103     _RL fVrDw(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
104 adcroft 1.1 C afFacMom - Tracer parameters for turning terms
105     C vfFacMom on and off.
106     C pfFacMom afFacMom - Advective terms
107     C cfFacMom vfFacMom - Eddy viscosity terms
108     C mTFacMom pfFacMom - Pressure terms
109     C cfFacMom - Coriolis terms
110     C foFacMom - Forcing
111     C mTFacMom - Metric term
112 jmc 1.23 C uDudxFac, AhDudxFac, etc ... individual term parameters for switching terms off
113 adcroft 1.1 _RS hFacZ(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
114     _RS r_hFacZ(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
115     _RS xA(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
116     _RS yA(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
117     _RL uTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
118     _RL vTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
119     _RL uFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
120     _RL vFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
121 jmc 1.8 _RL rTransU(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
122     _RL rTransV(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
123 adcroft 1.18 _RL KE(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
124 jmc 1.21 c _RL viscAh_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
125     c _RL viscAh_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
126     c _RL viscA4_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
127     c _RL viscA4_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
128     c _RL vort3(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
129     c _RL hDiv(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
130 adcroft 1.18 _RL strain(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
131     _RL tension(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
132 adcroft 1.1 _RL uDudxFac
133     _RL AhDudxFac
134     _RL vDudyFac
135     _RL AhDudyFac
136     _RL rVelDudrFac
137     _RL ArDudrFac
138     _RL fuFac
139     _RL mtFacU
140     _RL uDvdxFac
141     _RL AhDvdxFac
142     _RL vDvdyFac
143     _RL AhDvdyFac
144     _RL rVelDvdrFac
145     _RL ArDvdrFac
146     _RL fvFac
147     _RL mtFacV
148     LOGICAL bottomDragTerms
149 adcroft 1.3 CEOP
150 adcroft 1.1
151     C Initialise intermediate terms
152 jmc 1.23 DO j=1-OLy,sNy+OLy
153     DO i=1-OLx,sNx+OLx
154 adcroft 1.1 vF(i,j) = 0.
155     v4F(i,j) = 0.
156     cF(i,j) = 0.
157     mT(i,j) = 0.
158     fZon(i,j) = 0.
159     fMer(i,j) = 0.
160 jmc 1.23 fVrUp(i,j)= 0.
161     fVrDw(i,j)= 0.
162     rTransU(i,j)= 0.
163     rTransV(i,j)= 0.
164 adcroft 1.18 strain(i,j) = 0.
165 jmc 1.23 tension(i,j)= 0.
166     guDiss(i,j) = 0.
167     gvDiss(i,j) = 0.
168 adcroft 1.1 ENDDO
169     ENDDO
170    
171     C-- Term by term tracer parmeters
172     C o U momentum equation
173     uDudxFac = afFacMom*1.
174     AhDudxFac = vfFacMom*1.
175     vDudyFac = afFacMom*1.
176     AhDudyFac = vfFacMom*1.
177     rVelDudrFac = afFacMom*1.
178     ArDudrFac = vfFacMom*1.
179     mTFacU = mtFacMom*1.
180     fuFac = cfFacMom*1.
181     C o V momentum equation
182     uDvdxFac = afFacMom*1.
183     AhDvdxFac = vfFacMom*1.
184     vDvdyFac = afFacMom*1.
185     AhDvdyFac = vfFacMom*1.
186     rVelDvdrFac = afFacMom*1.
187     ArDvdrFac = vfFacMom*1.
188     mTFacV = mtFacMom*1.
189     fvFac = cfFacMom*1.
190 jmc 1.23
191     IF (implicitViscosity) THEN
192     ArDudrFac = 0.
193     ArDvdrFac = 0.
194     ENDIF
195 adcroft 1.1
196     IF ( no_slip_bottom
197     & .OR. bottomDragQuadratic.NE.0.
198     & .OR. bottomDragLinear.NE.0.) THEN
199     bottomDragTerms=.TRUE.
200     ELSE
201     bottomDragTerms=.FALSE.
202     ENDIF
203    
204     C-- Calculate open water fraction at vorticity points
205     CALL MOM_CALC_HFACZ(bi,bj,k,hFacZ,r_hFacZ,myThid)
206    
207     C---- Calculate common quantities used in both U and V equations
208     C Calculate tracer cell face open areas
209     DO j=1-OLy,sNy+OLy
210     DO i=1-OLx,sNx+OLx
211     xA(i,j) = _dyG(i,j,bi,bj)
212     & *drF(k)*_hFacW(i,j,k,bi,bj)
213     yA(i,j) = _dxG(i,j,bi,bj)
214     & *drF(k)*_hFacS(i,j,k,bi,bj)
215     ENDDO
216     ENDDO
217    
218     C Make local copies of horizontal flow field
219     DO j=1-OLy,sNy+OLy
220     DO i=1-OLx,sNx+OLx
221     uFld(i,j) = uVel(i,j,k,bi,bj)
222     vFld(i,j) = vVel(i,j,k,bi,bj)
223     ENDDO
224     ENDDO
225    
226     C Calculate velocity field "volume transports" through tracer cell faces.
227     DO j=1-OLy,sNy+OLy
228     DO i=1-OLx,sNx+OLx
229     uTrans(i,j) = uFld(i,j)*xA(i,j)
230     vTrans(i,j) = vFld(i,j)*yA(i,j)
231     ENDDO
232     ENDDO
233    
234 jmc 1.23 IF (bottomDragTerms) THEN
235     CALL MOM_CALC_KE(bi,bj,k,3,uFld,vFld,KE,myThid)
236     ENDIF
237 adcroft 1.1
238 jmc 1.21 IF (viscAstrain.NE.0. .OR. viscAtension.NE.0.) THEN
239 adcroft 1.18 CALL MOM_CALC_TENSION(bi,bj,k,uFld,vFld,
240     O tension,
241     I myThid)
242     CALL MOM_CALC_STRAIN(bi,bj,k,uFld,vFld,hFacZ,
243     O strain,
244     I myThid)
245 jmc 1.21 ENDIF
246 adcroft 1.18
247 jmc 1.8 C--- First call (k=1): compute vertical adv. flux fVerU(kUp) & fVerV(kUp)
248     IF (momAdvection.AND.k.EQ.1) THEN
249    
250     C- Calculate vertical transports above U & V points (West & South face):
251 jmc 1.23 CALL MOM_CALC_RTRANS( k, bi, bj,
252     O rTransU, rTransV,
253     I myTime, myIter, myThid)
254 jmc 1.8
255     C- Free surface correction term (flux at k=1)
256 jmc 1.23 CALL MOM_U_ADV_WU( bi,bj,k,uVel,wVel,rTransU,
257     O fVerU(1-OLx,1-OLy,kUp), myThid )
258 jmc 1.8
259 jmc 1.23 CALL MOM_V_ADV_WV( bi,bj,k,vVel,wVel,rTransV,
260     O fVerV(1-OLx,1-OLy,kUp), myThid )
261 jmc 1.8
262     C--- endif momAdvection & k=1
263     ENDIF
264    
265    
266     C--- Calculate vertical transports (at k+1) below U & V points :
267     IF (momAdvection) THEN
268 jmc 1.23 CALL MOM_CALC_RTRANS( k+1, bi, bj,
269     O rTransU, rTransV,
270     I myTime, myIter, myThid)
271 jmc 1.8 ENDIF
272    
273 adcroft 1.18 c IF (momViscosity) THEN
274     c & CALL MOM_CALC_VISCOSITY(bi,bj,k,
275     c I uFld,vFld,
276 jmc 1.21 c O viscAh_D,viscAh_Z,myThid)
277 jmc 1.8
278 jmc 1.23 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
279    
280 adcroft 1.1 C---- Zonal momentum equation starts here
281    
282 jmc 1.23 IF (momAdvection) THEN
283     C--- Calculate mean fluxes (advection) between cells for zonal flow.
284 adcroft 1.1
285     C-- Zonal flux (fZon is at east face of "u" cell)
286 jmc 1.23 C Mean flow component of zonal flux -> fZon
287     CALL MOM_U_ADV_UU(bi,bj,k,uTrans,uFld,fZon,myThid)
288 adcroft 1.1
289     C-- Meridional flux (fMer is at south face of "u" cell)
290 jmc 1.23 C Mean flow component of meridional flux -> fMer
291     CALL MOM_U_ADV_VU(bi,bj,k,vTrans,uFld,fMer,myThid)
292 adcroft 1.1
293     C-- Vertical flux (fVer is at upper face of "u" cell)
294 jmc 1.23 C Mean flow component of vertical flux (at k+1) -> fVer
295     CALL MOM_U_ADV_WU(
296     I bi,bj,k+1,uVel,wVel,rTransU,
297     O fVerU(1-OLx,1-OLy,kDown), myThid )
298 adcroft 1.1
299     C-- Tendency is minus divergence of the fluxes + coriolis + pressure term
300 jmc 1.23 DO j=jMin,jMax
301     DO i=iMin,iMax
302     gU(i,j,k,bi,bj) =
303 adcroft 1.1 #ifdef OLD_UV_GEOM
304 jmc 1.23 & -_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)/
305     & ( 0.5 _d 0*(rA(i,j,bi,bj)+rA(i-1,j,bi,bj)) )
306 adcroft 1.1 #else
307 jmc 1.23 & -_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)
308     & *recip_rAw(i,j,bi,bj)
309 adcroft 1.1 #endif
310 jmc 1.23 & *( ( fZon(i,j ) - fZon(i-1,j) )*uDudxFac
311     & +( fMer(i,j+1) - fMer(i, j) )*vDudyFac
312     & +(fVerU(i,j,kDown) - fVerU(i,j,kUp))*rkSign*rVelDudrFac
313     & )
314     ENDDO
315     ENDDO
316 adcroft 1.1
317 jmc 1.24 #ifdef ALLOW_DIAGNOSTICS
318     IF ( useDiagnostics ) THEN
319     CALL DIAGNOSTICS_FILL(fZon,'ADVx_Um ',k,1,2,bi,bj,myThid)
320     CALL DIAGNOSTICS_FILL(fMer,'ADVy_Um ',k,1,2,bi,bj,myThid)
321     CALL DIAGNOSTICS_FILL(fVerU(1-Olx,1-Oly,kUp),
322     & 'ADVrE_Um',k,1,2,bi,bj,myThid)
323     ENDIF
324     #endif
325    
326 jmc 1.8 #ifdef NONLIN_FRSURF
327     C-- account for 3.D divergence of the flow in rStar coordinate:
328 jmc 1.23 IF ( select_rStar.GT.0 ) THEN
329     DO j=jMin,jMax
330     DO i=iMin,iMax
331     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)
332 jmc 1.8 & - (rStarExpW(i,j,bi,bj) - 1. _d 0)/deltaTfreesurf
333     & *uVel(i,j,k,bi,bj)
334 jmc 1.23 ENDDO
335     ENDDO
336     ENDIF
337     IF ( select_rStar.LT.0 ) THEN
338     DO j=jMin,jMax
339     DO i=iMin,iMax
340     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)
341     & - rStarDhWDt(i,j,bi,bj)*uVel(i,j,k,bi,bj)
342     ENDDO
343     ENDDO
344     ENDIF
345     #endif /* NONLIN_FRSURF */
346    
347     ELSE
348     C- if momAdvection / else
349     DO j=1-OLy,sNy+OLy
350     DO i=1-OLx,sNx+OLx
351     gU(i,j,k,bi,bj) = 0. _d 0
352     ENDDO
353 jmc 1.8 ENDDO
354 jmc 1.23
355     C- endif momAdvection.
356 jmc 1.8 ENDIF
357 jmc 1.23
358     IF (momViscosity) THEN
359     C--- Calculate eddy fluxes (dissipation) between cells for zonal flow.
360    
361     C Bi-harmonic term del^2 U -> v4F
362     IF ( viscA4.NE.0. )
363     & CALL MOM_U_DEL2U(bi,bj,k,uFld,hFacZ,v4f,myThid)
364    
365     C Laplacian and bi-harmonic terms, Zonal Fluxes -> fZon
366     CALL MOM_U_XVISCFLUX(bi,bj,k,uFld,v4F,fZon,myThid)
367    
368     C Laplacian and bi-harmonic termis, Merid Fluxes -> fMer
369     CALL MOM_U_YVISCFLUX(bi,bj,k,uFld,v4F,hFacZ,fMer,myThid)
370    
371     C Eddy component of vertical flux (interior component only) -> fVrUp & fVrDw
372     IF (.NOT.implicitViscosity) THEN
373     CALL MOM_U_RVISCFLUX(bi,bj, k, uVel,KappaRU,fVrUp,myThid)
374     CALL MOM_U_RVISCFLUX(bi,bj,k+1,uVel,KappaRU,fVrDw,myThid)
375     ENDIF
376    
377     C-- Tendency is minus divergence of the fluxes
378     DO j=jMin,jMax
379     DO i=iMin,iMax
380     guDiss(i,j) =
381     #ifdef OLD_UV_GEOM
382     & -_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)/
383     & ( 0.5 _d 0*(rA(i,j,bi,bj)+rA(i-1,j,bi,bj)) )
384     #else
385     & -_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)
386     & *recip_rAw(i,j,bi,bj)
387     #endif
388     & *( ( fZon(i,j ) - fZon(i-1,j) )*AhDudxFac
389     & +( fMer(i,j+1) - fMer(i, j) )*AhDudyFac
390     & +( fVrDw(i,j) - fVrUp(i,j) )*rkSign*ArDudrFac
391     & )
392     ENDDO
393 jmc 1.8 ENDDO
394    
395 jmc 1.24 #ifdef ALLOW_DIAGNOSTICS
396     IF ( useDiagnostics ) THEN
397     CALL DIAGNOSTICS_FILL(fZon, 'VISCx_Um',k,1,2,bi,bj,myThid)
398     CALL DIAGNOSTICS_FILL(fMer, 'VISCy_Um',k,1,2,bi,bj,myThid)
399     IF (.NOT.implicitViscosity)
400     & CALL DIAGNOSTICS_FILL(fVrUp,'VISrE_Um',k,1,2,bi,bj,myThid)
401     ENDIF
402     #endif
403    
404 adcroft 1.1 C-- No-slip and drag BCs appear as body forces in cell abutting topography
405 jmc 1.23 IF (no_slip_sides) THEN
406 adcroft 1.1 C- No-slip BCs impose a drag at walls...
407 jmc 1.23 CALL MOM_U_SIDEDRAG(bi,bj,k,uFld,v4F,hFacZ,vF,myThid)
408     DO j=jMin,jMax
409     DO i=iMin,iMax
410     gUdiss(i,j) = gUdiss(i,j) + vF(i,j)
411     ENDDO
412     ENDDO
413     ENDIF
414 adcroft 1.1 C- No-slip BCs impose a drag at bottom
415 jmc 1.23 IF (bottomDragTerms) THEN
416     CALL MOM_U_BOTTOMDRAG(bi,bj,k,uFld,KE,KappaRU,vF,myThid)
417     DO j=jMin,jMax
418     DO i=iMin,iMax
419     gUdiss(i,j) = gUdiss(i,j) + vF(i,j)
420     ENDDO
421     ENDDO
422     ENDIF
423    
424     C- endif momViscosity
425 adcroft 1.1 ENDIF
426    
427 jmc 1.12 C-- Forcing term (moved to timestep.F)
428     c IF (momForcing)
429     c & CALL EXTERNAL_FORCING_U(
430     c I iMin,iMax,jMin,jMax,bi,bj,k,
431     c I myTime,myThid)
432 adcroft 1.1
433     C-- Metric terms for curvilinear grid systems
434 adcroft 1.5 IF (useNHMTerms) THEN
435     C o Non-hydrosatic metric terms
436 adcroft 1.1 CALL MOM_U_METRIC_NH(bi,bj,k,uFld,wVel,mT,myThid)
437     DO j=jMin,jMax
438     DO i=iMin,iMax
439     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+mTFacU*mT(i,j)
440     ENDDO
441     ENDDO
442 adcroft 1.5 ENDIF
443     IF (usingSphericalPolarMTerms) THEN
444 adcroft 1.1 CALL MOM_U_METRIC_SPHERE(bi,bj,k,uFld,vFld,mT,myThid)
445     DO j=jMin,jMax
446     DO i=iMin,iMax
447     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+mTFacU*mT(i,j)
448     ENDDO
449     ENDDO
450 afe 1.20 ENDIF
451 afe 1.19 IF (usingCylindricalGrid) THEN
452     CALL MOM_U_METRIC_CYLINDER(bi,bj,k,uFld,vFld,mT,myThid)
453     DO j=jMin,jMax
454     DO i=iMin,iMax
455     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+mTFacU*mT(i,j)
456     ENDDO
457     ENDDO
458 adcroft 1.1 ENDIF
459    
460 jmc 1.23 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
461 adcroft 1.1
462     C---- Meridional momentum equation starts here
463    
464 jmc 1.23 IF (momAdvection) THEN
465     C--- Calculate mean fluxes (advection) between cells for meridional flow.
466     C Mean flow component of zonal flux -> fZon
467     CALL MOM_V_ADV_UV(bi,bj,k,uTrans,vFld,fZon,myThid)
468 adcroft 1.1
469     C-- Meridional flux (fMer is at north face of "v" cell)
470 jmc 1.23 C Mean flow component of meridional flux -> fMer
471     CALL MOM_V_ADV_VV(bi,bj,k,vTrans,vFld,fMer,myThid)
472 adcroft 1.1
473     C-- Vertical flux (fVer is at upper face of "v" cell)
474 jmc 1.23 C Mean flow component of vertical flux (at k+1) -> fVerV
475     CALL MOM_V_ADV_WV(
476     I bi,bj,k+1,vVel,wVel,rTransV,
477     O fVerV(1-OLx,1-OLy,kDown), myThid )
478 adcroft 1.1
479     C-- Tendency is minus divergence of the fluxes + coriolis + pressure term
480 jmc 1.23 DO j=jMin,jMax
481     DO i=iMin,iMax
482     gV(i,j,k,bi,bj) =
483 adcroft 1.1 #ifdef OLD_UV_GEOM
484 jmc 1.23 & -_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)/
485     & ( 0.5 _d 0*(_rA(i,j,bi,bj)+_rA(i,j-1,bi,bj)) )
486 adcroft 1.1 #else
487 jmc 1.23 & -_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
488     & *recip_rAs(i,j,bi,bj)
489 adcroft 1.1 #endif
490 jmc 1.23 & *( ( fZon(i+1,j) - fZon(i,j ) )*uDvdxFac
491     & +( fMer(i, j) - fMer(i,j-1) )*vDvdyFac
492     & +(fVerV(i,j,kDown) - fVerV(i,j,kUp))*rkSign*rVelDvdrFac
493     & )
494 jmc 1.24 ENDDO
495     ENDDO
496    
497     #ifdef ALLOW_DIAGNOSTICS
498     IF ( useDiagnostics ) THEN
499     CALL DIAGNOSTICS_FILL(fZon,'ADVx_Vm ',k,1,2,bi,bj,myThid)
500     CALL DIAGNOSTICS_FILL(fMer,'ADVy_Vm ',k,1,2,bi,bj,myThid)
501     CALL DIAGNOSTICS_FILL(fVerV(1-Olx,1-Oly,kUp),
502     & 'ADVrE_Vm',k,1,2,bi,bj,myThid)
503     ENDIF
504     #endif
505 adcroft 1.1
506 jmc 1.8 #ifdef NONLIN_FRSURF
507     C-- account for 3.D divergence of the flow in rStar coordinate:
508 jmc 1.23 IF ( select_rStar.GT.0 ) THEN
509     DO j=jMin,jMax
510     DO i=iMin,iMax
511     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)
512 jmc 1.8 & - (rStarExpS(i,j,bi,bj) - 1. _d 0)/deltaTfreesurf
513     & *vVel(i,j,k,bi,bj)
514 jmc 1.23 ENDDO
515     ENDDO
516     ENDIF
517     IF ( select_rStar.LT.0 ) THEN
518     DO j=jMin,jMax
519     DO i=iMin,iMax
520     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)
521     & - rStarDhSDt(i,j,bi,bj)*vVel(i,j,k,bi,bj)
522     ENDDO
523     ENDDO
524     ENDIF
525     #endif /* NONLIN_FRSURF */
526    
527     ELSE
528     C- if momAdvection / else
529     DO j=1-OLy,sNy+OLy
530     DO i=1-OLx,sNx+OLx
531     gV(i,j,k,bi,bj) = 0. _d 0
532     ENDDO
533 jmc 1.8 ENDDO
534 jmc 1.23
535     C- endif momAdvection.
536 jmc 1.8 ENDIF
537 jmc 1.23
538     IF (momViscosity) THEN
539     C--- Calculate eddy fluxes (dissipation) between cells for meridional flow.
540     C Bi-harmonic term del^2 V -> v4F
541     IF ( viscA4.NE.0. )
542     & CALL MOM_V_DEL2V(bi,bj,k,vFld,hFacZ,v4f,myThid)
543    
544     C Laplacian and bi-harmonic terms, Zonal Fluxes -> fZon
545     CALL MOM_V_XVISCFLUX(bi,bj,k,vFld,v4f,hFacZ,fZon,myThid)
546    
547     C Laplacian and bi-harmonic termis, Merid Fluxes -> fMer
548     CALL MOM_V_YVISCFLUX(bi,bj,k,vFld,v4f,fMer,myThid)
549    
550     C Eddy component of vertical flux (interior component only) -> fVrUp & fVrDw
551     IF (.NOT.implicitViscosity) THEN
552     CALL MOM_V_RVISCFLUX(bi,bj, k, vVel,KappaRV,fVrUp,myThid)
553     CALL MOM_V_RVISCFLUX(bi,bj,k+1,vVel,KappaRV,fVrDw,myThid)
554     ENDIF
555    
556     C-- Tendency is minus divergence of the fluxes + coriolis + pressure term
557     DO j=jMin,jMax
558     DO i=iMin,iMax
559     gvDiss(i,j) =
560     #ifdef OLD_UV_GEOM
561     & -_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)/
562     & ( 0.5 _d 0*(_rA(i,j,bi,bj)+_rA(i,j-1,bi,bj)) )
563     #else
564     & -_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
565     & *recip_rAs(i,j,bi,bj)
566     #endif
567     & *( ( fZon(i+1,j) - fZon(i,j ) )*AhDvdxFac
568     & +( fMer(i, j) - fMer(i,j-1) )*AhDvdyFac
569     & +( fVrDw(i,j) - fVrUp(i,j) )*rkSign*ArDvdrFac
570     & )
571     ENDDO
572 jmc 1.8 ENDDO
573    
574 jmc 1.24 #ifdef ALLOW_DIAGNOSTICS
575     IF ( useDiagnostics ) THEN
576     CALL DIAGNOSTICS_FILL(fZon, 'VISCx_Vm',k,1,2,bi,bj,myThid)
577     CALL DIAGNOSTICS_FILL(fMer, 'VISCy_Vm',k,1,2,bi,bj,myThid)
578     IF (.NOT.implicitViscosity)
579     & CALL DIAGNOSTICS_FILL(fVrUp,'VISrE_Vm',k,1,2,bi,bj,myThid)
580     ENDIF
581     #endif
582    
583 adcroft 1.1 C-- No-slip and drag BCs appear as body forces in cell abutting topography
584 jmc 1.23 IF (no_slip_sides) THEN
585 adcroft 1.1 C- No-slip BCs impose a drag at walls...
586 jmc 1.23 CALL MOM_V_SIDEDRAG(bi,bj,k,vFld,v4F,hFacZ,vF,myThid)
587     DO j=jMin,jMax
588     DO i=iMin,iMax
589     gvDiss(i,j) = gvDiss(i,j) + vF(i,j)
590     ENDDO
591     ENDDO
592     ENDIF
593 adcroft 1.1 C- No-slip BCs impose a drag at bottom
594 jmc 1.23 IF (bottomDragTerms) THEN
595     CALL MOM_V_BOTTOMDRAG(bi,bj,k,vFld,KE,KappaRV,vF,myThid)
596     DO j=jMin,jMax
597     DO i=iMin,iMax
598     gvDiss(i,j) = gvDiss(i,j) + vF(i,j)
599     ENDDO
600     ENDDO
601     ENDIF
602    
603     C- endif momViscosity
604 adcroft 1.1 ENDIF
605    
606 jmc 1.12 C-- Forcing term (moved to timestep.F)
607     c IF (momForcing)
608     c & CALL EXTERNAL_FORCING_V(
609     c I iMin,iMax,jMin,jMax,bi,bj,k,
610     c I myTime,myThid)
611 adcroft 1.1
612     C-- Metric terms for curvilinear grid systems
613 adcroft 1.5 IF (useNHMTerms) THEN
614 adcroft 1.1 C o Spherical polar grid metric terms
615     CALL MOM_V_METRIC_NH(bi,bj,k,vFld,wVel,mT,myThid)
616     DO j=jMin,jMax
617     DO i=iMin,iMax
618     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+mTFacV*mT(i,j)
619     ENDDO
620     ENDDO
621 adcroft 1.5 ENDIF
622     IF (usingSphericalPolarMTerms) THEN
623 adcroft 1.1 CALL MOM_V_METRIC_SPHERE(bi,bj,k,uFld,mT,myThid)
624     DO j=jMin,jMax
625     DO i=iMin,iMax
626     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+mTFacV*mT(i,j)
627     ENDDO
628     ENDDO
629     ENDIF
630 afe 1.19 IF (usingCylindricalGrid) THEN
631     CALL MOM_V_METRIC_CYLINDER(bi,bj,k,uFld,vFld,mT,myThid)
632     DO j=jMin,jMax
633     DO i=iMin,iMax
634     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+mTFacV*mT(i,j)
635     ENDDO
636     ENDDO
637     ENDIF
638 adcroft 1.1
639 jmc 1.23 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
640 adcroft 1.1
641     C-- Coriolis term
642     C Note. As coded here, coriolis will not work with "thin walls"
643 jmc 1.12 c IF (useCDscheme) THEN
644     c CALL MOM_CDSCHEME(bi,bj,k,dPhiHydX,dPhiHydY,myThid)
645     c ELSE
646     IF (.NOT.useCDscheme) THEN
647     CALL MOM_U_CORIOLIS(bi,bj,k,vFld,cf,myThid)
648     DO j=jMin,jMax
649     DO i=iMin,iMax
650     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+fuFac*cf(i,j)
651     ENDDO
652     ENDDO
653 jmc 1.24 #ifdef ALLOW_DIAGNOSTICS
654     IF ( useDiagnostics )
655     & CALL DIAGNOSTICS_FILL(cf,'Um_Cori ',k,1,2,bi,bj,myThid)
656     #endif
657 jmc 1.12 CALL MOM_V_CORIOLIS(bi,bj,k,uFld,cf,myThid)
658     DO j=jMin,jMax
659     DO i=iMin,iMax
660     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+fvFac*cf(i,j)
661     ENDDO
662     ENDDO
663 jmc 1.24 #ifdef ALLOW_DIAGNOSTICS
664     IF ( useDiagnostics )
665     & CALL DIAGNOSTICS_FILL(cf,'Vm_Cori ',k,1,2,bi,bj,myThid)
666     #endif
667 jmc 1.12 ENDIF
668    
669 adcroft 1.7 IF (nonHydrostatic.OR.quasiHydrostatic) THEN
670 adcroft 1.6 CALL MOM_U_CORIOLIS_NH(bi,bj,k,wVel,cf,myThid)
671     DO j=jMin,jMax
672     DO i=iMin,iMax
673     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+fuFac*cf(i,j)
674     ENDDO
675     ENDDO
676     ENDIF
677 adcroft 1.1
678 jmc 1.23 C-- Set du/dt & dv/dt on boundaries to zero
679     DO j=jMin,jMax
680     DO i=iMin,iMax
681     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)*_maskW(i,j,k,bi,bj)
682     guDiss(i,j) = guDiss(i,j) *_maskW(i,j,k,bi,bj)
683     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)*_maskS(i,j,k,bi,bj)
684     gvDiss(i,j) = gvDiss(i,j) *_maskS(i,j,k,bi,bj)
685     ENDDO
686     ENDDO
687    
688 jmc 1.24 #ifdef ALLOW_DIAGNOSTICS
689     IF ( useDiagnostics ) THEN
690     IF (bottomDragTerms)
691     & CALL DIAGNOSTICS_FILL(KE, 'momKE ',k,1,2,bi,bj,myThid)
692     CALL DIAGNOSTICS_FILL(gU(1-Olx,1-Oly,k,bi,bj),
693     & 'Um_Advec',k,1,2,bi,bj,myThid)
694     CALL DIAGNOSTICS_FILL(gV(1-Olx,1-Oly,k,bi,bj),
695     & 'Vm_Advec',k,1,2,bi,bj,myThid)
696     IF (momViscosity) THEN
697     CALL DIAGNOSTICS_FILL(guDiss,'Um_Diss ',k,1,2,bi,bj,myThid)
698     CALL DIAGNOSTICS_FILL(gvDiss,'Vm_Diss ',k,1,2,bi,bj,myThid)
699     ENDIF
700     ENDIF
701     #endif /* ALLOW_DIAGNOSTICS */
702    
703 adcroft 1.1 RETURN
704     END
  ViewVC Help
Powered by ViewVC 1.1.22