Back to home page

MITgcm
 
 

    

File indexing completed on 2024-01-13 06:10:33 UTC

view on githubraw file Latest commit 005af54e on 2024-01-12 20:10:27 UTC
e9cb9b9110 Jean*0001 #include "OBCS_OPTIONS.h"
                0002 
                0003 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0004 CBOP
                0005 C     !ROUTINE: OBCS_DIAG_BALANCE
                0006 
                0007 C     !INTERFACE:
                0008       SUBROUTINE OBCS_DIAG_BALANCE(
26fc470c45 Jean*0009      U                              div2d,
                0010      I                              uTrans, vTrans, k,
e9cb9b9110 Jean*0011      I                              myTime, myIter, myThid )
                0012 
                0013 C     !DESCRIPTION:
                0014 C     *==========================================================*
                0015 C     | SUBROUTINE OBCS_DIAG_BALANCE
26fc470c45 Jean*0016 C     | o For diagnostics purpose, modify horizontal divergence
                0017 C     |   next (but outside) OB to ensure zero net inflow
e9cb9b9110 Jean*0018 C     *==========================================================*
                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 #include "OBCS_PARAMS.h"
                0029 #include "OBCS_GRID.h"
                0030 #include "OBCS_FIELDS.h"
                0031 
                0032 C     !INPUT/OUTPUT PARAMETERS:
26fc470c45 Jean*0033 C     div2d   :: horizontal divergence x grid-cell volume [m^3/s]
e9cb9b9110 Jean*0034 C     uTrans  :: horizontal transport to balance [m^3/s]
                0035 C     vTrans  :: horizontal transport to balance [m^3/s]
26fc470c45 Jean*0036 C     k       :: current level index
e9cb9b9110 Jean*0037 C     myTime  :: current time of simulation (s)
                0038 C     myIter  :: current iteration number
                0039 C     myThid  :: my Thread Id number
26fc470c45 Jean*0040       _RL div2d (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
e9cb9b9110 Jean*0041       _RL uTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0042       _RL vTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
26fc470c45 Jean*0043       INTEGER k
e9cb9b9110 Jean*0044       _RL myTime
                0045       INTEGER myIter
                0046       INTEGER myThid
                0047 CEOP
                0048 
                0049 #ifdef ALLOW_OBCS
                0050 #ifdef ALLOW_DIAGNOSTICS
                0051 
                0052 C     !FUNCTIONS:
                0053 
                0054 C     !LOCAL VARIABLES:
                0055 C     bi, bj       :: tile indices
                0056 C     i,j,k        :: loop indices
                0057 C     iB, jB       :: local index of open boundary
                0058 C     msgBuf       :: Informational/error message buffer
9317497c8c Jean*0059       INTEGER bi, bj, n
e9cb9b9110 Jean*0060       CHARACTER*(MAX_LEN_MBUF) msgBuf
9317497c8c Jean*0061       _RL areaOB(OBCS_maxConnect), tmpA
                0062       _RL inFlow(OBCS_maxConnect)
                0063       _RL tileAreaOB(nSx,nSy,OBCS_maxConnect)
                0064       _RL tileInFlow(nSx,nSy,OBCS_maxConnect)
005af54e38 Jean*0065 #if (defined ALLOW_OBCS_EAST ) || (defined ALLOW_OBCS_WEST )
                0066       INTEGER j, iB
                0067 #endif
                0068 #ifdef ALLOW_OBCS_WEST
                0069       INTEGER iBt
                0070 #endif
                0071 #if (defined ALLOW_OBCS_NORTH) || (defined ALLOW_OBCS_SOUTH)
                0072       INTEGER i, jB
                0073 #endif
                0074 #ifdef ALLOW_OBCS_SOUTH
                0075       INTEGER jBt
                0076 #endif
e9cb9b9110 Jean*0077 
                0078 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0079 
                0080 #ifdef ALLOW_DEBUG
                0081       IF (debugMode) CALL DEBUG_ENTER('OBCS_DIAG_BALANCE',myThid)
                0082 #endif
                0083 
                0084 C--   Integrate the transport through each OB
9317497c8c Jean*0085       DO n=1,OB_connectNumber(k)
                0086         areaOB(n)= 0. _d 0
                0087         inFlow(n)= 0. _d 0
90ee351aa6 Jean*0088         DO bj=myByLo(myThid),myByHi(myThid)
                0089          DO bi=myBxLo(myThid),myBxHi(myThid)
                0090           tileAreaOB(bi,bj,n) = 0.
                0091           tileInFlow(bi,bj,n) = 0.
                0092          ENDDO
9317497c8c Jean*0093         ENDDO
e9cb9b9110 Jean*0094       ENDDO
                0095 
                0096 #ifdef ALLOW_OBCS_EAST
                0097       DO bj=myByLo(myThid),myByHi(myThid)
                0098        DO bi=myBxLo(myThid),myBxHi(myThid)
                0099         IF ( tileHasOBE(bi,bj) ) THEN
                0100          DO j=1,sNy
                0101            iB = OB_Ie(j,bi,bj)
586e4775df Jean*0102            IF ( iB.NE.OB_indexNone .AND. iB.GT.1 ) THEN
9317497c8c Jean*0103             n = OBE_connect(j,k,bi,bj)
                0104             IF ( n.GE.1 ) THEN
                0105               tmpA = drF(k)*hFacW(iB,j,k,bi,bj)*dyG(iB,j,bi,bj)
                0106      &                     *maskInW(iB,j,bi,bj)
                0107               tileAreaOB(bi,bj,n) = tileAreaOB(bi,bj,n) + tmpA
                0108               tileInFlow(bi,bj,n) = tileInFlow(bi,bj,n)
                0109      &                            - uTrans(iB,j,bi,bj)
                0110             ENDIF
e9cb9b9110 Jean*0111            ENDIF
                0112          ENDDO
                0113         ENDIF
                0114        ENDDO
                0115       ENDDO
                0116 #endif /* ALLOW_OBCS_EAST */
                0117 
                0118 #ifdef ALLOW_OBCS_WEST
                0119       DO bj=myByLo(myThid),myByHi(myThid)
                0120        DO bi=myBxLo(myThid),myBxHi(myThid)
                0121         IF ( tileHasOBW(bi,bj) ) THEN
                0122          DO j=1,sNy
586e4775df Jean*0123            iB = OB_Iw(j,bi,bj)
                0124            IF ( iB.NE.OB_indexNone .AND. iB.LT.sNx ) THEN
9317497c8c Jean*0125             n = OBW_connect(j,k,bi,bj)
                0126             IF ( n.GE.1 ) THEN
                0127               iB = 1+iB
                0128               tmpA = drF(k)*hFacW(iB,j,k,bi,bj)*dyG(iB,j,bi,bj)
                0129      &                     *maskInW(iB,j,bi,bj)
                0130               tileAreaOB(bi,bj,n) = tileAreaOB(bi,bj,n) + tmpA
                0131               tileInFlow(bi,bj,n) = tileInFlow(bi,bj,n)
                0132      &                            + uTrans(iB,j,bi,bj)
                0133             ENDIF
e9cb9b9110 Jean*0134            ENDIF
                0135          ENDDO
                0136         ENDIF
                0137        ENDDO
                0138       ENDDO
                0139 #endif /* ALLOW_OBCS_WEST */
                0140 
                0141 #ifdef ALLOW_OBCS_NORTH
                0142       DO bj=myByLo(myThid),myByHi(myThid)
                0143        DO bi=myBxLo(myThid),myBxHi(myThid)
                0144         IF ( tileHasOBN(bi,bj) ) THEN
                0145          DO i=1,sNx
                0146            jB = OB_Jn(i,bi,bj)
586e4775df Jean*0147            IF ( jB.NE.OB_indexNone .AND. jB.GT.1 ) THEN
9317497c8c Jean*0148             n = OBN_connect(i,k,bi,bj)
                0149             IF ( n.GE.1 ) THEN
                0150               tmpA = drF(k)*hFacS(i,jB,k,bi,bj)*dxG(i,jB,bi,bj)
                0151      &                     *maskInS(i,jB,bi,bj)
                0152               tileAreaOB(bi,bj,n) = tileAreaOB(bi,bj,n) + tmpA
                0153               tileInFlow(bi,bj,n) = tileInFlow(bi,bj,n)
                0154      &                            - vTrans(i,jB,bi,bj)
                0155             ENDIF
e9cb9b9110 Jean*0156            ENDIF
                0157          ENDDO
                0158         ENDIF
                0159        ENDDO
                0160       ENDDO
                0161 #endif /* ALLOW_OBCS_NORTH */
                0162 
                0163 #ifdef ALLOW_OBCS_SOUTH
                0164       DO bj=myByLo(myThid),myByHi(myThid)
                0165        DO bi=myBxLo(myThid),myBxHi(myThid)
                0166         IF ( tileHasOBS(bi,bj) ) THEN
                0167          DO i=1,sNx
586e4775df Jean*0168            jB = OB_Js(i,bi,bj)
                0169            IF ( jB.NE.OB_indexNone .AND. jB.LT.sNy ) THEN
9317497c8c Jean*0170             n = OBS_connect(i,k,bi,bj)
                0171             IF ( n.GE.1 ) THEN
                0172               jB = 1+jB
                0173               tmpA = drF(k)*hFacS(i,jB,k,bi,bj)*dxG(i,jB,bi,bj)
                0174      &                     *maskInS(i,jB,bi,bj)
                0175               tileAreaOB(bi,bj,n) = tileAreaOB(bi,bj,n) + tmpA
                0176               tileInFlow(bi,bj,n) = tileInFlow(bi,bj,n)
                0177      &                            + vTrans(i,jB,bi,bj)
                0178             ENDIF
e9cb9b9110 Jean*0179            ENDIF
                0180          ENDDO
                0181         ENDIF
                0182        ENDDO
                0183       ENDDO
                0184 #endif /* ALLOW_OBCS_SOUTH */
                0185 
                0186 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
9317497c8c Jean*0187 C--   For each set of OB connected points, calculate a unique velocity
                0188 C     correction and correct to the corresponding OB point divergence
e9cb9b9110 Jean*0189 
9317497c8c Jean*0190       DO n=1,OB_connectNumber(k)
                0191        CALL GLOBAL_SUM_TILE_RL( tileAreaOB(1,1,n), areaOB(n), myThid )
                0192        IF ( areaOB(n).GT.0. ) THEN
                0193         CALL GLOBAL_SUM_TILE_RL( tileInFlow(1,1,n), inFlow(n), myThid )
e9cb9b9110 Jean*0194         IF ( debugLevel.GE.debLevC ) THEN
9317497c8c Jean*0195           WRITE(msgBuf,'(A,I9,2I4,A,1P2E16.8)')
                0196      &     'OBCS_DIAG_BALANCE (it,k,n=', myIter, k, n,
                0197      &       ' ) inFlow:',inFlow(n),inFlow(n)/areaOB(n)
e9cb9b9110 Jean*0198           CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
                0199      &                        SQUEEZE_RIGHT, myThid )
                0200         ENDIF
9317497c8c Jean*0201         inFlow(n) = inFlow(n) / areaOB(n)
                0202        ENDIF
                0203       ENDDO
e9cb9b9110 Jean*0204 
                0205 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0206 C--   Add correction:
                0207 
26fc470c45 Jean*0208       DO bj=myByLo(myThid),myByHi(myThid)
                0209        DO bi=myBxLo(myThid),myBxHi(myThid)
e9cb9b9110 Jean*0210 #ifdef ALLOW_OBCS_EAST
                0211          IF ( tileHasOBE(bi,bj) ) THEN
9317497c8c Jean*0212 c          DO j=1-OLy,sNy+OLy
                0213            DO j=1,sNy
586e4775df Jean*0214             IF ( OB_Ie(j,bi,bj).NE.OB_indexNone ) THEN
                0215              iB = OB_Ie(j,bi,bj)
9317497c8c Jean*0216              n = OBE_connect(j,k,bi,bj)
                0217              IF ( n.EQ.0 ) THEN
                0218                div2d(iB ,j,bi,bj) = div2d(iB ,j,bi,bj)
                0219      &                            - uTrans(iB,j,bi,bj)
                0220              ELSE
                0221                div2d(iB ,j,bi,bj) = div2d(iB,j,bi,bj)
                0222      &                 + inFlow(n)*drF(k)*hFacW(iB,j,k,bi,bj)
                0223      &                            *dyG(iB,j,bi,bj)*maskInW(iB,j,bi,bj)
                0224              ENDIF
e9cb9b9110 Jean*0225             ENDIF
                0226            ENDDO
                0227          ENDIF
                0228 #endif /* ALLOW_OBCS_EAST */
                0229 
                0230 #ifdef ALLOW_OBCS_WEST
                0231          IF ( tileHasOBW(bi,bj) ) THEN
9317497c8c Jean*0232 c          DO j=1-OLy,sNy+OLy
                0233            DO j=1,sNy
586e4775df Jean*0234             IF ( OB_Iw(j,bi,bj).NE.OB_indexNone ) THEN
9317497c8c Jean*0235              iBt = OB_Iw(j,bi,bj)
                0236              iB = 1+iBt
                0237              n = OBW_connect(j,k,bi,bj)
                0238              IF ( n.EQ.0 ) THEN
                0239                div2d(iBt,j,bi,bj) = div2d(iBt,j,bi,bj)
                0240      &                            + uTrans(iB,j,bi,bj)
                0241              ELSE
                0242                div2d(iBt,j,bi,bj) = div2d(iBt,j,bi,bj)
                0243      &                 + inFlow(n)*drF(k)*hFacW(iB,j,k,bi,bj)
                0244      &                            *dyG(iB,j,bi,bj)*maskInW(iB,j,bi,bj)
                0245              ENDIF
e9cb9b9110 Jean*0246             ENDIF
                0247            ENDDO
                0248          ENDIF
                0249 #endif /* ALLOW_OBCS_WEST */
                0250 
                0251 #ifdef ALLOW_OBCS_NORTH
                0252          IF ( tileHasOBN(bi,bj) ) THEN
9317497c8c Jean*0253 c          DO i=1-OLx,sNx+OLx
                0254            DO i=1,sNx
586e4775df Jean*0255             IF ( OB_Jn(i,bi,bj).NE.OB_indexNone ) THEN
                0256              jB = OB_Jn(i,bi,bj)
9317497c8c Jean*0257              n = OBN_connect(i,k,bi,bj)
                0258              IF ( n.EQ.0 ) THEN
                0259                div2d(i,jB ,bi,bj) = div2d(i,jB ,bi,bj)
                0260      &                            - vTrans(i,jB,bi,bj)
                0261              ELSE
                0262                div2d(i,jB ,bi,bj) = div2d(i,jB ,bi,bj)
                0263      &                 + inFlow(n)*drF(k)*hFacS(i,jB,k,bi,bj)
                0264      &                            *dxG(i,jB,bi,bj)*maskInS(i,jB,bi,bj)
                0265              ENDIF
e9cb9b9110 Jean*0266             ENDIF
                0267            ENDDO
                0268          ENDIF
                0269 #endif /* ALLOW_OBCS_NORTH */
                0270 
                0271 #ifdef ALLOW_OBCS_SOUTH
                0272          IF ( tileHasOBS(bi,bj) ) THEN
9317497c8c Jean*0273 c          DO i=1-OLx,sNx+OLx
                0274            DO i=1,sNx
586e4775df Jean*0275             IF ( OB_Js(i,bi,bj).NE.OB_indexNone ) THEN
9317497c8c Jean*0276              jBt = OB_Js(i,bi,bj)
                0277              jB = 1+jBt
90ee351aa6 Jean*0278              n = OBS_connect(i,k,bi,bj)
9317497c8c Jean*0279              IF ( n.EQ.0 ) THEN
                0280                div2d(i,jBt,bi,bj) = div2d(i,jBt,bi,bj)
                0281      &                            + vTrans(i,jB,bi,bj)
                0282              ELSE
                0283                div2d(i,jBt,bi,bj) = div2d(i,jBt,bi,bj)
                0284      &                 + inFlow(n)*drF(k)*hFacS(i,jB,k,bi,bj)
                0285      &                            *dxG(i,jB,bi,bj)*maskInS(i,jB,bi,bj)
                0286              ENDIF
e9cb9b9110 Jean*0287             ENDIF
                0288            ENDDO
                0289          ENDIF
                0290 #endif /* ALLOW_OBCS_SOUTH */
                0291 
26fc470c45 Jean*0292        ENDDO
                0293       ENDDO
                0294 
e9cb9b9110 Jean*0295 #ifdef ALLOW_DEBUG
                0296       IF (debugMode) CALL DEBUG_LEAVE('OBCS_DIAG_BALANCE',myThid)
                0297 #endif
                0298 
                0299 #endif /* ALLOW_DIAGNOSTICS */
                0300 #endif /* ALLOW_OBCS */
                0301 
                0302       RETURN
                0303       END
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated from https://github.com/MITgcm/MITgcm by the 2.2.1-MITgcm-0.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!