Back to home page

MITgcm

 
 

    


File indexing completed on 2024-08-14 05:10:50 UTC

view on githubraw file Latest commit a340904e on 2024-08-13 12:35:11 UTC
cec2469d72 Alis*0001 #include "MOM_VECINV_OPTIONS.h"
a340904e5a Ou W*0002 #include "MOM_COMMON_OPTIONS.h"
88b07ffa37 Jean*0003 #ifdef ALLOW_AUTODIFF
                0004 # include "AUTODIFF_OPTIONS.h"
                0005 #endif
9293d3c672 Hajo*0006 #ifdef ALLOW_GGL90
                0007 # include "GGL90_OPTIONS.h"
                0008 #endif
aea29c8517 Alis*0009 
039fe61d35 Jean*0010       SUBROUTINE MOM_VECINV(
07e17fa184 Jean*0011      I        bi,bj,k,iMin,iMax,jMin,jMax,
1833b564cb Jean*0012      I        kappaRU, kappaRV,
07e17fa184 Jean*0013      I        fVerUkm, fVerVkm,
                0014      O        fVerUkp, fVerVkp,
4667e4066d Jean*0015      O        guDiss, gvDiss,
07e17fa184 Jean*0016      I        myTime, myIter, myThid )
                0017 C     *==========================================================*
66d5e1666c Alis*0018 C     | S/R MOM_VECINV                                           |
aea29c8517 Alis*0019 C     | o Form the right hand-side of the momentum equation.     |
07e17fa184 Jean*0020 C     *==========================================================*
aea29c8517 Alis*0021 C     | Terms are evaluated one layer at a time working from     |
                0022 C     | the bottom to the top. The vertically integrated         |
                0023 C     | barotropic flow tendency term is evluated by summing the |
                0024 C     | tendencies.                                              |
                0025 C     | Notes:                                                   |
                0026 C     | We have not sorted out an entirely satisfactory formula  |
                0027 C     | for the diffusion equation bc with lopping. The present  |
                0028 C     | form produces a diffusive flux that does not scale with  |
                0029 C     | open-area. Need to do something to solidfy this and to   |
                0030 C     | deal "properly" with thin walls.                         |
07e17fa184 Jean*0031 C     *==========================================================*
aea29c8517 Alis*0032       IMPLICIT NONE
                0033 
                0034 C     == Global variables ==
                0035 #include "SIZE.h"
                0036 #include "EEPARAMS.h"
                0037 #include "PARAMS.h"
                0038 #include "GRID.h"
01f860d49e Jean*0039 #include "SURFACE.h"
f0501c53d1 Jean*0040 #include "DYNVARS.h"
79074ef66b Jean*0041 #include "FFIELDS.h"
a340904e5a Ou W*0042 #include "MOM_VISC.h"
9293d3c672 Hajo*0043 #if ( defined ALLOW_GGL90 && defined ALLOW_GGL90_LANGMUIR )
                0044 # include "GGL90.h"
                0045 #endif
cab1a69b8a Jean*0046 #ifdef ALLOW_TIMEAVE
f0501c53d1 Jean*0047 # include "TIMEAVE_STATV.h"
                0048 #endif
                0049 #ifdef ALLOW_MNC
                0050 # include "MNC_PARAMS.h"
cab1a69b8a Jean*0051 #endif
cd3c16aeda Patr*0052 #ifdef ALLOW_AUTODIFF_TAMC
                0053 # include "tamc.h"
                0054 #endif
aea29c8517 Alis*0055 
                0056 C     == Routine arguments ==
07e17fa184 Jean*0057 C     bi,bj   :: current tile indices
                0058 C     k       :: current vertical level
                0059 C     iMin,iMax,jMin,jMax :: loop ranges
                0060 C     fVerU   :: Flux of momentum in the vertical direction, out of the upper
f1a4cec01a Jean*0061 C     fVerV   :: face of a cell k ( flux into the cell above ).
07e17fa184 Jean*0062 C     fVerUkm :: vertical viscous flux of U, interface above (k-1/2)
                0063 C     fVerVkm :: vertical viscous flux of V, interface above (k-1/2)
                0064 C     fVerUkp :: vertical viscous flux of U, interface below (k+1/2)
                0065 C     fVerVkp :: vertical viscous flux of V, interface below (k+1/2)
                0066 
                0067 C     guDiss  :: dissipation tendency (all explicit terms), u component
                0068 C     gvDiss  :: dissipation tendency (all explicit terms), v component
                0069 C     myTime  :: current time
                0070 C     myIter  :: current time-step number
                0071 C     myThid  :: my Thread Id number
                0072       INTEGER bi,bj,k
                0073       INTEGER iMin,iMax,jMin,jMax
1833b564cb Jean*0074       _RL kappaRU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1)
                0075       _RL kappaRV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1)
07e17fa184 Jean*0076       _RL fVerUkm(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0077       _RL fVerVkm(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0078       _RL fVerUkp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0079       _RL fVerVkp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
4667e4066d Jean*0080       _RL guDiss(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0081       _RL gvDiss(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
f7d48db11c Jean*0082       _RL     myTime
3a279374db Alis*0083       INTEGER myIter
                0084       INTEGER myThid
aea29c8517 Alis*0085 
99bc607d7a Ed H*0086 #ifdef ALLOW_MOM_VECINV
3a279374db Alis*0087 C     == Functions ==
94a46dfe0d Jean*0088       LOGICAL  DIFFERENT_MULTIPLE
                0089       EXTERNAL DIFFERENT_MULTIPLE
a340904e5a Ou W*0090 #ifdef ALLOW_DIAGNOSTICS
                0091       LOGICAL  DIAGNOSTICS_IS_ON
                0092       EXTERNAL DIAGNOSTICS_IS_ON
                0093 #endif
3a279374db Alis*0094 
aea29c8517 Alis*0095 C     == Local variables ==
ed2dbbe83d Jean*0096 C     strainBC :: same as strain but account for no-slip BC
                0097 C     vort3BC  :: same as vort3  but account for no-slip BC
a340904e5a Ou W*0098 C     i,j      :: Loop counters
                0099       INTEGER i,j
aea29c8517 Alis*0100       _RL      vF (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
df1ec59c8a Jean*0101       _RL      vrF(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0102       _RL      uCf(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0103       _RL      vCf(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0104       _RS hFacZ   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
01f860d49e Jean*0105       _RS h0FacZ  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
df1ec59c8a Jean*0106       _RS r_hFacZ (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0107       _RL uFld    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0108       _RL vFld    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
9293d3c672 Hajo*0109 #if ( defined ALLOW_GGL90 && defined ALLOW_GGL90_LANGMUIR )
                0110       _RL uRes    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0111       _RL vRes    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0112 #endif
df1ec59c8a Jean*0113       _RL del2u   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0114       _RL del2v   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0115       _RL dStar   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0116       _RL zStar   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0117       _RL tension (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0118       _RL strain  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
ed2dbbe83d Jean*0119       _RL strainBC(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
f59d76b0dd Ed D*0120       _RL stretching(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
06244a5e4f Jean*0121 #ifdef ALLOW_LEITH_QG
f59d76b0dd Ed D*0122       _RL Nsquare (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
06244a5e4f Jean*0123 #endif
79074ef66b Jean*0124       _RL cDrag   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
df1ec59c8a Jean*0125       _RL KE      (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0126       _RL omega3  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0127       _RL vort3   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
ed2dbbe83d Jean*0128       _RL vort3BC (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
df1ec59c8a Jean*0129       _RL hDiv    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0130       _RL viscAh_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0131       _RL viscAh_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0132       _RL viscA4_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0133       _RL viscA4_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
07e17fa184 Jean*0134 C     xxxFac :: On-off tracer parameters used for switching terms off.
aea29c8517 Alis*0135       _RL  ArDudrFac
                0136       _RL  ArDvdrFac
df1ec59c8a Jean*0137       _RL  sideMaskFac
aea29c8517 Alis*0138       LOGICAL bottomDragTerms
f7d48db11c Jean*0139       LOGICAL writeDiag
a340904e5a Ou W*0140 #if ( defined ALLOW_MOM_TEND_EXTRA_DIAGS && defined ALLOW_DIAGNOSTICS )
                0141       LOGICAL botDragDiagIsOn
                0142       LOGICAL shelfDragDiagIsOn
                0143 #endif /* ALLOW_MOM_TEND_EXTRA_DIAGS && ALLOW_DIAGNOSTICS */
cd3c16aeda Patr*0144 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0145 C     kkey :: tape key (depends on level and tile indices)
                0146       INTEGER kkey
cd3c16aeda Patr*0147 #endif
cb356b4c5f Ed H*0148 #ifdef ALLOW_MNC
                0149       INTEGER offsets(9)
b22b541fe9 Ed H*0150       CHARACTER*(1) pf
cb356b4c5f Ed H*0151 #endif
                0152 
88b07ffa37 Jean*0153 #ifdef ALLOW_AUTODIFF
7d3b758ca2 Patr*0154 C--   only the kDown part of fverU/V is set in this subroutine
                0155 C--   the kUp is still required
f1a4cec01a Jean*0156 C--   In the case of mom_fluxform kUp is set as well
7d3b758ca2 Patr*0157 C--   (at least in part)
07e17fa184 Jean*0158       fVerUkm(1,1) = fVerUkm(1,1)
                0159       fVerVkm(1,1) = fVerVkm(1,1)
7d3b758ca2 Patr*0160 #endif
                0161 
cd3c16aeda Patr*0162 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0163       kkey = bi + (bj-1)*nSx + (ikey_dynamics-1)*nSx*nSy
                0164       kkey = k  + (kkey-1)*Nr
cd3c16aeda Patr*0165 #endif /* ALLOW_AUTODIFF_TAMC */
                0166 
94a46dfe0d Jean*0167       writeDiag = DIFFERENT_MULTIPLE(diagFreq, myTime, deltaTClock)
aea29c8517 Alis*0168 
ef1e764710 Ed H*0169 #ifdef ALLOW_MNC
                0170       IF (useMNC .AND. snapshot_mnc .AND. writeDiag) THEN
b22b541fe9 Ed H*0171         IF ( writeBinaryPrec .EQ. precFloat64 ) THEN
                0172           pf(1:1) = 'D'
                0173         ELSE
                0174           pf(1:1) = 'R'
                0175         ENDIF
cb356b4c5f Ed H*0176         IF ((bi .EQ. 1).AND.(bj .EQ. 1).AND.(k .EQ. 1)) THEN
                0177           CALL MNC_CW_SET_UDIM('mom_vi', -1, myThid)
987ff12cb6 Ed H*0178           CALL MNC_CW_RL_W_S('D','mom_vi',0,0,'T',myTime,myThid)
cb356b4c5f Ed H*0179           CALL MNC_CW_SET_UDIM('mom_vi', 0, myThid)
987ff12cb6 Ed H*0180           CALL MNC_CW_I_W_S('I','mom_vi',0,0,'iter',myIter,myThid)
cb356b4c5f Ed H*0181         ENDIF
                0182         DO i = 1,9
                0183           offsets(i) = 0
                0184         ENDDO
                0185         offsets(3) = k
9c98e82bbb Jean*0186 c       write(*,*) 'offsets = ',(offsets(i),i=1,9)
ef1e764710 Ed H*0187       ENDIF
                0188 #endif /*  ALLOW_MNC  */
                0189 
9c98e82bbb Jean*0190 C--   Initialise intermediate terms
                0191       DO j=1-OLy,sNy+OLy
                0192        DO i=1-OLx,sNx+OLx
4667e4066d Jean*0193         vF(i,j)    = 0.
                0194         vrF(i,j)   = 0.
aea29c8517 Alis*0195         uCf(i,j)   = 0.
                0196         vCf(i,j)   = 0.
                0197         del2u(i,j) = 0.
                0198         del2v(i,j) = 0.
                0199         dStar(i,j) = 0.
                0200         zStar(i,j) = 0.
4667e4066d Jean*0201         guDiss(i,j)= 0.
                0202         gvDiss(i,j)= 0.
9293d3c672 Hajo*0203 #if ( defined ALLOW_GGL90 && defined ALLOW_GGL90_LANGMUIR )
                0204 c       uRes(i,j)  = 0.
                0205 c       vRes(i,j)  = 0.
                0206 #endif
aea29c8517 Alis*0207         vort3(i,j) = 0.
4667e4066d Jean*0208         omega3(i,j)= 0.
df1ec59c8a Jean*0209         KE(i,j)    = 0.
9c98e82bbb Jean*0210 C-    need to initialise hDiv for MOM_VI_DEL2UV(call FILL_CS_CORNER_TR_RL)
                0211         hDiv(i,j)  = 0.
34e137f064 Jean*0212 c       viscAh_Z(i,j) = 0.
                0213 c       viscAh_D(i,j) = 0.
                0214 c       viscA4_Z(i,j) = 0.
                0215 c       viscA4_D(i,j) = 0.
629d440dd4 Patr*0216         strain(i,j)  = 0. _d 0
ed2dbbe83d Jean*0217         strainBC(i,j)= 0. _d 0
f59d76b0dd Ed D*0218         stretching(i,j) = 0. _d 0
06244a5e4f Jean*0219 #ifdef ALLOW_LEITH_QG
f59d76b0dd Ed D*0220         Nsquare(i,j) = 0. _d 0
06244a5e4f Jean*0221 #endif
629d440dd4 Patr*0222         tension(i,j) = 0. _d 0
88b07ffa37 Jean*0223 #ifdef ALLOW_AUTODIFF
cdc9f269ae Patr*0224         hFacZ(i,j)   = 0. _d 0
629d440dd4 Patr*0225 #endif
aea29c8517 Alis*0226        ENDDO
                0227       ENDDO
                0228 
                0229 C--   Term by term tracer parmeters
                0230 C     o U momentum equation
                0231       ArDudrFac    = vfFacMom*1.
                0232 C     o V momentum equation
                0233       ArDvdrFac    = vfFacMom*1.
                0234 
df1ec59c8a Jean*0235 C note: using standard stencil (no mask) results in under-estimating
                0236 C       vorticity at a no-slip boundary by a factor of 2 = sideDragFactor
                0237       IF ( no_slip_sides ) THEN
                0238         sideMaskFac = sideDragFactor
                0239       ELSE
                0240         sideMaskFac = 0. _d 0
                0241       ENDIF
                0242 
99731c717f Jean*0243       IF ( selectImplicitDrag.EQ.0 .AND.
                0244      &      (  no_slip_bottom
932b38363b Jean*0245      &    .OR. selectBotDragQuadr.GE.0
99731c717f Jean*0246      &    .OR. bottomDragLinear.NE.0. ) ) THEN
aea29c8517 Alis*0247        bottomDragTerms=.TRUE.
                0248       ELSE
                0249        bottomDragTerms=.FALSE.
                0250       ENDIF
                0251 
a340904e5a Ou W*0252 #if ( defined ALLOW_MOM_TEND_EXTRA_DIAGS && defined ALLOW_DIAGNOSTICS )
                0253       botDragDiagIsOn = .FALSE.
                0254       shelfDragDiagIsOn = .FALSE.
                0255       IF ( useDiagnostics ) THEN
                0256        IF ( bottomDragTerms ) botDragDiagIsOn =
                0257      &            DIAGNOSTICS_IS_ON( 'UBotDrag', myThid )
                0258      &       .OR. DIAGNOSTICS_IS_ON( 'VBotDrag', myThid )
                0259        IF ( useShelfIce ) shelfDragDiagIsOn =
                0260      &            DIAGNOSTICS_IS_ON( 'UShIDrag', myThid )
                0261      &       .OR. DIAGNOSTICS_IS_ON( 'VShIDrag', myThid )
                0262       ENDIF
                0263 #endif /* ALLOW_MOM_TEND_EXTRA_DIAGS && ALLOW_DIAGNOSTICS */
                0264 
aea29c8517 Alis*0265 C--   Calculate open water fraction at vorticity points
                0266       CALL MOM_CALC_HFACZ(bi,bj,k,hFacZ,r_hFacZ,myThid)
                0267 
                0268 C     Make local copies of horizontal flow field
                0269       DO j=1-OLy,sNy+OLy
                0270        DO i=1-OLx,sNx+OLx
                0271         uFld(i,j) = uVel(i,j,k,bi,bj)
                0272         vFld(i,j) = vVel(i,j,k,bi,bj)
                0273        ENDDO
                0274       ENDDO
                0275 
d09cd10c60 Gael*0276 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0277 CADJ STORE hFacZ(:,:)   = comlev1_bibj_k, key = kkey, byte = isbyte
                0278 CADJ STORE r_hFacZ(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
                0279 CADJ STORE fverukm(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
                0280 CADJ STORE fvervkm(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
d09cd10c60 Gael*0281 #endif
                0282 
cab1a69b8a Jean*0283 C note (jmc) : Dissipation and Vort3 advection do not necesary
                0284 C              use the same maskZ (and hFacZ)  => needs 2 call(s)
                0285 c     CALL MOM_VI_HFACZ_DISS(bi,bj,k,hFacZ,r_hFacZ,myThid)
                0286 
b8082fc644 Jean*0287       CALL MOM_CALC_KE(bi,bj,k,selectKEscheme,uFld,vFld,KE,myThid)
aea29c8517 Alis*0288 
7c7b0b4a46 Alis*0289       CALL MOM_CALC_RELVORT3(bi,bj,k,uFld,vFld,hFacZ,vort3,myThid)
aea29c8517 Alis*0290 
ed2dbbe83d Jean*0291 C-    mask vort3 and account for no-slip / free-slip BC in vort3BC:
                0292       DO j=1-OLy,sNy+OLy
                0293        DO i=1-OLx,sNx+OLx
                0294          vort3BC(i,j) = vort3(i,j)
                0295          IF ( hFacZ(i,j).EQ.zeroRS ) THEN
                0296            vort3BC(i,j) = sideMaskFac*vort3BC(i,j)
                0297            vort3(i,j)   = 0.
                0298          ENDIF
                0299        ENDDO
                0300       ENDDO
                0301 
d09cd10c60 Gael*0302 #ifdef ALLOW_AUTODIFF_TAMC
7448700841 Mart*0303 CADJ STORE KE(:,:)      = comlev1_bibj_k, key = kkey, byte = isbyte
edb6656069 Mart*0304 CADJ STORE vort3(:,:)   = comlev1_bibj_k, key = kkey, byte = isbyte
                0305 CADJ STORE vort3bc(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
d09cd10c60 Gael*0306 #endif
                0307 
aea29c8517 Alis*0308       IF (momViscosity) THEN
039fe61d35 Jean*0309 C--    For viscous term, compute horizontal divergence, tension & strain
df1ec59c8a Jean*0310 C      and mask relative vorticity (free-slip case):
                0311 
01f860d49e Jean*0312        DO j=1-OLy,sNy+OLy
                0313         DO i=1-OLx,sNx+OLx
                0314           h0FacZ(i,j) = hFacZ(i,j)
                0315         ENDDO
                0316        ENDDO
                0317 #ifdef NONLIN_FRSURF
                0318        IF ( no_slip_sides .AND. nonlinFreeSurf.GT.0 ) THEN
                0319         DO j=2-OLy,sNy+OLy
                0320          DO i=2-OLx,sNx+OLx
                0321           h0FacZ(i,j) = MIN(
                0322      &       MIN( h0FacW(i,j,k,bi,bj), h0FacW(i,j-1,k,bi,bj) ),
                0323      &       MIN( h0FacS(i,j,k,bi,bj), h0FacS(i-1,j,k,bi,bj) ) )
                0324          ENDDO
                0325         ENDDO
                0326        ENDIF
7448700841 Mart*0327 # ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0328 CADJ STORE h0FacZ(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
7448700841 Mart*0329 # endif
                0330 #endif /* NONLIN_FRSURF */
d09cd10c60 Gael*0331 
df1ec59c8a Jean*0332        CALL MOM_CALC_HDIV(bi,bj,k,2,uFld,vFld,hDiv,myThid)
                0333 
ed2dbbe83d Jean*0334        IF ( useVariableVisc .OR. useStrainTensionVisc ) THEN
                0335         CALL MOM_CALC_TENSION( bi,bj,k,uFld,vFld,tension,myThid )
                0336         CALL MOM_CALC_STRAIN( bi,bj,k,uFld,vFld,hFacZ,strain,myThid )
                0337 C-    mask strain and account for no-slip / free-slip BC in strainBC:
                0338         DO j=1-OLy,sNy+OLy
                0339          DO i=1-OLx,sNx+OLx
                0340            strainBC(i,j) = strain(i,j)
                0341            IF ( hFacZ(i,j).EQ.zeroRS ) THEN
                0342              strainBC(i,j) = sideMaskFac*strainBC(i,j)
                0343              strain(i,j)   = 0.
                0344            ENDIF
                0345          ENDDO
df1ec59c8a Jean*0346         ENDDO
f59d76b0dd Ed D*0347 #ifdef ALLOW_LEITH_QG
a125ab7be7 Jean*0348         IF ( viscC2LeithQG.NE.zeroRL ) THEN
f59d76b0dd Ed D*0349           CALL MOM_VISC_QGL_STRETCH(bi,bj,k,
                0350      O                            stretching, Nsquare,
                0351      I                            myTime, myIter, myThid )
                0352           CALL MOM_VISC_QGL_LIMIT(bi,bj,k,
                0353      O                          stretching,
                0354      I                          Nsquare, uFld,vFld,vort3,
                0355      I                          myTime, myIter, myThid )
                0356         ENDIF
                0357 #endif /* ALLOW_LEITH_QG */
ed2dbbe83d Jean*0358        ENDIF
df1ec59c8a Jean*0359 
d09cd10c60 Gael*0360 #ifdef ALLOW_AUTODIFF_TAMC
7448700841 Mart*0361 CADJ STORE hDiv(:,:)     = comlev1_bibj_k, key = kkey, byte = isbyte
edb6656069 Mart*0362 CADJ STORE tension(:,:)  = comlev1_bibj_k, key = kkey, byte = isbyte
                0363 CADJ STORE strain(:,:)   = comlev1_bibj_k, key = kkey, byte = isbyte
7448700841 Mart*0364 CADJ STORE strainBC(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
d09cd10c60 Gael*0365 #endif
                0366 
34e137f064 Jean*0367 C--    Calculate Lateral Viscosities
                0368        DO j=1-OLy,sNy+OLy
                0369         DO i=1-OLx,sNx+OLx
                0370          viscAh_D(i,j) = viscAhD
                0371          viscAh_Z(i,j) = viscAhZ
                0372          viscA4_D(i,j) = viscA4D
                0373          viscA4_Z(i,j) = viscA4Z
                0374         ENDDO
                0375        ENDDO
                0376        IF ( useVariableVisc ) THEN
ed2dbbe83d Jean*0377 C-     uses vort3BC & strainBC which account for no-slip / free-slip BC
79074ef66b Jean*0378 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0379 CADJ STORE KE(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
79074ef66b Jean*0380 #endif
34e137f064 Jean*0381          CALL MOM_CALC_VISC( bi, bj, k,
                0382      O            viscAh_Z, viscAh_D, viscA4_Z, viscA4_D,
f59d76b0dd Ed D*0383      I            hDiv, vort3BC, tension, strainBC, stretching,
                0384      I            KE, hfacZ, myThid )
d09cd10c60 Gael*0385 #ifdef ALLOW_AUTODIFF_TAMC
4240547d2d Mart*0386 # ifndef AUTODIFF_ALLOW_VISCFACADJ
                0387 C     These store directive must not be here, if you want to recompute
                0388 C     these viscosity coefficients with a modified viscFacAdj = viscFacInAd
                0389 C     because the store directives intentionally prevent the recomputation.
edb6656069 Mart*0390 CADJ STORE viscAh_Z(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
                0391 CADJ STORE viscAh_D(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
                0392 CADJ STORE viscA4_Z(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
                0393 CADJ STORE viscA4_D(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
4240547d2d Mart*0394 # endif /* AUTODIFF_ALLOW_VISCFACADJ */
                0395 #endif /* ALLOW_AUTODIFF_TAMC */
                0396        ENDIF
d09cd10c60 Gael*0397 
                0398 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0399 CADJ STORE hDiv(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
d09cd10c60 Gael*0400 #endif
                0401 
aea29c8517 Alis*0402 C      Calculate del^2 u and del^2 v for bi-harmonic term
f0501c53d1 Jean*0403        IF (useBiharmonicVisc) THEN
3a279374db Alis*0404          CALL MOM_VI_DEL2UV(bi,bj,k,hDiv,vort3,hFacZ,
                0405      O                      del2u,del2v,
ed2dbbe83d Jean*0406      I                      myThid)
88e5e58941 Jean*0407          CALL MOM_CALC_HDIV(bi,bj,k,2,del2u,del2v,dStar,myThid)
                0408          CALL MOM_CALC_RELVORT3(bi,bj,k,
                0409      &                          del2u,del2v,hFacZ,zStar,myThid)
df1ec59c8a Jean*0410        ENDIF
                0411 
d09cd10c60 Gael*0412 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0413 CADJ STORE del2u(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
                0414 CADJ STORE del2v(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
                0415 CADJ STORE dStar(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
                0416 CADJ STORE zStar(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
d09cd10c60 Gael*0417 #endif
                0418 
df1ec59c8a Jean*0419 C---   Calculate dissipation terms for U and V equations
b0c3bd7ab0 Bayl*0420 
ed2dbbe83d Jean*0421 C-     in terms of tension and strain
b0c3bd7ab0 Bayl*0422        IF (useStrainTensionVisc) THEN
ed2dbbe83d Jean*0423 C      use masked strain as if free-slip since side-drag is computed separately
f0501c53d1 Jean*0424          CALL MOM_HDISSIP( bi, bj, k,
ed2dbbe83d Jean*0425      I            tension, strain, hFacZ,
f0501c53d1 Jean*0426      I            viscAh_Z, viscAh_D, viscA4_Z, viscA4_D,
                0427      I            useHarmonicVisc, useBiharmonicVisc, useVariableVisc,
                0428      O            guDiss, gvDiss,
                0429      I            myThid )
b0c3bd7ab0 Bayl*0430        ELSE
ed2dbbe83d Jean*0431 C-     in terms of vorticity and divergence
f0501c53d1 Jean*0432          CALL MOM_VI_HDISSIP( bi, bj, k,
ed2dbbe83d Jean*0433      I            hDiv, vort3, dStar, zStar, hFacZ,
f0501c53d1 Jean*0434      I            viscAh_Z, viscAh_D, viscA4_Z, viscA4_D,
                0435      I            useHarmonicVisc, useBiharmonicVisc, useVariableVisc,
                0436      O            guDiss, gvDiss,
ed2dbbe83d Jean*0437      I            myThid )
07cc642809 Alis*0438        ENDIF
cab1a69b8a Jean*0439 
df1ec59c8a Jean*0440 C---  Other dissipation terms in Zonal momentum equation
aea29c8517 Alis*0441 
                0442 C--   Vertical flux (fVer is at upper face of "u" cell)
                0443 C     Eddy component of vertical flux (interior component only) -> vrF
79074ef66b Jean*0444        IF ( .NOT.implicitViscosity ) THEN
                0445         CALL MOM_U_RVISCFLUX(bi,bj,k+1,uVel,kappaRU,vrF,myThid)
aea29c8517 Alis*0446 C     Combine fluxes
79074ef66b Jean*0447         DO j=jMin,jMax
                0448          DO i=iMin,iMax
                0449           fVerUkp(i,j) = ArDudrFac*vrF(i,j)
                0450          ENDDO
4667e4066d Jean*0451         ENDDO
d09cd10c60 Gael*0452 
                0453 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0454 CADJ STORE fVerUkp(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
d09cd10c60 Gael*0455 #endif
                0456 
4667e4066d Jean*0457 C--   Tendency is minus divergence of the fluxes
f1a4cec01a Jean*0458 C     vert.visc.flx is scaled by deepFac2F (deep-atmos) and rhoFac (anelastic)
79074ef66b Jean*0459         DO j=jMin,jMax
                0460          DO i=iMin,iMax
                0461           guDiss(i,j) = guDiss(i,j)
                0462      &    -_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)
                0463      &    *recip_rAw(i,j,bi,bj)
                0464      &    *( fVerUkp(i,j) - fVerUkm(i,j) )*rkSign
                0465      &    *recip_deepFac2C(k)*recip_rhoFacC(k)
                0466          ENDDO
4667e4066d Jean*0467         ENDDO
79074ef66b Jean*0468        ENDIF
aea29c8517 Alis*0469 
039fe61d35 Jean*0470 C-- No-slip and drag BCs appear as body forces in cell abutting topography
79074ef66b Jean*0471        IF ( no_slip_sides ) THEN
aea29c8517 Alis*0472 C-     No-slip BCs impose a drag at walls...
79074ef66b Jean*0473         CALL MOM_U_SIDEDRAG( bi, bj, k,
                0474      I           uFld, del2u, h0FacZ,
                0475      I           viscAh_Z, viscA4_Z,
                0476      I           useHarmonicVisc, useBiharmonicVisc, useVariableVisc,
                0477      O           vF,
                0478      I           myThid )
                0479         DO j=jMin,jMax
                0480          DO i=iMin,iMax
                0481           guDiss(i,j) = guDiss(i,j)+vF(i,j)
                0482          ENDDO
aea29c8517 Alis*0483         ENDDO
79074ef66b Jean*0484        ENDIF
34e137f064 Jean*0485 
aea29c8517 Alis*0486 C-    No-slip BCs impose a drag at bottom
79074ef66b Jean*0487        IF ( bottomDragTerms ) THEN
                0488 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0489 CADJ STORE KE(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
79074ef66b Jean*0490 #endif
a125ab7be7 Jean*0491         CALL MOM_U_BOTDRAG_COEFF( bi, bj, k, .TRUE.,
79074ef66b Jean*0492      I                  uFld, vFld, kappaRU, KE,
                0493      O                  cDrag,
                0494      I                  myIter, myThid )
7448700841 Mart*0495 #ifdef ALLOW_AUTODIFF_TAMC
                0496 CADJ STORE cDrag(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
                0497 #endif
79074ef66b Jean*0498         DO j=jMin,jMax
                0499          DO i=iMin,iMax
a340904e5a Ou W*0500            vF(i,j) = -cDrag(i,j)*uFld(i,j)
                0501      &             *_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)
                0502            gUdiss(i,j) = gUdiss(i,j) + vF(i,j)
79074ef66b Jean*0503          ENDDO
aea29c8517 Alis*0504         ENDDO
79074ef66b Jean*0505         IF ( useDiagnostics ) THEN
                0506          DO j=jMin,jMax
                0507           DO i=iMin,iMax
                0508            botDragU(i,j,bi,bj) = botDragU(i,j,bi,bj)
                0509      &                         - cDrag(i,j)*uFld(i,j)*rUnit2mass
                0510           ENDDO
                0511          ENDDO
                0512         ENDIF
a340904e5a Ou W*0513 #if ( defined ALLOW_MOM_TEND_EXTRA_DIAGS && defined ALLOW_DIAGNOSTICS )
                0514         IF ( botDragDiagIsOn ) THEN
                0515          CALL DIAGNOSTICS_FILL( vF, 'UBotDrag', k,1,2,bi,bj, myThid )
                0516         ENDIF
                0517 #endif /* ALLOW_MOM_TEND_EXTRA_DIAGS && ALLOW_DIAGNOSTICS */
79074ef66b Jean*0518        ENDIF
dd49642a1d Mart*0519 #ifdef ALLOW_SHELFICE
e2d988bd46 Jean*0520        IF ( useShelfIce .AND. selectImplicitDrag.EQ.0 ) THEN
79074ef66b Jean*0521 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0522 CADJ STORE KE(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
79074ef66b Jean*0523 #endif
e2d988bd46 Jean*0524         CALL SHELFICE_U_DRAG_COEFF( bi, bj, k, .TRUE.,
                0525      I                  uFld, vFld, kappaRU, KE,
                0526      O                  cDrag,
                0527      I                  myIter, myThid )
7448700841 Mart*0528 #ifdef ALLOW_AUTODIFF_TAMC
                0529 CADJ STORE cDrag(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
                0530 #endif
79074ef66b Jean*0531         DO j=jMin,jMax
                0532          DO i=iMin,iMax
e2d988bd46 Jean*0533            gUdiss(i,j) = gUdiss(i,j)
                0534      &                 - cDrag(i,j)*uFld(i,j)
                0535      &                 *_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)
79074ef66b Jean*0536          ENDDO
dd49642a1d Mart*0537         ENDDO
a340904e5a Ou W*0538 #if ( defined ALLOW_MOM_TEND_EXTRA_DIAGS && defined ALLOW_DIAGNOSTICS )
                0539         IF ( shelfDragDiagIsOn ) THEN
                0540          DO j=jMin,jMax
                0541           DO i=iMin,iMax
                0542            vrF(i,j) = -cDrag(i,j)*uFld(i,j)
                0543      &              *_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)
                0544           ENDDO
                0545          ENDDO
                0546          CALL DIAGNOSTICS_FILL( vrF, 'UShIDrag', k,1,2,bi,bj, myThid )
                0547         ENDIF
                0548 #endif /* ALLOW_MOM_TEND_EXTRA_DIAGS && ALLOW_DIAGNOSTICS */
79074ef66b Jean*0549        ENDIF
dd49642a1d Mart*0550 #endif /* ALLOW_SHELFICE */
                0551 
df1ec59c8a Jean*0552 C---  Other dissipation terms in Meridional momentum equation
aea29c8517 Alis*0553 
                0554 C--   Vertical flux (fVer is at upper face of "v" cell)
                0555 C     Eddy component of vertical flux (interior component only) -> vrF
79074ef66b Jean*0556        IF ( .NOT.implicitViscosity ) THEN
                0557         CALL MOM_V_RVISCFLUX(bi,bj,k+1,vVel,kappaRV,vrF,myThid)
aea29c8517 Alis*0558 C     Combine fluxes -> fVerV
79074ef66b Jean*0559         DO j=jMin,jMax
                0560          DO i=iMin,iMax
                0561           fVerVkp(i,j) = ArDvdrFac*vrF(i,j)
                0562          ENDDO
4667e4066d Jean*0563         ENDDO
d09cd10c60 Gael*0564 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0565 CADJ STORE fVerVkp(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
d09cd10c60 Gael*0566 #endif
4667e4066d Jean*0567 C--   Tendency is minus divergence of the fluxes
f1a4cec01a Jean*0568 C     vert.visc.flx is scaled by deepFac2F (deep-atmos) and rhoFac (anelastic)
79074ef66b Jean*0569         DO j=jMin,jMax
                0570          DO i=iMin,iMax
                0571           gvDiss(i,j) = gvDiss(i,j)
                0572      &    -_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
                0573      &    *recip_rAs(i,j,bi,bj)
                0574      &    *( fVerVkp(i,j) - fVerVkm(i,j) )*rkSign
                0575      &    *recip_deepFac2C(k)*recip_rhoFacC(k)
                0576          ENDDO
4667e4066d Jean*0577         ENDDO
79074ef66b Jean*0578        ENDIF
aea29c8517 Alis*0579 
039fe61d35 Jean*0580 C-- No-slip and drag BCs appear as body forces in cell abutting topography
79074ef66b Jean*0581        IF ( no_slip_sides ) THEN
aea29c8517 Alis*0582 C-     No-slip BCs impose a drag at walls...
79074ef66b Jean*0583         CALL MOM_V_SIDEDRAG( bi, bj, k,
                0584      I           vFld, del2v, h0FacZ,
                0585      I           viscAh_Z, viscA4_Z,
                0586      I           useHarmonicVisc, useBiharmonicVisc, useVariableVisc,
                0587      O           vF,
                0588      I           myThid )
                0589         DO j=jMin,jMax
                0590          DO i=iMin,iMax
                0591           gvDiss(i,j) = gvDiss(i,j)+vF(i,j)
                0592          ENDDO
aea29c8517 Alis*0593         ENDDO
79074ef66b Jean*0594        ENDIF
34e137f064 Jean*0595 
aea29c8517 Alis*0596 C-    No-slip BCs impose a drag at bottom
79074ef66b Jean*0597        IF ( bottomDragTerms ) THEN
                0598 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0599 CADJ STORE KE(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
79074ef66b Jean*0600 #endif
a125ab7be7 Jean*0601         CALL MOM_V_BOTDRAG_COEFF( bi, bj, k, .TRUE.,
79074ef66b Jean*0602      I                  uFld, vFld, kappaRV, KE,
                0603      O                  cDrag,
                0604      I                  myIter, myThid )
7448700841 Mart*0605 #ifdef ALLOW_AUTODIFF_TAMC
                0606 CADJ STORE cDrag(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
                0607 #endif
79074ef66b Jean*0608         DO j=jMin,jMax
                0609          DO i=iMin,iMax
a340904e5a Ou W*0610            vF(i,j) = -cDrag(i,j)*vFld(i,j)
                0611      &             *_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
                0612            gvDiss(i,j) = gvDiss(i,j) + vF(i,j)
79074ef66b Jean*0613          ENDDO
aea29c8517 Alis*0614         ENDDO
79074ef66b Jean*0615         IF ( useDiagnostics ) THEN
                0616          DO j=jMin,jMax
                0617           DO i=iMin,iMax
                0618            botDragV(i,j,bi,bj) = botDragV(i,j,bi,bj)
                0619      &                         - cDrag(i,j)*vFld(i,j)*rUnit2mass
                0620           ENDDO
                0621          ENDDO
                0622         ENDIF
a340904e5a Ou W*0623 #if ( defined ALLOW_MOM_TEND_EXTRA_DIAGS && defined ALLOW_DIAGNOSTICS )
                0624         IF ( botDragDiagIsOn ) THEN
                0625          CALL DIAGNOSTICS_FILL( vF, 'VBotDrag', k,1,2,bi,bj, myThid )
                0626         ENDIF
                0627 #endif /* ALLOW_MOM_TEND_EXTRA_DIAGS && ALLOW_DIAGNOSTICS */
79074ef66b Jean*0628        ENDIF
dd49642a1d Mart*0629 #ifdef ALLOW_SHELFICE
e2d988bd46 Jean*0630        IF ( useShelfIce .AND. selectImplicitDrag.EQ.0 ) THEN
79074ef66b Jean*0631 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0632 CADJ STORE KE(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
79074ef66b Jean*0633 #endif
e2d988bd46 Jean*0634         CALL SHELFICE_V_DRAG_COEFF( bi, bj, k, .TRUE.,
                0635      I                  uFld, vFld, kappaRV, KE,
                0636      O                  cDrag,
                0637      I                  myIter, myThid )
7448700841 Mart*0638 #ifdef ALLOW_AUTODIFF_TAMC
                0639 CADJ STORE cDrag(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
                0640 #endif
79074ef66b Jean*0641         DO j=jMin,jMax
                0642          DO i=iMin,iMax
e2d988bd46 Jean*0643            gvDiss(i,j) = gvDiss(i,j)
                0644      &                 - cDrag(i,j)*vFld(i,j)
                0645      &                 *_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
79074ef66b Jean*0646          ENDDO
932b38363b Jean*0647         ENDDO
a340904e5a Ou W*0648 #if ( defined ALLOW_MOM_TEND_EXTRA_DIAGS && defined ALLOW_DIAGNOSTICS )
                0649         IF ( shelfDragDiagIsOn ) THEN
                0650          DO j=jMin,jMax
                0651           DO i=iMin,iMax
                0652            vrF(i,j) = -cDrag(i,j)*vFld(i,j)
                0653      &              *_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
                0654           ENDDO
                0655          ENDDO
                0656          CALL DIAGNOSTICS_FILL( vrF, 'VShIDrag', k,1,2,bi,bj, myThid )
                0657         ENDIF
                0658 #endif /* ALLOW_MOM_TEND_EXTRA_DIAGS && ALLOW_DIAGNOSTICS */
79074ef66b Jean*0659        ENDIF
dd49642a1d Mart*0660 #endif /* ALLOW_SHELFICE */
                0661 
34e137f064 Jean*0662 C--   if (momViscosity) end of block.
                0663       ENDIF
                0664 
                0665 C-    Return to standard hfacZ (min-4) and mask vort3 accordingly:
                0666 c     CALL MOM_VI_MASK_VORT3(bi,bj,k,hFacZ,r_hFacZ,vort3,myThid)
aea29c8517 Alis*0667 
df1ec59c8a Jean*0668 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0669 
                0670 C---  Prepare for Advection & Coriolis terms:
ed2dbbe83d Jean*0671 C-    calculate absolute vorticity
df1ec59c8a Jean*0672       IF (useAbsVorticity)
                0673      &  CALL MOM_CALC_ABSVORT3(bi,bj,k,vort3,omega3,myThid)
aea29c8517 Alis*0674 
d09cd10c60 Gael*0675 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0676 CADJ STORE omega3(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
d09cd10c60 Gael*0677 #endif
                0678 
aea29c8517 Alis*0679 C--   Horizontal Coriolis terms
a75214ff78 Jean*0680 c     IF (useCoriolis .AND. .NOT.useCDscheme
                0681 c    &    .AND. .NOT. useAbsVorticity) THEN
                0682 C- jmc: change it to keep the Coriolis terms when useAbsVorticity=T & momAdvection=F
711329b07b Jean*0683       IF ( useCoriolis .AND.
a75214ff78 Jean*0684      &     .NOT.( useCDscheme .OR. useAbsVorticity.AND.momAdvection )
                0685      &   ) THEN
                0686        IF (useAbsVorticity) THEN
7c3c2cec96 Jean*0687         CALL MOM_VI_U_CORIOLIS(bi,bj,k,selectVortScheme,useJamartMomAdv,
                0688      &                         vFld,omega3,hFacZ,r_hFacZ,
a75214ff78 Jean*0689      &                         uCf,myThid)
7c3c2cec96 Jean*0690         CALL MOM_VI_V_CORIOLIS(bi,bj,k,selectVortScheme,useJamartMomAdv,
                0691      &                         uFld,omega3,hFacZ,r_hFacZ,
a75214ff78 Jean*0692      &                         vCf,myThid)
9293d3c672 Hajo*0693 #if ( defined ALLOW_GGL90 && defined ALLOW_GGL90_LANGMUIR )
                0694        ELSEIF ( useLANGMUIR ) THEN
                0695         CALL GGL90_ADD_STOKESDRIFT(
                0696      O                 uRes, vRes,
                0697      I                 uFld, vFld, k, bi, bj, myThid )
                0698         CALL MOM_VI_CORIOLIS(bi,bj,k,uRes,vRes,hFacZ,r_hFacZ,
                0699      &                       uCf,vCf,myThid)
                0700 #endif
a75214ff78 Jean*0701        ELSE
                0702         CALL MOM_VI_CORIOLIS(bi,bj,k,uFld,vFld,hFacZ,r_hFacZ,
                0703      &                       uCf,vCf,myThid)
                0704        ENDIF
481f592dad Jean*0705        DO j=jMin,jMax
                0706         DO i=iMin,iMax
1aff67ca82 Jean*0707          gU(i,j,k,bi,bj) = uCf(i,j)
                0708          gV(i,j,k,bi,bj) = vCf(i,j)
481f592dad Jean*0709         ENDDO
aea29c8517 Alis*0710        ENDDO
f7d48db11c Jean*0711        IF ( writeDiag ) THEN
ef1e764710 Ed H*0712          IF (snapshot_mdsio) THEN
                0713            CALL WRITE_LOCAL_RL('fV','I10',1,uCf,bi,bj,k,myIter,myThid)
                0714            CALL WRITE_LOCAL_RL('fU','I10',1,vCf,bi,bj,k,myIter,myThid)
                0715          ENDIF
                0716 #ifdef ALLOW_MNC
                0717          IF (useMNC .AND. snapshot_mnc) THEN
b22b541fe9 Ed H*0718            CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj, 'fV', uCf,
cb356b4c5f Ed H*0719      &          offsets, myThid)
b22b541fe9 Ed H*0720            CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj, 'fU', vCf,
cb356b4c5f Ed H*0721      &          offsets, myThid)
ef1e764710 Ed H*0722          ENDIF
                0723 #endif /*  ALLOW_MNC  */
f7d48db11c Jean*0724        ENDIF
711329b07b Jean*0725 #ifdef ALLOW_DIAGNOSTICS
                0726        IF ( useDiagnostics ) THEN
                0727          CALL DIAGNOSTICS_FILL(uCf,'Um_Cori ',k,1,2,bi,bj,myThid)
                0728          CALL DIAGNOSTICS_FILL(vCf,'Vm_Cori ',k,1,2,bi,bj,myThid)
                0729        ENDIF
                0730 #endif /* ALLOW_DIAGNOSTICS */
4667e4066d Jean*0731       ELSE
                0732        DO j=jMin,jMax
                0733         DO i=iMin,iMax
1aff67ca82 Jean*0734          gU(i,j,k,bi,bj) = 0. _d 0
                0735          gV(i,j,k,bi,bj) = 0. _d 0
4667e4066d Jean*0736         ENDDO
                0737        ENDDO
481f592dad Jean*0738       ENDIF
                0739 
d09cd10c60 Gael*0740 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0741 CADJ STORE ucf(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
                0742 CADJ STORE vcf(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
d09cd10c60 Gael*0743 #endif
                0744 
481f592dad Jean*0745       IF (momAdvection) THEN
3add9696e1 Jean*0746 C--   Horizontal advection of relative (or absolute) vorticity
7fe6343684 Jean*0747        IF ( (highOrderVorticity.OR.upwindVorticity)
                0748      &     .AND.useAbsVorticity ) THEN
79074ef66b Jean*0749         CALL MOM_VI_U_CORIOLIS_C4(bi,bj,k,selectVortScheme,
3370e71df9 Mart*0750      &                         highOrderVorticity,upwindVorticity,
                0751      &                         vFld,omega3,r_hFacZ,
5d5fff3a81 Alis*0752      &                         uCf,myThid)
7fe6343684 Jean*0753        ELSEIF ( (highOrderVorticity.OR.upwindVorticity) ) THEN
79074ef66b Jean*0754         CALL MOM_VI_U_CORIOLIS_C4(bi,bj,k,selectVortScheme,
3370e71df9 Mart*0755      &                         highOrderVorticity, upwindVorticity,
                0756      &                         vFld,vort3, r_hFacZ,
3add9696e1 Jean*0757      &                         uCf,myThid)
7fe6343684 Jean*0758        ELSEIF ( useAbsVorticity ) THEN
3370e71df9 Mart*0759         CALL MOM_VI_U_CORIOLIS(bi,bj,k,selectVortScheme,useJamartMomAdv,
                0760      &                         vFld,omega3,hFacZ,r_hFacZ,
d4c99033f5 Jean*0761      &                         uCf,myThid)
5d5fff3a81 Alis*0762        ELSE
3370e71df9 Mart*0763         CALL MOM_VI_U_CORIOLIS(bi,bj,k,selectVortScheme,useJamartMomAdv,
                0764      &                         vFld,vort3, hFacZ,r_hFacZ,
5d5fff3a81 Alis*0765      &                         uCf,myThid)
                0766        ENDIF
481f592dad Jean*0767        DO j=jMin,jMax
                0768         DO i=iMin,iMax
                0769          gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
                0770         ENDDO
aea29c8517 Alis*0771        ENDDO
7fe6343684 Jean*0772        IF ( (highOrderVorticity.OR.upwindVorticity)
                0773      &     .AND.useAbsVorticity ) THEN
79074ef66b Jean*0774         CALL MOM_VI_V_CORIOLIS_C4(bi,bj,k,selectVortScheme,
3370e71df9 Mart*0775      &                         highOrderVorticity, upwindVorticity,
                0776      &                         uFld,omega3,r_hFacZ,
5d5fff3a81 Alis*0777      &                         vCf,myThid)
7fe6343684 Jean*0778        ELSEIF ( (highOrderVorticity.OR.upwindVorticity) ) THEN
79074ef66b Jean*0779         CALL MOM_VI_V_CORIOLIS_C4(bi,bj,k,selectVortScheme,
3370e71df9 Mart*0780      &                         highOrderVorticity, upwindVorticity,
                0781      &                         uFld,vort3, r_hFacZ,
3add9696e1 Jean*0782      &                         vCf,myThid)
7fe6343684 Jean*0783        ELSEIF ( useAbsVorticity ) THEN
3370e71df9 Mart*0784         CALL MOM_VI_V_CORIOLIS(bi,bj,k,selectVortScheme,useJamartMomAdv,
                0785      &                         uFld,omega3,hFacZ,r_hFacZ,
d4c99033f5 Jean*0786      &                         vCf,myThid)
5d5fff3a81 Alis*0787        ELSE
3370e71df9 Mart*0788         CALL MOM_VI_V_CORIOLIS(bi,bj,k,selectVortScheme,useJamartMomAdv,
                0789      &                         uFld,vort3, hFacZ,r_hFacZ,
5d5fff3a81 Alis*0790      &                         vCf,myThid)
                0791        ENDIF
481f592dad Jean*0792        DO j=jMin,jMax
                0793         DO i=iMin,iMax
                0794          gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
                0795         ENDDO
aea29c8517 Alis*0796        ENDDO
                0797 
d09cd10c60 Gael*0798 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0799 CADJ STORE ucf(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
                0800 CADJ STORE vcf(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte
d09cd10c60 Gael*0801 #endif
                0802 
f7d48db11c Jean*0803        IF ( writeDiag ) THEN
ef1e764710 Ed H*0804          IF (snapshot_mdsio) THEN
                0805            CALL WRITE_LOCAL_RL('zV','I10',1,uCf,bi,bj,k,myIter,myThid)
                0806            CALL WRITE_LOCAL_RL('zU','I10',1,vCf,bi,bj,k,myIter,myThid)
                0807          ENDIF
                0808 #ifdef ALLOW_MNC
                0809          IF (useMNC .AND. snapshot_mnc) THEN
b22b541fe9 Ed H*0810            CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj, 'zV', uCf,
cb356b4c5f Ed H*0811      &          offsets, myThid)
b22b541fe9 Ed H*0812            CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj, 'zU', vCf,
cb356b4c5f Ed H*0813      &          offsets, myThid)
ef1e764710 Ed H*0814          ENDIF
                0815 #endif /*  ALLOW_MNC  */
f7d48db11c Jean*0816        ENDIF
ef1e764710 Ed H*0817 
cab1a69b8a Jean*0818 #ifdef ALLOW_TIMEAVE
                0819        IF (taveFreq.GT.0.) THEN
                0820          CALL TIMEAVE_CUMUL_1K1T(uZetatave,vCf,deltaTClock,
                0821      &                           Nr, k, bi, bj, myThid)
                0822          CALL TIMEAVE_CUMUL_1K1T(vZetatave,uCf,deltaTClock,
                0823      &                           Nr, k, bi, bj, myThid)
                0824        ENDIF
5751fc37e0 Jean*0825 #endif /* ALLOW_TIMEAVE */
711329b07b Jean*0826 #ifdef ALLOW_DIAGNOSTICS
                0827        IF ( useDiagnostics ) THEN
                0828          CALL DIAGNOSTICS_FILL(uCf,'Um_AdvZ3',k,1,2,bi,bj,myThid)
                0829          CALL DIAGNOSTICS_FILL(vCf,'Vm_AdvZ3',k,1,2,bi,bj,myThid)
                0830        ENDIF
                0831 #endif /* ALLOW_DIAGNOSTICS */
cab1a69b8a Jean*0832 
481f592dad Jean*0833 C--   Vertical shear terms (-w*du/dr & -w*dv/dr)
6d4b75df68 Jean*0834        IF ( .NOT. momImplVertAdv ) THEN
f1a4cec01a Jean*0835         CALL MOM_VI_U_VERTSHEAR(bi,bj,k,uVel,wVel,uCf,myThid)
6d4b75df68 Jean*0836         DO j=jMin,jMax
                0837          DO i=iMin,iMax
                0838           gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
                0839          ENDDO
481f592dad Jean*0840         ENDDO
f1a4cec01a Jean*0841         CALL MOM_VI_V_VERTSHEAR(bi,bj,k,vVel,wVel,vCf,myThid)
6d4b75df68 Jean*0842         DO j=jMin,jMax
                0843          DO i=iMin,iMax
                0844           gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
                0845          ENDDO
481f592dad Jean*0846         ENDDO
711329b07b Jean*0847 #ifdef ALLOW_DIAGNOSTICS
                0848         IF ( useDiagnostics ) THEN
                0849          CALL DIAGNOSTICS_FILL(uCf,'Um_AdvRe',k,1,2,bi,bj,myThid)
                0850          CALL DIAGNOSTICS_FILL(vCf,'Vm_AdvRe',k,1,2,bi,bj,myThid)
                0851         ENDIF
                0852 #endif /* ALLOW_DIAGNOSTICS */
6d4b75df68 Jean*0853        ENDIF
aea29c8517 Alis*0854 
                0855 C--   Bernoulli term
f1a4cec01a Jean*0856        CALL MOM_VI_U_GRAD_KE(bi,bj,k,KE,uCf,myThid)
481f592dad Jean*0857        DO j=jMin,jMax
                0858         DO i=iMin,iMax
                0859          gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
                0860         ENDDO
aea29c8517 Alis*0861        ENDDO
f1a4cec01a Jean*0862        CALL MOM_VI_V_GRAD_KE(bi,bj,k,KE,vCf,myThid)
481f592dad Jean*0863        DO j=jMin,jMax
                0864         DO i=iMin,iMax
                0865          gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
                0866         ENDDO
                0867        ENDDO
f7d48db11c Jean*0868        IF ( writeDiag ) THEN
ef1e764710 Ed H*0869          IF (snapshot_mdsio) THEN
                0870            CALL WRITE_LOCAL_RL('KEx','I10',1,uCf,bi,bj,k,myIter,myThid)
                0871            CALL WRITE_LOCAL_RL('KEy','I10',1,vCf,bi,bj,k,myIter,myThid)
                0872          ENDIF
                0873 #ifdef ALLOW_MNC
                0874          IF (useMNC .AND. snapshot_mnc) THEN
b22b541fe9 Ed H*0875            CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj, 'KEx', uCf,
cb356b4c5f Ed H*0876      &          offsets, myThid)
b22b541fe9 Ed H*0877            CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj, 'KEy', vCf,
cb356b4c5f Ed H*0878      &          offsets, myThid)
df1ec59c8a Jean*0879          ENDIF
ef1e764710 Ed H*0880 #endif /*  ALLOW_MNC  */
f7d48db11c Jean*0881        ENDIF
                0882 
481f592dad Jean*0883 C--   end if momAdvection
                0884       ENDIF
                0885 
3daafce20b Jean*0886 C--   3.D Coriolis term (horizontal momentum, Eastward component: -fprime*w)
427e24e121 Jean*0887       IF ( select3dCoriScheme.GE.1 ) THEN
039fe61d35 Jean*0888         CALL MOM_U_CORIOLIS_NH(bi,bj,k,wVel,uCf,myThid)
                0889         DO j=jMin,jMax
                0890          DO i=iMin,iMax
                0891           gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
                0892          ENDDO
                0893         ENDDO
427e24e121 Jean*0894        IF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
039fe61d35 Jean*0895 C-     presently, non zero angleSinC array only supported with Curvilinear-Grid
                0896         CALL MOM_V_CORIOLIS_NH(bi,bj,k,wVel,vCf,myThid)
                0897         DO j=jMin,jMax
                0898          DO i=iMin,iMax
                0899           gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
                0900          ENDDO
                0901         ENDDO
                0902        ENDIF
                0903       ENDIF
                0904 
                0905 C--   Non-Hydrostatic (spherical) metric terms
                0906       IF ( useNHMTerms ) THEN
                0907        CALL MOM_U_METRIC_NH(bi,bj,k,uFld,wVel,uCf,myThid)
                0908        DO j=jMin,jMax
                0909         DO i=iMin,iMax
                0910          gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+uCf(i,j)
                0911         ENDDO
                0912        ENDDO
                0913        CALL MOM_V_METRIC_NH(bi,bj,k,vFld,wVel,vCf,myThid)
                0914        DO j=jMin,jMax
                0915         DO i=iMin,iMax
                0916          gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+vCf(i,j)
                0917         ENDDO
                0918        ENDDO
                0919       ENDIF
df1ec59c8a Jean*0920 
481f592dad Jean*0921 C--   Set du/dt & dv/dt on boundaries to zero
aea29c8517 Alis*0922       DO j=jMin,jMax
                0923        DO i=iMin,iMax
481f592dad Jean*0924         gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)*_maskW(i,j,k,bi,bj)
                0925         gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)*_maskS(i,j,k,bi,bj)
aea29c8517 Alis*0926        ENDDO
                0927       ENDDO
481f592dad Jean*0928 
5751fc37e0 Jean*0929 #ifdef ALLOW_DEBUG
8830b8f970 Jean*0930       IF ( debugLevel .GE. debLevC
5751fc37e0 Jean*0931      &   .AND. k.EQ.4 .AND. myIter.EQ.nIter0
                0932      &   .AND. nPx.EQ.1 .AND. nPy.EQ.1
                0933      &   .AND. useCubedSphereExchange ) THEN
e85db1faec Jean*0934         CALL DEBUG_CS_CORNER_UV( ' uDiss,vDiss from MOM_VECINV',
4667e4066d Jean*0935      &             guDiss,gvDiss, k, standardMessageUnit,bi,bj,myThid )
5751fc37e0 Jean*0936       ENDIF
                0937 #endif /* ALLOW_DEBUG */
aea29c8517 Alis*0938 
f7d48db11c Jean*0939       IF ( writeDiag ) THEN
ed2dbbe83d Jean*0940         IF (useBiharmonicVisc) THEN
                0941          CALL WRITE_LOCAL_RL( 'del2u', 'I10', 1, del2u,
                0942      &                        bi,bj,k, myIter, myThid )
                0943          CALL WRITE_LOCAL_RL( 'del2v', 'I10', 1, del2v,
                0944      &                        bi,bj,k, myIter, myThid )
                0945          CALL WRITE_LOCAL_RL( 'dStar', 'I10', 1, dStar,
                0946      &                        bi,bj,k, myIter, myThid )
                0947          CALL WRITE_LOCAL_RL( 'zStar', 'I10', 1, zStar,
                0948      &                        bi,bj,k, myIter, myThid )
                0949         ENDIF
ef1e764710 Ed H*0950         IF (snapshot_mdsio) THEN
df1ec59c8a Jean*0951          CALL WRITE_LOCAL_RL('W3','I10',1,omega3, bi,bj,k,myIter,myThid)
ed2dbbe83d Jean*0952          CALL WRITE_LOCAL_RL('Z3','I10',1,vort3BC,bi,bj,k,myIter,myThid)
df1ec59c8a Jean*0953          CALL WRITE_LOCAL_RL('KE','I10',1,KE,     bi,bj,k,myIter,myThid)
                0954          CALL WRITE_LOCAL_RL('D', 'I10',1,hDiv,   bi,bj,k,myIter,myThid)
                0955          CALL WRITE_LOCAL_RL('Dt','I10',1,tension,bi,bj,k,myIter,myThid)
ed2dbbe83d Jean*0956          CALL WRITE_LOCAL_RL( 'Ds', 'I10', 1, strainBC,
                0957      &                        bi,bj,k, myIter, myThid )
df1ec59c8a Jean*0958          CALL WRITE_LOCAL_RL('Du','I10',1,guDiss, bi,bj,k,myIter,myThid)
                0959          CALL WRITE_LOCAL_RL('Dv','I10',1,gvDiss, bi,bj,k,myIter,myThid)
ef1e764710 Ed H*0960         ENDIF
                0961 #ifdef ALLOW_MNC
                0962         IF (useMNC .AND. snapshot_mnc) THEN
df1ec59c8a Jean*0963           CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'W3',omega3,
                0964      &          offsets, myThid)
ed2dbbe83d Jean*0965           CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'Z3',vort3BC,
                0966      &          offsets, myThid)
df1ec59c8a Jean*0967           CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'KE',KE,
                0968      &          offsets, myThid)
                0969           CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'D', hDiv,
cb356b4c5f Ed H*0970      &          offsets, myThid)
b22b541fe9 Ed H*0971           CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'Dt',tension,
cb356b4c5f Ed H*0972      &          offsets, myThid)
ed2dbbe83d Jean*0973           CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'Ds',strainBC,
                0974      &          offsets, myThid)
b22b541fe9 Ed H*0975           CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'Du',guDiss,
cb356b4c5f Ed H*0976      &          offsets, myThid)
b22b541fe9 Ed H*0977           CALL MNC_CW_RL_W_OFFSET(pf,'mom_vi',bi,bj,'Dv',gvDiss,
cb356b4c5f Ed H*0978      &          offsets, myThid)
ef1e764710 Ed H*0979         ENDIF
                0980 #endif /*  ALLOW_MNC  */
3a279374db Alis*0981       ENDIF
3add9696e1 Jean*0982 
711329b07b Jean*0983 #ifdef ALLOW_DIAGNOSTICS
                0984       IF ( useDiagnostics ) THEN
ed2dbbe83d Jean*0985         CALL DIAGNOSTICS_FILL(vort3BC,'momVort3',k,1,2,bi,bj,myThid)
df1ec59c8a Jean*0986         CALL DIAGNOSTICS_FILL(KE,     'momKE   ',k,1,2,bi,bj,myThid)
711329b07b Jean*0987        IF (momViscosity) THEN
df1ec59c8a Jean*0988         CALL DIAGNOSTICS_FILL(hDiv,   'momHDiv ',k,1,2,bi,bj,myThid)
711329b07b Jean*0989        ENDIF
ed2dbbe83d Jean*0990        IF ( useVariableVisc .OR. useStrainTensionVisc ) THEN
                0991         CALL DIAGNOSTICS_FILL(tension, 'Tension ',k,1,2,bi,bj,myThid)
                0992         CALL DIAGNOSTICS_FILL(strainBC,'Strain  ',k,1,2,bi,bj,myThid)
f59d76b0dd Ed D*0993 C         stretching will be zero, unless using QG Leith
a125ab7be7 Jean*0994         IF ( viscC2LeithQG.NE.zeroRL ) THEN
f59d76b0dd Ed D*0995           CALL DIAGNOSTICS_FILL(stretching,
                0996      I                          'Stretch ',k,1,2,bi,bj,myThid)
                0997         ENDIF
ed2dbbe83d Jean*0998        ENDIF
07e17fa184 Jean*0999         CALL DIAGNOSTICS_FILL(gU(1-OLx,1-OLy,k,bi,bj),
df1ec59c8a Jean*1000      &                                'Um_Advec',k,1,2,bi,bj,myThid)
07e17fa184 Jean*1001         CALL DIAGNOSTICS_FILL(gV(1-OLx,1-OLy,k,bi,bj),
df1ec59c8a Jean*1002      &                                'Vm_Advec',k,1,2,bi,bj,myThid)
711329b07b Jean*1003       ENDIF
                1004 #endif /* ALLOW_DIAGNOSTICS */
                1005 
99bc607d7a Ed H*1006 #endif /* ALLOW_MOM_VECINV */
cab1a69b8a Jean*1007 
aea29c8517 Alis*1008       RETURN
                1009       END