/[MITgcm]/MITgcm/model/src/update_masks_etc.F
ViewVC logotype

Diff of /MITgcm/model/src/update_masks_etc.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.7 by mlosch, Mon Mar 19 14:32:47 2012 UTC revision 1.8 by jmc, Sun Jun 17 14:17:26 2012 UTC
# Line 3  C $Name$ Line 3  C $Name$
3    
4  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
5    
6    C--  File update_masks_etc.F:
7    C--   Contents
8    C--   o S/R UPDATE_MASKS_ETC
9    C--   o FCT SMOOTHMIN_RS( a, b )
10    C--   o FCT SMOOTHMIN_RL( a, b )
11    C--   o FCT SMOOTHABS_RS( x )
12    C--   o FCT SMOOTHABS_RL( x )
13    Cml   o S/R LIMIT_HFACC_TO_ONE
14    Cml   o S/R ADLIMIT_HFACC_TO_ONE
15    
16    C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
17  CBOP  CBOP
18  C     !ROUTINE: UPDATE_MASKS_ETC  C     !ROUTINE: UPDATE_MASKS_ETC
19  C     !INTERFACE:  C     !INTERFACE:
# Line 42  C     myThid -  Number of this instance Line 53  C     myThid -  Number of this instance
53        INTEGER myThid        INTEGER myThid
54    
55  #ifdef ALLOW_DEPTH_CONTROL  #ifdef ALLOW_DEPTH_CONTROL
56    C     !FUNCTIONS:
57          _RS SMOOTHMIN_RS
58          EXTERNAL SMOOTHMIN_RS
59    
60  C     !LOCAL VARIABLES:  C     !LOCAL VARIABLES:
61  C     == Local variables ==  C     == Local variables ==
62  C     bi,bj   :: Loop counters  C     bi,bj   :: Loop counters
# Line 55  Cml( Line 70  Cml(
70        INTEGER Im1, Jm1        INTEGER Im1, Jm1
71        _RL hFacCtmp, hFacCtmp2        _RL hFacCtmp, hFacCtmp2
72        _RL hFacMnSz        _RL hFacMnSz
       _RS smoothMin_R4  
       EXTERNAL smoothMin_R4  
73  Cml)  Cml)
74  CEOP  CEOP
75    
# Line 66  C    taking into account the lower_R Bou Line 79  C    taking into account the lower_R Bou
79         DO bi=myBxLo(myThid), myBxHi(myThid)         DO bi=myBxLo(myThid), myBxHi(myThid)
80          DO K=1, Nr          DO K=1, Nr
81           hFacMnSz=max( hFacMin, min(hFacMinDr*recip_drF(k),1. _d 0) )           hFacMnSz=max( hFacMin, min(hFacMinDr*recip_drF(k),1. _d 0) )
82           DO J=1-Oly,sNy+Oly           DO J=1-OLy,sNy+OLy
83            DO I=1-Olx,sNx+Olx            DO I=1-OLx,sNx+OLx
84  C      o Non-dimensional distance between grid bound. and domain lower_R bound.  C      o Non-dimensional distance between grid bound. and domain lower_R bound.
85  #ifdef ALLOW_DEPTH_CONTROL  #ifdef ALLOW_DEPTH_CONTROL
86             hFacCtmp = (rF(K)-xx_r_low(I,J,bi,bj))*recip_drF(K)             hFacCtmp = (rF(K)-xx_r_low(I,J,bi,bj))*recip_drF(K)
87  #else  #else
88             hFacCtmp = (rF(K)-R_low(I,J,bi,bj))*recip_drF(K)             hFacCtmp = (rF(K)-R_low(I,J,bi,bj))*recip_drF(K)
89  #endif /* ALLOW_DEPTH_CONTROL */  #endif /* ALLOW_DEPTH_CONTROL */
90  Cml           IF ( hFacCtmp .le. 0. _d 0 ) THEN  Cml           IF ( hFacCtmp .LE. 0. _d 0 ) THEN
91  CmlC           IF ( hFacCtmp .lt. 0.5*hfacMnSz ) THEN  CmlC           IF ( hFacCtmp .LT. 0.5*hfacMnSz ) THEN
92  Cml            hFacCtmp2 = 0. _d 0  Cml            hFacCtmp2 = 0. _d 0
93  Cml           ELSE  Cml           ELSE
94  Cml            hFacCtmp2 = hFacCtmp + hFacMnSz*(  Cml            hFacCtmp2 = hFacCtmp + hFacMnSz*(
95  Cml     &           EXP(-hFacCtmp/hFacMnSz)-EXP(-1./hFacMnSz) )  Cml     &           EXP(-hFacCtmp/hFacMnSz)-EXP(-1./hFacMnSz) )
96  Cml           ENDIF  Cml           ENDIF
97  Cml           call limit_hfacc_to_one( hFacCtmp2 )  Cml           CALL limit_hfacc_to_one( hFacCtmp2 )
98  Cml           hFacC(I,J,K,bi,bj) = hFacCtmp2  Cml           hFacC(I,J,K,bi,bj) = hFacCtmp2
99             IF ( hFacCtmp .le. 0. _d 0 ) THEN             IF ( hFacCtmp .LE. 0. _d 0 ) THEN
100  C           IF ( hFacCtmp .lt. 0.5*hfacMnSz ) THEN  C           IF ( hFacCtmp .LT. 0.5*hfacMnSz ) THEN
101              hFacC(I,J,K,bi,bj) = 0. _d 0              hFacC(I,J,K,bi,bj) = 0. _d 0
102             ELSEIF ( hFacCtmp .gt. 1. _d 0 ) THEN             ELSEIF ( hFacCtmp .GT. 1. _d 0 ) THEN
103              hFacC(I,J,K,bi,bj) = 1. _d 0              hFacC(I,J,K,bi,bj) = 1. _d 0
104             ELSE             ELSE
105              hFacC(I,J,K,bi,bj) = hFacCtmp + hFacMnSz*(              hFacC(I,J,K,bi,bj) = hFacCtmp + hFacMnSz*(
# Line 113  Cml           print '(A,F15.4,F20.16)', Line 126  Cml           print '(A,F15.4,F20.16)',
126  C - end bi,bj loops.  C - end bi,bj loops.
127         ENDDO         ENDDO
128        ENDDO        ENDDO
129  C  
130  C      _EXCH_XYZ_RS(hFacC,myThid)  C      _EXCH_XYZ_RS(hFacC,myThid)
131  C  
132  C-  Re-calculate lower-R Boundary position, taking into account hFacC  C-  Re-calculate lower-R Boundary position, taking into account hFacC
133        DO bj=myByLo(myThid), myByHi(myThid)        DO bj=myByLo(myThid), myByHi(myThid)
134         DO bi=myBxLo(myThid), myBxHi(myThid)         DO bi=myBxLo(myThid), myBxHi(myThid)
135          DO J=1-Oly,sNy+Oly          DO J=1-OLy,sNy+OLy
136           DO I=1-Olx,sNx+Olx           DO I=1-OLx,sNx+OLx
137            R_low(i,j,bi,bj) = rF(1)            R_low(i,j,bi,bj) = rF(1)
138           ENDDO           ENDDO
139          ENDDO          ENDDO
140          DO K=Nr,1,-1          DO K=Nr,1,-1
141           DO J=1-Oly,sNy+Oly           DO J=1-OLy,sNy+OLy
142            DO I=1-Olx,sNx+Olx            DO I=1-OLx,sNx+OLx
143             R_low(I,J,bi,bj) = R_low(I,J,bi,bj)             R_low(I,J,bi,bj) = R_low(I,J,bi,bj)
144       &                      - drF(K)*hFacC(I,J,K,bi,bj)       &                      - drF(K)*hFacC(I,J,K,bi,bj)
145            ENDDO            ENDDO
# Line 135  C-  Re-calculate lower-R Boundary positi Line 148  C-  Re-calculate lower-R Boundary positi
148  C - end bi,bj loops.  C - end bi,bj loops.
149         ENDDO         ENDDO
150        ENDDO        ENDDO
 C  
151    
152  Cml      DO bj=myByLo(myThid), myByHi(myThid)  Cml      DO bj=myByLo(myThid), myByHi(myThid)
153  Cml       DO bi=myBxLo(myThid), myBxHi(myThid)  Cml       DO bi=myBxLo(myThid), myBxHi(myThid)
154  CmlC-  Re-calculate Reference surface position, taking into account hFacC  CmlC-  Re-calculate Reference surface position, taking into account hFacC
155  CmlC   initialize Total column fluid thickness and surface k index  Cml        DO J=1-OLy,sNy+OLy
156  CmlC       Note: if no fluid (continent) ==> ksurf = Nr+1  Cml         DO I=1-OLx,sNx+OLx
 Cml        DO J=1-Oly,sNy+Oly  
 Cml         DO I=1-Olx,sNx+Olx  
 Cml          tmpfld(I,J,bi,bj) = 0.  
 Cml          ksurfC(I,J,bi,bj) = Nr+1  
157  Cml          Ro_surf(I,J,bi,bj) = R_low(I,J,bi,bj)  Cml          Ro_surf(I,J,bi,bj) = R_low(I,J,bi,bj)
158  Cml          DO K=Nr,1,-1  Cml          DO K=Nr,1,-1
159  Cml           Ro_surf(I,J,bi,bj) = Ro_surf(I,J,bi,bj)  Cml           Ro_surf(I,J,bi,bj) = Ro_surf(I,J,bi,bj)
160  Cml     &                        + drF(k)*hFacC(I,J,K,bi,bj)  Cml     &                        + drF(k)*hFacC(I,J,K,bi,bj)
 Cml           IF (maskC(I,J,K,bi,bj).NE.0.) THEN  
 Cml            ksurfC(I,J,bi,bj) = k  
 Cml            tmpfld(i,j,bi,bj) = tmpfld(i,j,bi,bj) + 1.  
 Cml           ENDIF  
161  Cml          ENDDO  Cml          ENDDO
162  Cml         ENDDO  Cml         ENDDO
163  Cml        ENDDO  Cml        ENDDO
# Line 175  CML  &         'Model Ro_surf (update_ma Line 179  CML  &         'Model Ro_surf (update_ma
179  C     Calculate quantities derived from XY depth map  C     Calculate quantities derived from XY depth map
180        DO bj = myByLo(myThid), myByHi(myThid)        DO bj = myByLo(myThid), myByHi(myThid)
181         DO bi = myBxLo(myThid), myBxHi(myThid)         DO bi = myBxLo(myThid), myBxHi(myThid)
182          DO j=1-Oly,sNy+Oly          DO j=1-OLy,sNy+OLy
183           DO i=1-Olx,sNx+Olx           DO i=1-OLx,sNx+OLx
184  C         Total fluid column thickness (r_unit) :  C         Total fluid column thickness (r_unit) :
185            tmpfld(i,j,bi,bj) = Ro_surf(i,j,bi,bj) - R_low(i,j,bi,bj)            tmpfld(i,j,bi,bj) = Ro_surf(i,j,bi,bj) - R_low(i,j,bi,bj)
186  C         Inverse of fluid column thickness (1/r_unit)  C         Inverse of fluid column thickness (1/r_unit)
# Line 202  CML   configurations. Line 206  CML   configurations.
206        DO bj=myByLo(myThid), myByHi(myThid)        DO bj=myByLo(myThid), myByHi(myThid)
207         DO bi=myBxLo(myThid), myBxHi(myThid)         DO bi=myBxLo(myThid), myBxHi(myThid)
208          DO K=1, Nr          DO K=1, Nr
209  CML         DO J=1-Oly+1,sNy+Oly           DO J=1-OLy,sNy+OLy
210  CML          DO I=1-Olx+1,sNx+Olx            DO I=1-OLx,sNx+OLx
 CML         DO J=1,sNy+1  
 CML          DO I=1,sNx+1  
          DO J=1-Oly,sNy+Oly  
           DO I=1-Olx,sNx+Olx  
211             Im1=MAX(I-1,1-OLx)             Im1=MAX(I-1,1-OLx)
212             Jm1=MAX(J-1,1-OLy)             Jm1=MAX(J-1,1-OLy)
213             IF (DYG(I,J,bi,bj).EQ.0.) THEN             IF (DYG(I,J,bi,bj).EQ.0.) THEN
# Line 216  C     boundaries such as happen on the l Line 216  C     boundaries such as happen on the l
216  C     We should really supply a flag for doing this.  C     We should really supply a flag for doing this.
217                hFacW(I,J,K,bi,bj)=0.                hFacW(I,J,K,bi,bj)=0.
218             ELSE             ELSE
 Cml              hFacW(I,J,K,bi,bj)=  
219                hFacW(I,J,K,bi,bj)=maskW(I,J,K,bi,bj)*                hFacW(I,J,K,bi,bj)=maskW(I,J,K,bi,bj)*
220  #ifdef USE_SMOOTH_MIN  #ifdef USE_SMOOTH_MIN
221       &           smoothMin_R4(hFacC(I,J,K,bi,bj),hFacC(Im1,J,K,bi,bj))       &           SMOOTHMIN_RS(hFacC(I,J,K,bi,bj),hFacC(Im1,J,K,bi,bj))
222  #else  #else
223       &                    MIN(hFacC(I,J,K,bi,bj),hFacC(Im1,J,K,bi,bj))       &                    MIN(hFacC(I,J,K,bi,bj),hFacC(Im1,J,K,bi,bj))
224  #endif /* USE_SMOOTH_MIN */  #endif /* USE_SMOOTH_MIN */
# Line 227  Cml              hFacW(I,J,K,bi,bj)= Line 226  Cml              hFacW(I,J,K,bi,bj)=
226             IF (DXG(I,J,bi,bj).EQ.0.) THEN             IF (DXG(I,J,bi,bj).EQ.0.) THEN
227                hFacS(I,J,K,bi,bj)=0.                hFacS(I,J,K,bi,bj)=0.
228             ELSE             ELSE
 Cml              hFacS(I,J,K,bi,bj)=  
229                hFacS(I,J,K,bi,bj)=maskS(I,J,K,bi,bj)*                hFacS(I,J,K,bi,bj)=maskS(I,J,K,bi,bj)*
230  #ifdef USE_SMOOTH_MIN  #ifdef USE_SMOOTH_MIN
231       &           smoothMin_R4(hFacC(I,J,K,bi,bj),hFacC(I,Jm1,K,bi,bj))       &           SMOOTHMIN_RS(hFacC(I,J,K,bi,bj),hFacC(I,Jm1,K,bi,bj))
232  #else  #else
233       &                    MIN(hFacC(I,J,K,bi,bj),hFacC(I,Jm1,K,bi,bj))       &                    MIN(hFacC(I,J,K,bi,bj),hFacC(I,Jm1,K,bi,bj))
234  #endif /* USE_SMOOTH_MIN */  #endif /* USE_SMOOTH_MIN */
# Line 240  Cml              hFacS(I,J,K,bi,bj)= Line 238  Cml              hFacS(I,J,K,bi,bj)=
238          ENDDO          ENDDO
239         ENDDO         ENDDO
240        ENDDO        ENDDO
241  #if (defined (ALLOW_AUTODIFF_TAMC) && \  #if ( defined (ALLOW_AUTODIFF_TAMC) && \
242       defined (ALLOW_AUTODIFF_MONITOR) && \        defined (ALLOW_AUTODIFF_MONITOR) && \
243       defined (ALLOW_DEPTH_CONTROL))        defined (ALLOW_DEPTH_CONTROL) )
244  C     Include call to a dummy routine. Its adjoint will be  C     Include call to a dummy routine. Its adjoint will be called at the proper
245  C     called at the proper place in the adjoint code.  C     place in the adjoint code. The adjoint routine will print out adjoint
246  C     The adjoint routine will print out adjoint values  C     values if requested. The location of the call is important, it has to be
247  C     if requested. The location of the call is important,  C     after the adjoint of the exchanges (DO_GTERM_BLOCKING_EXCHANGES).
 C     it has to be after the adjoint of the exchanges  
 C     (DO_GTERM_BLOCKING_EXCHANGES).  
248  Cml      CALL DUMMY_IN_HFAC( 'W', 0, myThid )  Cml      CALL DUMMY_IN_HFAC( 'W', 0, myThid )
249  Cml      CALL DUMMY_IN_HFAC( 'S', 0, myThid )  Cml      CALL DUMMY_IN_HFAC( 'S', 0, myThid )
250  #endif  #endif
 Cml      CALL EXCH_UV_XYZ_RL(hFacW,hFacS,.FALSE.,myThid)  
251        CALL EXCH_UV_XYZ_RS(hFacW,hFacS,.FALSE.,myThid)        CALL EXCH_UV_XYZ_RS(hFacW,hFacS,.FALSE.,myThid)
252  #if (defined (ALLOW_AUTODIFF_TAMC) && \  #if ( defined (ALLOW_AUTODIFF_TAMC) && \
253       defined (ALLOW_AUTODIFF_MONITOR) && \        defined (ALLOW_AUTODIFF_MONITOR) && \
254       defined (ALLOW_DEPTH_CONTROL))        defined (ALLOW_DEPTH_CONTROL) )
255  C     Include call to a dummy routine. Its adjoint will be  C     Include call to a dummy routine. Its adjoint will be called at the proper
256  C     called at the proper place in the adjoint code.  C     place in the adjoint code. The adjoint routine will print out adjoint
257  C     The adjoint routine will print out adjoint values  C     values if requested. The location of the call is important, it has to be
258  C     if requested. The location of the call is important,  C     after the adjoint of the exchanges (DO_GTERM_BLOCKING_EXCHANGES).
 C     it has to be after the adjoint of the exchanges  
 C     (DO_GTERM_BLOCKING_EXCHANGES).  
259  Cml      CALL DUMMY_IN_HFAC( 'W', 1, myThid )  Cml      CALL DUMMY_IN_HFAC( 'W', 1, myThid )
260  Cml      CALL DUMMY_IN_HFAC( 'S', 1, myThid )  Cml      CALL DUMMY_IN_HFAC( 'S', 1, myThid )
261  #endif  #endif
# Line 295  Cml   optimization! Line 288  Cml   optimization!
288        DO bj = myByLo(myThid), myByHi(myThid)        DO bj = myByLo(myThid), myByHi(myThid)
289         DO bi = myBxLo(myThid), myBxHi(myThid)         DO bi = myBxLo(myThid), myBxHi(myThid)
290          DO K=1,Nr          DO K=1,Nr
291           DO J=1-Oly,sNy+Oly           DO J=1-OLy,sNy+OLy
292            DO I=1-Olx,sNx+Olx            DO I=1-OLx,sNx+OLx
293             IF (hFacC(I,J,K,bi,bj) .NE. 0. ) THEN             IF (hFacC(I,J,K,bi,bj) .NE. 0. ) THEN
294  Cml           IF (maskC(I,J,K,bi,bj) .NE. 0. ) THEN  Cml           IF (maskC(I,J,K,bi,bj) .NE. 0. ) THEN
295              recip_hFacC(I,J,K,bi,bj) = 1. _d 0 / hFacC(I,J,K,bi,bj)              recip_hFacC(I,J,K,bi,bj) = 1. _d 0 / hFacC(I,J,K,bi,bj)
# Line 336  Cml      _EXCH_XYZ_RS(maskS    , myThid Line 329  Cml      _EXCH_XYZ_RS(maskS    , myThid
329  Cml      DO bj = myByLo(myThid), myByHi(myThid)  Cml      DO bj = myByLo(myThid), myByHi(myThid)
330  Cml       DO bi = myBxLo(myThid), myBxHi(myThid)  Cml       DO bi = myBxLo(myThid), myBxHi(myThid)
331  CmlCml)  CmlCml)
332  C-    Calculate surface k index for interface W & S (U & V points)  #ifdef NONLIN_FRSURF
333          DO J=1-Oly,sNy+Oly  C--   Save initial geometrical hFac factor into h0Fac (fixed in time):
334           DO I=1-Olx,sNx+Olx  C     Note: In case 1 pkg modifies hFac (from packages_init_fixed, called
335            ksurfW(I,J,bi,bj) = Nr+1  C     later in sequence of calls) this pkg would need also to update h0Fac.
336            ksurfS(I,J,bi,bj) = Nr+1          DO k=1,Nr
337            DO k=Nr,1,-1           DO j=1-OLy,sNy+OLy
338  Cml           IF (hFacW(I,J,K,bi,bj).NE.0.) THEN            DO i=1-OLx,sNx+OLx
339             IF (maskW(I,J,K,bi,bj).NE.0.) THEN             h0FacC(i,j,k,bi,bj) = _hFacC(i,j,k,bi,bj)
340                ksurfW(I,J,bi,bj) = k             h0FacW(i,j,k,bi,bj) = _hFacW(i,j,k,bi,bj)
341             ENDIF             h0FacS(i,j,k,bi,bj) = _hFacS(i,j,k,bi,bj)
 Cml           IF (hFacS(I,J,K,bi,bj).NE.0.) THEN  
            IF (maskS(I,J,K,bi,bj).NE.0.) THEN  
               ksurfS(I,J,bi,bj) = k  
   
            ENDIF  
342            ENDDO            ENDDO
343           ENDDO           ENDDO
344          ENDDO          ENDDO
345    #endif /* NONLIN_FRSURF */
346  C - end bi,bj loops.  C - end bi,bj loops.
347         ENDDO         ENDDO
348        ENDDO        ENDDO
349    
 c #ifdef ALLOW_NONHYDROSTATIC  
 C--   Calculate "recip_hFacU" = reciprocal hfac distance/volume for W cells  
 C     not used ; computed locally in CALC_GW  
 c #endif  
   
350  #endif /* ALLOW_DEPTH_CONTROL */  #endif /* ALLOW_DEPTH_CONTROL */
351        RETURN        RETURN
352        END        END
353    
354  #ifdef USE_SMOOTH_MIN  #ifdef USE_SMOOTH_MIN
355        _RS function smoothMin_R4( a, b )  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
356    
357          _RS FUNCTION SMOOTHMIN_RS( a, b )
358    
359        implicit none        IMPLICIT NONE
360    
361        _RS a, b        _RS a, b
362    
363        _RS smoothAbs_R4        _RS SMOOTHABS_RS
364        external smoothAbs_R4        EXTERNAL SMOOTHABS_RS
365    
366  Cml      smoothMin_R4 = .5*(a+b)  Cml      smoothMin_R4 = .5*(a+b)
367        smoothMin_R4 = .5*( a+b - smoothAbs_R4(a-b) )        SMOOTHMIN_RS = .5*( a+b - SMOOTHABS_RS(a-b) )
368  CML      smoothMin_R4 = MIN(a,b)  CML      smoothMin_R4 = MIN(a,b)
369    
370        return        RETURN
371        end        END
372    
373        _RL function smoothMin_R8( a, b )        _RL FUNCTION SMOOTHMIN_RL( a, b )
374    
375        implicit none        IMPLICIT NONE
376    
377        _RL a, b        _RL a, b
378    
379        _RL smoothAbs_R8        _RL SMOOTHABS_RL
380        external smoothAbs_R8        EXTERNAL SMOOTHABS_RL
381    
382  Cml      smoothMin_R8 = .5*(a+b)  Cml      smoothMin_R8 = .5*(a+b)
383        smoothMin_R8 = .5*( a+b - smoothAbs_R8(a-b) )        SMOOTHMIN_RL = .5*( a+b - SMOOTHABS_RL(a-b) )
384  Cml      smoothMin_R8 = MIN(a,b)  Cml      smoothMin_R8 = MIN(a,b)
385    
386        return        RETURN
387        end        END
388    
389        _RS function smoothAbs_R4( x )        _RS FUNCTION SMOOTHABS_RS( x )
390    
391        implicit none        IMPLICIT NONE
392  C     === Global variables ===  C     === Global variables ===
393  #include "SIZE.h"  #include "SIZE.h"
394  #include "EEPARAMS.h"  #include "EEPARAMS.h"
# Line 412  C     input parameter Line 398  C     input parameter
398  c     local variable  c     local variable
399        _RS sf, rsf        _RS sf, rsf
400    
401        if ( smoothAbsFuncRange .lt. 0.0 ) then        IF ( smoothAbsFuncRange .LT. 0.0 ) THEN
402  c     limit of smoothMin(a,b) = .5*(a+b)  c     limit of smoothMin(a,b) = .5*(a+b)
403           smoothAbs_R4 = 0.           SMOOTHABS_RS = 0.
404        else        ELSE
405           if ( smoothAbsFuncRange .ne. 0.0 ) then           IF ( smoothAbsFuncRange .NE. 0.0 ) THEN
406              sf  = 10.0/smoothAbsFuncRange              sf  = 10.0/smoothAbsFuncRange
407              rsf = 1./sf              rsf = 1./sf
408           else           ELSE
409  c     limit of smoothMin(a,b) = min(a,b)  c     limit of smoothMin(a,b) = min(a,b)
410              sf  = 0.              sf  = 0.
411              rsf = 0.              rsf = 0.
412           end if           ENDIF
413  c  c
414           if ( x .gt. smoothAbsFuncRange ) then           IF ( x .GT. smoothAbsFuncRange ) THEN
415              smoothAbs_R4 = x              SMOOTHABS_RS = x
416           else if ( x .lt. -smoothAbsFuncRange ) then           ELSEIF ( x .LT. -smoothAbsFuncRange ) THEN
417              smoothAbs_R4 = -x              SMOOTHABS_RS = -x
418           else           ELSE
419              smoothAbs_R4 = log(.5*(exp(x*sf)+exp(-x*sf)))*rsf              SMOOTHABS_RS = log(.5*(exp(x*sf)+exp(-x*sf)))*rsf
420           end if           ENDIF
421        end if        ENDIF
422    
423        return        RETURN
424        end        END
425    
426        _RL function smoothAbs_R8( x )        _RL FUNCTION SMOOTHABS_RL( x )
427    
428        implicit none        IMPLICIT NONE
429  C     === Global variables ===  C     === Global variables ===
430  #include "SIZE.h"  #include "SIZE.h"
431  #include "EEPARAMS.h"  #include "EEPARAMS.h"
# Line 449  C     input parameter Line 435  C     input parameter
435  c     local variable  c     local variable
436        _RL sf, rsf        _RL sf, rsf
437    
438        if ( smoothAbsFuncRange .lt. 0.0 ) then        IF ( smoothAbsFuncRange .LT. 0.0 ) THEN
439  c     limit of smoothMin(a,b) = .5*(a+b)  c     limit of smoothMin(a,b) = .5*(a+b)
440           smoothAbs_R8 = 0.           SMOOTHABS_RL = 0.
441        else        ELSE
442           if ( smoothAbsFuncRange .ne. 0.0 ) then           IF ( smoothAbsFuncRange .NE. 0.0 ) THEN
443              sf  = 10.0D0/smoothAbsFuncRange              sf  = 10.0D0/smoothAbsFuncRange
444              rsf = 1.D0/sf              rsf = 1.D0/sf
445           else           ELSE
446  c     limit of smoothMin(a,b) = min(a,b)  c     limit of smoothMin(a,b) = min(a,b)
447              sf  = 0.D0              sf  = 0.D0
448              rsf = 0.D0              rsf = 0.D0
449           end if           ENDIF
450  c  c
451           if ( x .ge. smoothAbsFuncRange ) then           IF ( x .GE. smoothAbsFuncRange ) THEN
452              smoothAbs_R8 = x              SMOOTHABS_RL = x
453           else if ( x .le. -smoothAbsFuncRange ) then           ELSEIF ( x .LE. -smoothAbsFuncRange ) THEN
454              smoothAbs_R8 = -x              SMOOTHABS_RL = -x
455           else           ELSE
456              smoothAbs_R8 = log(.5*(exp(x*sf)+exp(-x*sf)))*rsf              SMOOTHABS_RL = log(.5*(exp(x*sf)+exp(-x*sf)))*rsf
457           end if           ENDIF
458        end if        ENDIF
459    
460        return        RETURN
461        end        END
462  #endif /* USE_SMOOTH_MIN */  #endif /* USE_SMOOTH_MIN */
463    
464  Cml#ifdef ALLOW_DEPTH_CONTROL  Cml#ifdef ALLOW_DEPTH_CONTROL
# Line 483  Cmlcadj SUBROUTINE limit_hfacc_to_one DE Line 469  Cmlcadj SUBROUTINE limit_hfacc_to_one DE
469  Cmlcadj SUBROUTINE limit_hfacc_to_one REQUIRED  Cmlcadj SUBROUTINE limit_hfacc_to_one REQUIRED
470  Cmlcadj SUBROUTINE limit_hfacc_to_one ADNAME  = adlimit_hfacc_to_one  Cmlcadj SUBROUTINE limit_hfacc_to_one ADNAME  = adlimit_hfacc_to_one
471  Cml#endif /* ALLOW_DEPTH_CONTROL */  Cml#endif /* ALLOW_DEPTH_CONTROL */
472  Cml      subroutine limit_hfacc_to_one( hf )  Cml      SUBROUTINE LIMIT_HFACC_TO_ONE( hf )
473  Cml  Cml
474  Cml      _RL hf  Cml      _RL hf
475  Cml  Cml
476  Cml      if ( hf .gt. 1. _d 0 ) then  Cml      IF ( hf .GT. 1. _d 0 ) THEN
477  Cml       hf = 1. _d 0  Cml       hf = 1. _d 0
478  Cml      endif  Cml      ENDIF
479  Cml  Cml
480  Cml      return  Cml      RETURN
481  Cml      end  Cml      END
482  Cml  Cml
483  Cml      subroutine adlimit_hfacc_to_one( hf, adhf )  Cml      SUBROUTINE ADLIMIT_HFACC_TO_ONE( hf, adhf )
484  Cml  Cml
485  Cml      _RL hf, adhf  Cml      _RL hf, adhf
486  Cml  Cml
487  Cml      return  Cml      RETURN
488  Cml      end  Cml      END
489    
490  #ifdef ALLOW_DEPTH_CONTROL  #ifdef ALLOW_DEPTH_CONTROL
491  cadj SUBROUTINE dummy_in_hfac INPUT   = 1, 2, 3  cadj SUBROUTINE dummy_in_hfac INPUT   = 1, 2, 3
# Line 511  cadj SUBROUTINE dummy_in_hfac INFLUENCED Line 497  cadj SUBROUTINE dummy_in_hfac INFLUENCED
497  cadj SUBROUTINE dummy_in_hfac ADNAME  = addummy_in_hfac  cadj SUBROUTINE dummy_in_hfac ADNAME  = addummy_in_hfac
498  cadj SUBROUTINE dummy_in_hfac FTLNAME = g_dummy_in_hfac  cadj SUBROUTINE dummy_in_hfac FTLNAME = g_dummy_in_hfac
499  #endif /* ALLOW_DEPTH_CONTROL */  #endif /* ALLOW_DEPTH_CONTROL */
   
Legend:
Removed from v.1.7  
changed lines
  Added in v.1.8
  ViewVC Help
Powered by ViewVC 1.1.22