Back to home page

MITgcm

 
 

    


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

view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
bd0c662004 Jean*0001 #include "SHAP_FILT_OPTIONS.h"
                0002 
                0003 CBOP
                0004 C     !ROUTINE: SHAP_FILT_COMPUTVORT
                0005 C     !INTERFACE:
                0006       SUBROUTINE SHAP_FILT_COMPUTVORT(
                0007      I           uFld, vFld,
                0008      O           vort,
                0009      I           k, bi,bj, myThid )
                0010 C     !DESCRIPTION: \bv
                0011 C     *==========================================================*
                0012 C     | S/R SHAP_FILT_COMPUTVORT
                0013 C     | o Calculate delta_i[vFld]-delta_j[uFld]
                0014 C     *==========================================================*
                0015 C     | o used in computational-mode filter to replace relative
                0016 C     |   vorticity
                0017 C     *==========================================================*
                0018 C     \ev
                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"
                0028 #ifdef ALLOW_EXCH2
f9f661930b Jean*0029 #include "W2_EXCH2_SIZE.h"
bd0c662004 Jean*0030 #include "W2_EXCH2_TOPOLOGY.h"
                0031 #endif /* ALLOW_EXCH2 */
                0032 
                0033 C     !INPUT/OUTPUT PARAMETERS:
                0034 C     == Routine arguments
                0035 C     uFld :: velocity field (U component) on which filter applies
                0036 C     vFld :: velocity field (V component) on which filter applies
                0037 C     myThid :: Thread number for this instance of SHAP_FILT_UV_S2
                0038       _RL uFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0039       _RL vFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0040       _RL vort(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0041       INTEGER k, bi,bj
                0042       INTEGER myThid
                0043 
                0044 #ifdef ALLOW_SHAP_FILT
                0045 
                0046 C     !LOCAL VARIABLES:
                0047 C     == Local variables ==
                0048       INTEGER i,j
                0049       _RS maskZ
                0050       LOGICAL northWestCorner, northEastCorner,
                0051      &        southWestCorner, southEastCorner
                0052       INTEGER myFace
                0053 #ifdef ALLOW_EXCH2
                0054       INTEGER myTile
                0055 #endif /* ALLOW_EXCH2 */
                0056 CEOP
                0057 
41b8eb6e45 Jean*0058 #ifdef ALLOW_AUTODIFF
bd0c662004 Jean*0059 C-    Initialisation :
41b8eb6e45 Jean*0060       DO j=1-OLy,sNy+OLy
                0061         DO i=1-OLx,sNx+OLx
bd0c662004 Jean*0062           vort(i,j)= 0.
                0063         ENDDO
                0064       ENDDO
3cff26a467 Jean*0065 #endif
                0066 
bd0c662004 Jean*0067 C-    replace Physical calc Div & Vort by computational one :
41b8eb6e45 Jean*0068       DO j=2-OLy,sNy+OLy
                0069         DO i=2-OLx,sNx+OLx
bd0c662004 Jean*0070           vort(i,j) = ( vFld(i,j)-vFld(i-1,j) )
                0071      &              - ( uFld(i,j)-uFld(i,j-1) )
                0072           maskZ = (maskW(i,j,k,bi,bj)+maskW(i,j-1,k,bi,bj))
                0073      &           *(maskS(i,j,k,bi,bj)+maskS(i-1,j,k,bi,bj))
                0074           IF (maskZ.LT.1.) vort(i,j)=0.
                0075         ENDDO
                0076       ENDDO
                0077 
                0078 C-    Special stuff for Cubed Sphere
                0079       IF (useCubedSphereExchange) THEN
                0080 #ifdef ALLOW_EXCH2
c424ee7cc7 Jean*0081          myTile = W2_myTileList(bi,bj)
bd0c662004 Jean*0082          myFace = exch2_myFace(myTile)
                0083          southWestCorner = exch2_isWedge(myTile).EQ.1
                0084      &               .AND. exch2_isSedge(myTile).EQ.1
                0085          southEastCorner = exch2_isEedge(myTile).EQ.1
                0086      &               .AND. exch2_isSedge(myTile).EQ.1
                0087          northEastCorner = exch2_isEedge(myTile).EQ.1
                0088      &               .AND. exch2_isNedge(myTile).EQ.1
                0089          northWestCorner = exch2_isWedge(myTile).EQ.1
                0090      &               .AND. exch2_isNedge(myTile).EQ.1
                0091 #else
                0092          myFace = bi
                0093          southWestCorner = .TRUE.
                0094          southEastCorner = .TRUE.
                0095          northWestCorner = .TRUE.
                0096          northEastCorner = .TRUE.
                0097 #endif /* ALLOW_EXCH2 */
                0098 C---
                0099          IF ( southWestCorner ) THEN
                0100            i=1
                0101            j=1
                0102            maskZ = maskW(i,j,k,bi,bj)+maskW(i,j-1,k,bi,bj)
                0103      &            +maskS(i,j,k,bi,bj)
                0104            IF (maskZ.GE.2.) THEN
                0105              vort(i,j)=
                0106      &          (+vFld(i,j) -uFld(i,j) ) +uFld(i,j-1)
                0107              vort(i,j)=vort(i,j)*4. _d 0 / 3. _d 0
                0108            ELSE
                0109              vort(i,j)=0.
                0110            ENDIF
                0111          ENDIF
                0112 C---
                0113          IF ( southEastCorner ) THEN
                0114            i=sNx+1
                0115            j=1
                0116            maskZ = maskW(i,j,k,bi,bj)+maskW(i,j-1,k,bi,bj)
                0117      &                               +maskS(i-1,j,k,bi,bj)
                0118            IF (maskZ.GE.2.) THEN
                0119              IF ( myFace.EQ.2 ) THEN
                0120               vort(i,j)=
                0121      &          (-uFld(i,j) -vFld(i-1,j) ) +uFld(i,j-1)
                0122              ELSEIF ( myFace.EQ.4 ) THEN
                0123               vort(i,j)=
                0124      &          (-vFld(i-1,j) +uFld(i,j-1) ) -uFld(i,j)
                0125              ELSE
                0126               vort(i,j)=
                0127      &          (+uFld(i,j-1) -uFld(i,j) ) -vFld(i-1,j)
                0128              ENDIF
                0129              vort(i,j)=vort(i,j)*4. _d 0 / 3. _d 0
                0130            ELSE
                0131              vort(i,j)=0.
                0132            ENDIF
                0133          ENDIF
                0134 C---
                0135          IF ( northWestCorner ) THEN
                0136            i=1
                0137            j=sNy+1
                0138            maskZ = maskW(i,j,k,bi,bj)+maskW(i,j-1,k,bi,bj)
                0139      &            +maskS(i,j,k,bi,bj)
                0140            IF (maskZ.GE.2.) THEN
                0141              IF ( myFace.EQ.1 ) THEN
                0142               vort(i,j)=
                0143      &          (+uFld(i,j-1) +vFld(i,j) ) -uFld(i,j)
                0144              ELSEIF ( myFace.EQ.3 ) THEN
                0145               vort(i,j)=
                0146      &          (-uFld(i,j) +uFld(i,j-1) ) +vFld(i,j)
                0147              ELSE
                0148               vort(i,j)=
                0149      &          (+vFld(i,j) -uFld(i,j) ) +uFld(i,j-1)
                0150              ENDIF
                0151              vort(i,j)=vort(i,j)*4. _d 0 / 3. _d 0
                0152            ELSE
                0153              vort(i,j)=0.
                0154            ENDIF
                0155          ENDIF
                0156 C---
                0157          IF ( northEastCorner ) THEN
                0158            i=sNx+1
                0159            j=sNy+1
                0160            maskZ = maskW(i,j,k,bi,bj)+maskW(i,j-1,k,bi,bj)
                0161      &                               +maskS(i-1,j,k,bi,bj)
                0162            IF (maskZ.GE.2.) THEN
                0163              IF ( MOD(myFace,2).EQ.1 ) THEN
                0164               vort(i,j)=
                0165      &          (-uFld(i,j) -vFld(i-1,j) ) +uFld(i,j-1)
                0166              ELSE
                0167               vort(i,j)=
                0168      &          (+uFld(i,j-1) -uFld(i,j) ) -vFld(i-1,j)
                0169              ENDIF
                0170              vort(i,j)=vort(i,j)*4. _d 0 / 3. _d 0
                0171            ELSE
                0172              vort(i,j)=0.
                0173            ENDIF
                0174          ENDIF
                0175 C---  end if useCubedSphereExchange:
                0176       ENDIF
                0177 
                0178 c---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0179 
                0180 #endif /* ALLOW_SHAP_FILT */
                0181 
                0182       RETURN
                0183       END