Back to home page

MITgcm

 
 

    


File indexing completed on 2018-03-02 18:42:17 UTC

view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
cec2469d72 Alis*0001 #include "MOM_VECINV_OPTIONS.h"
aea29c8517 Alis*0002 
355718db32 Jean*0003 CBOP
                0004 C     !ROUTINE: MOM_VI_U_CORIOLIS
                0005 C     !INTERFACE:
5d7e0a8948 Jean*0006       SUBROUTINE MOM_VI_U_CORIOLIS(
3370e71df9 Mart*0007      I                     bi, bj, k, 
                0008      I                     selectVortScheme, useJamartMomAdv,
355718db32 Jean*0009      I                     vFld, omega3, hFacZ, r_hFacZ,
                0010      O                     uCoriolisTerm,
                0011      I                     myThid )
                0012 C     !DESCRIPTION: \bv
cab1a69b8a Jean*0013 C     *==========================================================*
                0014 C     | S/R MOM_VI_U_CORIOLIS
355718db32 Jean*0015 C     |==========================================================*
                0016 C     | o Calculate flux (in Y-dir.) of vorticity at U point
                0017 C     |   using 2nd order interpolation
cab1a69b8a Jean*0018 C     *==========================================================*
355718db32 Jean*0019 C     \ev
                0020 
                0021 C     !USES:
                0022       IMPLICIT NONE
aea29c8517 Alis*0023 
                0024 C     == Global variables ==
                0025 #include "SIZE.h"
                0026 #include "EEPARAMS.h"
                0027 #include "GRID.h"
                0028 
355718db32 Jean*0029 C     !INPUT/OUTPUT PARAMETERS:
aea29c8517 Alis*0030 C     == Routine arguments ==
355718db32 Jean*0031       INTEGER bi, bj, k
                0032       _RL     vFld     (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0033       _RL     omega3   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0034       _RS     hFacZ    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0035       _RS     r_hFacZ  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
aea29c8517 Alis*0036       _RL uCoriolisTerm(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
3370e71df9 Mart*0037       INTEGER selectVortScheme
                0038       LOGICAL useJamartMomAdv
aea29c8517 Alis*0039       INTEGER myThid
355718db32 Jean*0040 CEOP
aea29c8517 Alis*0041 
                0042 C     == Local variables ==
355718db32 Jean*0043 C     msgBuf :: Informational/error meesage buffer
                0044       CHARACTER*(MAX_LEN_MBUF) msgBuf
                0045       LOGICAL upwindVort3
                0046       INTEGER i, j
                0047       _RL     vBarXY, vBarXm, vBarXp
                0048       _RL     vort3u
4e5f84a272 Jean*0049       _RL     vort3mj, vort3ij, vort3mp, vort3ip
                0050       _RL     oneThird, tmpFac
355718db32 Jean*0051       _RS     epsil
                0052       PARAMETER( upwindVort3 = .FALSE. )
aea29c8517 Alis*0053 
5d7e0a8948 Jean*0054       epsil = 1. _d -9
4e5f84a272 Jean*0055       tmpFac = 1. _d 0
                0056 c     oneThird = 1. _d 0 / ( 1. _d 0 + 2.*tmpFac )
                0057       oneThird = 1. _d 0 / 3. _d 0
aea29c8517 Alis*0058 
355718db32 Jean*0059       IF ( selectVortScheme.EQ.0 ) THEN
                0060 C--   using enstrophy conserving scheme (Shallow-Water Eq.) by Sadourny, JAS 75
                0061 
                0062        DO j=1-Oly,sNy+Oly-1
                0063         DO i=2-Olx,sNx+Olx
aea29c8517 Alis*0064          vBarXY=0.25*(
616600b8d2 Patr*0065      &      (vFld( i , j )*dxG( i , j ,bi,bj)*_hFacS( i , j ,k,bi,bj)
                0066      &      +vFld(i-1, j )*dxG(i-1, j ,bi,bj)*_hFacS(i-1, j ,k,bi,bj))
                0067      &     +(vFld( i ,j+1)*dxG( i ,j+1,bi,bj)*_hFacS( i ,j+1,k,bi,bj)
355718db32 Jean*0068      &      +vFld(i-1,j+1)*dxG(i-1,j+1,bi,bj)*_hFacS(i-1,j+1,k,bi,bj))
                0069      &               )
                0070          IF (upwindVort3) THEN
aea29c8517 Alis*0071           IF (vBarXY.GT.0.) THEN
355718db32 Jean*0072            vort3u=omega3(i,j)*r_hFacZ(i,j)
aea29c8517 Alis*0073           ELSE
355718db32 Jean*0074            vort3u=omega3(i,j+1)*r_hFacZ(i,j+1)
aea29c8517 Alis*0075           ENDIF
                0076          ELSE
cab1a69b8a Jean*0077            vort3u=0.5*(omega3(i,j)*r_hFacZ(i,j)
                0078      &                +omega3(i,j+1)*r_hFacZ(i,j+1))
aea29c8517 Alis*0079          ENDIF
355718db32 Jean*0080          uCoriolisTerm(i,j)= +vort3u*vBarXY*recip_dxC(i,j,bi,bj)
                0081      &                              * _maskW(i,j,k,bi,bj)
                0082         ENDDO
                0083        ENDDO
                0084 
                0085       ELSEIF ( selectVortScheme.EQ.1 ) THEN
                0086 C--   same as above, with different formulation (relatively to hFac)
                0087 
4e5f84a272 Jean*0088        DO j=1-Oly,sNy+Oly-1
355718db32 Jean*0089         DO i=2-Olx,sNx+Olx
                0090          vBarXY= 0.5*(
                0091      &      (vFld( i , j )*dxG( i , j ,bi,bj)*hFacZ(i, j )
                0092      &      +vFld(i-1, j )*dxG(i-1, j ,bi,bj)*hFacZ(i, j ))
5d7e0a8948 Jean*0093      &     +(vFld( i ,j+1)*dxG( i ,j+1,bi,bj)*hFacZ(i,j+1)
                0094      &      +vFld(i-1,j+1)*dxG(i-1,j+1,bi,bj)*hFacZ(i,j+1))
355718db32 Jean*0095      &               )/MAX( epsil, hFacZ(i,j)+hFacZ(i,j+1) )
                0096          IF (upwindVort3) THEN
cab1a69b8a Jean*0097           IF (vBarXY.GT.0.) THEN
                0098            vort3u=omega3(i,j)
                0099           ELSE
                0100            vort3u=omega3(i,j+1)
                0101           ENDIF
                0102          ELSE
5d7e0a8948 Jean*0103            vort3u=0.5*(omega3(i,j)+omega3(i,j+1))
cab1a69b8a Jean*0104          ENDIF
355718db32 Jean*0105          uCoriolisTerm(i,j)= +vort3u*vBarXY*recip_dxC(i,j,bi,bj)
                0106      &                              * _maskW(i,j,k,bi,bj)
                0107         ENDDO
                0108        ENDDO
                0109 
                0110       ELSEIF ( selectVortScheme.EQ.2 ) THEN
                0111 C--   using energy conserving scheme (used by Sadourny in JAS 75 paper)
                0112 
                0113        DO j=1-Oly,sNy+Oly-1
                0114         DO i=2-Olx,sNx+Olx
                0115          vBarXm=0.5*(
                0116      &       vFld( i , j )*dxG( i , j ,bi,bj)*_hFacS( i , j ,k,bi,bj)
                0117      &      +vFld(i-1, j )*dxG(i-1, j ,bi,bj)*_hFacS(i-1, j ,k,bi,bj) )
                0118          vBarXp=0.5*(
                0119      &       vFld( i ,j+1)*dxG( i ,j+1,bi,bj)*_hFacS( i ,j+1,k,bi,bj)
                0120      &      +vFld(i-1,j+1)*dxG(i-1,j+1,bi,bj)*_hFacS(i-1,j+1,k,bi,bj) )
                0121          IF (upwindVort3) THEN
                0122           IF ( (vBarXm+vBarXp) .GT.0.) THEN
                0123            vort3u=vBarXm*r_hFacZ(i, j )*omega3(i, j )
                0124           ELSE
                0125            vort3u=vBarXp*r_hFacZ(i,j+1)*omega3(i,j+1)
                0126           ENDIF
                0127          ELSE
                0128            vort3u = ( vBarXm*r_hFacZ(i, j )*omega3(i, j )
                0129      &               +vBarXp*r_hFacZ(i,j+1)*omega3(i,j+1)
                0130      &              )*0.5 _d 0
                0131          ENDIF
                0132          uCoriolisTerm(i,j)= +vort3u*recip_dxC(i,j,bi,bj)
                0133      &                              * _maskW(i,j,k,bi,bj)
                0134         ENDDO
                0135        ENDDO
                0136 
4e5f84a272 Jean*0137       ELSEIF ( selectVortScheme.EQ.3 ) THEN
                0138 C--   using energy & enstrophy conserving scheme
                0139 C     (from Sadourny, described by Burridge & Haseler, ECMWF Rep.4, 1977)
                0140 
                0141 C     domain where uCoriolisTerm is valid :
                0142 C     [ 3-Olx : sNx+Olx-1 ] x [ 2-Oly : sNy+Oly-1 ]
                0143 C     (=> might need overlap of 3 if using CD-scheme)
                0144        DO j=1-Oly,sNy+Oly-1
                0145         DO i=2-Olx,sNx+Olx-1
                0146          vort3mj= ( r_hFacZ(i, j )*omega3(i, j )
                0147      &            +(r_hFacZ(i,j+1)*omega3(i,j+1)
                0148      &             +r_hFacZ(i-1,j)*omega3(i-1,j)
                0149      &            ))*oneThird
                0150 c    &            )*tmpFac)*oneThird
                0151      &      *vFld(i-1, j )*dxG(i-1, j ,bi,bj)*_hFacS(i-1, j ,k,bi,bj)
                0152          vort3ij= ( r_hFacZ(i, j )*omega3(i, j )
                0153      &            +(r_hFacZ(i,j+1)*omega3(i,j+1)
                0154      &             +r_hFacZ(i+1,j)*omega3(i+1,j)
                0155      &            ))*oneThird
                0156 c    &            )*tmpFac)*oneThird
                0157      &      *vFld( i , j )*dxG( i , j ,bi,bj)*_hFacS( i , j ,k,bi,bj)
                0158          vort3mp= ( r_hFacZ(i,j+1)*omega3(i,j+1)
                0159      &            +(r_hFacZ(i, j )*omega3(i, j )
                0160      &             +r_hFacZ(i-1,j+1)*omega3(i-1,j+1)
                0161      &            ))*oneThird
                0162 c    &            )*tmpFac)*oneThird
                0163      &      *vFld(i-1,j+1)*dxG(i-1,j+1,bi,bj)*_hFacS(i-1,j+1,k,bi,bj)
                0164          vort3ip= ( r_hFacZ(i,j+1)*omega3(i,j+1)
                0165      &            +(r_hFacZ(i, j )*omega3(i, j )
                0166      &             +r_hFacZ(i+1,j+1)*omega3(i+1,j+1)
                0167      &            ))*oneThird
                0168 c    &            )*tmpFac)*oneThird
                0169      &      *vFld( i ,j+1)*dxG( i ,j+1,bi,bj)*_hFacS( i ,j+1,k,bi,bj)
                0170 C---
                0171          uCoriolisTerm(i,j)= +( (vort3mj+vort3ij)+(vort3mp+vort3ip) )
                0172      &                     *0.25 _d 0 *recip_dxC(i,j,bi,bj)
                0173      &                                * _maskW(i,j,k,bi,bj)
                0174         ENDDO
                0175        ENDDO
                0176 
355718db32 Jean*0177       ELSE
                0178         WRITE(msgBuf,'(A,I5,A)')
                0179      &   'MOM_VI_U_CORIOLIS: selectVortScheme=', selectVortScheme,
                0180      &   ' not implemented'
                0181         CALL PRINT_ERROR( msgBuf, myThid )
                0182         STOP 'ABNORMAL END: S/R MOM_VI_U_CORIOLIS'
                0183 
                0184       ENDIF
                0185 
                0186       IF ( useJamartMomAdv ) THEN
                0187        DO j=1-Oly,sNy+Oly-1
4e5f84a272 Jean*0188         DO i=2-Olx,sNx+Olx-1
355718db32 Jean*0189          uCoriolisTerm(i,j) = uCoriolisTerm(i,j)
                0190      &           * 4. _d 0 * _hFacW(i,j,k,bi,bj)
                0191      &           / MAX( epsil,
                0192      &                  (_hFacS(i, j ,k,bi,bj)+_hFacS(i-1, j ,k,bi,bj))
                0193      &                 +(_hFacS(i,j+1,k,bi,bj)+_hFacS(i-1,j+1,k,bi,bj))
                0194      &                )
                0195         ENDDO
aea29c8517 Alis*0196        ENDDO
355718db32 Jean*0197       ENDIF
aea29c8517 Alis*0198 
                0199       RETURN
                0200       END