Back to home page

MITgcm
 
 

    

File indexing completed on 2023-11-05 05:09:59 UTC

view on githubraw file Latest commit 65754df4 on 2023-11-04 17:55:24 UTC
1eadaea85b Jean*0001 #include "PACKAGES_CONFIG.h"
05b6d6742b Patr*0002 #include "CPP_OPTIONS.h"
517dbdc414 Jean*0003 #ifdef ALLOW_AUTODIFF
                0004 # include "AUTODIFF_OPTIONS.h"
                0005 #endif
aecc8b0f47 Mart*0006 #ifdef ALLOW_CTRL
                0007 # include "CTRL_OPTIONS.h"
                0008 #endif
05b6d6742b Patr*0009 
169b2214c5 Jean*0010 C--  File update_masks_etc.F:
                0011 C--   Contents
                0012 C--   o S/R UPDATE_MASKS_ETC
                0013 C--   o FCT SMOOTHMIN_RS( a, b )
                0014 C--   o FCT SMOOTHMIN_RL( a, b )
                0015 C--   o FCT SMOOTHABS_RS( x )
                0016 C--   o FCT SMOOTHABS_RL( x )
                0017 Cml   o S/R LIMIT_HFACC_TO_ONE
                0018 Cml   o S/R ADLIMIT_HFACC_TO_ONE
                0019 
                0020 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
05b6d6742b Patr*0021 CBOP
                0022 C     !ROUTINE: UPDATE_MASKS_ETC
                0023 C     !INTERFACE:
                0024       SUBROUTINE UPDATE_MASKS_ETC( myThid )
                0025 C     !DESCRIPTION: \bv
                0026 C     *==========================================================*
ab5a98a4ed Jean*0027 C     | SUBROUTINE UPDATE_MASKS_ETC
                0028 C     | o Re-initialise masks and topography factors after a new
                0029 C     |   hFacC has been calculated by the minimizer
05b6d6742b Patr*0030 C     *==========================================================*
65e083f882 Jean*0031 C     | These arrays are used throughout the code and describe
                0032 C     | the topography of the domain through masks (0s and 1s)
                0033 C     | and fractional height factors (0<hFac<1). The latter
                0034 C     | distinguish between the lopped-cell and full-step
                0035 C     | topographic representations.
05b6d6742b Patr*0036 C     *==========================================================*
                0037 C     | code taken from ini_masks_etc.F
                0038 C     *==========================================================*
                0039 C     \ev
                0040 
                0041 C     !USES:
                0042       IMPLICIT NONE
                0043 C     === Global variables ===
                0044 #include "SIZE.h"
                0045 #include "EEPARAMS.h"
                0046 #include "PARAMS.h"
                0047 #include "GRID.h"
                0048 #include "SURFACE.h"
                0049 Cml we need optimcycle for storing the new hFaC(C/W/S) and depth
517dbdc414 Jean*0050 #ifdef ALLOW_AUTODIFF
65754df434 Mart*0051 # include "OPTIMCYCLE.h"
65e083f882 Jean*0052 #endif
05b6d6742b Patr*0053 
                0054 C     !INPUT/OUTPUT PARAMETERS:
                0055 C     == Routine arguments ==
                0056 C     myThid -  Number of this instance of INI_MASKS_ETC
                0057       INTEGER myThid
                0058 
                0059 #ifdef ALLOW_DEPTH_CONTROL
169b2214c5 Jean*0060 C     !FUNCTIONS:
                0061       _RS SMOOTHMIN_RS
                0062       EXTERNAL SMOOTHMIN_RS
                0063 
05b6d6742b Patr*0064 C     !LOCAL VARIABLES:
                0065 C     == Local variables ==
ab5a98a4ed Jean*0066 C     bi,bj   :: Loop counters
05b6d6742b Patr*0067 C     I,J,K
ab5a98a4ed Jean*0068 C     tmpfld  :: Temporary array used to compute & write Total Depth
05b6d6742b Patr*0069       INTEGER bi, bj
ab5a98a4ed Jean*0070       INTEGER I, J, K
                0071       _RS tmpfld(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
05b6d6742b Patr*0072       CHARACTER*(MAX_LEN_MBUF) suff
                0073 Cml(
                0074       INTEGER Im1, Jm1
aecc8b0f47 Mart*0075       _RL hFacCtmp
                0076 Cml   _RL hFacCtmp2
05b6d6742b Patr*0077       _RL hFacMnSz
                0078 Cml)
                0079 CEOP
                0080 
                0081 C- Calculate lopping factor hFacC : over-estimate the part inside of the domain
                0082 C    taking into account the lower_R Boundary (Bathymetrie / Top of Atmos)
                0083       DO bj=myByLo(myThid), myByHi(myThid)
                0084        DO bi=myBxLo(myThid), myBxHi(myThid)
                0085         DO K=1, Nr
                0086          hFacMnSz=max( hFacMin, min(hFacMinDr*recip_drF(k),1. _d 0) )
169b2214c5 Jean*0087          DO J=1-OLy,sNy+OLy
                0088           DO I=1-OLx,sNx+OLx
05b6d6742b Patr*0089 C      o Non-dimensional distance between grid bound. and domain lower_R bound.
11c3150c71 Mart*0090            hFacCtmp = ( rF(K) - xx_r_low(I,J,bi,bj) )*recip_drF(K)
169b2214c5 Jean*0091 Cml           IF ( hFacCtmp .LE. 0. _d 0 ) THEN
                0092 CmlC           IF ( hFacCtmp .LT. 0.5*hfacMnSz ) THEN
05b6d6742b Patr*0093 Cml            hFacCtmp2 = 0. _d 0
                0094 Cml           ELSE
                0095 Cml            hFacCtmp2 = hFacCtmp + hFacMnSz*(
                0096 Cml     &           EXP(-hFacCtmp/hFacMnSz)-EXP(-1./hFacMnSz) )
                0097 Cml           ENDIF
169b2214c5 Jean*0098 Cml           CALL limit_hfacc_to_one( hFacCtmp2 )
05b6d6742b Patr*0099 Cml           hFacC(I,J,K,bi,bj) = hFacCtmp2
169b2214c5 Jean*0100            IF ( hFacCtmp .LE. 0. _d 0 ) THEN
                0101 C           IF ( hFacCtmp .LT. 0.5*hfacMnSz ) THEN
05b6d6742b Patr*0102             hFacC(I,J,K,bi,bj) = 0. _d 0
169b2214c5 Jean*0103            ELSEIF ( hFacCtmp .GT. 1. _d 0 ) THEN
05b6d6742b Patr*0104             hFacC(I,J,K,bi,bj) = 1. _d 0
                0105            ELSE
                0106             hFacC(I,J,K,bi,bj) = hFacCtmp + hFacMnSz*(
                0107      &           EXP(-hFacCtmp/hFacMnSz)-EXP(-1./hFacMnSz) )
                0108            ENDIF
                0109 Cml           print '(A,3I5,F20.16)', 'ml-hfac:', I,J,K,hFacC(I,J,K,bi,bj)
                0110 CmlC      o Select between, closed, open or partial (0,1,0-1)
                0111 Cml            hFacCtmp=min( max( hFacCtmp, 0. _d 0) , 1. _d 0)
                0112 CmlC      o Impose minimum fraction and/or size (dimensional)
                0113 Cml            IF (hFacCtmp.LT.hFacMnSz) THEN
                0114 Cml             IF (hFacCtmp.LT.hFacMnSz*0.5) THEN
                0115 Cml              hFacC(I,J,K,bi,bj)=0.
                0116 Cml             ELSE
                0117 Cml              hFacC(I,J,K,bi,bj)=hFacMnSz
                0118 Cml             ENDIF
                0119 Cml            ELSE
                0120 Cml             hFacC(I,J,K,bi,bj)=hFacCtmp
                0121 Cml            ENDIF
                0122 Cml           ENDIF
                0123 Cml           print '(A,F15.4,F20.16)', 'ml-hfac:', R_low(i,j,bi,bj),hFacC(I,J,K,bi,bj)
                0124           ENDDO
                0125          ENDDO
                0126         ENDDO
                0127 C - end bi,bj loops.
                0128        ENDDO
                0129       ENDDO
169b2214c5 Jean*0130 
12c8b75709 Jean*0131 C      _EXCH_XYZ_RS(hFacC,myThid)
169b2214c5 Jean*0132 
05b6d6742b Patr*0133 C-  Re-calculate lower-R Boundary position, taking into account hFacC
                0134       DO bj=myByLo(myThid), myByHi(myThid)
                0135        DO bi=myBxLo(myThid), myBxHi(myThid)
169b2214c5 Jean*0136         DO J=1-OLy,sNy+OLy
                0137          DO I=1-OLx,sNx+OLx
bbebcf4c4e Mart*0138           R_low(i,j,bi,bj) = rF(1)
                0139          ENDDO
                0140         ENDDO
                0141         DO K=Nr,1,-1
169b2214c5 Jean*0142          DO J=1-OLy,sNy+OLy
                0143           DO I=1-OLx,sNx+OLx
05b6d6742b Patr*0144            R_low(I,J,bi,bj) = R_low(I,J,bi,bj)
bbebcf4c4e Mart*0145      &                      - drF(K)*hFacC(I,J,K,bi,bj)
05b6d6742b Patr*0146           ENDDO
                0147          ENDDO
                0148         ENDDO
                0149 C - end bi,bj loops.
                0150        ENDDO
                0151       ENDDO
                0152 
                0153 Cml      DO bj=myByLo(myThid), myByHi(myThid)
                0154 Cml       DO bi=myBxLo(myThid), myBxHi(myThid)
                0155 CmlC-  Re-calculate Reference surface position, taking into account hFacC
169b2214c5 Jean*0156 Cml        DO J=1-OLy,sNy+OLy
                0157 Cml         DO I=1-OLx,sNx+OLx
05b6d6742b Patr*0158 Cml          Ro_surf(I,J,bi,bj) = R_low(I,J,bi,bj)
                0159 Cml          DO K=Nr,1,-1
                0160 Cml           Ro_surf(I,J,bi,bj) = Ro_surf(I,J,bi,bj)
                0161 Cml     &                        + drF(k)*hFacC(I,J,K,bi,bj)
                0162 Cml          ENDDO
                0163 Cml         ENDDO
                0164 Cml        ENDDO
                0165 CmlC     - end bi,bj loops.
                0166 Cml       ENDDO
                0167 Cml      ENDDO
                0168 
522c728681 Jean*0169       IF ( plotLevel.GE.debLevC ) THEN
65e083f882 Jean*0170         _BARRIER
                0171         CALL PLOT_FIELD_XYRS( R_low,
                0172      &         'Model R_low (update_masks_etc)', 1, myThid )
                0173 CML I assume that Ro_surf is not changed anywhere else in the code
                0174 CML and since it is not changed in this routine, we do not need to
05b6d6742b Patr*0175 CML print it again.
65e083f882 Jean*0176 CML     CALL PLOT_FIELD_XYRS( Ro_surf,
                0177 CML  &         'Model Ro_surf (update_masks_etc)', 1, myThid )
                0178       ENDIF
05b6d6742b Patr*0179 
                0180 C     Calculate quantities derived from XY depth map
                0181       DO bj = myByLo(myThid), myByHi(myThid)
                0182        DO bi = myBxLo(myThid), myBxHi(myThid)
169b2214c5 Jean*0183         DO j=1-OLy,sNy+OLy
                0184          DO i=1-OLx,sNx+OLx
05b6d6742b Patr*0185 C         Total fluid column thickness (r_unit) :
                0186           tmpfld(i,j,bi,bj) = Ro_surf(i,j,bi,bj) - R_low(i,j,bi,bj)
                0187 C         Inverse of fluid column thickness (1/r_unit)
                0188           IF ( tmpfld(i,j,bi,bj) .LE. 0. ) THEN
                0189            recip_Rcol(i,j,bi,bj) = 0.
                0190           ELSE
ab5a98a4ed Jean*0191            recip_Rcol(i,j,bi,bj) = 1. _d 0 / tmpfld(i,j,bi,bj)
05b6d6742b Patr*0192           ENDIF
                0193          ENDDO
                0194         ENDDO
                0195        ENDDO
                0196       ENDDO
12c8b75709 Jean*0197 C     _EXCH_XY_RS(   recip_Rcol, myThid )
05b6d6742b Patr*0198 
                0199 C     hFacW and hFacS (at U and V points)
                0200 CML   This will be the crucial part of the code, because here the minimum
                0201 CML   function MIN is involved which does not have a continuous derivative
                0202 CML   for MIN(x,y) at y=x.
                0203 CML   The thin walls representation has been moved into this loop, that is
                0204 CML   before the call to EXCH_UV_XVY_RS, because TAMC will prefer it this
65e083f882 Jean*0205 CML   way. On the other hand, this might cause difficulties in some
05b6d6742b Patr*0206 CML   configurations.
                0207       DO bj=myByLo(myThid), myByHi(myThid)
                0208        DO bi=myBxLo(myThid), myBxHi(myThid)
                0209         DO K=1, Nr
169b2214c5 Jean*0210          DO J=1-OLy,sNy+OLy
                0211           DO I=1-OLx,sNx+OLx
05b6d6742b Patr*0212            Im1=MAX(I-1,1-OLx)
                0213            Jm1=MAX(J-1,1-OLy)
                0214            IF (DYG(I,J,bi,bj).EQ.0.) THEN
                0215 C     thin walls representation of non-periodic
                0216 C     boundaries such as happen on the lat-lon grid at the N/S poles.
                0217 C     We should really supply a flag for doing this.
                0218               hFacW(I,J,K,bi,bj)=0.
                0219            ELSE
                0220               hFacW(I,J,K,bi,bj)=maskW(I,J,K,bi,bj)*
                0221 #ifdef USE_SMOOTH_MIN
169b2214c5 Jean*0222      &           SMOOTHMIN_RS(hFacC(I,J,K,bi,bj),hFacC(Im1,J,K,bi,bj))
05b6d6742b Patr*0223 #else
                0224      &                    MIN(hFacC(I,J,K,bi,bj),hFacC(Im1,J,K,bi,bj))
                0225 #endif /* USE_SMOOTH_MIN */
                0226            ENDIF
                0227            IF (DXG(I,J,bi,bj).EQ.0.) THEN
                0228               hFacS(I,J,K,bi,bj)=0.
                0229            ELSE
                0230               hFacS(I,J,K,bi,bj)=maskS(I,J,K,bi,bj)*
                0231 #ifdef USE_SMOOTH_MIN
169b2214c5 Jean*0232      &           SMOOTHMIN_RS(hFacC(I,J,K,bi,bj),hFacC(I,Jm1,K,bi,bj))
ab5a98a4ed Jean*0233 #else
05b6d6742b Patr*0234      &                    MIN(hFacC(I,J,K,bi,bj),hFacC(I,Jm1,K,bi,bj))
                0235 #endif /* USE_SMOOTH_MIN */
ab5a98a4ed Jean*0236            ENDIF
05b6d6742b Patr*0237           ENDDO
                0238          ENDDO
                0239         ENDDO
                0240        ENDDO
                0241       ENDDO
11c3150c71 Mart*0242 #if ( defined ALLOW_AUTODIFF && defined ALLOW_AUTODIFF_MONITOR )
169b2214c5 Jean*0243 C     Include call to a dummy routine. Its adjoint will be called at the proper
                0244 C     place in the adjoint code. The adjoint routine will print out adjoint
                0245 C     values if requested. The location of the call is important, it has to be
                0246 C     after the adjoint of the exchanges (DO_GTERM_BLOCKING_EXCHANGES).
05b6d6742b Patr*0247 Cml      CALL DUMMY_IN_HFAC( 'W', 0, myThid )
                0248 Cml      CALL DUMMY_IN_HFAC( 'S', 0, myThid )
                0249 #endif
                0250       CALL EXCH_UV_XYZ_RS(hFacW,hFacS,.FALSE.,myThid)
11c3150c71 Mart*0251 #if ( defined ALLOW_AUTODIFF && defined ALLOW_AUTODIFF_MONITOR )
169b2214c5 Jean*0252 C     Include call to a dummy routine. Its adjoint will be called at the proper
                0253 C     place in the adjoint code. The adjoint routine will print out adjoint
                0254 C     values if requested. The location of the call is important, it has to be
                0255 C     after the adjoint of the exchanges (DO_GTERM_BLOCKING_EXCHANGES).
05b6d6742b Patr*0256 Cml      CALL DUMMY_IN_HFAC( 'W', 1, myThid )
                0257 Cml      CALL DUMMY_IN_HFAC( 'S', 1, myThid )
                0258 #endif
                0259 
                0260 C-    Write to disk: Total Column Thickness & hFac(C,W,S):
                0261       WRITE(suff,'(I10.10)') optimcycle
                0262       CALL WRITE_FLD_XY_RS( 'Depth.',suff,tmpfld,optimcycle,myThid)
                0263       CALL WRITE_FLD_XYZ_RS( 'hFacC.',suff,hFacC,optimcycle,myThid)
                0264       CALL WRITE_FLD_XYZ_RS( 'hFacW.',suff,hFacW,optimcycle,myThid)
                0265       CALL WRITE_FLD_XYZ_RS( 'hFacS.',suff,hFacS,optimcycle,myThid)
                0266 
522c728681 Jean*0267       IF ( plotLevel.GE.debLevC ) THEN
65e083f882 Jean*0268         _BARRIER
05b6d6742b Patr*0269 C--   Write to monitor file (standard output)
65e083f882 Jean*0270         CALL PLOT_FIELD_XYZRS( hFacC,'hFacC (update_masks_etc)',
                0271      &                                          Nr, 1, myThid )
                0272         CALL PLOT_FIELD_XYZRS( hFacW,'hFacW (update_masks_etc)',
                0273      &                                          Nr, 1, myThid )
                0274         CALL PLOT_FIELD_XYZRS( hFacS,'hFacS (update_masks_etc)',
                0275      &                                          Nr, 1, myThid )
                0276       ENDIF
05b6d6742b Patr*0277 
                0278 C     Masks and reciprocals of hFac[CWS]
                0279 Cml   The masks should stay constant, so they are not recomputed at this time
                0280 Cml   implicitly implying that no cell that is wet in the begin will ever dry
65e083f882 Jean*0281 Cml   up! This is a strong constraint and should be implementent as a hard
05b6d6742b Patr*0282 Cml   inequality contraint when performing optimization (m1qn3 cannot do that)
ff02675122 Jean*0283 Cml   Also, I am assuming here that the new hFac(s) never become zero during
05b6d6742b Patr*0284 Cml   optimization!
                0285       DO bj = myByLo(myThid), myByHi(myThid)
                0286        DO bi = myBxLo(myThid), myBxHi(myThid)
                0287         DO K=1,Nr
169b2214c5 Jean*0288          DO J=1-OLy,sNy+OLy
                0289           DO I=1-OLx,sNx+OLx
05b6d6742b Patr*0290            IF (hFacC(I,J,K,bi,bj) .NE. 0. ) THEN
                0291 Cml           IF (maskC(I,J,K,bi,bj) .NE. 0. ) THEN
ab5a98a4ed Jean*0292             recip_hFacC(I,J,K,bi,bj) = 1. _d 0 / hFacC(I,J,K,bi,bj)
05b6d6742b Patr*0293 Cml            maskC(I,J,K,bi,bj) = 1.
                0294            ELSE
                0295             recip_hFacC(I,J,K,bi,bj) = 0.
                0296 Cml            maskC(I,J,K,bi,bj) = 0.
                0297            ENDIF
                0298            IF (hFacW(I,J,K,bi,bj) .NE. 0. ) THEN
                0299 Cml           IF (maskW(I,J,K,bi,bj) .NE. 0. ) THEN
ab5a98a4ed Jean*0300             recip_hFacW(I,J,K,bi,bj) = 1. _d 0 / hFacw(I,J,K,bi,bj)
05b6d6742b Patr*0301 Cml            maskW(I,J,K,bi,bj) = 1.
                0302            ELSE
                0303             recip_hFacW(I,J,K,bi,bj) = 0.
                0304 Cml            maskW(I,J,K,bi,bj) = 0.
                0305            ENDIF
                0306            IF (hFacS(I,J,K,bi,bj) .NE. 0. ) THEN
                0307 Cml           IF (maskS(I,J,K,bi,bj) .NE. 0. ) THEN
ab5a98a4ed Jean*0308             recip_hFacS(I,J,K,bi,bj) = 1. _d 0 / hFacS(I,J,K,bi,bj)
05b6d6742b Patr*0309 Cml            maskS(I,J,K,bi,bj) = 1.
                0310            ELSE
                0311             recip_hFacS(I,J,K,bi,bj) = 0.
                0312 Cml            maskS(I,J,K,bi,bj) = 0.
                0313            ENDIF
                0314           ENDDO
                0315          ENDDO
                0316         ENDDO
                0317 CmlCml(
                0318 Cml       ENDDO
                0319 Cml      ENDDO
12c8b75709 Jean*0320 Cml      _EXCH_XYZ_RS(recip_hFacC    , myThid )
                0321 Cml      _EXCH_XYZ_RS(recip_hFacW    , myThid )
                0322 Cml      _EXCH_XYZ_RS(recip_hFacS    , myThid )
                0323 Cml      _EXCH_XYZ_RS(maskC    , myThid )
                0324 Cml      _EXCH_XYZ_RS(maskW    , myThid )
                0325 Cml      _EXCH_XYZ_RS(maskS    , myThid )
05b6d6742b Patr*0326 Cml      DO bj = myByLo(myThid), myByHi(myThid)
                0327 Cml       DO bi = myBxLo(myThid), myBxHi(myThid)
                0328 CmlCml)
169b2214c5 Jean*0329 #ifdef NONLIN_FRSURF
                0330 C--   Save initial geometrical hFac factor into h0Fac (fixed in time):
                0331 C     Note: In case 1 pkg modifies hFac (from packages_init_fixed, called
                0332 C     later in sequence of calls) this pkg would need also to update h0Fac.
                0333         DO k=1,Nr
                0334          DO j=1-OLy,sNy+OLy
                0335           DO i=1-OLx,sNx+OLx
                0336            h0FacC(i,j,k,bi,bj) = _hFacC(i,j,k,bi,bj)
                0337            h0FacW(i,j,k,bi,bj) = _hFacW(i,j,k,bi,bj)
                0338            h0FacS(i,j,k,bi,bj) = _hFacS(i,j,k,bi,bj)
05b6d6742b Patr*0339           ENDDO
                0340          ENDDO
                0341         ENDDO
169b2214c5 Jean*0342 #endif /* NONLIN_FRSURF */
05b6d6742b Patr*0343 C - end bi,bj loops.
                0344        ENDDO
                0345       ENDDO
ab5a98a4ed Jean*0346 
05b6d6742b Patr*0347 #endif /* ALLOW_DEPTH_CONTROL */
                0348       RETURN
                0349       END
                0350 
                0351 #ifdef USE_SMOOTH_MIN
169b2214c5 Jean*0352 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0353 
                0354       _RS FUNCTION SMOOTHMIN_RS( a, b )
05b6d6742b Patr*0355 
169b2214c5 Jean*0356       IMPLICIT NONE
05b6d6742b Patr*0357 
                0358       _RS a, b
                0359 
169b2214c5 Jean*0360       _RS SMOOTHABS_RS
                0361       EXTERNAL SMOOTHABS_RS
05b6d6742b Patr*0362 
                0363 Cml      smoothMin_R4 = .5*(a+b)
169b2214c5 Jean*0364       SMOOTHMIN_RS = .5*( a+b - SMOOTHABS_RS(a-b) )
05b6d6742b Patr*0365 CML      smoothMin_R4 = MIN(a,b)
                0366 
169b2214c5 Jean*0367       RETURN
                0368       END
05b6d6742b Patr*0369 
169b2214c5 Jean*0370       _RL FUNCTION SMOOTHMIN_RL( a, b )
05b6d6742b Patr*0371 
169b2214c5 Jean*0372       IMPLICIT NONE
05b6d6742b Patr*0373 
                0374       _RL a, b
                0375 
169b2214c5 Jean*0376       _RL SMOOTHABS_RL
                0377       EXTERNAL SMOOTHABS_RL
05b6d6742b Patr*0378 
                0379 Cml      smoothMin_R8 = .5*(a+b)
169b2214c5 Jean*0380       SMOOTHMIN_RL = .5*( a+b - SMOOTHABS_RL(a-b) )
05b6d6742b Patr*0381 Cml      smoothMin_R8 = MIN(a,b)
                0382 
169b2214c5 Jean*0383       RETURN
                0384       END
05b6d6742b Patr*0385 
169b2214c5 Jean*0386       _RS FUNCTION SMOOTHABS_RS( x )
ab5a98a4ed Jean*0387 
169b2214c5 Jean*0388       IMPLICIT NONE
05b6d6742b Patr*0389 C     === Global variables ===
                0390 #include "SIZE.h"
                0391 #include "EEPARAMS.h"
                0392 #include "PARAMS.h"
                0393 C     input parameter
                0394       _RS x
                0395 c     local variable
                0396       _RS sf, rsf
                0397 
169b2214c5 Jean*0398       IF ( smoothAbsFuncRange .LT. 0.0 ) THEN
05b6d6742b Patr*0399 c     limit of smoothMin(a,b) = .5*(a+b)
169b2214c5 Jean*0400          SMOOTHABS_RS = 0.
                0401       ELSE
                0402          IF ( smoothAbsFuncRange .NE. 0.0 ) THEN
05b6d6742b Patr*0403             sf  = 10.0/smoothAbsFuncRange
                0404             rsf = 1./sf
169b2214c5 Jean*0405          ELSE
05b6d6742b Patr*0406 c     limit of smoothMin(a,b) = min(a,b)
                0407             sf  = 0.
                0408             rsf = 0.
169b2214c5 Jean*0409          ENDIF
05b6d6742b Patr*0410 c
169b2214c5 Jean*0411          IF ( x .GT. smoothAbsFuncRange ) THEN
                0412             SMOOTHABS_RS = x
                0413          ELSEIF ( x .LT. -smoothAbsFuncRange ) THEN
                0414             SMOOTHABS_RS = -x
                0415          ELSE
                0416             SMOOTHABS_RS = log(.5*(exp(x*sf)+exp(-x*sf)))*rsf
                0417          ENDIF
                0418       ENDIF
05b6d6742b Patr*0419 
169b2214c5 Jean*0420       RETURN
                0421       END
05b6d6742b Patr*0422 
169b2214c5 Jean*0423       _RL FUNCTION SMOOTHABS_RL( x )
ab5a98a4ed Jean*0424 
169b2214c5 Jean*0425       IMPLICIT NONE
05b6d6742b Patr*0426 C     === Global variables ===
                0427 #include "SIZE.h"
                0428 #include "EEPARAMS.h"
                0429 #include "PARAMS.h"
                0430 C     input parameter
                0431       _RL x
                0432 c     local variable
                0433       _RL sf, rsf
                0434 
169b2214c5 Jean*0435       IF ( smoothAbsFuncRange .LT. 0.0 ) THEN
05b6d6742b Patr*0436 c     limit of smoothMin(a,b) = .5*(a+b)
169b2214c5 Jean*0437          SMOOTHABS_RL = 0.
                0438       ELSE
                0439          IF ( smoothAbsFuncRange .NE. 0.0 ) THEN
05b6d6742b Patr*0440             sf  = 10.0D0/smoothAbsFuncRange
                0441             rsf = 1.D0/sf
169b2214c5 Jean*0442          ELSE
05b6d6742b Patr*0443 c     limit of smoothMin(a,b) = min(a,b)
                0444             sf  = 0.D0
                0445             rsf = 0.D0
169b2214c5 Jean*0446          ENDIF
ab5a98a4ed Jean*0447 c
169b2214c5 Jean*0448          IF ( x .GE. smoothAbsFuncRange ) THEN
                0449             SMOOTHABS_RL = x
                0450          ELSEIF ( x .LE. -smoothAbsFuncRange ) THEN
                0451             SMOOTHABS_RL = -x
                0452          ELSE
                0453             SMOOTHABS_RL = log(.5*(exp(x*sf)+exp(-x*sf)))*rsf
                0454          ENDIF
                0455       ENDIF
                0456 
                0457       RETURN
                0458       END
05b6d6742b Patr*0459 #endif /* USE_SMOOTH_MIN */
                0460 
                0461 Cml#ifdef ALLOW_DEPTH_CONTROL
                0462 Cmlcadj SUBROUTINE limit_hfacc_to_one INPUT   = 1
                0463 Cmlcadj SUBROUTINE limit_hfacc_to_one OUTPUT  = 1
                0464 Cmlcadj SUBROUTINE limit_hfacc_to_one ACTIVE  = 1
                0465 Cmlcadj SUBROUTINE limit_hfacc_to_one DEPEND  = 1
                0466 Cmlcadj SUBROUTINE limit_hfacc_to_one REQUIRED
                0467 Cmlcadj SUBROUTINE limit_hfacc_to_one ADNAME  = adlimit_hfacc_to_one
                0468 Cml#endif /* ALLOW_DEPTH_CONTROL */
169b2214c5 Jean*0469 Cml      SUBROUTINE LIMIT_HFACC_TO_ONE( hf )
05b6d6742b Patr*0470 Cml
                0471 Cml      _RL hf
ab5a98a4ed Jean*0472 Cml
169b2214c5 Jean*0473 Cml      IF ( hf .GT. 1. _d 0 ) THEN
05b6d6742b Patr*0474 Cml       hf = 1. _d 0
169b2214c5 Jean*0475 Cml      ENDIF
05b6d6742b Patr*0476 Cml
169b2214c5 Jean*0477 Cml      RETURN
                0478 Cml      END
05b6d6742b Patr*0479 Cml
169b2214c5 Jean*0480 Cml      SUBROUTINE ADLIMIT_HFACC_TO_ONE( hf, adhf )
05b6d6742b Patr*0481 Cml
                0482 Cml      _RL hf, adhf
ab5a98a4ed Jean*0483 Cml
169b2214c5 Jean*0484 Cml      RETURN
                0485 Cml      END
05b6d6742b Patr*0486 
                0487 #ifdef ALLOW_DEPTH_CONTROL
                0488 cadj SUBROUTINE dummy_in_hfac INPUT   = 1, 2, 3
65e083f882 Jean*0489 cadj SUBROUTINE dummy_in_hfac OUTPUT  =
                0490 cadj SUBROUTINE dummy_in_hfac ACTIVE  =
05b6d6742b Patr*0491 cadj SUBROUTINE dummy_in_hfac DEPEND  = 1, 2, 3
                0492 cadj SUBROUTINE dummy_in_hfac REQUIRED
                0493 cadj SUBROUTINE dummy_in_hfac INFLUENCED
                0494 cadj SUBROUTINE dummy_in_hfac ADNAME  = addummy_in_hfac
                0495 cadj SUBROUTINE dummy_in_hfac FTLNAME = g_dummy_in_hfac
                0496 #endif /* ALLOW_DEPTH_CONTROL */
[ 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!