/[MITgcm]/MITgcm/pkg/gmredi/gmredi_slope_psi.F
ViewVC logotype

Diff of /MITgcm/pkg/gmredi/gmredi_slope_psi.F

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

revision 1.2 by heimbach, Fri Feb 15 21:28:07 2002 UTC revision 1.3 by heimbach, Thu Nov 14 22:43:49 2002 UTC
# Line 6  C $Name$ Line 6  C $Name$
6  CStartOfInterface  CStartOfInterface
7        SUBROUTINE GMREDI_SLOPE_PSI_B(        SUBROUTINE GMREDI_SLOPE_PSI_B(
8       I             dSigmaDrW,dSigmaDrS,       I             dSigmaDrW,dSigmaDrS,
9       I             depthZ,       I             depthZ,K,
10       U             SlopeX, SlopeY,       U             SlopeX, SlopeY,
11       O             taperX, taperY,       O             taperX, taperY,
12       I             bi,bj, myThid )       I             bi,bj, myThid )
# Line 15  C     | SUBROUTINE GMREDI_SLOPE_PSI_B Line 15  C     | SUBROUTINE GMREDI_SLOPE_PSI_B
15  C     | o Calculate slopes for use in GM/Redi tensor             |  C     | o Calculate slopes for use in GM/Redi tensor             |
16  C     |==========================================================|  C     |==========================================================|
17  C     | On entry:                                                |  C     | On entry:                                                |
18  C     |            dSigmaDrW conatins the d/dz Sigma          |  C     |            dSigmaDrW conatins the d/dz Sigma             |
19  C     |            SlopeX/Y     contains X/Y gradients of sigma  |  C     |            SlopeX/Y     contains X/Y gradients of sigma  |
20  C     |            depthZ       conatins the height (m) of level |  C     |            depthZ       conatins the height (m) of level |
21  C     | On exit:                                                 |  C     | On exit:                                                 |
22  C     |            dSigmaDrW conatins the effective dSig/dz   |  C     |            dSigmaDrW conatins the effective dSig/dz      |
23  C     |            SlopeX/Y     contains X/Y slopes              |  C     |            SlopeX/Y     contains X/Y slopes              |
24  C     |            taperFct     contains tapering funct. value ; |  C     |            taperFct     contains tapering funct. value ; |
25  C     |                         = 1 when using no tapering       |  C     |                         = 1 when using no tapering       |
# Line 32  C     == Global variables == Line 32  C     == Global variables ==
32  #include "GMREDI.h"  #include "GMREDI.h"
33  #include "PARAMS.h"  #include "PARAMS.h"
34    
35    #ifdef ALLOW_AUTODIFF_TAMC
36    #include "tamc.h"
37    #include "tamc_keys.h"
38    #endif /* ALLOW_AUTODIFF_TAMC */
39    
40  C     == Routine arguments ==  C     == Routine arguments ==
41  C  C
42        _RL SlopeX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL SlopeX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
# Line 41  C Line 46  C
46        _RL taperX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL taperX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
47        _RL taperY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL taperY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
48        _RL depthZ        _RL depthZ
49        INTEGER bi,bj,myThid        INTEGER K,bi,bj,myThid
50  CEndOfInterface  CEndOfInterface
51    
52  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
53    #ifdef GM_BOLUS_ADVEC
54    
55  C     == Local variables ==  C     == Local variables ==
56        _RL Small_Number        _RL gradSmod(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
57        PARAMETER(Small_Number=1.D-12)        _RL dSigmaDrLtd(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
58        _RL gradSmod,f1,Smod,f2,Rnondim,Cspd,Lrho        _RL drdsigmaltd(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
59        _RL dSigmaDrLtd, dRdSigmaLtd        _RL f1,Smod,f2,Rnondim,Cspd,Lrho
60        _RL maxSlopeSqr        _RL maxSlopeSqr
61        _RL fpi        _RL fpi
62        PARAMETER(fpi=3.141592653589793047592d0)        PARAMETER(fpi=3.141592653589793047592d0)
63        INTEGER i,j        INTEGER i,j
64    
65    #ifdef ALLOW_AUTODIFF_TAMC
66          act1 = bi - myBxLo(myThid)
67          max1 = myBxHi(myThid) - myBxLo(myThid) + 1
68          act2 = bj - myByLo(myThid)
69          max2 = myByHi(myThid) - myByLo(myThid) + 1
70          act3 = myThid - 1
71          max3 = nTx*nTy
72          act4 = ikey_dynamics - 1
73          ikey = (act1 + 1) + act2*max1
74         &                  + act3*max1*max2
75         &                  + act4*max1*max2*max3
76          kkey = (ikey-1)*Nr + k
77    #endif /* ALLOW_AUTODIFF_TAMC */
78    
79    #ifdef ALLOW_AUTODIFF_TAMC
80    CADJ STORE slopeX(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
81    CADJ STORE slopeY(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
82    #endif
83    
84        IF (GM_taper_scheme.EQ.'orig' .OR.        IF (GM_taper_scheme.EQ.'orig' .OR.
85       &    GM_taper_scheme.EQ.'clipping') THEN       &    GM_taper_scheme.EQ.'clipping') THEN
86    
87  C-      Original implementation in mitgcmuv  C-      Original implementation in mitgcmuv
88  C       (this turns out to be the same as Cox slope clipping)  C       (this turns out to be the same as Cox slope clipping)
89    
90    C-- X-comp
91    
92    #ifdef ALLOW_AUTODIFF_TAMC
93          DO j=1-Oly+1,sNy+Oly-1
94           DO i=1-Olx+1,sNx+Olx-1
95            dSigmaDrLtd(i,j) = 0. _d 0
96           ENDDO
97          ENDDO
98    #endif /* ALLOW_AUTODIFF_TAMC */
99    
100  C-      Cox 1987 "Slope clipping"  C-      Cox 1987 "Slope clipping"
101          DO j=1-Oly+1,sNy+Oly-1          DO j=1-Oly+1,sNy+Oly-1
102           DO i=1-Olx+1,sNx+Olx-1           DO i=1-Olx+1,sNx+Olx-1
103              dsigmadrltd(i,j) = -(GM_Small_Number+
104  c         gradSmod=SlopeX(i,j)*SlopeX(i,j)       &     abs(SlopeX(i,j))*GM_rMaxSlope)
105  c    &            +SlopeY(i,j)*SlopeY(i,j)           ENDDO
106  c         if (gradSmod .NE. 0.) gradSmod=sqrt(gradSmod)          ENDDO
107            gradSmod=abs(SlopeX(i,j))  #ifdef ALLOW_AUTODIFF_TAMC
108    CADJ STORE dSigmaDrltd(:,:)  = comlev1_bibj_k, key=kkey, byte=isbyte
109            dSigmaDrLtd = -(Small_Number+gradSmod*GM_rMaxSlope)  CADJ STORE dSigmaDrW(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
110            IF (dSigmaDrW(i,j).GE.dSigmaDrLtd)  #endif
111       &        dSigmaDrW(i,j) = dSigmaDrLtd          DO j=1-Oly+1,sNy+Oly-1
112             DO i=1-Olx+1,sNx+Olx-1
113              IF (dSigmaDrW(i,j).GE.dsigmadrltd(i,j))
114         &        dSigmaDrW(i,j) = dsigmadrltd(i,j)
115             ENDDO
116            ENDDO
117    #ifdef ALLOW_AUTODIFF_TAMC
118    CADJ STORE dSigmaDrW(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
119    #endif
120            DO j=1-Oly+1,sNy+Oly-1
121             DO i=1-Olx+1,sNx+Olx-1
122            SlopeX(i,j) = -SlopeX(i,j)/dSigmaDrW(i,j)            SlopeX(i,j) = -SlopeX(i,j)/dSigmaDrW(i,j)
123              taperX(i,j)=1. _d 0
124             ENDDO
125            ENDDO
126    
127            gradSmod=abs(SlopeY(i,j))  C-- Y-comp
           dSigmaDrLtd = -(Small_Number+gradSmod*GM_rMaxSlope)  
           IF (dSigmaDrS(i,j).GE.dSigmaDrLtd)  
      &        dSigmaDrS(i,j) = dSigmaDrLtd  
           SlopeY(i,j) = -SlopeY(i,j)/dSigmaDrS(i,j)  
128    
129            taperX(i,j)=1. _d 0  #ifdef ALLOW_AUTODIFF_TAMC
130          DO j=1-Oly+1,sNy+Oly-1
131           DO i=1-Olx+1,sNx+Olx-1
132            dSigmaDrLtd(i,j) = 0. _d 0
133           ENDDO
134          ENDDO
135    #endif /* ALLOW_AUTODIFF_TAMC */
136            DO j=1-Oly+1,sNy+Oly-1
137             DO i=1-Olx+1,sNx+Olx-1
138              dsigmadrltd(i,j) = -(GM_Small_Number+
139         &     abs(SlopeY(i,j))*GM_rMaxSlope)
140             ENDDO
141            ENDDO
142    #ifdef ALLOW_AUTODIFF_TAMC
143    CADJ STORE dSigmaDrltd(:,:)  = comlev1_bibj_k, key=kkey, byte=isbyte
144    CADJ STORE dSigmaDrS(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
145    #endif
146            DO j=1-Oly+1,sNy+Oly-1
147             DO i=1-Olx+1,sNx+Olx-1
148              IF (dSigmaDrS(i,j).GE.dsigmadrltd(i,j))
149         &        dSigmaDrS(i,j) = dsigmadrltd(i,j)
150             ENDDO
151            ENDDO
152    #ifdef ALLOW_AUTODIFF_TAMC
153    CADJ STORE dSigmaDrS(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
154    #endif
155            DO j=1-Oly+1,sNy+Oly-1
156             DO i=1-Olx+1,sNx+Olx-1
157              SlopeY(i,j) = -SlopeY(i,j)/dSigmaDrS(i,j)
158            taperY(i,j)=1. _d 0            taperY(i,j)=1. _d 0
   
159           ENDDO           ENDDO
160          ENDDO          ENDDO
161    
162        ELSE        ELSE
163    
164    #ifdef ALLOW_AUTODIFF_TAMC
165    CADJ STORE slopeX(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
166    CADJ STORE dSigmaDrW(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
167    #endif
168    
169  C- Compute the slope, no clipping, but avoid reverse slope in negatively  C- Compute the slope, no clipping, but avoid reverse slope in negatively
170  C                                  stratified (Sigma_Z > 0) region :  C                                  stratified (Sigma_Z > 0) region :
171          DO j=1-Oly+1,sNy+Oly-1          DO j=1-Oly+1,sNy+Oly-1
172           DO i=1-Olx+1,sNx+Olx-1           DO i=1-Olx+1,sNx+Olx-1
173              dsigmadrltd(i,j) = -GM_Small_Number
174            dSigmaDrLtd = -Small_Number            IF (dSigmaDrW(i,j).GE.dsigmadrltd(i,j))
175            IF (dSigmaDrW(i,j).GE.dSigmaDrLtd)       &        dSigmaDrW(i,j) = dsigmadrltd(i,j)
176       &        dSigmaDrW(i,j) = dSigmaDrLtd           ENDDO
177            dRdSigmaLtd = 1./dSigmaDrW(i,j)          ENDDO
178    #ifdef ALLOW_AUTODIFF_TAMC
179    CADJ STORE dsigmadrW(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
180    #endif
181            DO j=1-Oly+1,sNy+Oly-1
182             DO i=1-Olx+1,sNx+Olx-1
183              drdsigmaltd(i,j) = 1./dSigmaDrW(i,j)
184            SlopeX(i,j) = -SlopeX(i,j)/dSigmaDrW(i,j)            SlopeX(i,j) = -SlopeX(i,j)/dSigmaDrW(i,j)
185              taperX(i,j)=1. _d 0
186             ENDDO
187            ENDDO
188    
189            dSigmaDrLtd = -Small_Number  #ifdef ALLOW_AUTODIFF_TAMC
190            IF (dSigmaDrS(i,j).GE.dSigmaDrLtd)  CADJ STORE slopeY(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
191       &        dSigmaDrS(i,j) = dSigmaDrLtd  CADJ STORE dSigmaDrS(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
192            SlopeY(i,j) = -SlopeY(i,j)/dSigmaDrS(i,j)  #endif
   
 c         SlopeSqr(i,j)=SlopeX(i,j)*SlopeX(i,j)  
 c    &                 +SlopeY(i,j)*SlopeY(i,j)  
193    
194            taperX(i,j)=1. _d 0          DO j=1-Oly+1,sNy+Oly-1
195             DO i=1-Olx+1,sNx+Olx-1
196              dsigmadrltd(i,j) = -GM_Small_Number
197              IF (dSigmaDrS(i,j).GE.dsigmadrltd(i,j))
198         &        dSigmaDrS(i,j) = dsigmadrltd(i,j)
199             ENDDO
200            ENDDO
201    #ifdef ALLOW_AUTODIFF_TAMC
202    CADJ STORE dsigmadrS(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
203    #endif
204            DO j=1-Oly+1,sNy+Oly-1
205             DO i=1-Olx+1,sNx+Olx-1
206              SlopeY(i,j) = -SlopeY(i,j)/dSigmaDrS(i,j)
207            taperY(i,j)=1. _d 0            taperY(i,j)=1. _d 0
   
208           ENDDO           ENDDO
209          ENDDO          ENDDO
210    
211  C- Compute the tapering function for the GM+Redi tensor :  C- Compute the tapering function for the GM+Redi tensor :
212    
213    #ifdef ALLOW_AUTODIFF_TAMC
214    CADJ STORE slopeX(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
215    CADJ STORE slopeY(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
216    #endif
217    
218         IF (GM_taper_scheme.EQ.'linear') THEN         IF (GM_taper_scheme.EQ.'linear') THEN
219    
220  C-      Simplest adiabatic tapering = Smax/Slope (linear)  C-      Simplest adiabatic tapering = Smax/Slope (linear)
# Line 202  C-      Large, Danabasoglu and Doney, JP Line 300  C-      Large, Danabasoglu and Doney, JP
300        ENDIF        ENDIF
301    
302    
303    #endif /* BOLUS_ADVEC */
304  #endif /* ALLOW_GMREDI */  #endif /* ALLOW_GMREDI */
305    
306        RETURN        RETURN
Legend:
Removed from v.1.2  
changed lines
  Added in v.1.3
  ViewVC Help
Powered by ViewVC 1.1.22