Back to home page

MITgcm

 
 

    


File indexing completed on 2020-07-16 05:11:23 UTC

view on githubraw file Latest commit 5fe78992 on 2020-07-16 01:49:33 UTC
925df4f4b9 Mart*0001 #include "SEAICE_OPTIONS.h"
                0002 #ifdef ALLOW_AUTODIFF
                0003 # include "AUTODIFF_OPTIONS.h"
                0004 #endif
                0005 
                0006 C--  File seaice_calc_stressdiv.F
                0007 C--   Contents
                0008 C--   o SEAICE_CALC_STRESSDIV
                0009 C--   o SEAICE_CALC_STRESS
                0010 
                0011 CBOP
                0012 C !ROUTINE: SEAICE_CALC_STRESSDIV
                0013 C !INTERFACE: ==========================================================
cbf501ab81 Jean*0014       SUBROUTINE SEAICE_CALC_STRESSDIV(
925df4f4b9 Mart*0015      I     e11, e22, e12, press, zeta, eta, etaZ,
                0016      O     stressDivergenceX, stressDivergenceY,
                0017      I     bi, bj, myTime, myIter, myThid )
                0018 
                0019 C !DESCRIPTION: \bv
                0020 C     *===========================================================*
                0021 C     | SUBROUTINE SEAICE_CALC_STRESSDIV
                0022 C     | o compute ice internal stress divergence
                0023 C     |
                0024 C     | Martin Losch, May 2017, Martin.Losch@awi.de
                0025 C     *===========================================================*
                0026 C \ev
                0027 
                0028 C !USES: ===============================================================
                0029       IMPLICIT NONE
                0030 
                0031 #include "SIZE.h"
                0032 #include "EEPARAMS.h"
                0033 #include "GRID.h"
                0034 
                0035 C     !INPUT/OUTPUT PARAMETERS:
                0036 C     === Routine arguments ===
                0037 C     e11/e22/e12 :: strain rate tensor components
                0038 C     press  :: maximal compressive strength
                0039 C     zeta   :: bulk viscosity
cbf501ab81 Jean*0040 C     eta    :: shear viscosity
925df4f4b9 Mart*0041 C     etaZ   :: shear viscosity at vorticity points
                0042 C     stressDivergenceX/Y :: x/y-component of stress divergence
cbf501ab81 Jean*0043 C     myTime :: Simulation time
                0044 C     myIter :: Simulation timestep number
                0045 C     myThid :: my Thread Id. number
925df4f4b9 Mart*0046       _RL e11  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0047       _RL e22  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0048       _RL e12  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0049       _RL press(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0050       _RL zeta (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0051       _RL eta  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0052       _RL etaZ (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0053       _RL stressDivergenceX(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0054       _RL stressDivergenceY(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
cbf501ab81 Jean*0055       INTEGER bi, bj
                0056       _RL     myTime
                0057       INTEGER myIter
                0058       INTEGER myThid
925df4f4b9 Mart*0059 CEOP
                0060 
                0061 #ifdef SEAICE_CGRID
                0062 C !LOCAL VARIABLES: ====================================================
                0063 C     === Local variables ===
                0064 C     i,j       :: inner loop counters
cbf501ab81 Jean*0065       _RL sig11(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0066       _RL sig22(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0067       _RL sig12(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
925df4f4b9 Mart*0068       INTEGER i, j
                0069 
cbf501ab81 Jean*0070       CALL SEAICE_CALC_STRESS(
925df4f4b9 Mart*0071      I     e11, e22, e12, press, zeta, eta, etaZ,
                0072      O     sig11, sig22, sig12,
                0073      I     bi, bj, myTime, myIter, myThid )
                0074 
                0075 C     compute divergence of stress tensor
cbf501ab81 Jean*0076 
5fe78992ba Mart*0077       DO j=1,sNy+1
                0078        DO i=1,sNx+1
925df4f4b9 Mart*0079         stressDivergenceX(i,j,bi,bj) =
                0080      &       ( sig11(i  ,j  ) * _dyF(i  ,j,bi,bj)
5fe78992ba Mart*0081      &       - sig11(i-1,j  ) * _dyF(i-1,j,bi,bj)
925df4f4b9 Mart*0082      &       + sig12(i  ,j+1) * _dxV(i,j+1,bi,bj)
                0083      &       - sig12(i  ,j  ) * _dxV(i,j  ,bi,bj)
                0084      &       ) * recip_rAw(i,j,bi,bj)
                0085         stressDivergenceY(i,j,bi,bj) =
                0086      &       ( sig22(i  ,j  ) * _dxF(i,j  ,bi,bj)
                0087      &       - sig22(i  ,j-1) * _dxF(i,j-1,bi,bj)
5fe78992ba Mart*0088      &       + sig12(i+1,j  ) * _dyU(i+1,j,bi,bj)
925df4f4b9 Mart*0089      &       - sig12(i  ,j  ) * _dyU(i  ,j,bi,bj)
                0090      &       ) * recip_rAs(i,j,bi,bj)
                0091        ENDDO
                0092       ENDDO
                0093 
                0094       RETURN
                0095       END
                0096 
                0097 CBOP
                0098 C !ROUTINE: SEAICE_CALC_STRESS
                0099 C !INTERFACE: ==========================================================
cbf501ab81 Jean*0100       SUBROUTINE SEAICE_CALC_STRESS(
925df4f4b9 Mart*0101      I     e11, e22, e12, press, zeta, eta, etaZ,
                0102      O     sig11, sig22, sig12,
                0103      I     bi, bj, myTime, myIter, myThid )
                0104 
                0105 C !DESCRIPTION: \bv
                0106 C     *===========================================================*
                0107 C     | SUBROUTINE SEAICE_CALC_STRESS
                0108 C     | o compute ice internal stress according to rheology
                0109 C     |
                0110 C     | Martin Losch, May 2017, Martin.Losch@awi.de
                0111 C     *===========================================================*
                0112 C \ev
                0113 
                0114 C !USES: ===============================================================
                0115       IMPLICIT NONE
                0116 
                0117 #include "SIZE.h"
                0118 #include "EEPARAMS.h"
                0119 #include "GRID.h"
                0120 
                0121 C     !INPUT/OUTPUT PARAMETERS:
                0122 C     === Routine arguments ===
                0123 C     myTime :: Simulation time
                0124 C     myIter :: Simulation timestep number
                0125 C     myThid :: my Thread Id. number
                0126 C     e11/e22/e12 :: strain rate tensor components
                0127 C     press  :: maximal compressive strength
                0128 C     zeta   :: bulk viscosity
cbf501ab81 Jean*0129 C     eta    :: shear viscosity
925df4f4b9 Mart*0130 C     etaZ   :: shear viscosity at vorticity points
                0131 C     sig11/sig22/sig12   :: stress tensor components
                0132       _RL     myTime
                0133       INTEGER myIter
                0134       INTEGER myThid
                0135       INTEGER bi, bj
                0136       _RL e11  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0137       _RL e22  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0138       _RL e12  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0139       _RL press(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0140       _RL zeta (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0141       _RL eta  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0142       _RL etaZ (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0143       _RL sig11(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0144       _RL sig22(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0145       _RL sig12(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0146 CEOP
                0147 
                0148 C !LOCAL VARIABLES: ====================================================
                0149 C     === Local variables ===
                0150 C     i,j           :: inner loop counters
                0151 C     eplus, eminus :: convenient abbreviations for e11+e22, e11-e22
                0152       INTEGER i, j
                0153       _RL     eplus, eminus
                0154 
                0155 C     compute components of stress tensor from current strain rate tensor
                0156       DO j=1-OLy,sNy+OLy
                0157        DO i=1-OLx,sNx+OLx
                0158         sig11(i,j) = 0. _d 0
                0159         sig22(i,j) = 0. _d 0
                0160         sig12(i,j) = 0. _d 0
                0161        ENDDO
                0162       ENDDO
                0163 
5fe78992ba Mart*0164       DO j=0,sNy+1
                0165        DO i=0,sNx+1
925df4f4b9 Mart*0166         eplus = e11(i,j,bi,bj) + e22(i,j,bi,bj)
                0167         eminus= e11(i,j,bi,bj) - e22(i,j,bi,bj)
                0168         sig11(i,j) = zeta(i,j,bi,bj)*eplus + eta(i,j,bi,bj)*eminus
                0169      &       - 0.5 _d 0 * press(i,j,bi,bj)
                0170         sig22(i,j) = zeta(i,j,bi,bj)*eplus - eta(i,j,bi,bj)*eminus
                0171      &       - 0.5 _d 0 * press(i,j,bi,bj)
                0172        ENDDO
                0173       ENDDO
cbf501ab81 Jean*0174 
5fe78992ba Mart*0175       DO j=1,sNy+2
                0176        DO i=1,sNx+2
925df4f4b9 Mart*0177         sig12(i,j) = 2. _d 0 * e12(i,j,bi,bj) * etaZ(i,j,bi,bj)
                0178        ENDDO
                0179       ENDDO
                0180 #endif /* SEAICE_CGRID */
                0181       RETURN
                0182       END