/[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.29 - (hide annotations) (download)
Wed Oct 12 01:10:10 2005 UTC (18 years, 8 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint57v_post, checkpoint57w_post
Changes since 1.28: +32 -9 lines
apply free-slip / no-slip BC on vorticity & strain.

1 jmc 1.29 C $Header: /u/gcmpack/MITgcm/pkg/mom_fluxform/mom_fluxform.F,v 1.28 2005/09/27 13:38:21 baylor 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 baylor 1.25 _RL viscAh_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
125     _RL viscAh_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
126     _RL viscA4_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
127     _RL viscA4_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
128     _RL vort3(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
129     _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 jmc 1.29 _RL sideMaskFac
149 baylor 1.25 LOGICAL bottomDragTerms,harmonic,biharmonic,useVariableViscosity
150 adcroft 1.3 CEOP
151 adcroft 1.1
152     C Initialise intermediate terms
153 jmc 1.23 DO j=1-OLy,sNy+OLy
154     DO i=1-OLx,sNx+OLx
155 adcroft 1.1 vF(i,j) = 0.
156     v4F(i,j) = 0.
157     cF(i,j) = 0.
158     mT(i,j) = 0.
159     fZon(i,j) = 0.
160     fMer(i,j) = 0.
161 jmc 1.23 fVrUp(i,j)= 0.
162     fVrDw(i,j)= 0.
163     rTransU(i,j)= 0.
164     rTransV(i,j)= 0.
165 adcroft 1.18 strain(i,j) = 0.
166 jmc 1.23 tension(i,j)= 0.
167     guDiss(i,j) = 0.
168     gvDiss(i,j) = 0.
169 adcroft 1.1 ENDDO
170     ENDDO
171    
172     C-- Term by term tracer parmeters
173     C o U momentum equation
174     uDudxFac = afFacMom*1.
175     AhDudxFac = vfFacMom*1.
176     vDudyFac = afFacMom*1.
177     AhDudyFac = vfFacMom*1.
178     rVelDudrFac = afFacMom*1.
179     ArDudrFac = vfFacMom*1.
180     mTFacU = mtFacMom*1.
181     fuFac = cfFacMom*1.
182     C o V momentum equation
183     uDvdxFac = afFacMom*1.
184     AhDvdxFac = vfFacMom*1.
185     vDvdyFac = afFacMom*1.
186     AhDvdyFac = vfFacMom*1.
187     rVelDvdrFac = afFacMom*1.
188     ArDvdrFac = vfFacMom*1.
189     mTFacV = mtFacMom*1.
190     fvFac = cfFacMom*1.
191 jmc 1.23
192     IF (implicitViscosity) THEN
193     ArDudrFac = 0.
194     ArDvdrFac = 0.
195     ENDIF
196 adcroft 1.1
197 jmc 1.29 C note: using standard stencil (no mask) results in under-estimating
198     C vorticity at a no-slip boundary by a factor of 2 = sideDragFactor
199     IF ( no_slip_sides ) THEN
200     sideMaskFac = sideDragFactor
201     ELSE
202     sideMaskFac = 0. _d 0
203     ENDIF
204    
205 adcroft 1.1 IF ( no_slip_bottom
206     & .OR. bottomDragQuadratic.NE.0.
207     & .OR. bottomDragLinear.NE.0.) THEN
208     bottomDragTerms=.TRUE.
209     ELSE
210     bottomDragTerms=.FALSE.
211     ENDIF
212    
213     C-- Calculate open water fraction at vorticity points
214     CALL MOM_CALC_HFACZ(bi,bj,k,hFacZ,r_hFacZ,myThid)
215    
216     C---- Calculate common quantities used in both U and V equations
217     C Calculate tracer cell face open areas
218     DO j=1-OLy,sNy+OLy
219     DO i=1-OLx,sNx+OLx
220     xA(i,j) = _dyG(i,j,bi,bj)
221     & *drF(k)*_hFacW(i,j,k,bi,bj)
222     yA(i,j) = _dxG(i,j,bi,bj)
223     & *drF(k)*_hFacS(i,j,k,bi,bj)
224     ENDDO
225     ENDDO
226    
227     C Make local copies of horizontal flow field
228     DO j=1-OLy,sNy+OLy
229     DO i=1-OLx,sNx+OLx
230     uFld(i,j) = uVel(i,j,k,bi,bj)
231     vFld(i,j) = vVel(i,j,k,bi,bj)
232     ENDDO
233     ENDDO
234    
235     C Calculate velocity field "volume transports" through tracer cell faces.
236     DO j=1-OLy,sNy+OLy
237     DO i=1-OLx,sNx+OLx
238     uTrans(i,j) = uFld(i,j)*xA(i,j)
239     vTrans(i,j) = vFld(i,j)*yA(i,j)
240     ENDDO
241     ENDDO
242    
243 baylor 1.25 CALL MOM_CALC_KE(bi,bj,k,2,uFld,vFld,KE,myThid)
244 jmc 1.29 IF ( momViscosity) THEN
245     CALL MOM_CALC_HDIV(bi,bj,k,2,uFld,vFld,hDiv,myThid)
246     CALL MOM_CALC_RELVORT3(bi,bj,k,uFld,vFld,hFacZ,vort3,myThid)
247     CALL MOM_CALC_TENSION(bi,bj,k,uFld,vFld,tension,myThid)
248     CALL MOM_CALC_STRAIN(bi,bj,k,uFld,vFld,hFacZ,strain,myThid)
249     DO j=1-Oly,sNy+Oly
250     DO i=1-Olx,sNx+Olx
251     IF ( hFacZ(i,j).EQ.0. ) THEN
252     vort3(i,j) = sideMaskFac*vort3(i,j)
253     strain(i,j) = sideMaskFac*strain(i,j)
254     ENDIF
255     ENDDO
256     ENDDO
257     #ifdef ALLOW_DIAGNOSTICS
258     IF ( useDiagnostics ) THEN
259     CALL DIAGNOSTICS_FILL(hDiv, 'momHDiv ',k,1,2,bi,bj,myThid)
260     CALL DIAGNOSTICS_FILL(vort3, 'momVort3',k,1,2,bi,bj,myThid)
261     CALL DIAGNOSTICS_FILL(tension,'Tension ',k,1,2,bi,bj,myThid)
262     CALL DIAGNOSTICS_FILL(strain, 'Strain ',k,1,2,bi,bj,myThid)
263     ENDIF
264     #endif
265     ENDIF
266 adcroft 1.18
267 jmc 1.8 C--- First call (k=1): compute vertical adv. flux fVerU(kUp) & fVerV(kUp)
268     IF (momAdvection.AND.k.EQ.1) THEN
269    
270     C- Calculate vertical transports above U & V points (West & South face):
271 jmc 1.23 CALL MOM_CALC_RTRANS( k, bi, bj,
272     O rTransU, rTransV,
273     I myTime, myIter, myThid)
274 jmc 1.8
275     C- Free surface correction term (flux at k=1)
276 jmc 1.23 CALL MOM_U_ADV_WU( bi,bj,k,uVel,wVel,rTransU,
277     O fVerU(1-OLx,1-OLy,kUp), myThid )
278 jmc 1.8
279 jmc 1.23 CALL MOM_V_ADV_WV( bi,bj,k,vVel,wVel,rTransV,
280     O fVerV(1-OLx,1-OLy,kUp), myThid )
281 jmc 1.8
282     C--- endif momAdvection & k=1
283     ENDIF
284    
285    
286     C--- Calculate vertical transports (at k+1) below U & V points :
287     IF (momAdvection) THEN
288 jmc 1.23 CALL MOM_CALC_RTRANS( k+1, bi, bj,
289     O rTransU, rTransV,
290     I myTime, myIter, myThid)
291 jmc 1.8 ENDIF
292    
293 baylor 1.25 IF (momViscosity) THEN
294     CALL MOM_CALC_VISC(
295     I bi,bj,k,
296     O viscAh_Z,viscAh_D,viscA4_Z,viscA4_D,
297     O harmonic,biharmonic,useVariableViscosity,
298 jmc 1.26 I hDiv,vort3,tension,strain,KE,hFacZ,
299 baylor 1.25 I myThid)
300     ENDIF
301 jmc 1.8
302 jmc 1.23 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
303    
304 adcroft 1.1 C---- Zonal momentum equation starts here
305    
306 jmc 1.23 IF (momAdvection) THEN
307     C--- Calculate mean fluxes (advection) between cells for zonal flow.
308 adcroft 1.1
309     C-- Zonal flux (fZon is at east face of "u" cell)
310 jmc 1.23 C Mean flow component of zonal flux -> fZon
311     CALL MOM_U_ADV_UU(bi,bj,k,uTrans,uFld,fZon,myThid)
312 adcroft 1.1
313     C-- Meridional flux (fMer is at south face of "u" cell)
314 jmc 1.23 C Mean flow component of meridional flux -> fMer
315     CALL MOM_U_ADV_VU(bi,bj,k,vTrans,uFld,fMer,myThid)
316 adcroft 1.1
317     C-- Vertical flux (fVer is at upper face of "u" cell)
318 jmc 1.23 C Mean flow component of vertical flux (at k+1) -> fVer
319     CALL MOM_U_ADV_WU(
320     I bi,bj,k+1,uVel,wVel,rTransU,
321     O fVerU(1-OLx,1-OLy,kDown), myThid )
322 adcroft 1.1
323     C-- Tendency is minus divergence of the fluxes + coriolis + pressure term
324 jmc 1.23 DO j=jMin,jMax
325     DO i=iMin,iMax
326     gU(i,j,k,bi,bj) =
327 adcroft 1.1 #ifdef OLD_UV_GEOM
328 jmc 1.23 & -_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)/
329     & ( 0.5 _d 0*(rA(i,j,bi,bj)+rA(i-1,j,bi,bj)) )
330 adcroft 1.1 #else
331 jmc 1.23 & -_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)
332     & *recip_rAw(i,j,bi,bj)
333 adcroft 1.1 #endif
334 jmc 1.23 & *( ( fZon(i,j ) - fZon(i-1,j) )*uDudxFac
335     & +( fMer(i,j+1) - fMer(i, j) )*vDudyFac
336     & +(fVerU(i,j,kDown) - fVerU(i,j,kUp))*rkSign*rVelDudrFac
337     & )
338     ENDDO
339     ENDDO
340 adcroft 1.1
341 jmc 1.24 #ifdef ALLOW_DIAGNOSTICS
342     IF ( useDiagnostics ) THEN
343     CALL DIAGNOSTICS_FILL(fZon,'ADVx_Um ',k,1,2,bi,bj,myThid)
344     CALL DIAGNOSTICS_FILL(fMer,'ADVy_Um ',k,1,2,bi,bj,myThid)
345     CALL DIAGNOSTICS_FILL(fVerU(1-Olx,1-Oly,kUp),
346     & 'ADVrE_Um',k,1,2,bi,bj,myThid)
347     ENDIF
348     #endif
349    
350 jmc 1.8 #ifdef NONLIN_FRSURF
351     C-- account for 3.D divergence of the flow in rStar coordinate:
352 jmc 1.23 IF ( select_rStar.GT.0 ) THEN
353     DO j=jMin,jMax
354     DO i=iMin,iMax
355     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)
356 jmc 1.8 & - (rStarExpW(i,j,bi,bj) - 1. _d 0)/deltaTfreesurf
357     & *uVel(i,j,k,bi,bj)
358 jmc 1.23 ENDDO
359     ENDDO
360     ENDIF
361     IF ( select_rStar.LT.0 ) THEN
362     DO j=jMin,jMax
363     DO i=iMin,iMax
364     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)
365     & - rStarDhWDt(i,j,bi,bj)*uVel(i,j,k,bi,bj)
366     ENDDO
367     ENDDO
368     ENDIF
369     #endif /* NONLIN_FRSURF */
370    
371     ELSE
372     C- if momAdvection / else
373     DO j=1-OLy,sNy+OLy
374     DO i=1-OLx,sNx+OLx
375     gU(i,j,k,bi,bj) = 0. _d 0
376     ENDDO
377 jmc 1.8 ENDDO
378 jmc 1.23
379     C- endif momAdvection.
380 jmc 1.8 ENDIF
381 jmc 1.23
382     IF (momViscosity) THEN
383     C--- Calculate eddy fluxes (dissipation) between cells for zonal flow.
384    
385     C Bi-harmonic term del^2 U -> v4F
386 baylor 1.25 IF (biharmonic)
387 jmc 1.23 & CALL MOM_U_DEL2U(bi,bj,k,uFld,hFacZ,v4f,myThid)
388    
389     C Laplacian and bi-harmonic terms, Zonal Fluxes -> fZon
390 baylor 1.25 CALL MOM_U_XVISCFLUX(bi,bj,k,uFld,v4F,fZon,
391 baylor 1.27 I viscAh_D,viscA4_D,myThid)
392 jmc 1.23
393     C Laplacian and bi-harmonic termis, Merid Fluxes -> fMer
394 baylor 1.25 CALL MOM_U_YVISCFLUX(bi,bj,k,uFld,v4F,hFacZ,fMer,
395 baylor 1.27 I viscAh_Z,viscA4_Z,myThid)
396 jmc 1.23
397     C Eddy component of vertical flux (interior component only) -> fVrUp & fVrDw
398     IF (.NOT.implicitViscosity) THEN
399     CALL MOM_U_RVISCFLUX(bi,bj, k, uVel,KappaRU,fVrUp,myThid)
400     CALL MOM_U_RVISCFLUX(bi,bj,k+1,uVel,KappaRU,fVrDw,myThid)
401     ENDIF
402    
403     C-- Tendency is minus divergence of the fluxes
404     DO j=jMin,jMax
405     DO i=iMin,iMax
406     guDiss(i,j) =
407     #ifdef OLD_UV_GEOM
408     & -_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)/
409     & ( 0.5 _d 0*(rA(i,j,bi,bj)+rA(i-1,j,bi,bj)) )
410     #else
411     & -_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)
412     & *recip_rAw(i,j,bi,bj)
413     #endif
414     & *( ( fZon(i,j ) - fZon(i-1,j) )*AhDudxFac
415     & +( fMer(i,j+1) - fMer(i, j) )*AhDudyFac
416     & +( fVrDw(i,j) - fVrUp(i,j) )*rkSign*ArDudrFac
417     & )
418     ENDDO
419 jmc 1.8 ENDDO
420    
421 jmc 1.24 #ifdef ALLOW_DIAGNOSTICS
422     IF ( useDiagnostics ) THEN
423     CALL DIAGNOSTICS_FILL(fZon, 'VISCx_Um',k,1,2,bi,bj,myThid)
424     CALL DIAGNOSTICS_FILL(fMer, 'VISCy_Um',k,1,2,bi,bj,myThid)
425     IF (.NOT.implicitViscosity)
426     & CALL DIAGNOSTICS_FILL(fVrUp,'VISrE_Um',k,1,2,bi,bj,myThid)
427     ENDIF
428     #endif
429    
430 adcroft 1.1 C-- No-slip and drag BCs appear as body forces in cell abutting topography
431 jmc 1.23 IF (no_slip_sides) THEN
432 adcroft 1.1 C- No-slip BCs impose a drag at walls...
433 baylor 1.27 CALL MOM_U_SIDEDRAG(
434     I bi,bj,k,
435     I uFld, v4f, hFacZ,
436     I viscAh_Z,viscA4_Z,
437     I harmonic,biharmonic,useVariableViscosity,
438     O vF,
439     I myThid)
440 jmc 1.23 DO j=jMin,jMax
441     DO i=iMin,iMax
442     gUdiss(i,j) = gUdiss(i,j) + vF(i,j)
443     ENDDO
444     ENDDO
445     ENDIF
446 adcroft 1.1 C- No-slip BCs impose a drag at bottom
447 jmc 1.23 IF (bottomDragTerms) THEN
448     CALL MOM_U_BOTTOMDRAG(bi,bj,k,uFld,KE,KappaRU,vF,myThid)
449     DO j=jMin,jMax
450     DO i=iMin,iMax
451     gUdiss(i,j) = gUdiss(i,j) + vF(i,j)
452     ENDDO
453     ENDDO
454     ENDIF
455    
456     C- endif momViscosity
457 adcroft 1.1 ENDIF
458    
459 jmc 1.12 C-- Forcing term (moved to timestep.F)
460     c IF (momForcing)
461     c & CALL EXTERNAL_FORCING_U(
462     c I iMin,iMax,jMin,jMax,bi,bj,k,
463     c I myTime,myThid)
464 adcroft 1.1
465     C-- Metric terms for curvilinear grid systems
466 adcroft 1.5 IF (useNHMTerms) THEN
467     C o Non-hydrosatic metric terms
468 adcroft 1.1 CALL MOM_U_METRIC_NH(bi,bj,k,uFld,wVel,mT,myThid)
469     DO j=jMin,jMax
470     DO i=iMin,iMax
471     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+mTFacU*mT(i,j)
472     ENDDO
473     ENDDO
474 adcroft 1.5 ENDIF
475     IF (usingSphericalPolarMTerms) THEN
476 adcroft 1.1 CALL MOM_U_METRIC_SPHERE(bi,bj,k,uFld,vFld,mT,myThid)
477     DO j=jMin,jMax
478     DO i=iMin,iMax
479     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+mTFacU*mT(i,j)
480     ENDDO
481     ENDDO
482 afe 1.20 ENDIF
483 afe 1.19 IF (usingCylindricalGrid) THEN
484     CALL MOM_U_METRIC_CYLINDER(bi,bj,k,uFld,vFld,mT,myThid)
485     DO j=jMin,jMax
486     DO i=iMin,iMax
487     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+mTFacU*mT(i,j)
488     ENDDO
489     ENDDO
490 adcroft 1.1 ENDIF
491    
492 jmc 1.23 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
493 adcroft 1.1
494     C---- Meridional momentum equation starts here
495    
496 jmc 1.23 IF (momAdvection) THEN
497     C--- Calculate mean fluxes (advection) between cells for meridional flow.
498     C Mean flow component of zonal flux -> fZon
499     CALL MOM_V_ADV_UV(bi,bj,k,uTrans,vFld,fZon,myThid)
500 adcroft 1.1
501     C-- Meridional flux (fMer is at north face of "v" cell)
502 jmc 1.23 C Mean flow component of meridional flux -> fMer
503     CALL MOM_V_ADV_VV(bi,bj,k,vTrans,vFld,fMer,myThid)
504 adcroft 1.1
505     C-- Vertical flux (fVer is at upper face of "v" cell)
506 jmc 1.23 C Mean flow component of vertical flux (at k+1) -> fVerV
507     CALL MOM_V_ADV_WV(
508     I bi,bj,k+1,vVel,wVel,rTransV,
509     O fVerV(1-OLx,1-OLy,kDown), myThid )
510 adcroft 1.1
511     C-- Tendency is minus divergence of the fluxes + coriolis + pressure term
512 jmc 1.23 DO j=jMin,jMax
513     DO i=iMin,iMax
514     gV(i,j,k,bi,bj) =
515 adcroft 1.1 #ifdef OLD_UV_GEOM
516 jmc 1.23 & -_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)/
517     & ( 0.5 _d 0*(_rA(i,j,bi,bj)+_rA(i,j-1,bi,bj)) )
518 adcroft 1.1 #else
519 jmc 1.23 & -_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
520     & *recip_rAs(i,j,bi,bj)
521 adcroft 1.1 #endif
522 jmc 1.23 & *( ( fZon(i+1,j) - fZon(i,j ) )*uDvdxFac
523     & +( fMer(i, j) - fMer(i,j-1) )*vDvdyFac
524     & +(fVerV(i,j,kDown) - fVerV(i,j,kUp))*rkSign*rVelDvdrFac
525     & )
526 jmc 1.24 ENDDO
527     ENDDO
528    
529     #ifdef ALLOW_DIAGNOSTICS
530     IF ( useDiagnostics ) THEN
531     CALL DIAGNOSTICS_FILL(fZon,'ADVx_Vm ',k,1,2,bi,bj,myThid)
532     CALL DIAGNOSTICS_FILL(fMer,'ADVy_Vm ',k,1,2,bi,bj,myThid)
533     CALL DIAGNOSTICS_FILL(fVerV(1-Olx,1-Oly,kUp),
534     & 'ADVrE_Vm',k,1,2,bi,bj,myThid)
535     ENDIF
536     #endif
537 adcroft 1.1
538 jmc 1.8 #ifdef NONLIN_FRSURF
539     C-- account for 3.D divergence of the flow in rStar coordinate:
540 jmc 1.23 IF ( select_rStar.GT.0 ) THEN
541     DO j=jMin,jMax
542     DO i=iMin,iMax
543     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)
544 jmc 1.8 & - (rStarExpS(i,j,bi,bj) - 1. _d 0)/deltaTfreesurf
545     & *vVel(i,j,k,bi,bj)
546 jmc 1.23 ENDDO
547     ENDDO
548     ENDIF
549     IF ( select_rStar.LT.0 ) THEN
550     DO j=jMin,jMax
551     DO i=iMin,iMax
552     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)
553     & - rStarDhSDt(i,j,bi,bj)*vVel(i,j,k,bi,bj)
554     ENDDO
555     ENDDO
556     ENDIF
557     #endif /* NONLIN_FRSURF */
558    
559     ELSE
560     C- if momAdvection / else
561     DO j=1-OLy,sNy+OLy
562     DO i=1-OLx,sNx+OLx
563     gV(i,j,k,bi,bj) = 0. _d 0
564     ENDDO
565 jmc 1.8 ENDDO
566 jmc 1.23
567     C- endif momAdvection.
568 jmc 1.8 ENDIF
569 jmc 1.23
570     IF (momViscosity) THEN
571     C--- Calculate eddy fluxes (dissipation) between cells for meridional flow.
572     C Bi-harmonic term del^2 V -> v4F
573 baylor 1.25 IF (biharmonic)
574 jmc 1.23 & CALL MOM_V_DEL2V(bi,bj,k,vFld,hFacZ,v4f,myThid)
575    
576     C Laplacian and bi-harmonic terms, Zonal Fluxes -> fZon
577 baylor 1.25 CALL MOM_V_XVISCFLUX(bi,bj,k,vFld,v4f,hFacZ,fZon,
578 baylor 1.27 I viscAh_Z,viscA4_Z,myThid)
579 jmc 1.23
580     C Laplacian and bi-harmonic termis, Merid Fluxes -> fMer
581 baylor 1.25 CALL MOM_V_YVISCFLUX(bi,bj,k,vFld,v4f,fMer,
582 baylor 1.27 I viscAh_D,viscA4_D,myThid)
583 jmc 1.23
584     C Eddy component of vertical flux (interior component only) -> fVrUp & fVrDw
585     IF (.NOT.implicitViscosity) THEN
586     CALL MOM_V_RVISCFLUX(bi,bj, k, vVel,KappaRV,fVrUp,myThid)
587     CALL MOM_V_RVISCFLUX(bi,bj,k+1,vVel,KappaRV,fVrDw,myThid)
588     ENDIF
589    
590     C-- Tendency is minus divergence of the fluxes + coriolis + pressure term
591     DO j=jMin,jMax
592     DO i=iMin,iMax
593     gvDiss(i,j) =
594     #ifdef OLD_UV_GEOM
595     & -_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)/
596     & ( 0.5 _d 0*(_rA(i,j,bi,bj)+_rA(i,j-1,bi,bj)) )
597     #else
598     & -_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
599     & *recip_rAs(i,j,bi,bj)
600     #endif
601     & *( ( fZon(i+1,j) - fZon(i,j ) )*AhDvdxFac
602     & +( fMer(i, j) - fMer(i,j-1) )*AhDvdyFac
603     & +( fVrDw(i,j) - fVrUp(i,j) )*rkSign*ArDvdrFac
604     & )
605     ENDDO
606 jmc 1.8 ENDDO
607    
608 jmc 1.24 #ifdef ALLOW_DIAGNOSTICS
609     IF ( useDiagnostics ) THEN
610     CALL DIAGNOSTICS_FILL(fZon, 'VISCx_Vm',k,1,2,bi,bj,myThid)
611     CALL DIAGNOSTICS_FILL(fMer, 'VISCy_Vm',k,1,2,bi,bj,myThid)
612     IF (.NOT.implicitViscosity)
613     & CALL DIAGNOSTICS_FILL(fVrUp,'VISrE_Vm',k,1,2,bi,bj,myThid)
614     ENDIF
615     #endif
616    
617 adcroft 1.1 C-- No-slip and drag BCs appear as body forces in cell abutting topography
618 jmc 1.23 IF (no_slip_sides) THEN
619 adcroft 1.1 C- No-slip BCs impose a drag at walls...
620 baylor 1.27 CALL MOM_V_SIDEDRAG(
621     I bi,bj,k,
622     I vFld, v4f, hFacZ,
623     I viscAh_Z,viscA4_Z,
624     I harmonic,biharmonic,useVariableViscosity,
625     O vF,
626     I myThid)
627 jmc 1.23 DO j=jMin,jMax
628     DO i=iMin,iMax
629     gvDiss(i,j) = gvDiss(i,j) + vF(i,j)
630     ENDDO
631     ENDDO
632     ENDIF
633 adcroft 1.1 C- No-slip BCs impose a drag at bottom
634 jmc 1.23 IF (bottomDragTerms) THEN
635     CALL MOM_V_BOTTOMDRAG(bi,bj,k,vFld,KE,KappaRV,vF,myThid)
636     DO j=jMin,jMax
637     DO i=iMin,iMax
638     gvDiss(i,j) = gvDiss(i,j) + vF(i,j)
639     ENDDO
640     ENDDO
641     ENDIF
642    
643     C- endif momViscosity
644 adcroft 1.1 ENDIF
645    
646 jmc 1.12 C-- Forcing term (moved to timestep.F)
647     c IF (momForcing)
648     c & CALL EXTERNAL_FORCING_V(
649     c I iMin,iMax,jMin,jMax,bi,bj,k,
650     c I myTime,myThid)
651 adcroft 1.1
652     C-- Metric terms for curvilinear grid systems
653 adcroft 1.5 IF (useNHMTerms) THEN
654 adcroft 1.1 C o Spherical polar grid metric terms
655     CALL MOM_V_METRIC_NH(bi,bj,k,vFld,wVel,mT,myThid)
656     DO j=jMin,jMax
657     DO i=iMin,iMax
658     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+mTFacV*mT(i,j)
659     ENDDO
660     ENDDO
661 adcroft 1.5 ENDIF
662     IF (usingSphericalPolarMTerms) THEN
663 adcroft 1.1 CALL MOM_V_METRIC_SPHERE(bi,bj,k,uFld,mT,myThid)
664     DO j=jMin,jMax
665     DO i=iMin,iMax
666     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+mTFacV*mT(i,j)
667     ENDDO
668     ENDDO
669     ENDIF
670 afe 1.19 IF (usingCylindricalGrid) THEN
671     CALL MOM_V_METRIC_CYLINDER(bi,bj,k,uFld,vFld,mT,myThid)
672     DO j=jMin,jMax
673     DO i=iMin,iMax
674     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+mTFacV*mT(i,j)
675     ENDDO
676     ENDDO
677     ENDIF
678 adcroft 1.1
679 jmc 1.23 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
680 adcroft 1.1
681     C-- Coriolis term
682     C Note. As coded here, coriolis will not work with "thin walls"
683 jmc 1.12 c IF (useCDscheme) THEN
684     c CALL MOM_CDSCHEME(bi,bj,k,dPhiHydX,dPhiHydY,myThid)
685     c ELSE
686     IF (.NOT.useCDscheme) THEN
687     CALL MOM_U_CORIOLIS(bi,bj,k,vFld,cf,myThid)
688     DO j=jMin,jMax
689     DO i=iMin,iMax
690     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+fuFac*cf(i,j)
691     ENDDO
692     ENDDO
693 jmc 1.24 #ifdef ALLOW_DIAGNOSTICS
694     IF ( useDiagnostics )
695     & CALL DIAGNOSTICS_FILL(cf,'Um_Cori ',k,1,2,bi,bj,myThid)
696     #endif
697 jmc 1.12 CALL MOM_V_CORIOLIS(bi,bj,k,uFld,cf,myThid)
698     DO j=jMin,jMax
699     DO i=iMin,iMax
700     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+fvFac*cf(i,j)
701     ENDDO
702     ENDDO
703 jmc 1.24 #ifdef ALLOW_DIAGNOSTICS
704     IF ( useDiagnostics )
705     & CALL DIAGNOSTICS_FILL(cf,'Vm_Cori ',k,1,2,bi,bj,myThid)
706     #endif
707 jmc 1.12 ENDIF
708    
709 adcroft 1.7 IF (nonHydrostatic.OR.quasiHydrostatic) THEN
710 adcroft 1.6 CALL MOM_U_CORIOLIS_NH(bi,bj,k,wVel,cf,myThid)
711     DO j=jMin,jMax
712     DO i=iMin,iMax
713     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+fuFac*cf(i,j)
714     ENDDO
715     ENDDO
716     ENDIF
717 adcroft 1.1
718 jmc 1.23 C-- Set du/dt & dv/dt on boundaries to zero
719     DO j=jMin,jMax
720     DO i=iMin,iMax
721     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)*_maskW(i,j,k,bi,bj)
722     guDiss(i,j) = guDiss(i,j) *_maskW(i,j,k,bi,bj)
723     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)*_maskS(i,j,k,bi,bj)
724     gvDiss(i,j) = gvDiss(i,j) *_maskS(i,j,k,bi,bj)
725     ENDDO
726     ENDDO
727    
728 jmc 1.24 #ifdef ALLOW_DIAGNOSTICS
729     IF ( useDiagnostics ) THEN
730 baylor 1.28 CALL DIAGNOSTICS_FILL(KE, 'momKE ',k,1,2,bi,bj,myThid)
731 jmc 1.24 CALL DIAGNOSTICS_FILL(gU(1-Olx,1-Oly,k,bi,bj),
732     & 'Um_Advec',k,1,2,bi,bj,myThid)
733     CALL DIAGNOSTICS_FILL(gV(1-Olx,1-Oly,k,bi,bj),
734     & 'Vm_Advec',k,1,2,bi,bj,myThid)
735     IF (momViscosity) THEN
736     CALL DIAGNOSTICS_FILL(guDiss,'Um_Diss ',k,1,2,bi,bj,myThid)
737     CALL DIAGNOSTICS_FILL(gvDiss,'Vm_Diss ',k,1,2,bi,bj,myThid)
738     ENDIF
739     ENDIF
740     #endif /* ALLOW_DIAGNOSTICS */
741    
742 adcroft 1.1 RETURN
743     END

  ViewVC Help
Powered by ViewVC 1.1.22