Back to home page

MITgcm

 
 

    


File indexing completed on 2021-02-17 06:11:44 UTC

view on githubraw file Latest commit 2d5bb917 on 2020-09-13 20:55:43 UTC
af960ebfb4 Jean*0001 #include "MOM_COMMON_OPTIONS.h"
                0002 
                0003 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0004 CBOP
                0005 C     !ROUTINE: MOM_CALC_SMAG_3D
                0006 
                0007 C     !INTERFACE:
                0008       SUBROUTINE MOM_CALC_SMAG_3D(
                0009      I        str11, str22, str33, str12, str13, str23,
                0010      O        viscAh3d_00, viscAh3d_12,
                0011      O        viscAh3d_13, viscAh3d_23,
                0012      I        smag3D_hLsC, smag3D_hLsW, smag3D_hLsS, smag3D_hLsZ,
                0013      I        k, bi, bj, myThid )
                0014 
                0015 C     !DESCRIPTION:
                0016 C     Calculate Smagorinsky 3-D (harmonic) viscosities
                0017 C      at current level k (for viscAh3d_00 & viscAh3d_12)
                0018 C      and at level k+1   (for viscAh3d_13 & viscAh3d_23)
                0019 
                0020 C     !USES:
                0021       IMPLICIT NONE
                0022 
                0023 C     == Global variables ==
                0024 #include "SIZE.h"
                0025 #include "EEPARAMS.h"
                0026 #include "PARAMS.h"
                0027 #include "GRID.h"
2d5bb917cc Jean*0028 #ifdef ALLOW_SMAG_3D_DIFFUSIVITY
                0029 # include "DYNVARS.h"
                0030 #endif
af960ebfb4 Jean*0031 c#include "MOM_VISC.h"
                0032 
                0033 C     !INPUT PARAMETERS:
                0034 C     str11       :: strain component Vxx @ grid-cell center
                0035 C     str22       :: strain component Vyy @ grid-cell center
                0036 C     str33       :: strain component Vzz @ grid-cell center
                0037 C     str12       :: strain component Vxy @ grid-cell corner
                0038 C     str13       :: strain component Vxz @ above uVel
                0039 C     str23       :: strain component Vyz @ above vVel
                0040 C     smag3D_hLsC :: horiz. grid length scale (power 2/3) at grid cell center
                0041 C     smag3D_hLsW :: horiz. grid length scale (power 2/3) at western  edge
                0042 C     smag3D_hLsS :: horiz. grid length scale (power 2/3) at southern egde
                0043 C     smag3D_hLsZ :: horiz. grid length scale (power 2/3) at grid cell corner
                0044 C     k           :: current level index
                0045 C     bi, bj      :: tile indices
                0046 C     myThid      :: my Thread Id number
                0047       _RL str11(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0048       _RL str22(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0049       _RL str33(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0050       _RL str12(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0051       _RL str13(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1)
                0052       _RL str23(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1)
                0053       _RS smag3D_hLsC(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0054       _RS smag3D_hLsW(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0055       _RS smag3D_hLsS(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0056       _RS smag3D_hLsZ(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0057       INTEGER k, bi, bj
                0058       INTEGER myThid
                0059 
                0060 C     !OUTPUT PARAMETERS:
                0061 C     viscAh3d_00 :: Smagorinsky viscosity @ grid-cell center
                0062 C     viscAh3d_12 :: Smagorinsky viscosity @ grid-cell corner
                0063 C     viscAh3d_13 :: Smagorinsky viscosity @ above uVel
                0064 C     viscAh3d_23 :: Smagorinsky viscosity @ above vVel
                0065       _RL viscAh3d_00(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0066       _RL viscAh3d_12(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0067       _RL viscAh3d_13(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1)
                0068       _RL viscAh3d_23(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1)
                0069 CEOP
                0070 
245cf5c24c Jean*0071 #ifdef ALLOW_SMAG_3D
af960ebfb4 Jean*0072 C     !LOCAL VARIABLES:
6a6adb74e6 Jean*0073 C     S66   :: Sum of squared of the 3 strain component @ grid-cell center
                0074 C     S12   :: Squared of strain component Vxy @ grid-cell corner
                0075 C     S13   :: Squared of strain component Vxz @ above uVel
                0076 C     S23   :: Squared of strain component Vyz @ above vVel
af960ebfb4 Jean*0077       INTEGER i, j
6a6adb74e6 Jean*0078       INTEGER kl, n
2d5bb917cc Jean*0079       _RL twoThird, tmpFac, locVisc
6a6adb74e6 Jean*0080       _RL S66(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
                0081       _RL S12(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
                0082       _RL S13(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
                0083       _RL S23(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
af960ebfb4 Jean*0084 
                0085 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0086 
                0087       twoThird = 2. _d 0 / 3. _d 0
                0088 
6a6adb74e6 Jean*0089       DO n=1,2
                0090         kl = k + n - 1
                0091         IF ( kl.LE.Nr ) THEN
af960ebfb4 Jean*0092           DO j=1-OLy,sNy+OLy
6a6adb74e6 Jean*0093            DO i=1-OLx,sNx+OLx
                0094              S66(i,j,n) = str11(i,j,kl)**2
                0095      &                  + str22(i,j,kl)**2
                0096      &                  + str33(i,j,kl)**2
                0097              S12(i,j,n) = str12(i,j,kl)**2
                0098            ENDDO
af960ebfb4 Jean*0099           ENDDO
                0100         ELSE
                0101           DO j=1-OLy,sNy+OLy
6a6adb74e6 Jean*0102            DO i=1-OLx,sNx+OLx
                0103              S66(i,j,n) = 0. _d 0
                0104              S12(i,j,n) = 0. _d 0
                0105            ENDDO
af960ebfb4 Jean*0106           ENDDO
                0107         ENDIF
6a6adb74e6 Jean*0108           DO j=1-OLy,sNy+OLy
                0109            DO i=1-OLx,sNx+OLx
                0110              S13(i,j,n) = str13(i,j,kl)**2
                0111              S23(i,j,n) = str23(i,j,kl)**2
                0112            ENDDO
                0113           ENDDO
af960ebfb4 Jean*0114       ENDDO
                0115 
                0116 C--  ------------------------------------------------------------------
                0117 C--  calculate current level Smag viscosity coeff
                0118 C--  ------------------------------------------------------------------
                0119 
6a6adb74e6 Jean*0120 C     Current level k --> n=1
                0121       kl = k
                0122       n = 1
2d5bb917cc Jean*0123       tmpFac = twoRL*SQRT(twoRL)*drF(kl)**twoThird
af960ebfb4 Jean*0124 
                0125 C     viscAh3d_00 = sqrt( S11+S22+S33+2*(S12+S13+S23) ) @ grid-cell center
                0126 
                0127       DO j=1-OLy,sNy+OLy-1
                0128        DO i=1-OLx,sNx+OLx-1
2d5bb917cc Jean*0129          locVisc =
af960ebfb4 Jean*0130      &    smag3D_hLsC(i,j,bi,bj)*tmpFac*SQRT(
6a6adb74e6 Jean*0131      &            S66( i , j , n )
af960ebfb4 Jean*0132      &   +  0.5*( S12( i ,j+1, n )+S12(i+1,j+1, n )
                0133      &           +S12( i , j , n )+S12(i+1, j , n ) )
                0134      &   +  0.5*( S13( i , j , n )+S13(i+1, j , n )
                0135      &           +S13( i , j ,n+1)+S13(i+1, j ,n+1) )
                0136      &   +  0.5*( S23( i , j , n )+S23( i ,j+1, n )
                0137      &           +S23( i , j ,n+1)+S23( i ,j+1,n+1) )
                0138      &                                      )
2d5bb917cc Jean*0139          viscAh3d_00(i,j,kl) = locVisc*smag3D_coeff
                0140 #ifdef ALLOW_SMAG_3D_DIFFUSIVITY
                0141          smag3D_diffK(i,j,kl,bi,bj) = locVisc*smag3D_diffCoeff
                0142 #endif
af960ebfb4 Jean*0143        ENDDO
                0144       ENDDO
                0145 
                0146 C     viscAh3d_12 = sqrt( S11+S22+S33+2*(S12+S13+S23) ) @ grid-cell corner
                0147 
2d5bb917cc Jean*0148       tmpFac = smag3D_coeff*tmpFac
af960ebfb4 Jean*0149       DO j=2-OLy,sNy+OLy
                0150        DO i=2-OLx,sNx+OLx
6a6adb74e6 Jean*0151          viscAh3d_12(i,j,kl) =
af960ebfb4 Jean*0152      &    smag3D_hLsZ(i,j,bi,bj)*tmpFac*SQRT(
6a6adb74e6 Jean*0153      &     0.25*( S66(i-1, j , n )+S66( i , j , n )
                0154      &           +S66(i-1,j-1, n )+S66( i ,j-1, n ) )
af960ebfb4 Jean*0155      &   + 2.0 *  S12( i , j , n )
                0156      &   + 0.5 *( S13( i ,j-1, n )+S13( i , j , n )
                0157      &           +S13( i ,j-1,n+1)+S13( i , j ,n+1) )
                0158      &   + 0.5 *( S23(i-1, j , n )+S23( i , j , n )
                0159      &           +S23(i-1, j ,n+1)+S23( i , j ,n+1) )
                0160      &                                      )
                0161        ENDDO
                0162       ENDDO
                0163 
                0164 C--  ------------------------------------------------------------------
                0165 C--  calculate  next level (k+1) viscosity coeff (uz,vz)
                0166 C--  ------------------------------------------------------------------
                0167 
                0168       IF ( k.EQ.1 ) THEN
                0169         DO j=1-OLy,sNy+OLy
                0170          DO i=1-OLx,sNx+OLx
                0171            viscAh3d_13(i,j,k) = 0. _d 0
                0172            viscAh3d_23(i,j,k) = 0. _d 0
                0173          ENDDO
                0174         ENDDO
                0175       ENDIF
                0176 
6a6adb74e6 Jean*0177 C     Next level k+1 --> n=2
                0178       kl = k+1
                0179       n = 2
                0180       tmpFac = smag3D_coeff*twoRL*SQRT(twoRL)*drC(kl)**twoThird
af960ebfb4 Jean*0181 
                0182 C     viscAh3d_13 = sqrt( S11+S22+S33+2*(S12+S13+S23) ) @ above uVel
                0183 
                0184       DO j=1-OLy,sNy+OLy-1
                0185        DO i=2-OLx,sNx+OLx
6a6adb74e6 Jean*0186          viscAh3d_13(i,j,kl) =
af960ebfb4 Jean*0187      &    smag3D_hLsW(i,j,bi,bj)*tmpFac*SQRT(
6a6adb74e6 Jean*0188      &     0.25*( S66(i-1, j ,n-1)+S66( i , j ,n-1)
                0189      &           +S66(i-1, j , n )+S66( i , j , n ) )
af960ebfb4 Jean*0190      &   + 0.5 *( S12( i , j ,n-1)+S12( i ,j+1,n-1)
                0191      &           +S12( i , j , n )+S12( i ,j+1, n ) )
                0192      &   + 2.0 *  S13( i , j , n )
                0193      &   + 0.5 *( S23(i-1,j+1, n )+S23( i ,j+1, n )
                0194      &           +S23(i-1, j , n )+S23( i , j , n ) )
                0195      &                                      )
                0196        ENDDO
                0197       ENDDO
                0198 
                0199 C     viscAh3d_23 = sqrt( S11+S22+S33+2*(S12+S13+S23) ) @ above vVel
                0200 
                0201       DO j=2-OLy,sNy+OLy
                0202        DO i=1-OLx,sNx+OLx-1
6a6adb74e6 Jean*0203          viscAh3d_23(i,j,kl) =
af960ebfb4 Jean*0204      &    smag3D_hLsS(i,j,bi,bj)*tmpFac*SQRT(
6a6adb74e6 Jean*0205      &     0.25*( S66( i ,j-1,n-1)+S66( i , j ,n-1)
                0206      &           +S66( i ,j-1, n )+S66( i , j , n ) )
af960ebfb4 Jean*0207      &   + 0.5 *( S12( i , j ,n-1)+S12(i+1, j ,n-1)
                0208      &           +S12( i , j , n )+S12(i+1, j , n ) )
                0209      &   + 0.5 *( S13( i , j , n )+S13(i+1, j , n )
                0210      &           +S13( i ,j-1, n )+S13(i+1,j-1, n ) )
                0211      &   + 2.0 *  S23( i , j , n )
                0212      &                                      )
                0213        ENDDO
                0214       ENDDO
                0215 
                0216 #ifdef ALLOW_DIAGNOSTICS
6a6adb74e6 Jean*0217 c     IF ( useDiagnostics.AND. k.EQ.Nr ) THEN
                0218 c       CALL DIAGNOSTICS_FILL( viscAh3d_00, 'Smag3D_C',
                0219 c    &                         0, Nr, 2, bi, bj, myThid )
af960ebfb4 Jean*0220 c     ENDIF
                0221 #endif
                0222 
245cf5c24c Jean*0223 #endif /* ALLOW_SMAG_3D */
af960ebfb4 Jean*0224       RETURN
                0225       END