Back to home page

MITgcm

 
 

    


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

view on githubraw file Latest commit 005af54e on 2024-01-12 20:10:27 UTC
96bbd4e2a5 Patr*0001 #include "OBCS_OPTIONS.h"
                0002 
eb5718bf71 Jean*0003 C--   File obcs_sponge.F:
                0004 C--    Contents:
                0005 C--    o OBCS_SPONGE_U
                0006 C--    o OBCS_SPONGE_V
                0007 C--    o OBCS_SPONGE_T
                0008 C--    o OBCS_SPONGE_S
                0009 
                0010 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0011 
96bbd4e2a5 Patr*0012 CStartOfInterface
                0013       SUBROUTINE OBCS_SPONGE_U(
73b1dccda0 Jean*0014      U                gU_arr,
                0015      I                iMin,iMax,jMin,jMax, k, bi, bj,
                0016      I                myTime, myIter, myThid )
eb5718bf71 Jean*0017 C     *==========================================================*
                0018 C     | S/R OBCS_SPONGE_U
                0019 C     | o Contains problem specific forcing for zonal velocity.
                0020 C     *==========================================================*
                0021 C     | Adds a relaxation term to gU near Open-Boundaries
                0022 C     *==========================================================*
96bbd4e2a5 Patr*0023       IMPLICIT NONE
                0024 
                0025 C     == Global data ==
                0026 #include "SIZE.h"
                0027 #include "EEPARAMS.h"
                0028 #include "PARAMS.h"
                0029 #include "GRID.h"
                0030 #include "DYNVARS.h"
9b4f2a04e2 Jean*0031 #include "OBCS_PARAMS.h"
                0032 #include "OBCS_GRID.h"
                0033 #include "OBCS_FIELDS.h"
96bbd4e2a5 Patr*0034 
                0035 C     == Routine arguments ==
73b1dccda0 Jean*0036 C     gU_arr    :: the tendency array
                0037 C     iMin,iMax :: Working range of x-index for applying forcing.
                0038 C     jMin,jMax :: Working range of y-index for applying forcing.
                0039 C     k         :: Current vertical level index
                0040 C     bi,bj     :: Current tile indices
                0041       _RL gU_arr(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0042       INTEGER iMin, iMax, jMin, jMax
                0043       INTEGER k, bi, bj
eb5718bf71 Jean*0044       _RL myTime
73b1dccda0 Jean*0045       INTEGER myIter
96bbd4e2a5 Patr*0046       INTEGER myThid
                0047 CEndOfInterface
                0048 
b2a6135872 Dimi*0049 #if defined(ALLOW_OBCS) && defined(ALLOW_OBCS_SPONGE)
96bbd4e2a5 Patr*0050 C     == Local variables ==
                0051 C     Loop counters
005af54e38 Jean*0052       INTEGER i, j
                0053 #if (defined ALLOW_OBCS_EAST ) || (defined ALLOW_OBCS_WEST )
                0054       INTEGER isl
                0055 #endif
                0056 #if (defined ALLOW_OBCS_NORTH) || (defined ALLOW_OBCS_SOUTH)
                0057       INTEGER jsl
                0058 #endif
96bbd4e2a5 Patr*0059       _RL urelax, lambda_obcs_u
                0060 
74019f026d Jean*0061       IF ( useOBCSsponge .AND. spongeThickness.NE.0 ) THEN
96bbd4e2a5 Patr*0062 
                0063 C Northern Open Boundary
432b2a4585 Dimi*0064 #ifdef ALLOW_OBCS_NORTH
                0065        IF ( OBCSsponge_N .AND. OBCSsponge_UatNS ) THEN
                0066         DO i=iMin,iMax
                0067          IF ( OB_Jn(i,bi,bj).NE.OB_indexNone ) THEN
                0068           DO jsl= 1,spongeThickness
                0069            j=OB_Jn(i,bi,bj)-jsl
                0070            IF ((j.ge.jmin).and.(j.le.jmax)) THEN
                0071             IF (useLinearSponge) THEN
                0072              urelax = OBNu(i,k,bi,bj)
                0073             ELSE
                0074              urelax=(
                0075      &            float(spongeThickness-jsl)*OBNu(i,k,bi,bj)
                0076      &            + float(jsl)*uVel(i,j,k,bi,bj) )
                0077      &            / float(spongeThickness)
                0078             ENDIF
                0079             lambda_obcs_u = (
                0080      &           float(spongeThickness-jsl)*Vrelaxobcsbound
                0081      &           + float(jsl)*Vrelaxobcsinner)
                0082      &           / float(spongeThickness)
                0083             IF (lambda_obcs_u.ne.0.) THEN
                0084              lambda_obcs_u = 1. _d 0 / lambda_obcs_u
                0085             ELSE
                0086              lambda_obcs_u = 0. _d 0
                0087             ENDIF
                0088             gU_arr(i,j) = gU_arr(i,j)
                0089      &           - _maskW(i,j,k,bi,bj) * lambda_obcs_u
                0090      &           * ( uVel(i,j,k,bi,bj) - urelax )
                0091            ENDIF
                0092           ENDDO
96bbd4e2a5 Patr*0093          ENDIF
                0094         ENDDO
                0095        ENDIF
432b2a4585 Dimi*0096 #endif /* ALLOW_OBCS_NORTH */
96bbd4e2a5 Patr*0097 
                0098 C Southern Open Boundary
432b2a4585 Dimi*0099 #ifdef ALLOW_OBCS_SOUTH
                0100        IF ( OBCSsponge_S .AND. OBCSsponge_UatNS ) THEN
                0101         DO i=iMin,iMax
                0102          IF ( OB_Js(i,bi,bj).NE.OB_indexNone ) THEN
                0103           DO jsl= 1,spongeThickness
                0104            j=OB_Js(i,bi,bj)+jsl
                0105            IF ((j.ge.jmin).and.(j.le.jmax)) THEN
                0106             IF (useLinearSponge) THEN
                0107              urelax= OBSu(i,k,bi,bj)
                0108             ELSE
                0109              urelax=(
                0110      &            float(spongeThickness-jsl)*OBSu(i,k,bi,bj)
                0111      &            + float(jsl)*uVel(i,j,k,bi,bj) )
                0112      &            / float(spongeThickness)
                0113             ENDIF
                0114             lambda_obcs_u = (
                0115      &           float(spongeThickness-jsl)*Vrelaxobcsbound
                0116      &           + float(jsl)*Vrelaxobcsinner)
                0117      &           / float(spongeThickness)
                0118             IF (lambda_obcs_u.ne.0.) THEN
                0119              lambda_obcs_u = 1. _d 0 / lambda_obcs_u
                0120             ELSE
                0121              lambda_obcs_u = 0. _d 0
                0122             ENDIF
                0123             gU_arr(i,j) = gU_arr(i,j)
                0124      &           - _maskW(i,j,k,bi,bj) * lambda_obcs_u
                0125      &           * ( uVel(i,j,k,bi,bj) - urelax )
                0126            ENDIF
                0127           ENDDO
96bbd4e2a5 Patr*0128          ENDIF
                0129         ENDDO
                0130        ENDIF
432b2a4585 Dimi*0131 #endif /* ALLOW_OBCS_SOUTH */
96bbd4e2a5 Patr*0132 
eb5718bf71 Jean*0133 C Eastern Open Boundary
432b2a4585 Dimi*0134 #ifdef ALLOW_OBCS_EAST
                0135        IF ( OBCSsponge_E .AND. OBCSsponge_UatEW ) THEN
                0136         DO j=jMin,jMax
                0137          IF ( OB_Ie(j,bi,bj).NE.OB_indexNone ) THEN
                0138           DO isl= 1,spongeThickness
                0139            i=OB_Ie(j,bi,bj)-isl
                0140            IF ((i.ge.imin).and.(i.le.imax)) THEN
                0141             IF (useLinearSponge) THEN
                0142              urelax=OBEu(j,k,bi,bj)
                0143             ELSE
                0144              urelax=(
                0145      &            float(spongeThickness-isl)*OBEu(j,k,bi,bj)
                0146      &            + float(isl)*uVel(i,j,k,bi,bj) )
                0147      &            / float(spongeThickness)
                0148             ENDIF
                0149             lambda_obcs_u = (
                0150      &           float(spongeThickness-isl)*Urelaxobcsbound
                0151      &           + float(isl)*Urelaxobcsinner)
                0152      &           / float(spongeThickness)
                0153             IF (lambda_obcs_u.ne.0.) THEN
                0154              lambda_obcs_u = 1. _d 0 / lambda_obcs_u
                0155             ELSE
                0156              lambda_obcs_u = 0. _d 0
                0157             ENDIF
                0158             gU_arr(i,j) = gU_arr(i,j)
                0159      &           - _maskW(i,j,k,bi,bj) * lambda_obcs_u
                0160      &           * ( uVel(i,j,k,bi,bj) - urelax )
                0161            ENDIF
                0162           ENDDO
96bbd4e2a5 Patr*0163          ENDIF
                0164         ENDDO
                0165        ENDIF
432b2a4585 Dimi*0166 #endif /* ALLOW_OBCS_EAST */
96bbd4e2a5 Patr*0167 
                0168 C Western Open Boundary
432b2a4585 Dimi*0169 #ifdef ALLOW_OBCS_WEST
                0170        IF ( OBCSsponge_W .AND. OBCSsponge_UatEW ) THEN
                0171         DO j=jMin,jMax
                0172          IF ( OB_Iw(j,bi,bj).NE.OB_indexNone ) THEN
                0173           DO isl= 1,spongeThickness
                0174            i=OB_Iw(j,bi,bj)+isl+1
                0175            IF ((i.ge.imin).and.(i.le.imax)) THEN
                0176             IF (useLinearSponge) THEN
                0177              urelax= OBWu(j,k,bi,bj)
                0178             ELSE
                0179              urelax=(
                0180      &            float(spongeThickness-isl)*OBWu(j,k,bi,bj)
                0181      &            + float(isl)*uVel(i,j,k,bi,bj) )
                0182      &            / float(spongeThickness)
                0183             ENDIF
                0184             lambda_obcs_u= (
                0185      &           float(spongeThickness-isl)*Urelaxobcsbound
                0186      &           + float(isl)*Urelaxobcsinner)
                0187      &           / float(spongeThickness)
                0188             IF (lambda_obcs_u.ne.0.) THEN
                0189              lambda_obcs_u = 1. _d 0 / lambda_obcs_u
                0190             ELSE
                0191              lambda_obcs_u = 0. _d 0
                0192             ENDIF
                0193             gU_arr(i,j) = gU_arr(i,j)
                0194      &           - _maskW(i,j,k,bi,bj) * lambda_obcs_u
                0195      &           * ( uVel(i,j,k,bi,bj) - urelax )
                0196            ENDIF
                0197           ENDDO
96bbd4e2a5 Patr*0198          ENDIF
                0199         ENDDO
                0200        ENDIF
432b2a4585 Dimi*0201 #endif /* ALLOW_OBCS_WEST */
96bbd4e2a5 Patr*0202 
                0203       ENDIF
                0204 
f7f6458a80 Jean*0205 #endif /* ALLOW_OBCS & ALLOW_OBCS_SPONGE */
96bbd4e2a5 Patr*0206 
                0207       RETURN
                0208       END
                0209 
eb5718bf71 Jean*0210 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0211 
96bbd4e2a5 Patr*0212 CStartOfInterface
                0213       SUBROUTINE OBCS_SPONGE_V(
73b1dccda0 Jean*0214      U                gV_arr,
                0215      I                iMin,iMax,jMin,jMax, k, bi, bj,
                0216      I                myTime, myIter, myThid )
eb5718bf71 Jean*0217 C     *==========================================================*
                0218 C     | S/R OBCS_SPONGE_V
                0219 C     | o Contains problem specific forcing for merid velocity.
                0220 C     *==========================================================*
                0221 C     | Adds a relaxation term to gV near Open-Boundaries
                0222 C     *==========================================================*
96bbd4e2a5 Patr*0223       IMPLICIT NONE
                0224 
                0225 C     == Global data ==
                0226 #include "SIZE.h"
                0227 #include "EEPARAMS.h"
                0228 #include "PARAMS.h"
                0229 #include "GRID.h"
                0230 #include "DYNVARS.h"
9b4f2a04e2 Jean*0231 #include "OBCS_PARAMS.h"
                0232 #include "OBCS_GRID.h"
                0233 #include "OBCS_FIELDS.h"
96bbd4e2a5 Patr*0234 
                0235 C     == Routine arguments ==
73b1dccda0 Jean*0236 C     gV_arr    :: the tendency array
                0237 C     iMin,iMax :: Working range of x-index for applying forcing.
                0238 C     jMin,jMax :: Working range of y-index for applying forcing.
                0239 C     k         :: Current vertical level index
                0240 C     bi,bj     :: Current tile indices
                0241       _RL gV_arr(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0242       INTEGER iMin, iMax, jMin, jMax
                0243       INTEGER k, bi, bj
eb5718bf71 Jean*0244       _RL myTime
73b1dccda0 Jean*0245       INTEGER myIter
96bbd4e2a5 Patr*0246       INTEGER myThid
                0247 CEndOfInterface
f7f6458a80 Jean*0248 
b2a6135872 Dimi*0249 #if defined(ALLOW_OBCS) && defined(ALLOW_OBCS_SPONGE)
96bbd4e2a5 Patr*0250 C     == Local variables ==
                0251 C     Loop counters
005af54e38 Jean*0252       INTEGER i, j
                0253 #if (defined ALLOW_OBCS_EAST ) || (defined ALLOW_OBCS_WEST )
                0254       INTEGER isl
                0255 #endif
                0256 #if (defined ALLOW_OBCS_NORTH) || (defined ALLOW_OBCS_SOUTH)
                0257       INTEGER jsl
                0258 #endif
96bbd4e2a5 Patr*0259       _RL vrelax,lambda_obcs_v
                0260 
74019f026d Jean*0261       IF ( useOBCSsponge .AND. spongeThickness.NE.0 ) THEN
96bbd4e2a5 Patr*0262 
                0263 C Northern Open Boundary
432b2a4585 Dimi*0264 #ifdef ALLOW_OBCS_NORTH
                0265        IF ( OBCSsponge_N .AND. OBCSsponge_VatNS ) THEN
                0266         DO i=iMin,iMax
                0267          IF ( OB_Jn(i,bi,bj).NE.OB_indexNone ) THEN
                0268           DO jsl= 1,spongeThickness
                0269            j=OB_Jn(i,bi,bj)-jsl
                0270            IF ((j.ge.jmin).and.(j.le.jmax)) THEN
                0271             IF (useLinearSponge) THEN
                0272              vrelax= OBNv(i,k,bi,bj)
                0273             ELSE
                0274              vrelax=(
                0275      &            float(spongeThickness-jsl)*OBNv(i,k,bi,bj)
                0276      &            + float(jsl)*vVel(i,j,k,bi,bj) )
                0277      &            / float(spongeThickness)
                0278             ENDIF
                0279             lambda_obcs_v = (
                0280      &           float(spongeThickness-jsl)*Vrelaxobcsbound
                0281      &           + float(jsl)*Vrelaxobcsinner)
                0282      &           / float(spongeThickness)
                0283             IF (lambda_obcs_v.ne.0.) THEN
                0284              lambda_obcs_v = 1. _d 0 / lambda_obcs_v
                0285             ELSE
                0286              lambda_obcs_v = 0. _d 0
                0287             ENDIF
                0288             gV_arr(i,j) = gV_arr(i,j)
                0289      &           - _maskS(i,j,k,bi,bj) * lambda_obcs_v
                0290      &           * ( vVel(i,j,k,bi,bj) - vrelax )
                0291            ENDIF
                0292           ENDDO
96bbd4e2a5 Patr*0293          ENDIF
                0294         ENDDO
                0295        ENDIF
432b2a4585 Dimi*0296 #endif /* ALLOW_OBCS_NORTH */
96bbd4e2a5 Patr*0297 
                0298 C Southern Open Boundary
432b2a4585 Dimi*0299 #ifdef ALLOW_OBCS_SOUTH
                0300        IF ( OBCSsponge_S .AND. OBCSsponge_VatNS ) THEN
                0301         DO i=iMin,iMax
                0302          IF ( OB_Js(i,bi,bj).NE.OB_indexNone ) THEN
                0303           DO jsl= 1,spongeThickness
                0304            j=OB_Js(i,bi,bj)+jsl+1
                0305            IF ((j.ge.jmin).and.(j.le.jmax)) THEN
                0306             IF (useLinearSponge) THEN
                0307              vrelax= OBSv(i,k,bi,bj)
                0308             ELSE
                0309              vrelax=(
                0310      &            float(spongeThickness-jsl)*OBSv(i,k,bi,bj)
                0311      &            + float(jsl)*vVel(i,j,k,bi,bj) )
                0312      &            / float(spongeThickness)
                0313             ENDIF
                0314             lambda_obcs_v = (
                0315      &           float(spongeThickness-jsl)*Vrelaxobcsbound
                0316      &           + float(jsl)*Vrelaxobcsinner)
                0317      &           / float(spongeThickness)
                0318             IF (lambda_obcs_v.ne.0.) THEN
                0319              lambda_obcs_v = 1. _d 0 / lambda_obcs_v
                0320             ELSE
                0321              lambda_obcs_v = 0. _d 0
                0322             ENDIF
                0323             gV_arr(i,j) = gV_arr(i,j)
                0324      &           - _maskS(i,j,k,bi,bj) * lambda_obcs_v
                0325      &           * ( vVel(i,j,k,bi,bj) - vrelax )
                0326            ENDIF
                0327           ENDDO
96bbd4e2a5 Patr*0328          ENDIF
                0329         ENDDO
                0330        ENDIF
432b2a4585 Dimi*0331 #endif /* ALLOW_OBCS_SOUTH */
96bbd4e2a5 Patr*0332 
eb5718bf71 Jean*0333 C Eastern Open Boundary
432b2a4585 Dimi*0334 #ifdef ALLOW_OBCS_EAST
                0335        IF ( OBCSsponge_E .AND. OBCSsponge_VatEW ) THEN
                0336         DO j=jMin,jMax
                0337          IF ( OB_Ie(j,bi,bj).NE.OB_indexNone ) THEN
                0338           DO isl= 1,spongeThickness
                0339            i=OB_Ie(j,bi,bj)-isl
                0340            IF ((i.ge.imin).and.(i.le.imax)) THEN
                0341             IF (useLinearSponge) THEN
                0342              vrelax= OBEv(j,k,bi,bj)
                0343             ELSE
                0344              vrelax=(
                0345      &            float(spongeThickness-isl)*OBEv(j,k,bi,bj)
                0346      &            + float(isl)*vVel(i,j,k,bi,bj) )
                0347      &            / float(spongeThickness)
                0348             ENDIF
                0349             lambda_obcs_v = (
                0350      &           float(spongeThickness-isl)*Urelaxobcsbound
                0351      &           + float(isl)*Urelaxobcsinner)
                0352      &           / float(spongeThickness)
                0353             If (lambda_obcs_v.ne.0.) THEN
                0354              lambda_obcs_v = 1. _d 0 / lambda_obcs_v
                0355             ELSE
                0356              lambda_obcs_v = 0. _d 0
                0357             ENDIF
                0358             gV_arr(i,j) = gV_arr(i,j)
                0359      &           - _maskS(i,j,k,bi,bj) * lambda_obcs_v
                0360      &           * ( vVel(i,j,k,bi,bj) - vrelax )
                0361            ENDIF
                0362           ENDDO
96bbd4e2a5 Patr*0363          ENDIF
                0364         ENDDO
                0365        ENDIF
432b2a4585 Dimi*0366 #endif /* ALLOW_OBCS_EAST */
96bbd4e2a5 Patr*0367 
                0368 C Western Open Boundary
432b2a4585 Dimi*0369 #ifdef ALLOW_OBCS_WEST
                0370        IF ( OBCSsponge_W .AND. OBCSsponge_VatEW ) THEN
                0371         DO j=jMin,jMax
                0372          IF ( OB_Iw(j,bi,bj).NE.OB_indexNone ) THEN
                0373           DO isl= 1,spongeThickness
                0374            i=OB_Iw(j,bi,bj)+isl
                0375            IF ((i.ge.imin).and.(i.le.imax)) THEN
                0376             IF (useLinearSponge) THEN
                0377              vrelax = OBWv(j,k,bi,bj)
                0378             ELSE
                0379              vrelax=(
                0380      &            float(spongeThickness-isl)*OBWv(j,k,bi,bj)
                0381      &            + float(isl)*vVel(i,j,k,bi,bj) )
                0382      &            / float(spongeThickness)
                0383             ENDIF
                0384             lambda_obcs_v = (
                0385      &           float(spongeThickness-isl)*Urelaxobcsbound
                0386      &           + float(isl)*Urelaxobcsinner)
                0387      &           / float(spongeThickness)
                0388             IF (lambda_obcs_v.ne.0.) THEN
                0389              lambda_obcs_v = 1. _d 0 / lambda_obcs_v
                0390             ELSE
                0391              lambda_obcs_v = 0. _d 0
                0392             ENDIF
                0393             gV_arr(i,j) = gV_arr(i,j)
                0394      &           - _maskS(i,j,k,bi,bj) * lambda_obcs_v
                0395      &           * ( vVel(i,j,k,bi,bj) - vrelax )
                0396            ENDIF
                0397           ENDDO
96bbd4e2a5 Patr*0398          ENDIF
                0399         ENDDO
                0400        ENDIF
432b2a4585 Dimi*0401 #endif /* ALLOW_OBCS_WEST */
eb5718bf71 Jean*0402 
96bbd4e2a5 Patr*0403       ENDIF
                0404 
f7f6458a80 Jean*0405 #endif /* ALLOW_OBCS & ALLOW_OBCS_SPONGE */
96bbd4e2a5 Patr*0406 
                0407       RETURN
                0408       END
                0409 
eb5718bf71 Jean*0410 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
96bbd4e2a5 Patr*0411 
                0412 CStartOfInterface
                0413       SUBROUTINE OBCS_SPONGE_T(
73b1dccda0 Jean*0414      U                gT_arr,
                0415      I                iMin,iMax,jMin,jMax, k, bi, bj,
                0416      I                myTime, myIter, myThid )
eb5718bf71 Jean*0417 C     *==========================================================*
                0418 C     | S/R OBCS_SPONGE_T
                0419 C     | o Contains problem specific forcing for temperature.
                0420 C     *==========================================================*
                0421 C     | Adds a relaxation term to gT near Open-Boundaries
                0422 C     *==========================================================*
96bbd4e2a5 Patr*0423       IMPLICIT NONE
                0424 
                0425 C     == Global data ==
                0426 #include "SIZE.h"
                0427 #include "EEPARAMS.h"
                0428 #include "PARAMS.h"
                0429 #include "GRID.h"
                0430 #include "DYNVARS.h"
9b4f2a04e2 Jean*0431 #include "OBCS_PARAMS.h"
                0432 #include "OBCS_GRID.h"
                0433 #include "OBCS_FIELDS.h"
96bbd4e2a5 Patr*0434 
                0435 C     == Routine arguments ==
73b1dccda0 Jean*0436 C     gT_arr    :: the tendency array
                0437 C     iMin,iMax :: Working range of x-index for applying forcing.
                0438 C     jMin,jMax :: Working range of y-index for applying forcing.
                0439 C     k         :: Current vertical level index
                0440 C     bi,bj     :: Current tile indices
                0441       _RL gT_arr(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0442       INTEGER iMin, iMax, jMin, jMax
                0443       INTEGER k, bi, bj
eb5718bf71 Jean*0444       _RL myTime
73b1dccda0 Jean*0445       INTEGER myIter
96bbd4e2a5 Patr*0446       INTEGER myThid
                0447 CEndOfInterface
                0448 
b2a6135872 Dimi*0449 #if defined(ALLOW_OBCS) && defined(ALLOW_OBCS_SPONGE)
96bbd4e2a5 Patr*0450 C     == Local variables ==
                0451 C     Loop counters
005af54e38 Jean*0452       INTEGER i, j
                0453 #if (defined ALLOW_OBCS_EAST ) || (defined ALLOW_OBCS_WEST )
                0454       INTEGER isl
                0455 #endif
                0456 #if (defined ALLOW_OBCS_NORTH) || (defined ALLOW_OBCS_SOUTH)
                0457       INTEGER jsl
                0458 #endif
96bbd4e2a5 Patr*0459       _RL trelax, lambda_obcs_t
                0460 
74019f026d Jean*0461       IF ( useOBCSsponge .AND. spongeThickness.NE.0 ) THEN
96bbd4e2a5 Patr*0462 
                0463 C Northern Open Boundary
432b2a4585 Dimi*0464 #ifdef ALLOW_OBCS_NORTH
                0465        IF ( OBCSsponge_N .AND. OBCSsponge_Theta ) THEN
                0466         DO i=iMin,iMax
                0467          IF ( OB_Jn(i,bi,bj).NE.OB_indexNone ) THEN
                0468           DO jsl= 1,spongeThickness
                0469            j=OB_Jn(i,bi,bj)-jsl
                0470            IF ((j.ge.jmin).and.(j.le.jmax)) THEN
                0471             IF  (OBNt(i,k,bi,bj).ne. 0.d0) then
                0472              IF (useLinearSponge) THEN
                0473               trelax = OBNt(i,k,bi,bj)
                0474              ELSE
                0475               trelax=(
                0476      &             float(spongeThickness-jsl)*OBNt(i,k,bi,bj)
                0477      &             + float(jsl)*theta(i,j,k,bi,bj) )
                0478      &             / float(spongeThickness)
                0479              ENDIF
                0480              lambda_obcs_t = (
                0481      &            float(spongeThickness-jsl)*Vrelaxobcsbound
                0482      &            + float(jsl)*Vrelaxobcsinner)
                0483      &            / float(spongeThickness)
                0484              IF (lambda_obcs_t.ne.0.) THEN
                0485               lambda_obcs_t = 1. _d 0 / lambda_obcs_t
                0486              ELSE
                0487               lambda_obcs_t = 0. _d 0
                0488              ENDIF
                0489              gT_arr(i,j) = gT_arr(i,j)
                0490      &            -  maskC(i,j,k,bi,bj) * lambda_obcs_t
                0491      &            * ( theta(i,j,k,bi,bj) - trelax )
                0492             ENDIF
                0493            ENDIF
                0494           ENDDO
96bbd4e2a5 Patr*0495          ENDIF
                0496         ENDDO
                0497        ENDIF
432b2a4585 Dimi*0498 #endif /* ALLOW_OBCS_NORTH */
96bbd4e2a5 Patr*0499 
                0500 C Southern Open Boundary
432b2a4585 Dimi*0501 #ifdef ALLOW_OBCS_SOUTH
                0502        IF ( OBCSsponge_S .AND. OBCSsponge_Theta ) THEN
                0503         DO i=iMin,iMax
                0504          IF ( OB_Js(i,bi,bj).NE.OB_indexNone ) THEN
                0505           DO jsl= 1,spongeThickness
                0506            j=OB_Js(i,bi,bj)+jsl
                0507            IF ((j.ge.jmin).and.(j.le.jmax)) THEN
                0508             IF  (OBSt(i,k,bi,bj).ne. 0.d0) then
                0509              IF (useLinearSponge) THEN
                0510               trelax= OBSt(i,k,bi,bj)
                0511              ELSE
                0512               trelax=(
                0513      &             float(spongeThickness-jsl)*OBSt(i,k,bi,bj)
                0514      &             + float(jsl)*theta(i,j,k,bi,bj) )
                0515      &             / float(spongeThickness)
                0516              ENDIF
                0517              lambda_obcs_t = (
                0518      &            float(spongeThickness-jsl)*Vrelaxobcsbound
                0519      &            + float(jsl)*Vrelaxobcsinner)
                0520      &            / float(spongeThickness)
                0521              IF (lambda_obcs_t.ne.0.) THEN
                0522               lambda_obcs_t = 1. _d 0 / lambda_obcs_t
                0523              ELSE
                0524               lambda_obcs_t = 0. _d 0
                0525              ENDIF
                0526              gT_arr(i,j) = gT_arr(i,j)
                0527      &            - maskC(i,j,k,bi,bj) * lambda_obcs_t
                0528      &            * ( theta(i,j,k,bi,bj) - trelax )
                0529             ENDIF
                0530            ENDIF
                0531           ENDDO
96bbd4e2a5 Patr*0532          ENDIF
                0533         ENDDO
                0534        ENDIF
432b2a4585 Dimi*0535 #endif /* ALLOW_OBCS_SOUTH */
96bbd4e2a5 Patr*0536 
eb5718bf71 Jean*0537 C Eastern Open Boundary
432b2a4585 Dimi*0538 #ifdef ALLOW_OBCS_EAST
                0539        IF ( OBCSsponge_E .AND. OBCSsponge_Theta ) THEN
                0540         DO j=jMin,jMax
                0541          IF ( OB_Ie(j,bi,bj).NE.OB_indexNone ) THEN
                0542           DO isl= 1,spongeThickness
                0543            i=OB_Ie(j,bi,bj)-isl
                0544            IF ((i.ge.imin).and.(i.le.imax)) THEN
                0545             IF  (OBEt(j,k,bi,bj).ne. 0.d0) then
                0546              IF (useLinearSponge) THEN
                0547               trelax = OBEt(j,k,bi,bj)
                0548              ELSE
                0549               trelax=(
                0550      &             float(spongeThickness-isl)*OBEt(j,k,bi,bj)
                0551      &             + float(isl)*theta(i,j,k,bi,bj) )
                0552      &             / float(spongeThickness)
                0553              ENDIF
                0554              lambda_obcs_t = (
                0555      &            float(spongeThickness-isl)*Urelaxobcsbound
                0556      &            + float(isl)*Urelaxobcsinner)
                0557      &            / float(spongeThickness)
                0558              IF (lambda_obcs_t.ne.0.) THEN
                0559               lambda_obcs_t = 1. _d 0 / lambda_obcs_t
                0560              ELSE
                0561               lambda_obcs_t = 0. _d 0
                0562              ENDIF
                0563              gT_arr(i,j) = gT_arr(i,j)
                0564      &            - maskC(i,j,k,bi,bj) * lambda_obcs_t
                0565      &            * ( theta(i,j,k,bi,bj) - trelax )
                0566             ENDIF
                0567            ENDIF
                0568           ENDDO
96bbd4e2a5 Patr*0569          ENDIF
                0570         ENDDO
                0571        ENDIF
432b2a4585 Dimi*0572 #endif /* ALLOW_OBCS_EAST */
96bbd4e2a5 Patr*0573 
                0574 C Western Open Boundary
432b2a4585 Dimi*0575 #ifdef ALLOW_OBCS_WEST
                0576        IF ( OBCSsponge_W .AND. OBCSsponge_Theta ) THEN
                0577         DO j=jMin,jMax
                0578          IF ( OB_Iw(j,bi,bj).NE.OB_indexNone ) THEN
                0579           DO isl= 1,spongeThickness
                0580            i=OB_Iw(j,bi,bj)+isl
                0581            IF ((i.ge.imin).and.(i.le.imax)) THEN
                0582             IF  (OBWt(j,k,bi,bj).ne. 0.d0) THEN
                0583              IF (useLinearSponge) THEN
                0584               trelax= OBWt(j,k,bi,bj)
                0585              ELSE
                0586               trelax=(
                0587      &             float(spongeThickness-isl)*OBWt(j,k,bi,bj)
                0588      &             + float(isl)*theta(i,j,k,bi,bj) )
                0589      &             / float(spongeThickness)
                0590              ENDIF
                0591              lambda_obcs_t= (
                0592      &            float(spongeThickness-isl)*Urelaxobcsbound
                0593      &            + float(isl)*Urelaxobcsinner)
                0594      &            / float(spongeThickness)
                0595              IF (lambda_obcs_t .ne. 0.) THEN
                0596               lambda_obcs_t = 1. _d 0 / lambda_obcs_t
                0597              ELSE
                0598               lambda_obcs_t = 0. _d 0
                0599              ENDIF
                0600              gT_arr(i,j) = gT_arr(i,j)
                0601      &            - maskC(i,j,k,bi,bj) * lambda_obcs_t
                0602      &            * ( theta(i,j,k,bi,bj) - trelax )
                0603             ENDIF
                0604            ENDIF
                0605           ENDDO
96bbd4e2a5 Patr*0606          ENDIF
                0607         ENDDO
                0608        ENDIF
432b2a4585 Dimi*0609 #endif /* ALLOW_OBCS_WEST */
96bbd4e2a5 Patr*0610 
                0611       ENDIF
                0612 
f7f6458a80 Jean*0613 #endif /* ALLOW_OBCS & ALLOW_OBCS_SPONGE */
96bbd4e2a5 Patr*0614 
                0615       RETURN
                0616       END
                0617 
eb5718bf71 Jean*0618 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
96bbd4e2a5 Patr*0619 
                0620 CStartOfInterface
                0621       SUBROUTINE OBCS_SPONGE_S(
73b1dccda0 Jean*0622      U                gS_arr,
                0623      I                iMin,iMax,jMin,jMax, k, bi, bj,
                0624      I                myTime, myIter, myThid )
eb5718bf71 Jean*0625 C     *==========================================================*
                0626 C     | S/R OBCS_SPONGE_S
                0627 C     | o Contains problem specific forcing for salinity.
                0628 C     *==========================================================*
                0629 C     | Adds a relaxation term to gS near Open-Boundaries
                0630 C     *==========================================================*
96bbd4e2a5 Patr*0631       IMPLICIT NONE
                0632 
                0633 C     == Global data ==
                0634 #include "SIZE.h"
                0635 #include "EEPARAMS.h"
                0636 #include "PARAMS.h"
                0637 #include "GRID.h"
                0638 #include "DYNVARS.h"
9b4f2a04e2 Jean*0639 #include "OBCS_PARAMS.h"
                0640 #include "OBCS_GRID.h"
                0641 #include "OBCS_FIELDS.h"
96bbd4e2a5 Patr*0642 
                0643 C     == Routine arguments ==
73b1dccda0 Jean*0644 C     gS_arr    :: the tendency array
                0645 C     iMin,iMax :: Working range of x-index for applying forcing.
                0646 C     jMin,jMax :: Working range of y-index for applying forcing.
                0647 C     k         :: Current vertical level index
                0648 C     bi,bj     :: Current tile indices
                0649       _RL gS_arr(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0650       INTEGER iMin, iMax, jMin, jMax
                0651       INTEGER k, bi, bj
eb5718bf71 Jean*0652       _RL myTime
73b1dccda0 Jean*0653       INTEGER myIter
96bbd4e2a5 Patr*0654       INTEGER myThid
                0655 CEndOfInterface
                0656 
b2a6135872 Dimi*0657 #if defined(ALLOW_OBCS) && defined(ALLOW_OBCS_SPONGE)
96bbd4e2a5 Patr*0658 C     == Local variables ==
                0659 C     Loop counters
005af54e38 Jean*0660       INTEGER i, j
                0661 #if (defined ALLOW_OBCS_EAST ) || (defined ALLOW_OBCS_WEST )
                0662       INTEGER isl
                0663 #endif
                0664 #if (defined ALLOW_OBCS_NORTH) || (defined ALLOW_OBCS_SOUTH)
                0665       INTEGER jsl
                0666 #endif
96bbd4e2a5 Patr*0667       _RL srelax, lambda_obcs_s
                0668 
74019f026d Jean*0669       IF ( useOBCSsponge .AND. spongeThickness.NE.0 ) THEN
96bbd4e2a5 Patr*0670 
                0671 C Northern Open Boundary
432b2a4585 Dimi*0672 #ifdef ALLOW_OBCS_NORTH
                0673        IF ( OBCSsponge_N .AND. OBCSsponge_Salt ) THEN
                0674         DO i=iMin,iMax
                0675          IF ( OB_Jn(i,bi,bj).NE.OB_indexNone ) THEN
                0676           DO jsl= 1,spongeThickness
                0677            j=OB_Jn(i,bi,bj)-jsl
                0678            IF ((j.ge.jmin).and.(j.le.jmax)) THEN
                0679             IF  (OBNs(i,k,bi,bj).ne. 0.d0) then
                0680              IF (useLinearSponge) THEN
                0681               srelax= OBNs(i,k,bi,bj)
                0682              ELSE
                0683               srelax=(
                0684      &             float(spongeThickness-jsl)*OBNs(i,k,bi,bj)
                0685      &             + float(jsl)*salt(i,j,k,bi,bj) )
                0686      &             / float(spongeThickness)
                0687              ENDIF
                0688              lambda_obcs_s = (
                0689      &            float(spongeThickness-jsl)*Vrelaxobcsbound
                0690      &            + float(jsl)*Vrelaxobcsinner)
                0691      &            / float(spongeThickness)
                0692              IF (lambda_obcs_s.ne.0.) THEN
                0693               lambda_obcs_s = 1. _d 0 / lambda_obcs_s
                0694              ELSE
                0695               lambda_obcs_s = 0. _d 0
                0696              ENDIF
                0697              gS_arr(i,j) = gS_arr(i,j)
                0698      &            - maskC(i,j,k,bi,bj) * lambda_obcs_s
                0699      &            * ( salt(i,j,k,bi,bj) - srelax )
                0700             ENDIF
                0701            ENDIF
                0702           ENDDO
96bbd4e2a5 Patr*0703          ENDIF
                0704         ENDDO
                0705        ENDIF
432b2a4585 Dimi*0706 #endif /* ALLOW_OBCS_NORTH */
96bbd4e2a5 Patr*0707 
                0708 C Southern Open Boundary
432b2a4585 Dimi*0709 #ifdef ALLOW_OBCS_SOUTH
                0710        IF ( OBCSsponge_S .AND. OBCSsponge_Salt ) THEN
                0711         DO i=iMin,iMax
                0712          IF ( OB_Js(i,bi,bj).NE.OB_indexNone ) THEN
                0713           DO jsl= 1,spongeThickness
                0714            j=OB_Js(i,bi,bj)+jsl
                0715            IF ((j.ge.jmin).and.(j.le.jmax)) THEN
                0716             IF  (OBSs(i,k,bi,bj).ne. 0.d0) THEN
                0717             IF (useLinearSponge) THEN
                0718              srelax= OBSs(i,k,bi,bj)
                0719             ELSE
                0720              srelax=(
                0721      &            float(spongeThickness-jsl)*OBSs(i,k,bi,bj)
                0722      &            + float(jsl)*salt(i,j,k,bi,bj) )
                0723      &            / float(spongeThickness)
                0724             ENDIF
                0725             lambda_obcs_s = (
                0726      &           float(spongeThickness)*Vrelaxobcsbound
                0727      &           + float(jsl)*Vrelaxobcsinner)
                0728      &           / float(spongeThickness)
                0729             IF (lambda_obcs_s.ne.0.) THEN
                0730              lambda_obcs_s = 1. _d 0 / lambda_obcs_s
                0731             ELSE
                0732              lambda_obcs_s = 0. _d 0
                0733             ENDIF
                0734             gS_arr(i,j) = gS_arr(i,j)
                0735      &           - maskC(i,j,k,bi,bj) * lambda_obcs_s
                0736      &           * ( salt(i,j,k,bi,bj) - srelax )
                0737            ENDIF
                0738           ENDIF
                0739          ENDDO
                0740         ENDIF
                0741        ENDDO
                0742       ENDIF
                0743 #endif /* ALLOW_OBCS_SOUTH */
96bbd4e2a5 Patr*0744 
eb5718bf71 Jean*0745 C Eastern Open Boundary
432b2a4585 Dimi*0746 #ifdef ALLOW_OBCS_EAST
                0747       IF ( OBCSsponge_E .AND. OBCSsponge_Salt ) THEN
                0748        DO j=jMin,jMax
                0749         IF ( OB_Ie(j,bi,bj).NE.OB_indexNone ) THEN
                0750          DO isl= 1,spongeThickness
                0751           i=OB_Ie(j,bi,bj)-isl
                0752           IF ((i.ge.imin).and.(i.le.imax)) THEN
                0753            IF  (OBEs(j,k,bi,bj).ne. 0.d0) THEN
                0754             IF (useLinearSponge) THEN
                0755              srelax=  OBEs(j,k,bi,bj)
                0756             ELSE
                0757              srelax=(
                0758      &            float(spongeThickness-isl)*OBEs(j,k,bi,bj)
                0759      &            + float(isl)*salt(i,j,k,bi,bj) )
                0760      &            / float(spongeThickness)
                0761             ENDIF
                0762             lambda_obcs_s = (
                0763      &           float(spongeThickness-isl)*Urelaxobcsbound
                0764      &           + float(isl)*Urelaxobcsinner)
                0765      &           / float(spongeThickness)
                0766             IF (lambda_obcs_s.ne.0.) THEN
                0767              lambda_obcs_s = 1. _d 0 / lambda_obcs_s
                0768             ELSE
                0769              lambda_obcs_s = 0. _d 0
                0770             ENDIF
                0771             gS_arr(i,j) = gS_arr(i,j)
                0772      &           - maskC(i,j,k,bi,bj) * lambda_obcs_s
                0773      &           * ( salt(i,j,k,bi,bj) - srelax )
                0774            ENDIF
                0775           ENDIF
                0776          ENDDO
                0777         ENDIF
                0778        ENDDO
                0779       ENDIF
                0780 #endif /* ALLOW_OBCS_EAST */
96bbd4e2a5 Patr*0781 
                0782 C Western Open Boundary
432b2a4585 Dimi*0783 #ifdef ALLOW_OBCS_WEST
                0784       IF ( OBCSsponge_W .AND. OBCSsponge_Salt ) THEN
                0785        DO j=jMin,jMax
                0786         IF ( OB_Iw(j,bi,bj).NE.OB_indexNone ) THEN
                0787          DO isl= 1,spongeThickness
                0788           i=OB_Iw(j,bi,bj)+isl
                0789           IF ((i.ge.imin).and.(i.le.imax)) THEN
73b1dccda0 Jean*0790            IF  (OBWs(j,k,bi,bj).ne. 0.d0) then
432b2a4585 Dimi*0791             IF (useLinearSponge) THEN
                0792              srelax= OBWs(j,k,bi,bj)
                0793             ELSE
                0794              srelax=(
                0795      &            float(spongeThickness-isl)*OBWs(j,k,bi,bj)
                0796      &            + float(isl)*salt(i,j,k,bi,bj) )
                0797      &            / float(spongeThickness)
                0798             ENDIF
                0799             lambda_obcs_s= (
                0800      &           float(spongeThickness-isl)*Urelaxobcsbound
                0801      &           + float(isl)*Urelaxobcsinner)
                0802      &           / float(spongeThickness)
                0803             IF (lambda_obcs_s.ne.0.) THEN
                0804              lambda_obcs_s = 1. _d 0 / lambda_obcs_s
                0805             ELSE
                0806              lambda_obcs_s = 0. _d 0
                0807             ENDIF
                0808             gS_arr(i,j) = gS_arr(i,j)
                0809      &           - maskC(i,j,k,bi,bj) * lambda_obcs_s
                0810      &           * ( salt(i,j,k,bi,bj) - srelax )
                0811            ENDIF
                0812           ENDIF
                0813          ENDDO
                0814         ENDIF
                0815        ENDDO
                0816       ENDIF
                0817 #endif /* ALLOW_OBCS_WEST */
96bbd4e2a5 Patr*0818 
                0819       ENDIF
                0820 
f7f6458a80 Jean*0821 #endif /* ALLOW_OBCS & ALLOW_OBCS_SPONGE */
96bbd4e2a5 Patr*0822 
                0823       RETURN
                0824       END