Back to home page

MITgcm

 
 

    


File indexing completed on 2018-03-02 18:43:53 UTC

view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
aea29c8517 Alis*0001 #include "SHAP_FILT_OPTIONS.h"
                0002 
15c70d7cd1 Jean*0003 CBOP
                0004 C     !ROUTINE: SHAP_FILT_UV_S1
                0005 C     !INTERFACE:
                0006       SUBROUTINE SHAP_FILT_UV_S1(
                0007      U           uFld, vFld, tmpFldU, tmpFldV,
                0008      I           kSize, myTime, myThid )
                0009 C     !DESCRIPTION: \bv
                0010 C     *==========================================================*
                0011 C     | S/R SHAP_FILT_UV_S1
                0012 C     | o Applies Shapiro filter to velocity field (u & v).
                0013 C     | o use filtering function "S1" = [1 - d_xx^n - d_yy^n]
                0014 C     |     with no grid spacing (computational Filter) ;
                0015 C     |     include No-Slip option.
                0016 C     *==========================================================*
                0017 C     \ev
                0018 
                0019 C     !USES:
aea29c8517 Alis*0020       IMPLICIT NONE
                0021 
                0022 C     == Global variables ===
                0023 #include "SIZE.h"
                0024 #include "EEPARAMS.h"
                0025 #include "PARAMS.h"
                0026 #include "GRID.h"
                0027 #include "SHAP_FILT.h"
                0028 
15c70d7cd1 Jean*0029 C     !INPUT/OUTPUT PARAMETERS:
aea29c8517 Alis*0030 C     == Routine arguments
15c70d7cd1 Jean*0031 C     uFld :: velocity field (U component) on which filter applies
                0032 C     vFld :: velocity field (V component) on which filter applies
                0033 C     tmpFldU :: working temporary array
                0034 C     tmpFldV :: working temporary array
                0035 C     kSize :: length of 3rd Dim : either =1 (2D field) or =Nr (3D field)
                0036 C     myTime :: Current time in simulation
                0037 C     myThid :: Thread number for this instance of SHAP_FILT_UV_S1
                0038       INTEGER kSize
                0039       _RL uFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy,kSize,nSx,nSy)
                0040       _RL vFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy,kSize,nSx,nSy)
                0041       _RL tmpFldU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,kSize,nSx,nSy)
                0042       _RL tmpFldV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,kSize,nSx,nSy)
aea29c8517 Alis*0043       _RL     myTime
                0044       INTEGER myThid
                0045 
                0046 #ifdef ALLOW_SHAP_FILT
                0047 
15c70d7cd1 Jean*0048 C     !LOCAL VARIABLES: 
aea29c8517 Alis*0049 C     == Local variables ==
15c70d7cd1 Jean*0050       INTEGER bi,bj,k,i,j,N
aea29c8517 Alis*0051       _RL tmpGrdU(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0052       _RL tmpGrdV(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0053       _RL maskZj,maskZp
15c70d7cd1 Jean*0054       _RL noSlipFact
aea29c8517 Alis*0055       _RL tmpScal
15c70d7cd1 Jean*0056 CEOP
                0057 
                0058       noSlipFact = Shap_noSlip*2. _d 0
aea29c8517 Alis*0059 
                0060       IF (nShapUV.gt.0) THEN
                0061 
                0062         DO bj=myByLo(myThid),myByHi(myThid)
                0063          DO bi=myBxLo(myThid),myBxHi(myThid)
15c70d7cd1 Jean*0064           DO K=1,kSize
                0065            DO J=1-OLy,sNy+OLy
                0066             DO I=1-OLx,sNx+OLx
aea29c8517 Alis*0067              tmpFldU(i,j,k,bi,bj)=uFld(i,j,k,bi,bj)
                0068      &                *_maskW(i,j,k,bi,bj)
                0069              tmpFldV(i,j,k,bi,bj)=vFld(i,j,k,bi,bj)
                0070      &                *_maskS(i,j,k,bi,bj)
                0071             ENDDO
                0072            ENDDO
                0073           ENDDO
                0074          ENDDO
                0075         ENDDO
                0076 
                0077 
                0078 C      d_xx^n tmpFld 
                0079 
                0080        DO N=1,nShapUV
                0081 
15c70d7cd1 Jean*0082         IF (kSize.EQ.Nr) THEN
                0083           CALL EXCH_UV_XYZ_RL(tmpFldU,tmpFldV,.TRUE.,myThid)
                0084         ELSE
                0085           CALL EXCH_UV_XY_RL(tmpFldU,tmpFldV,.TRUE.,myThid)
                0086         ENDIF
aea29c8517 Alis*0087 
                0088         DO bj=myByLo(myThid),myByHi(myThid)
                0089          DO bi=myBxLo(myThid),myBxHi(myThid)
15c70d7cd1 Jean*0090           DO K=1,kSize
aea29c8517 Alis*0091 
                0092 C          Uxx
                0093            DO J=1,sNy
                0094             DO I=1,sNx+1
                0095              tmpGrdU(i,j) = -0.25*(
                0096      &          tmpFldU(i-1,j,k,bi,bj) + tmpFldU(i+1,j,k,bi,bj)
                0097      &             - 2.*tmpFldU(i,j,k,bi,bj)
                0098      &            )*_maskW(i,j,k,bi,bj)
                0099             ENDDO
                0100            ENDDO
                0101 
                0102            DO J=1,sNy
                0103             DO I=1,sNx+1
                0104              tmpFldU(i,j,k,bi,bj) = tmpGrdU(i,j)
                0105             ENDDO
                0106            ENDDO
                0107 
                0108 C          Vyy
                0109            DO J=1,sNy+1
                0110             DO I=1,sNx
                0111              tmpGrdV(i,j) = -0.25*(
                0112      &          tmpFldV(i,j-1,k,bi,bj) + tmpFldV(i,j+1,k,bi,bj)
                0113      &             - 2.*tmpFldV(i,j,k,bi,bj)
                0114      &            )*_maskS(i,j,k,bi,bj)
                0115             ENDDO
                0116            ENDDO
                0117 
                0118            DO J=1,sNy+1
                0119             DO I=1,sNx
                0120              tmpFldV(i,j,k,bi,bj) = tmpGrdV(i,j)
                0121             ENDDO
                0122            ENDDO
                0123 
                0124           ENDDO
                0125          ENDDO
                0126         ENDDO
                0127 
                0128        ENDDO
                0129 
15c70d7cd1 Jean*0130 C      F <-  [1 - d_xx^n *deltaT/tau].F
aea29c8517 Alis*0131        DO bj=myByLo(myThid),myByHi(myThid)
                0132         DO bi=myBxLo(myThid),myBxHi(myThid)
15c70d7cd1 Jean*0133          DO K=1,kSize
aea29c8517 Alis*0134           DO J=1,sNy
                0135            DO I=1,sNx+1
15c70d7cd1 Jean*0136             tmpScal = uFld(i,j,k,bi,bj)
                0137             uFld(i,j,k,bi,bj) = uFld(i,j,k,bi,bj)
                0138      &         -0.5*tmpFldU(i,j,k,bi,bj)*deltaTmom/Shap_uvtau
                0139             tmpFldU(i,j,k,bi,bj) = tmpScal*_maskW(i,j,k,bi,bj)
aea29c8517 Alis*0140            ENDDO
                0141           ENDDO
                0142           DO J=1,sNy+1
                0143            DO I=1,sNx
15c70d7cd1 Jean*0144             tmpScal = vFld(i,j,k,bi,bj)
                0145             vFld(i,j,k,bi,bj) = vFld(i,j,k,bi,bj)
                0146      &         -0.5*tmpFldV(i,j,k,bi,bj)*deltaTmom/Shap_uvtau
                0147             tmpFldV(i,j,k,bi,bj) = tmpScal*_maskS(i,j,k,bi,bj)
aea29c8517 Alis*0148            ENDDO
                0149           ENDDO
                0150          ENDDO
                0151         ENDDO
                0152        ENDDO
                0153 
                0154 
                0155 C      d_yy^n tmpFld 
                0156 
                0157        DO N=1,nShapUV
                0158 
15c70d7cd1 Jean*0159         IF (kSize.EQ.Nr) THEN
                0160           CALL EXCH_UV_XYZ_RL(tmpFldU,tmpFldV,.TRUE.,myThid)
                0161         ELSE
                0162           CALL EXCH_UV_XY_RL(tmpFldU,tmpFldV,.TRUE.,myThid)
                0163         ENDIF
aea29c8517 Alis*0164 
                0165         DO bj=myByLo(myThid),myByHi(myThid)
                0166          DO bi=myBxLo(myThid),myBxHi(myThid)
15c70d7cd1 Jean*0167           DO K=1,kSize
aea29c8517 Alis*0168 
                0169 C          Uyy
                0170            DO J=1,sNy
                0171             DO I=1,sNx+1
                0172              maskZj=_maskS(i-1, j ,k,bi,bj)
                0173      &             *_maskS( i , j ,k,bi,bj)
                0174              maskZp=_maskS(i-1,j+1,k,bi,bj)
                0175      &             *_maskS( i ,j+1,k,bi,bj)
                0176              tmpGrdU(i,j) = -0.25*(
                0177      &        (tmpFldU(i,j+1,k,bi,bj)-tmpFldU(i, j ,k,bi,bj))*maskZp
                0178      &       -(tmpFldU(i, j ,k,bi,bj)-tmpFldU(i,j-1,k,bi,bj))*maskZj
15c70d7cd1 Jean*0179      &       -noSlipFact*(2.-maskZj-maskZp)*tmpFldU(i,j,k,bi,bj)
aea29c8517 Alis*0180      &             )*_maskW(i,j,k,bi,bj)
                0181             ENDDO
                0182            ENDDO
                0183 
                0184            IF (useCubedSphereExchange) THEN
                0185             DO I=1,sNx+1,sNx
15c70d7cd1 Jean*0186             J=1
aea29c8517 Alis*0187              maskZj=maskS(i-1, j ,k,bi,bj)*maskS( i , j ,k,bi,bj)
                0188              maskZp=maskS(i-1,j+1,k,bi,bj)*maskS( i ,j+1,k,bi,bj)
                0189              tmpGrdU(i,j) = -0.25*(
                0190      &        (tmpFldU(i,j+1,k,bi,bj)-tmpFldU(i, j ,k,bi,bj))*maskZp
                0191      &       -(tmpFldU(i, j ,k,bi,bj)-0*tmpFldU(i,j-1,k,bi,bj))*maskZj
15c70d7cd1 Jean*0192      &       -noSlipFact*(2.-maskZj-maskZp)*tmpFldU(i,j,k,bi,bj)
aea29c8517 Alis*0193      &             )*_maskW(i,j,k,bi,bj)
                0194             J=sNy
                0195              maskZj=maskS(i-1, j ,k,bi,bj)*maskS( i , j ,k,bi,bj)
                0196              maskZp=maskS(i-1,j+1,k,bi,bj)*maskS( i ,j+1,k,bi,bj)
                0197              tmpGrdU(i,j) = -0.25*(
                0198      &        (0*tmpFldU(i,j+1,k,bi,bj)-tmpFldU(i, j ,k,bi,bj))*maskZp
                0199      &       -(tmpFldU(i, j ,k,bi,bj)-tmpFldU(i,j-1,k,bi,bj))*maskZj
15c70d7cd1 Jean*0200      &       -noSlipFact*(2.-maskZj-maskZp)*tmpFldU(i,j,k,bi,bj)
aea29c8517 Alis*0201      &             )*_maskW(i,j,k,bi,bj)
                0202             ENDDO
                0203            ENDIF
                0204 
                0205            DO J=1,sNy
                0206             DO I=1,sNx+1
                0207              tmpFldU(i,j,k,bi,bj) = tmpGrdU(i,j)
                0208             ENDDO
                0209            ENDDO
                0210 
                0211 C          Vxx
                0212            DO J=1,sNy+1
                0213             DO I=1,sNx
                0214              maskZj=_maskW( i ,j-1,k,bi,bj)
                0215      &             *_maskW( i , j ,k,bi,bj)
                0216              maskZp=_maskW(i+1,j-1,k,bi,bj)
                0217      &             *_maskW(i+1, j ,k,bi,bj)
                0218              tmpGrdV(i,j) = -0.25*(
                0219      &        (tmpFldV(i+1,j,k,bi,bj)-tmpFldV( i ,j,k,bi,bj))*maskZp
                0220      &       -(tmpFldV( i ,j,k,bi,bj)-tmpFldV(i-1,j,k,bi,bj))*maskZj
15c70d7cd1 Jean*0221      &       -noSlipFact*(2.-maskZj-maskZp)*tmpFldV(i,j,k,bi,bj)
aea29c8517 Alis*0222      &             )*_maskS(i,j,k,bi,bj)
                0223             ENDDO
                0224            ENDDO
                0225 
                0226            IF (useCubedSphereExchange) THEN
                0227             DO J=1,sNy+1,sNy
                0228             I=1
                0229              maskZj=maskW( i ,j-1,k,bi,bj)*maskW( i , j ,k,bi,bj)
                0230              maskZp=maskW(i+1,j-1,k,bi,bj)*maskW(i+1, j ,k,bi,bj)
                0231              tmpGrdV(i,j) = -0.25*(
                0232      &        (tmpFldV(i+1,j,k,bi,bj)-tmpFldV( i ,j,k,bi,bj))*maskZp
                0233      &       -(tmpFldV( i ,j,k,bi,bj)-0*tmpFldV(i-1,j,k,bi,bj))*maskZj
15c70d7cd1 Jean*0234      &       -noSlipFact*(2.-maskZj-maskZp)*tmpFldV(i,j,k,bi,bj)
aea29c8517 Alis*0235      &             )*_maskS(i,j,k,bi,bj)
                0236             I=sNx
                0237              maskZj=maskW( i ,j-1,k,bi,bj)*maskW( i , j ,k,bi,bj)
                0238              maskZp=maskW(i+1,j-1,k,bi,bj)*maskW(i+1, j ,k,bi,bj)
                0239              tmpGrdV(i,j) = -0.25*(
                0240      &        (0*tmpFldV(i+1,j,k,bi,bj)-tmpFldV( i ,j,k,bi,bj))*maskZp
                0241      &       -(tmpFldV( i ,j,k,bi,bj)-tmpFldV(i-1,j,k,bi,bj))*maskZj
15c70d7cd1 Jean*0242      &       -noSlipFact*(2.-maskZj-maskZp)*tmpFldV(i,j,k,bi,bj)
aea29c8517 Alis*0243      &             )*_maskS(i,j,k,bi,bj)
                0244             ENDDO
                0245            ENDIF
                0246 
                0247            DO J=1,sNy+1
                0248             DO I=1,sNx
                0249              tmpFldV(i,j,k,bi,bj) = tmpGrdV(i,j)
                0250             ENDDO
                0251            ENDDO
                0252 
                0253           ENDDO
                0254          ENDDO
                0255         ENDDO
                0256 
                0257        ENDDO
                0258 
15c70d7cd1 Jean*0259 C      F <-  [1 - d_yy^n *deltaT/tau].F
aea29c8517 Alis*0260        DO bj=myByLo(myThid),myByHi(myThid)
                0261         DO bi=myBxLo(myThid),myBxHi(myThid)
15c70d7cd1 Jean*0262          DO K=1,kSize
aea29c8517 Alis*0263           DO J=1,sNy
                0264            DO I=1,sNx+1
15c70d7cd1 Jean*0265             uFld(i,j,k,bi,bj) = uFld(i,j,k,bi,bj)
                0266      &         -0.5*tmpFldU(i,j,k,bi,bj)*deltaTmom/Shap_uvtau
aea29c8517 Alis*0267            ENDDO
                0268           ENDDO
                0269           DO J=1,sNy+1
                0270            DO I=1,sNx
15c70d7cd1 Jean*0271             vFld(i,j,k,bi,bj) = vFld(i,j,k,bi,bj)
                0272      &         -0.5*tmpFldV(i,j,k,bi,bj)*deltaTmom/Shap_uvtau
aea29c8517 Alis*0273            ENDDO
                0274           ENDDO
                0275          ENDDO
                0276         ENDDO
                0277        ENDDO
                0278 
15c70d7cd1 Jean*0279         IF (kSize.EQ.Nr) THEN
                0280           CALL EXCH_UV_XYZ_RL(uFld,vFld,.TRUE.,myThid)
                0281         ELSEIF (kSize.EQ.1) THEN
                0282           CALL EXCH_UV_XY_RL(uFld,vFld,.TRUE.,myThid)
                0283         ELSE
                0284           STOP 'S/R SHAP_FILT_UV_S1: kSize is wrong'
                0285         ENDIF
                0286 
aea29c8517 Alis*0287       ENDIF
                0288 #endif /* ALLOW_SHAP_FILT */
                0289 
                0290       RETURN
                0291       END