pkg/mom_vecinv/mom_vi_hdissip.F

cec2469d72 Alis*0001 #include "MOM_VECINV_OPTIONS.h"
aea29c8517 Alis*0002
                0003       SUBROUTINE MOM_VI_HDISSIP(
2f384c201c Jean*0004      I        bi, bj, k,
                0005      I        hDiv, vort3, dStar, zStar, hFacZ,
                0006      I        viscAh_Z, viscAh_D, viscA4_Z, viscA4_D,
                0007      I        harmonic, biharmonic, useVariableViscosity,
                0008      O        uDissip, vDissip,
                0009      I        myThid )
460dc72355 Patr*0010
aea29c8517 Alis*0011       IMPLICIT NONE
2f384c201c Jean*0012
aea29c8517 Alis*0013 C     Calculate horizontal dissipation terms
                0014 C     [del^2 - del^4] (u,v)
                0015
                0016 C     == Global variables ==
                0017 #include "SIZE.h"
                0018 #include "EEPARAMS.h"
                0019 #include "PARAMS.h"
2f384c201c Jean*0020 #include "GRID.h"
aea29c8517 Alis*0021
                0022 C     == Routine arguments ==
2f384c201c Jean*0023       INTEGER bi, bj, k
                0024       _RL hDiv (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
aea29c8517 Alis*0025       _RL vort3(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0026       _RL dStar(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0027       _RL zStar(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
2f384c201c Jean*0028       _RS hFacZ(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
d0d29b39da Jean*0029       _RL viscAh_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0030       _RL viscAh_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0031       _RL viscA4_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0032       _RL viscA4_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
b0c3bd7ab0 Bayl*0033       LOGICAL harmonic, biharmonic, useVariableViscosity
2f384c201c Jean*0034       _RL uDissip(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0035       _RL vDissip(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
998681995e Bayl*0036       INTEGER myThid
e46ec53fc5 Alis*0037
998681995e Bayl*0038 C     == Local variables ==
2f384c201c Jean*0039       INTEGER i, j
                0040       _RL Zip, Zij, Zpj, Dim, Dij, Dmj, uD2, vD2, uD4, vD4
                0041       _RL Zip1, Zij1, Zpj1
aea29c8517 Alis*0042
86fd561884 Jean*0043 C     - Laplacian  terms
b0c3bd7ab0 Bayl*0044       IF (harmonic) THEN
396048594e Alis*0045 C This bit scales the harmonic dissipation operator to be proportional
                0046 C to the grid-cell area over the time-step. viscAh is then non-dimensional
2f384c201c Jean*0047 C and should be less than 1/8, for example viscAh=0.01
90ef0f383b Mart*0048        IF (useVariableViscosity) THEN
2f384c201c Jean*0049         DO j=2-OLy,sNy+OLy-1
                0050          DO i=2-OLx,sNx+OLx-1
90ef0f383b Mart*0051
bbf2faef13 Patr*0052           Dij=hDiv( i , j )*viscAh_D(i,j)
                0053           Dim=hDiv( i ,j-1)*viscAh_D(i,j-1)
                0054           Dmj=hDiv(i-1, j )*viscAh_D(i-1,j)
                0055           Zij=hFacZ( i , j )*vort3( i , j )*viscAh_Z(i,j)
                0056           Zip=hFacZ( i ,j+1)*vort3( i ,j+1)*viscAh_Z(i,j+1)
                0057           Zpj=hFacZ(i+1, j )*vort3(i+1, j )*viscAh_Z(i+1,j)
                0058
396048594e Alis*0059           uD2 = (
                0060      &               cosFacU(j,bi,bj)*( Dij-Dmj )*recip_DXC(i,j,bi,bj)
616600b8d2 Patr*0061      &     -_recip_hFacW(i,j,k,bi,bj)*( Zip-Zij )*recip_DYG(i,j,bi,bj) )
58f454ab40 Mart*0062 #ifdef ISOTROPIC_COS_SCALING
                0063      &                                           *cosFacU(j,bi,bj)
                0064 #endif /* ISOTROPIC_COS_SCALING */
396048594e Alis*0065           vD2 = (
616600b8d2 Patr*0066      &      _recip_hFacS(i,j,k,bi,bj)*( Zpj-Zij )*recip_DXG(i,j,bi,bj)
396048594e Alis*0067      &                                           *cosFacV(j,bi,bj)
                0068      &                               +( Dij-Dim )*recip_DYC(i,j,bi,bj) )
58f454ab40 Mart*0069 #ifdef ISOTROPIC_COS_SCALING
                0070      &                                           *cosFacV(j,bi,bj)
                0071 #endif /* ISOTROPIC_COS_SCALING */
90ef0f383b Mart*0072
                0073           uDissip(i,j) = uD2
                0074           vDissip(i,j) = vD2
                0075
                0076          ENDDO
                0077         ENDDO
                0078        ELSE
2f384c201c Jean*0079         DO j=2-OLy,sNy+OLy-1
                0080          DO i=2-OLx,sNx+OLx-1
90ef0f383b Mart*0081
bbf2faef13 Patr*0082           Dim=hDiv( i ,j-1)
                0083           Dij=hDiv( i , j )
                0084           Dmj=hDiv(i-1, j )
                0085           Zip=hFacZ( i ,j+1)*vort3( i ,j+1)
                0086           Zij=hFacZ( i , j )*vort3( i , j )
                0087           Zpj=hFacZ(i+1, j )*vort3(i+1, j )
                0088
86fd561884 Jean*0089           uD2 = viscAhD*
aea29c8517 Alis*0090      &               cosFacU(j,bi,bj)*( Dij-Dmj )*recip_DXC(i,j,bi,bj)
616600b8d2 Patr*0091      &        - viscAhZ*_recip_hFacW(i,j,k,bi,bj)*
86fd561884 Jean*0092      &                                ( Zip-Zij )*recip_DYG(i,j,bi,bj)
58f454ab40 Mart*0093 #ifdef ISOTROPIC_COS_SCALING
                0094      &                                           *cosFacU(j,bi,bj)
                0095 #endif /* ISOTROPIC_COS_SCALING */
616600b8d2 Patr*0096           vD2 = viscAhZ*_recip_hFacS(i,j,k,bi,bj)*
86fd561884 Jean*0097      &               cosFacV(j,bi,bj)*( Zpj-Zij )*recip_DXG(i,j,bi,bj)
                0098      &        + viscAhD*              ( Dij-Dim )*recip_DYC(i,j,bi,bj)
58f454ab40 Mart*0099 #ifdef ISOTROPIC_COS_SCALING
                0100      &                                           *cosFacV(j,bi,bj)
                0101 #endif /* ISOTROPIC_COS_SCALING */
86fd561884 Jean*0102
90ef0f383b Mart*0103           uDissip(i,j) = uD2
                0104           vDissip(i,j) = vD2
86fd561884 Jean*0105
90ef0f383b Mart*0106          ENDDO
86fd561884 Jean*0107         ENDDO
90ef0f383b Mart*0108        ENDIF
86fd561884 Jean*0109       ELSE
2f384c201c Jean*0110        DO j=2-OLy,sNy+OLy-1
                0111         DO i=2-OLx,sNx+OLx-1
86fd561884 Jean*0112          uDissip(i,j) = 0.
                0113          vDissip(i,j) = 0.
                0114         ENDDO
                0115        ENDDO
                0116       ENDIF
                0117
                0118 C     - Bi-harmonic terms
b0c3bd7ab0 Bayl*0119       IF (biharmonic) THEN
aea29c8517 Alis*0120
90ef0f383b Mart*0121 C This bit scales the harmonic dissipation operator to be proportional
                0122 C to the grid-cell area over the time-step. viscAh is then non-dimensional
2f384c201c Jean*0123 C and should be less than 1/8, for example viscAh=0.01
90ef0f383b Mart*0124        IF (useVariableViscosity) THEN
2f384c201c Jean*0125         DO j=2-OLy,sNy+OLy-1
                0126          DO i=2-OLx,sNx+OLx-1
aea29c8517 Alis*0127
90ef0f383b Mart*0128 #ifdef MOM_VI_ORIGINAL_VISCA4
                0129           Dim=dyF( i ,j-1,bi,bj)*dStar( i ,j-1)
                0130           Dij=dyF( i , j ,bi,bj)*dStar( i , j )
                0131           Dmj=dyF(i-1, j ,bi,bj)*dStar(i-1, j )
2f384c201c Jean*0132
90ef0f383b Mart*0133           Zip1=dxV( i ,j+1,bi,bj)*hFacZ( i ,j+1)*zStar( i ,j+1)
                0134           Zij1=dxV( i , j ,bi,bj)*hFacZ( i , j )*zStar( i , j )
                0135           Zpj1=dxV(i+1, j ,bi,bj)*hFacZ(i+1, j )*zStar(i+1, j )
fd0a408c49 Jean*0136 #else
90ef0f383b Mart*0137           Dim=dStar( i ,j-1)
                0138           Dij=dStar( i , j )
                0139           Dmj=dStar(i-1, j )
fd0a408c49 Jean*0140
90ef0f383b Mart*0141           Zip1=hFacZ( i ,j+1)*zStar( i ,j+1)
                0142           Zij1=hFacZ( i , j )*zStar( i , j )
                0143           Zpj1=hFacZ(i+1, j )*zStar(i+1, j )
fd0a408c49 Jean*0144 #endif
e46ec53fc5 Alis*0145           Dij=Dij*viscA4_D(i,j)
                0146           Dim=Dim*viscA4_D(i,j-1)
                0147           Dmj=Dmj*viscA4_D(i-1,j)
bbf2faef13 Patr*0148           Zij=Zij1*viscA4_Z(i,j)
                0149           Zip=Zip1*viscA4_Z(i,j+1)
                0150           Zpj=Zpj1*viscA4_Z(i+1,j)
fd0a408c49 Jean*0151
                0152 #ifdef MOM_VI_ORIGINAL_VISCA4
396048594e Alis*0153           uD4 = recip_rAw(i,j,bi,bj)*(
                0154      &                             ( (Dij-Dmj)*cosFacU(j,bi,bj) )
2f384c201c Jean*0155      &  -_recip_hFacW(i,j,k,bi,bj)*( Zip-Zij )
f8798ebccf Patr*0156 # ifdef ISOTROPIC_COS_SCALING
58f454ab40 Mart*0157      &                                        *cosFacU(j,bi,bj)
f8798ebccf Patr*0158 # endif /* ISOTROPIC_COS_SCALING */
58f454ab40 Mart*0159      &         )
396048594e Alis*0160           vD4 = recip_rAs(i,j,bi,bj)*(
616600b8d2 Patr*0161      &   _recip_hFacS(i,j,k,bi,bj)*( (Zpj-Zij)*cosFacV(j,bi,bj) )
2f384c201c Jean*0162      &   +                         ( Dij-Dim )
f8798ebccf Patr*0163 # ifdef ISOTROPIC_COS_SCALING
58f454ab40 Mart*0164      &                                        *cosFacV(j,bi,bj)
f8798ebccf Patr*0165 # endif /* ISOTROPIC_COS_SCALING */
58f454ab40 Mart*0166      &         )
90ef0f383b Mart*0167 #else /* MOM_VI_ORIGINAL_VISCA4 */
                0168           uD4 = (
                0169      &               cosFacU(j,bi,bj)*( Dij-Dmj )*recip_DXC(i,j,bi,bj)
                0170      &     -_recip_hFacW(i,j,k,bi,bj)*( Zip-Zij )*recip_DYG(i,j,bi,bj) )
                0171 # ifdef ISOTROPIC_COS_SCALING
                0172      &                                           *cosFacU(j,bi,bj)
                0173 # endif /* ISOTROPIC_COS_SCALING */
                0174           vD4 = (
                0175      &      _recip_hFacS(i,j,k,bi,bj)*( Zpj-Zij )*recip_DXG(i,j,bi,bj)
                0176      &                                           *cosFacV(j,bi,bj)
                0177      &                               +( Dij-Dim )*recip_DYC(i,j,bi,bj) )
                0178 # ifdef ISOTROPIC_COS_SCALING
                0179      &                                           *cosFacV(j,bi,bj)
                0180 # endif /* ISOTROPIC_COS_SCALING */
                0181 #endif  /* MOM_VI_ORIGINAL_VISCA4 */
                0182
                0183           uDissip(i,j) = uDissip(i,j) - uD4
                0184           vDissip(i,j) = vDissip(i,j) - vD4
                0185
                0186          ENDDO
                0187         ENDDO
                0188        ELSE
2f384c201c Jean*0189         DO j=2-OLy,sNy+OLy-1
                0190          DO i=2-OLx,sNx+OLx-1
90ef0f383b Mart*0191
                0192 #ifdef MOM_VI_ORIGINAL_VISCA4
                0193           Dim=dyF( i ,j-1,bi,bj)*dStar( i ,j-1)
                0194           Dij=dyF( i , j ,bi,bj)*dStar( i , j )
                0195           Dmj=dyF(i-1, j ,bi,bj)*dStar(i-1, j )
2f384c201c Jean*0196
90ef0f383b Mart*0197           Zip1=dxV( i ,j+1,bi,bj)*hFacZ( i ,j+1)*zStar( i ,j+1)
                0198           Zij1=dxV( i , j ,bi,bj)*hFacZ( i , j )*zStar( i , j )
                0199           Zpj1=dxV(i+1, j ,bi,bj)*hFacZ(i+1, j )*zStar(i+1, j )
                0200 #else
                0201           Dim=dStar( i ,j-1)
                0202           Dij=dStar( i , j )
                0203           Dmj=dStar(i-1, j )
                0204
                0205           Zip1=hFacZ( i ,j+1)*zStar( i ,j+1)
                0206           Zij1=hFacZ( i , j )*zStar( i , j )
                0207           Zpj1=hFacZ(i+1, j )*zStar(i+1, j )
                0208 #endif
bbf2faef13 Patr*0209           Zij=Zij1
                0210           Zip=Zip1
                0211           Zpj=Zpj1
90ef0f383b Mart*0212
                0213 #ifdef MOM_VI_ORIGINAL_VISCA4
396048594e Alis*0214           uD4 = recip_rAw(i,j,bi,bj)*(
2f384c201c Jean*0215      &                             viscA4D*( Dij-Dmj )*cosFacU(j,bi,bj)
                0216      &  -_recip_hFacW(i,j,k,bi,bj)*viscA4Z*( Zip-Zij )
f8798ebccf Patr*0217 # ifdef ISOTROPIC_COS_SCALING
2f384c201c Jean*0218      &                                                *cosFacU(j,bi,bj)
f8798ebccf Patr*0219 # endif /* ISOTROPIC_COS_SCALING */
2f384c201c Jean*0220      &                               )
396048594e Alis*0221           vD4 = recip_rAs(i,j,bi,bj)*(
2f384c201c Jean*0222      &   _recip_hFacS(i,j,k,bi,bj)*viscA4Z*( Zpj-Zij )*cosFacV(j,bi,bj)
                0223      &   +                         viscA4D*( Dij-Dim )
f8798ebccf Patr*0224 # ifdef ISOTROPIC_COS_SCALING
2f384c201c Jean*0225      &                                                *cosFacV(j,bi,bj)
f8798ebccf Patr*0226 # endif /* ISOTROPIC_COS_SCALING */
2f384c201c Jean*0227      &                               )
fd0a408c49 Jean*0228 #else /* MOM_VI_ORIGINAL_VISCA4 */
d0d29b39da Jean*0229           uD4 = viscA4D*
fd0a408c49 Jean*0230      &               cosFacU(j,bi,bj)*( Dij-Dmj )*recip_DXC(i,j,bi,bj)
616600b8d2 Patr*0231      &        - viscA4Z*_recip_hFacW(i,j,k,bi,bj)*
86fd561884 Jean*0232      &                                ( Zip-Zij )*recip_DYG(i,j,bi,bj)
f8798ebccf Patr*0233 # ifdef ISOTROPIC_COS_SCALING
58f454ab40 Mart*0234      &                                           *cosFacU(j,bi,bj)
f8798ebccf Patr*0235 # endif /* ISOTROPIC_COS_SCALING */
616600b8d2 Patr*0236           vD4 = viscA4Z*_recip_hFacS(i,j,k,bi,bj)*
86fd561884 Jean*0237      &               cosFacV(j,bi,bj)*( Zpj-Zij )*recip_DXG(i,j,bi,bj)
                0238      &        + viscA4D*              ( Dij-Dim )*recip_DYC(i,j,bi,bj)
f8798ebccf Patr*0239 # ifdef ISOTROPIC_COS_SCALING
58f454ab40 Mart*0240      &                                           *cosFacV(j,bi,bj)
f8798ebccf Patr*0241 # endif /* ISOTROPIC_COS_SCALING */
fd0a408c49 Jean*0242 #endif /* MOM_VI_ORIGINAL_VISCA4 */
aea29c8517 Alis*0243
90ef0f383b Mart*0244           uDissip(i,j) = uDissip(i,j) - uD4
                0245           vDissip(i,j) = vDissip(i,j) - vD4
2f384c201c Jean*0246
90ef0f383b Mart*0247          ENDDO
86fd561884 Jean*0248         ENDDO
90ef0f383b Mart*0249        ENDIF
86fd561884 Jean*0250       ENDIF
aea29c8517 Alis*0251
81e837e1e7 Jean*0252       IF ( harmonic .OR. biharmonic ) THEN
                0253        DO j=1-OLy,sNy+OLy-1
                0254         DO i=1-OLx,sNx+OLx-1
                0255          uDissip(i,j) = uDissip(i,j)*maskW(i,j,k,bi,bj)
f1a4cec01a Jean*0256      &                              *recip_deepFacC(k)
81e837e1e7 Jean*0257          vDissip(i,j) = vDissip(i,j)*maskS(i,j,k,bi,bj)
f1a4cec01a Jean*0258      &                              *recip_deepFacC(k)
81e837e1e7 Jean*0259         ENDDO
                0260        ENDDO
                0261       ENDIF
                0262
aea29c8517 Alis*0263       RETURN
                0264       END