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

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

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

revision 1.8 by jmc, Tue Dec 4 15:09:34 2001 UTC revision 1.13 by heimbach, Thu Nov 28 17:30:34 2002 UTC
# Line 5  C $Name$ Line 5  C $Name$
5    
6  CStartOfInterface  CStartOfInterface
7        SUBROUTINE GMREDI_CALC_TENSOR(        SUBROUTINE GMREDI_CALC_TENSOR(
8       I             bi, bj, iMin, iMax, jMin, jMax, K,       I             bi, bj, iMin, iMax, jMin, jMax,
9       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
10       I             myThid )       I             myThid )
11  C     /==========================================================\  C     /==========================================================\
# Line 24  C     == Global variables == Line 24  C     == Global variables ==
24  #include "GMREDI.h"  #include "GMREDI.h"
25  #include "GMREDI_DIAGS.h"  #include "GMREDI_DIAGS.h"
26    
27    #ifdef ALLOW_AUTODIFF_TAMC
28    #include "tamc.h"
29    #include "tamc_keys.h"
30    #endif /* ALLOW_AUTODIFF_TAMC */
31    
32  C     == Routine arguments ==  C     == Routine arguments ==
33  C  C
34        _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
35        _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
36        _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
37        INTEGER bi,bj,iMin,iMax,jMin,jMax,K        INTEGER bi,bj,iMin,iMax,jMin,jMax
38        INTEGER myThid        INTEGER myThid
39  CEndOfInterface  CEndOfInterface
40    
41  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
42    
43  C     == Local variables ==  C     == Local variables ==
44        INTEGER i,j,km1,kp1        INTEGER i,j,k,km1,kp1
45        _RL SlopeX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL SlopeX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
46        _RL SlopeY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL SlopeY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
47          _RL dSigmaDx(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
48          _RL dSigmaDy(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
49        _RL dSigmaDrReal(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL dSigmaDrReal(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
50        _RL SlopeSqr(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL SlopeSqr(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
51        _RL taperFct(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL taperFct(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
52          _RL maskp1, maskm1, Kgm_tmp
53    
54  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
55        _RS deltaH,zero_rs        _RS deltaH,zero_rs
56        PARAMETER(zero_rs=0.)        PARAMETER(zero_rs=0.)
57        _RL N2,SN        _RL N2,SN
58        _RL Ssq        _RL Ssq(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
59  #endif  #endif
60    
61    #ifdef ALLOW_AUTODIFF_TAMC
62              act1 = bi - myBxLo(myThid)
63              max1 = myBxHi(myThid) - myBxLo(myThid) + 1
64              act2 = bj - myByLo(myThid)
65              max2 = myByHi(myThid) - myByLo(myThid) + 1
66              act3 = myThid - 1
67              max3 = nTx*nTy
68              act4 = ikey_dynamics - 1
69              igmkey = (act1 + 1) + act2*max1
70         &                      + act3*max1*max2
71         &                      + act4*max1*max2*max3
72    #endif /* ALLOW_AUTODIFF_TAMC */
73    
74        km1=max(1,K-1)  #ifdef GM_VISBECK_VARIABLE_K
75        kp1=min(Nr,K)        DO j=1-Oly,sNy+Oly
76           DO i=1-Olx,sNx+Olx
77            VisbeckK(i,j,bi,bj) = 0. _d 0
78           ENDDO
79          ENDDO
80    #endif
81    
82          DO k=2,Nr
83    C-- 1rst loop on k : compute Tensor Coeff. at W points.
84           km1 = MAX(1,k-1)
85           maskm1 = 1. _d 0
86           IF (k.LE.1) maskm1 = 0. _d 0
87    
88  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
89  !HPF$ INDEPENDENT         kkey = (igmkey-1)*Nr + k
90           DO j=1-Oly,sNy+Oly
91            DO i=1-Olx,sNx+Olx
92             SlopeX(i,j)       = 0. _d 0
93             SlopeY(i,j)       = 0. _d 0
94             dSigmaDx(i,j)     = 0. _d 0
95             dSigmaDy(i,j)     = 0. _d 0
96             dSigmaDrReal(i,j) = 0. _d 0
97             SlopeSqr(i,j)     = 0. _d 0
98             taperFct(i,j)     = 0. _d 0
99             Kwx(i,j,k,bi,bj)  = 0. _d 0
100             Kwy(i,j,k,bi,bj)  = 0. _d 0
101             Kwz(i,j,k,bi,bj)  = 0. _d 0
102    # ifdef GM_NON_UNITY_DIAGONAL
103             Kux(i,j,k,bi,bj)  = 0. _d 0
104             Kvy(i,j,k,bi,bj)  = 0. _d 0
105    # endif
106    # ifdef GM_EXTRA_DIAGONAL
107             Kuz(i,j,k,bi,bj)  = 0. _d 0
108             Kvz(i,j,k,bi,bj)  = 0. _d 0
109    # endif
110    # ifdef GM_BOLUS_ADVEC
111             GM_PsiX(i,j,k,bi,bj)  = 0. _d 0
112             GM_PsiY(i,j,k,bi,bj)  = 0. _d 0
113    # endif
114            ENDDO
115           ENDDO
116  #endif  #endif
117    
118        DO j=1-Oly+1,sNy+Oly-1        DO j=1-Oly+1,sNy+Oly-1
 #ifdef ALLOW_AUTODIFF_TAMC  
 !HPF$ INDEPENDENT  
 #endif  
119         DO i=1-Olx+1,sNx+Olx-1         DO i=1-Olx+1,sNx+Olx-1
   
120  C      Gradient of Sigma at rVel points  C      Gradient of Sigma at rVel points
121          SlopeX(i,j)=0.25*( sigmaX(i+1, j ,km1) +sigmaX(i,j,km1)          dSigmaDx(i,j)=0.25*( sigmaX(i+1, j ,km1) +sigmaX(i,j,km1)
122       &                    +sigmaX(i+1, j , k ) +sigmaX(i,j, k ) )       &                    +sigmaX(i+1, j , k ) +sigmaX(i,j, k ) )
123       &                  *maskC(i,j,k,bi,bj)       &                  *maskC(i,j,k,bi,bj)*maskm1
124          SlopeY(i,j)=0.25*( sigmaY( i ,j+1,km1) +sigmaY(i,j,km1)          dSigmaDy(i,j)=0.25*( sigmaY( i ,j+1,km1) +sigmaY(i,j,km1)
125       &                    +sigmaY( i ,j+1, k ) +sigmaY(i,j, k ) )       &                    +sigmaY( i ,j+1, k ) +sigmaY(i,j, k ) )
126       &                  *maskC(i,j,k,bi,bj)       &                  *maskC(i,j,k,bi,bj)*maskm1
127          dSigmaDrReal(i,j)=sigmaR(i,j,k)          dSigmaDrReal(i,j)=sigmaR(i,j,k)*maskm1
   
128         ENDDO         ENDDO
129        ENDDO        ENDDO
130    
131    #ifdef ALLOW_AUTODIFF_TAMC
132    CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
133    CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
134    CADJ STORE dsigmadrreal(:,:)   = comlev1_bibj_k, key=kkey, byte=isbyte
135    #endif /* ALLOW_AUTODIFF_TAMC */
136    
137  C     Calculate slopes for use in tensor, taper and/or clip  C     Calculate slopes for use in tensor, taper and/or clip
138        CALL GMREDI_SLOPE_LIMIT(        CALL GMREDI_SLOPE_LIMIT(
139       I             dSigmadRReal,       U             dSigmadRReal,
140       I             rF(K),       I             rF(K),K,
141       U             SlopeX, SlopeY,       U             SlopeX, SlopeY,
142         U             dSigmaDx, dSigmaDy,
143       O             SlopeSqr, taperFct,       O             SlopeSqr, taperFct,
144       I             bi, bj, myThid )       I             bi, bj, myThid )
145    
# Line 88  C     Calculate slopes for use in tensor Line 147  C     Calculate slopes for use in tensor
147         DO i=1-Olx+1,sNx+Olx-1         DO i=1-Olx+1,sNx+Olx-1
148    
149  C       Mask Iso-neutral slopes  C       Mask Iso-neutral slopes
150          SlopeX(i,j)=SlopeX(i,j)*maskC(i,j,k,bi,bj)          SlopeX(i,j)=SlopeX(i,j)*maskC(i,j,k,bi,bj)*maskm1
151          SlopeY(i,j)=SlopeY(i,j)*maskC(i,j,k,bi,bj)          SlopeY(i,j)=SlopeY(i,j)*maskC(i,j,k,bi,bj)*maskm1
152          SlopeSqr(i,j)=SlopeSqr(i,j)*maskC(i,j,k,bi,bj)          SlopeSqr(i,j)=SlopeSqr(i,j)*maskC(i,j,k,bi,bj)*maskm1
153  c       Ssq=SlopeX(i,j)*SlopeX(i,j)+SlopeY(i,j)*SlopeY(i,j)  
154           ENDDO
155  C       Components of Redi/GM tensor        ENDDO
156          Kwx(i,j,k,bi,bj)=2.*SlopeX(i,j)*taperFct(i,j)  
157          Kwy(i,j,k,bi,bj)=2.*SlopeY(i,j)*taperFct(i,j)  #ifdef ALLOW_AUTODIFF_TAMC
158    CADJ STORE SlopeX(:,:)          = comlev1_bibj_k, key=kkey, byte=isbyte
159    CADJ STORE SlopeY(:,:)          = comlev1_bibj_k, key=kkey, byte=isbyte
160    CADJ STORE SlopeSqr(:,:)        = comlev1_bibj_k, key=kkey, byte=isbyte
161    CADJ STORE taperFct(:,:)        = comlev1_bibj_k, key=kkey, byte=isbyte
162    #ifdef GM_VISBECK_VARIABLE_K
163    CADJ STORE dSigmaDrReal(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
164    #endif
165    #endif /* ALLOW_AUTODIFF_TAMC */
166    
167          DO j=1-Oly+1,sNy+Oly-1
168           DO i=1-Olx+1,sNx+Olx-1
169    
170    C       Components of Redi/GM tensor
171            Kwx(i,j,k,bi,bj)= SlopeX(i,j)*taperFct(i,j)
172            Kwy(i,j,k,bi,bj)= SlopeY(i,j)*taperFct(i,j)
173          Kwz(i,j,k,bi,bj)= SlopeSqr(i,j)*taperFct(i,j)          Kwz(i,j,k,bi,bj)= SlopeSqr(i,j)*taperFct(i,j)
174    
175  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
176    
177  C- note (jmc) : moved here since only used in VISBECK_VARIABLE_K  C- note (jmc) : moved here since only used in VISBECK_VARIABLE_K
178  C           but don't know if *taperFct (or **2 ?) is necessary  C           but don't know if *taperFct (or **2 ?) is necessary
179          Ssq=SlopeSqr(i,j)*taperFct(i,j)          Ssq(i,j)=SlopeSqr(i,j)*taperFct(i,j)
180    
181  C--     Depth average of M^2/N^2 * N  C--     Depth average of M^2/N^2 * N
182    
# Line 118  C       Now we convert deltaH to a non-d Line 192  C       Now we convert deltaH to a non-d
192          deltaH=deltaH/GM_Visbeck_depth          deltaH=deltaH/GM_Visbeck_depth
193    
194          IF (K.eq.2) VisbeckK(i,j,bi,bj)=0.          IF (K.eq.2) VisbeckK(i,j,bi,bj)=0.
195          IF (Ssq.NE.0.) THEN          IF ( Ssq(i,j).NE.0. .AND. dSigmaDrReal(i,j).NE.0. ) THEN
196           N2= -Gravity*recip_Rhonil*dSigmaDrReal(i,j)           N2= -Gravity*recip_RhoConst*dSigmaDrReal(i,j)
197           SN=sqrt(Ssq*N2)           SN=sqrt(Ssq(i,j)*N2)
198           VisbeckK(i,j,bi,bj)=VisbeckK(i,j,bi,bj)+deltaH           VisbeckK(i,j,bi,bj)=VisbeckK(i,j,bi,bj)+deltaH
199       &      *GM_Visbeck_alpha*GM_Visbeck_length*GM_Visbeck_length*SN       &      *GM_Visbeck_alpha*GM_Visbeck_length*GM_Visbeck_length*SN
200          ENDIF          ENDIF
201    
 C       Limit range that KapGM can take  
         VisbeckK(i,j,bi,bj)=  
      &     min(VisbeckK(i,j,bi,bj),GM_Visbeck_maxval_K)  
   
202  #endif /* GM_VISBECK_VARIABLE_K */  #endif /* GM_VISBECK_VARIABLE_K */
203    
204           ENDDO
205          ENDDO
206    
207    C-- end 1rst loop on vertical level index k
208          ENDDO
209    
210    
211    #ifdef GM_VISBECK_VARIABLE_K
212    #ifdef ALLOW_AUTODIFF_TAMC
213    CADJ STORE VisbeckK(:,:,bi,bj) = comlev1_bibj, key=igmkey, byte=isbyte
214    #endif
215          IF ( GM_Visbeck_alpha.NE.0. ) THEN
216    C-     Limit range that KapGM can take
217           DO j=1-Oly+1,sNy+Oly-1
218            DO i=1-Olx+1,sNx+Olx-1
219             VisbeckK(i,j,bi,bj)=
220         &       MIN(VisbeckK(i,j,bi,bj),GM_Visbeck_maxval_K)
221  #ifdef ALLOW_TIMEAVE  #ifdef ALLOW_TIMEAVE
222  C--     Time-average           Visbeck_K_T(i,j,bi,bj)=Visbeck_K_T(i,j,bi,bj)
223          GM_Kwx_T(i,j,k,bi,bj)=GM_Kwx_T(i,j,k,bi,bj)       &                         +VisbeckK(i,j,bi,bj)*deltaTclock
224       &                       +Kwx(i,j,k,bi,bj)*deltaTclock  #endif
225          GM_Kwy_T(i,j,k,bi,bj)=GM_Kwy_T(i,j,k,bi,bj)          ENDDO
226       &                       +Kwy(i,j,k,bi,bj)*deltaTclock         ENDDO
227          GM_Kwz_T(i,j,k,bi,bj)=GM_Kwz_T(i,j,k,bi,bj)        ENDIF
228       &                       +Kwz(i,j,k,bi,bj)*deltaTclock  #endif /* GM_VISBECK_VARIABLE_K */
229    
230    
231    C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
232    
233    C-- 2nd loop on k : compute Tensor Coeff. at U,V levels.
234          DO k=1,Nr
235           kp1 = MIN(Nr,k+1)
236           maskp1 = 1. _d 0
237           IF (k.GE.Nr) maskp1 = 0. _d 0
238    
239    #ifdef ALLOW_AUTODIFF_TAMC
240           kkey = (igmkey-1)*Nr + k
241  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
242          IF (K.EQ.Nr)  CADJ STORE Kwx(:,:,k,bi,bj) = comlev1_bibj, key=kkey, byte=isbyte
243       &  Visbeck_K_T(i,j,bi,bj)=Visbeck_K_T(i,j,bi,bj)  CADJ STORE Kwy(:,:,k,bi,bj) = comlev1_bibj, key=kkey, byte=isbyte
244       &                       +VisbeckK(i,j,bi,bj)*deltaTclock  CADJ STORE Kwz(:,:,k,bi,bj) = comlev1_bibj, key=kkey, byte=isbyte
245  #endif  #endif
246  #endif /* ALLOW_TIMEAVE */  #endif
247    
248    C-    express the Tensor in term of Diffusivity (= m**2 / s )
249          DO j=1-Oly+1,sNy+Oly-1
250           DO i=1-Olx+1,sNx+Olx-1
251            Kgm_tmp = GM_isopycK + GM_skewflx*GM_background_K
252    #ifdef GM_VISBECK_VARIABLE_K
253         &          + VisbeckK(i,j,bi,bj)*(1.+GM_skewflx)    
254    #endif
255            Kwx(i,j,k,bi,bj)= Kgm_tmp*Kwx(i,j,k,bi,bj)
256            Kwy(i,j,k,bi,bj)= Kgm_tmp*Kwy(i,j,k,bi,bj)
257            Kwz(i,j,k,bi,bj)= ( GM_isopycK
258    #ifdef GM_VISBECK_VARIABLE_K
259         &                    + VisbeckK(i,j,bi,bj)
260    #endif
261         &                    )*Kwz(i,j,k,bi,bj)
262         ENDDO         ENDDO
263        ENDDO        ENDDO
264    #ifdef ALLOW_AUTODIFF_TAMC
265  #ifdef ALLOW_TIMEAVE  #ifdef GM_VISBECK_VARIABLE_K
266        GM_TimeAve(k,bi,bj)=GM_TimeAve(k,bi,bj)+deltaTclock  CADJ STORE VisbeckK(:,:,bi,bj) =
267    CADJ &     comlev1_bibj, key=kkey, byte=isbyte
268    #endif
269  #endif  #endif
270    
271    #if ( defined (GM_NON_UNITY_DIAGONAL) || defined (GM_EXTRA_DIAGONAL) )
272    
 #ifdef GM_NON_UNITY_DIAGONAL  
273  C     Gradient of Sigma at U points  C     Gradient of Sigma at U points
274        DO j=1-Oly+1,sNy+Oly-1        DO j=1-Oly+1,sNy+Oly-1
275         DO i=1-Olx+1,sNx+Olx-1         DO i=1-Olx+1,sNx+Olx-1
276          SlopeX(i,j)=sigmaX(i,j,km1)          dSigmaDx(i,j)=sigmaX(i,j,k)
277       &          *_maskW(i,j,k,bi,bj)       &          *_maskW(i,j,k,bi,bj)
278          SlopeY(i,j)=0.25*( sigmaY(i-1,j+1,k) +sigmaY(i,j+1,k)          dSigmaDy(i,j)=0.25*( sigmaY(i-1,j+1,k) +sigmaY(i,j+1,k)
279       &                    +sigmaY(i-1, j ,k) +sigmaY(i, j ,k) )       &                      +sigmaY(i-1, j ,k) +sigmaY(i, j ,k) )
280       &          *_maskW(i,j,k,bi,bj)       &          *_maskW(i,j,k,bi,bj)
281          dSigmaDrReal(i,j)=0.25*( sigmaR(i-1,j, k ) +sigmaR(i,j, k )          dSigmaDrReal(i,j)=0.25*( sigmaR(i-1,j, k ) +sigmaR(i,j, k )
282       &                          +sigmaR(i-1,j,kp1) +sigmaR(i,j,kp1) )       &                  +maskp1*(sigmaR(i-1,j,kp1) +sigmaR(i,j,kp1)) )
283       &          *_maskW(i,j,k,bi,bj)       &          *_maskW(i,j,k,bi,bj)
284         ENDDO         ENDDO
285        ENDDO        ENDDO
286    
287    #ifdef ALLOW_AUTODIFF_TAMC
288    CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
289    CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
290    CADJ STORE dsigmadrreal(:,:)   = comlev1_bibj_k, key=kkey, byte=isbyte
291    #endif /* ALLOW_AUTODIFF_TAMC */
292    
293  C     Calculate slopes for use in tensor, taper and/or clip  C     Calculate slopes for use in tensor, taper and/or clip
294        CALL GMREDI_SLOPE_LIMIT(        CALL GMREDI_SLOPE_LIMIT(
295       I             dSigmadRReal,       U             dSigmadRReal,
296       I             rF(K),       I             rF(K),K,
297       U             SlopeX, SlopeY,       U             SlopeX, SlopeY,
298         U             dSigmaDx, dSigmaDy,
299       O             SlopeSqr, taperFct,       O             SlopeSqr, taperFct,
300       I             bi, bj, myThid )       I             bi, bj, myThid )
301    
302        DO j=1-Oly+1,sNy+Oly-1  #ifdef GM_NON_UNITY_DIAGONAL
303         DO i=1-Olx+1,sNx+Olx-1          DO j=1-Oly+1,sNy+Oly-1
304          Kux(i,j,k,bi,bj)=taperFct(i,j)           DO i=1-Olx+1,sNx+Olx-1
305         ENDDO            Kux(i,j,k,bi,bj) =
306        ENDDO       &     ( GM_isopycK
307    #ifdef GM_VISBECK_VARIABLE_K
308         &     +0.5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))
309    #endif
310         &     )
311         &     *taperFct(i,j)
312             ENDDO
313            ENDDO
314    #ifdef ALLOW_AUTODIFF_TAMC
315    # ifndef GM_TAPER_ORIG_CLIPPING
316    CADJ STORE Kux(:,:,k,bi,bj)  = comlev1_bibj_k, key=kkey, byte=isbyte
317    # endif
318    #endif
319            DO j=1-Oly+1,sNy+Oly-1
320             DO i=1-Olx+1,sNx+Olx-1
321              Kux(i,j,k,bi,bj) = MAX( Kux(i,j,k,bi,bj), GM_Kmin_horiz )
322             ENDDO
323            ENDDO
324    #endif /* GM_NON_UNITY_DIAGONAL */
325    
326    #ifdef GM_EXTRA_DIAGONAL
327    
328    #ifdef ALLOW_AUTODIFF_TAMC
329    CADJ STORE SlopeX(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
330    CADJ STORE taperFct(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte
331    #endif /* ALLOW_AUTODIFF_TAMC */
332          IF (GM_ExtraDiag) THEN
333            DO j=1-Oly+1,sNy+Oly-1
334             DO i=1-Olx+1,sNx+Olx-1
335              Kuz(i,j,k,bi,bj) =
336         &     ( GM_isopycK - GM_skewflx*GM_background_K
337    #ifdef GM_VISBECK_VARIABLE_K
338         &     +0.5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))*GM_advect
339    #endif
340         &     )*SlopeX(i,j)*taperFct(i,j)
341             ENDDO
342            ENDDO
343          ENDIF
344    #endif /* GM_EXTRA_DIAGONAL */
345    
346  C     Gradient of Sigma at V points  C     Gradient of Sigma at V points
347        DO j=1-Oly+1,sNy+Oly-1        DO j=1-Oly+1,sNy+Oly-1
348         DO i=1-Olx+1,sNx+Olx-1         DO i=1-Olx+1,sNx+Olx-1
349          SlopeX(i,j)=0.25*( sigmaX(i, j ,k) +sigmaX(i+1, j ,k)          dSigmaDx(i,j)=0.25*( sigmaX(i, j ,k) +sigmaX(i+1, j ,k)
350       &                    +sigmaX(i,j-1,k) +sigmaX(i+1,j-1,k) )       &                    +sigmaX(i,j-1,k) +sigmaX(i+1,j-1,k) )
351       &          *_maskS(i,j,k,bi,bj)       &          *_maskS(i,j,k,bi,bj)
352          SlopeY(i,j)=sigmaY(i,j,km1)          dSigmaDy(i,j)=sigmaY(i,j,k)
353       &          *_maskS(i,j,k,bi,bj)       &          *_maskS(i,j,k,bi,bj)
354          dSigmaDrReal(i,j)=0.25*( sigmaR(i,j-1, k ) +sigmaR(i,j, k )          dSigmaDrReal(i,j)=0.25*( sigmaR(i,j-1, k ) +sigmaR(i,j, k )
355       &                          +sigmaR(i,j-1,kp1) +sigmaR(i,j,kp1) )       &                  +maskp1*(sigmaR(i,j-1,kp1) +sigmaR(i,j,kp1)) )
356       &          *_maskS(i,j,k,bi,bj)       &          *_maskS(i,j,k,bi,bj)
357         ENDDO         ENDDO
358        ENDDO        ENDDO
359    
360    #ifdef ALLOW_AUTODIFF_TAMC
361    CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
362    CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
363    CADJ STORE dsigmadrreal(:,:)   = comlev1_bibj_k, key=kkey, byte=isbyte
364    #endif /* ALLOW_AUTODIFF_TAMC */
365    
366  C     Calculate slopes for use in tensor, taper and/or clip  C     Calculate slopes for use in tensor, taper and/or clip
367        CALL GMREDI_SLOPE_LIMIT(        CALL GMREDI_SLOPE_LIMIT(
368       I             dSigmadRReal,       U             dSigmadRReal,
369       I             rF(K),       I             rF(K),K,
370       U             SlopeX, SlopeY,       U             SlopeX, SlopeY,
371         U             dSigmaDx, dSigmaDy,
372       O             SlopeSqr, taperFct,       O             SlopeSqr, taperFct,
373       I             bi, bj, myThid )       I             bi, bj, myThid )
374    
375    #ifdef GM_NON_UNITY_DIAGONAL
376            DO j=1-Oly+1,sNy+Oly-1
377             DO i=1-Olx+1,sNx+Olx-1
378              Kvy(i,j,k,bi,bj) =
379         &     ( GM_isopycK
380    #ifdef GM_VISBECK_VARIABLE_K
381         &     +0.5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))
382    #endif
383         &     )
384         &     *taperFct(i,j)
385             ENDDO
386            ENDDO
387    #ifdef ALLOW_AUTODIFF_TAMC
388    # ifndef GM_TAPER_ORIG_CLIPPING
389    CADJ STORE Kvy(:,:,k,bi,bj)  = comlev1_bibj_k, key=kkey, byte=isbyte
390    # endif
391    #endif
392            DO j=1-Oly+1,sNy+Oly-1
393             DO i=1-Olx+1,sNx+Olx-1
394              Kvy(i,j,k,bi,bj) = MAX( Kvy(i,j,k,bi,bj), GM_Kmin_horiz )
395             ENDDO
396            ENDDO
397    #endif /* GM_NON_UNITY_DIAGONAL */
398    
399    #ifdef GM_EXTRA_DIAGONAL
400    
401    #ifdef ALLOW_AUTODIFF_TAMC
402    CADJ STORE SlopeY(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
403    CADJ STORE taperFct(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte
404    #endif /* ALLOW_AUTODIFF_TAMC */
405          IF (GM_ExtraDiag) THEN
406            DO j=1-Oly+1,sNy+Oly-1
407             DO i=1-Olx+1,sNx+Olx-1
408              Kvz(i,j,k,bi,bj) =
409         &     ( GM_isopycK - GM_skewflx*GM_background_K
410    #ifdef GM_VISBECK_VARIABLE_K
411         &     +0.5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))*GM_advect
412    #endif
413         &     )*SlopeY(i,j)*taperFct(i,j)
414             ENDDO
415            ENDDO
416          ENDIF
417    #endif /* GM_EXTRA_DIAGONAL */
418    
419    #endif /* GM_NON_UNITY_DIAGONAL || GM_EXTRA_DIAGONAL */
420    
421    #ifdef ALLOW_TIMEAVE
422    C--   Time-average
423        DO j=1-Oly+1,sNy+Oly-1        DO j=1-Oly+1,sNy+Oly-1
424         DO i=1-Olx+1,sNx+Olx-1         DO i=1-Olx+1,sNx+Olx-1
425          Kvy(i,j,k,bi,bj)=taperFct(i,j)          GM_Kwx_T(i,j,k,bi,bj)=GM_Kwx_T(i,j,k,bi,bj)
426         &                       +Kwx(i,j,k,bi,bj)*deltaTclock
427            GM_Kwy_T(i,j,k,bi,bj)=GM_Kwy_T(i,j,k,bi,bj)
428         &                       +Kwy(i,j,k,bi,bj)*deltaTclock
429            GM_Kwz_T(i,j,k,bi,bj)=GM_Kwz_T(i,j,k,bi,bj)
430         &                       +Kwz(i,j,k,bi,bj)*deltaTclock
431         ENDDO         ENDDO
432        ENDDO        ENDDO
433          GM_TimeAve(k,bi,bj)=GM_TimeAve(k,bi,bj)+deltaTclock
434    #endif /* ALLOW_TIMEAVE */
435    
436  #endif /* GM_NON_UNITY_DIAGONAL */  C-- end 2nd loop on vertical level index k
437          ENDDO
438    
439    
440    #ifdef GM_BOLUS_ADVEC
441          IF (GM_AdvForm) THEN
442            CALL GMREDI_CALC_PSI_B(
443         I             bi, bj, iMin, iMax, jMin, jMax,
444         I             sigmaX, sigmaY, sigmaR,
445         I             myThid )
446          ENDIF
447    #endif
448    
449  #endif /* ALLOW_GMREDI */  #endif /* ALLOW_GMREDI */
450    
# Line 222  C     Calculate slopes for use in tensor Line 453  C     Calculate slopes for use in tensor
453    
454    
455        SUBROUTINE GMREDI_CALC_TENSOR_DUMMY(        SUBROUTINE GMREDI_CALC_TENSOR_DUMMY(
456       I             bi, bj, iMin, iMax, jMin, jMax, K,       I             bi, bj, iMin, iMax, jMin, jMax,
457       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
458       I             myThid )       I             myThid )
459  C     /==========================================================\  C     /==========================================================\
# Line 245  C Line 476  C
476        _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
477        _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
478        _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
479        INTEGER bi,bj,iMin,iMax,jMin,jMax,K        INTEGER bi,bj,iMin,iMax,jMin,jMax
480        INTEGER myThid        INTEGER myThid
481  CEndOfInterface  CEndOfInterface
482    
483        INTEGER i, j        INTEGER i, j, k
484    
485  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
486    
487        DO j=1-Oly+1,sNy+Oly-1        DO k=1,Nr
488         DO i=1-Olx+1,sNx+Olx-1         DO j=1-Oly+1,sNy+Oly-1
489          Kwx(i,j,k,bi,bj) = 0.0          DO i=1-Olx+1,sNx+Olx-1
490          Kwy(i,j,k,bi,bj) = 0.0           Kwx(i,j,k,bi,bj) = 0.0
491          Kwz(i,j,k,bi,bj) = 0.0           Kwy(i,j,k,bi,bj) = 0.0
492             Kwz(i,j,k,bi,bj) = 0.0
493            ENDDO
494         ENDDO         ENDDO
495        ENDDO        ENDDO
496  #endif /* ALLOW_GMREDI */  #endif /* ALLOW_GMREDI */
497    
498        end        RETURN
499          END
Legend:
Removed from v.1.8  
changed lines
  Added in v.1.13
  ViewVC Help
Powered by ViewVC 1.1.22