/[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.2 by heimbach, Mon Jan 8 20:11:04 2001 UTC revision 1.14 by heimbach, Fri Jan 10 00:48:39 2003 UTC
# Line 1  Line 1 
1  C $Header$  C $Header$
2    C $Name$
3    
4  #include "GMREDI_OPTIONS.h"  #include "GMREDI_OPTIONS.h"
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 23  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 dRdSigmaLtd(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL SlopeSqr(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
51        _RL Ssq        _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        _RL deltaH,zero_rs
56        PARAMETER(zero_rs=0.)        PARAMETER(zero_rs=0.D0)
57        _RL N2,SN        _RL N2,SN
58          _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)=op25*( 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          SlopeY(i,j)=0.25*( sigmaY( i ,j+1,km1) +sigmaY(i,j,km1)       &                  *maskC(i,j,k,bi,bj)*maskm1
124            dSigmaDy(i,j)=op25*( 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          dSigmaDrReal(i,j)=sigmaR(i,j,k)       &                  *maskC(i,j,k,bi,bj)*maskm1
127            dSigmaDrReal(i,j)=sigmaR(i,j,k)*maskm1
         if (hFacC(i,j,k,bi,bj).eq.0.) then  
          SlopeX(i,j)=0.  
          SlopeY(i,j)=0.  
         endif  
   
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       O             dRdSigmaLtd,       U             dSigmaDx, dSigmaDy,
143         O             SlopeSqr, taperFct,
144       I             bi, bj, myThid )       I             bi, bj, myThid )
145    
146        DO j=1-Oly+1,sNy+Oly-1        DO j=1-Oly+1,sNy+Oly-1
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          if (hFacC(i,j,k,bi,bj).eq.0.) then          SlopeX(i,j)=SlopeX(i,j)*maskC(i,j,k,bi,bj)*maskm1
151           SlopeX(i,j)=0.          SlopeY(i,j)=SlopeY(i,j)*maskC(i,j,k,bi,bj)*maskm1
152           SlopeY(i,j)=0.          SlopeSqr(i,j)=SlopeSqr(i,j)*maskC(i,j,k,bi,bj)*maskm1
153          endif  
154          Ssq=SlopeX(i,j)*SlopeX(i,j)+SlopeY(i,j)*SlopeY(i,j)         ENDDO
155          ENDDO
156  C       Components of Redi/GM tensor  
157          Kwx(i,j,k,myThid)=2.*SlopeX(i,j)  #ifdef ALLOW_AUTODIFF_TAMC
158          Kwy(i,j,k,myThid)=2.*SlopeY(i,j)  CADJ STORE SlopeSqr(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte
159          Kwz(i,j,k,myThid)=Ssq  #endif /* ALLOW_AUTODIFF_TAMC */
160    
161          DO j=1-Oly+1,sNy+Oly-1
162           DO i=1-Olx+1,sNx+Olx-1
163    
164    C       Components of Redi/GM tensor
165            Kwx(i,j,k,bi,bj)= SlopeX(i,j)*taperFct(i,j)
166            Kwy(i,j,k,bi,bj)= SlopeY(i,j)*taperFct(i,j)
167            Kwz(i,j,k,bi,bj)= SlopeSqr(i,j)*taperFct(i,j)
168    
169  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
170    
171    C- note (jmc) : moved here since only used in VISBECK_VARIABLE_K
172    C           but don't know if *taperFct (or **2 ?) is necessary
173            Ssq(i,j)=SlopeSqr(i,j)*taperFct(i,j)
174    
175  C--     Depth average of M^2/N^2 * N  C--     Depth average of M^2/N^2 * N
176    
177  C       Calculate terms for mean Richardson number  C       Calculate terms for mean Richardson number
# Line 114  C       If negative then we are below th Line 185  C       If negative then we are below th
185  C       Now we convert deltaH to a non-dimensional fraction  C       Now we convert deltaH to a non-dimensional fraction
186          deltaH=deltaH/GM_Visbeck_depth          deltaH=deltaH/GM_Visbeck_depth
187    
188          if (K.eq.2) VisbeckK(i,j,myThid)=0.          IF (K.eq.2) VisbeckK(i,j,bi,bj)=0.
189  Calt?   if (dSigmaDrReal(i,j).NE.0.) then          IF ( Ssq(i,j).NE.0. .AND. dSigmaDrReal(i,j).NE.0. ) THEN
190  Calt?    N2=(-Gravity*recip_Rhonil)*dSigmaDrReal(i,j)           N2= -Gravity*recip_RhoConst*dSigmaDrReal(i,j)
191          if (dRdSigmaLtd(i,j).NE.0.) then           SN=sqrt(Ssq(i,j)*N2)
192           N2=(-Gravity*recip_Rhonil)/dRdSigmaLtd(i,j)           VisbeckK(i,j,bi,bj)=VisbeckK(i,j,bi,bj)+deltaH
          SN=sqrt(Ssq*N2)  
          VisbeckK(i,j,myThid)=VisbeckK(i,j,myThid)+deltaH  
193       &      *GM_Visbeck_alpha*GM_Visbeck_length*GM_Visbeck_length*SN       &      *GM_Visbeck_alpha*GM_Visbeck_length*GM_Visbeck_length*SN
194          endif          ENDIF
   
 C       Limit range that KapGM can take  
         VisbeckK(i,j,myThid)=  
      &     min(VisbeckK(i,j,myThid),GM_Visbeck_maxval_K)  
195    
196  #endif /* GM_VISBECK_VARIABLE_K */  #endif /* GM_VISBECK_VARIABLE_K */
197    
198           ENDDO
199          ENDDO
200    
201    C-- end 1rst loop on vertical level index k
202          ENDDO
203    
204    
 #ifdef INCLUDE_DIAGNOSTICS_INTERFACE_CODE  
 C--     Time-average  
         GM_Kwx_T(i,j,k,bi,bj)=GM_Kwx_T(i,j,k,bi,bj)  
      &                       +Kwx(i,j,k,myThid)*deltaTclock  
         GM_Kwy_T(i,j,k,bi,bj)=GM_Kwy_T(i,j,k,bi,bj)  
      &                       +Kwy(i,j,k,myThid)*deltaTclock  
         GM_Kwz_T(i,j,k,bi,bj)=GM_Kwz_T(i,j,k,bi,bj)  
      &                       +Kwz(i,j,k,myThid)*deltaTclock  
205  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
206          IF (K.EQ.Nr)  #ifdef ALLOW_AUTODIFF_TAMC
207       &  Visbeck_K_T(i,j,bi,bj)=Visbeck_K_T(i,j,bi,bj)  CADJ STORE VisbeckK(:,:,bi,bj) = comlev1_bibj, key=igmkey, byte=isbyte
208       &                       +VisbeckK(i,j,myThid)*deltaTclock  #endif
209          IF ( GM_Visbeck_alpha.NE.0. ) THEN
210    C-     Limit range that KapGM can take
211           DO j=1-Oly+1,sNy+Oly-1
212            DO i=1-Olx+1,sNx+Olx-1
213             VisbeckK(i,j,bi,bj)=
214         &       MIN(VisbeckK(i,j,bi,bj),GM_Visbeck_maxval_K)
215    #ifdef ALLOW_TIMEAVE
216             Visbeck_K_T(i,j,bi,bj)=Visbeck_K_T(i,j,bi,bj)
217         &                         +VisbeckK(i,j,bi,bj)*deltaTclock
218  #endif  #endif
219            ENDDO
220         ENDDO         ENDDO
221        ENDDO        ENDIF
222        GM_TimeAve(k,bi,bj)=GM_TimeAve(k,bi,bj)+deltaTclock  #endif /* GM_VISBECK_VARIABLE_K */
223  #endif /* INCLUDE_DIAGNOSTICS_INTERFACE_CODE */  
224    
225    C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
226    
227    C-- 2nd loop on k : compute Tensor Coeff. at U,V levels.
228          DO k=1,Nr
229           kp1 = MIN(Nr,k+1)
230           maskp1 = 1. _d 0
231           IF (k.GE.Nr) maskp1 = 0. _d 0
232    
233    #ifdef ALLOW_AUTODIFF_TAMC
234           kkey = (igmkey-1)*Nr + k
235  #ifdef GM_NON_UNITY_DIAGONAL  #ifdef GM_NON_UNITY_DIAGONAL
236    CADJ STORE Kwx(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
237    CADJ STORE Kwy(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
238    CADJ STORE Kwz(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
239    #endif
240    #endif
241    
242    C-    express the Tensor in term of Diffusivity (= m**2 / s )
243          DO j=1-Oly+1,sNy+Oly-1
244           DO i=1-Olx+1,sNx+Olx-1
245            Kgm_tmp = GM_isopycK + GM_skewflx*GM_background_K
246    #ifdef GM_VISBECK_VARIABLE_K
247         &          + VisbeckK(i,j,bi,bj)*(1.+GM_skewflx)    
248    #endif
249            Kwx(i,j,k,bi,bj)= Kgm_tmp*Kwx(i,j,k,bi,bj)
250            Kwy(i,j,k,bi,bj)= Kgm_tmp*Kwy(i,j,k,bi,bj)
251            Kwz(i,j,k,bi,bj)= ( GM_isopycK
252    #ifdef GM_VISBECK_VARIABLE_K
253         &                    + VisbeckK(i,j,bi,bj)
254    #endif
255         &                    )*Kwz(i,j,k,bi,bj)
256           ENDDO
257          ENDDO
258    
259    #if ( defined (GM_NON_UNITY_DIAGONAL) || defined (GM_EXTRA_DIAGONAL) )
260    
261  C     Gradient of Sigma at U points  C     Gradient of Sigma at U points
262        DO j=1-Oly+1,sNy+Oly-1        DO j=1-Oly+1,sNy+Oly-1
263         DO i=1-Olx+1,sNx+Olx-1         DO i=1-Olx+1,sNx+Olx-1
264          SlopeX(i,j)=sigmaX(i,j,km1)          dSigmaDx(i,j)=sigmaX(i,j,k)
      &          *_maskW(i,j,k,bi,bj)  
         SlopeY(i,j)=0.25*( sigmaY(i-1,j+1,k) +sigmaY(i,j+1,k)  
      &                    +sigmaY(i-1, j ,k) +sigmaY(i, j ,k) )  
265       &          *_maskW(i,j,k,bi,bj)       &          *_maskW(i,j,k,bi,bj)
266          dSigmaDrReal(i,j)=0.25*( sigmaR(i-1,j, k ) +sigmaR(i,j, k )          dSigmaDy(i,j)=op25*( sigmaY(i-1,j+1,k) +sigmaY(i,j+1,k)
267       &                          +sigmaR(i-1,j,kp1) +sigmaR(i,j,kp1) )       &                      +sigmaY(i-1, j ,k) +sigmaY(i, j ,k) )
268       &          *_maskW(i,j,k,bi,bj)       &          *_maskW(i,j,k,bi,bj)
269            dSigmaDrReal(i,j)=op25*( sigmaR(i-1,j, k ) +sigmaR(i,j, k )
270         &                  +maskp1*(sigmaR(i-1,j,kp1) +sigmaR(i,j,kp1)) )
271         &          *_maskW(i,j,k,bi,bj)*maskp1
272         ENDDO         ENDDO
273        ENDDO        ENDDO
274    
275    #ifdef ALLOW_AUTODIFF_TAMC
276    CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
277    CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
278    CADJ STORE dsigmadrreal(:,:)   = comlev1_bibj_k, key=kkey, byte=isbyte
279    #endif /* ALLOW_AUTODIFF_TAMC */
280    
281  C     Calculate slopes for use in tensor, taper and/or clip  C     Calculate slopes for use in tensor, taper and/or clip
282        CALL GMREDI_SLOPE_LIMIT(        CALL GMREDI_SLOPE_LIMIT(
283       I             dSigmadRReal,       U             dSigmadRReal,
284       I             rF(K),       I             rF(K),K,
285       U             SlopeX, SlopeY,       U             SlopeX, SlopeY,
286       O             dRdSigmaLtd,       U             dSigmaDx, dSigmaDy,
287         O             SlopeSqr, taperFct,
288       I             bi, bj, myThid )       I             bi, bj, myThid )
289    
290        DO j=1-Oly+1,sNy+Oly-1  #ifdef GM_NON_UNITY_DIAGONAL
291         DO i=1-Olx+1,sNx+Olx-1          DO j=1-Oly+1,sNy+Oly-1
292          Kux(i,j,k,myThid)=(dSigmaDrReal(i,j)*dRdSigmaLtd(i,j))**2           DO i=1-Olx+1,sNx+Olx-1
293         ENDDO            Kux(i,j,k,bi,bj) =
294        ENDDO       &     ( GM_isopycK
295    #ifdef GM_VISBECK_VARIABLE_K
296         &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))
297    #endif
298         &     )
299         &     *taperFct(i,j)
300             ENDDO
301            ENDDO
302    #ifdef ALLOW_AUTODIFF_TAMC
303    # ifndef GM_TAPER_ORIG_CLIPPING
304    CADJ STORE Kux(:,:,k,bi,bj)  = comlev1_bibj_k, key=kkey, byte=isbyte
305    # endif
306    #endif
307            DO j=1-Oly+1,sNy+Oly-1
308             DO i=1-Olx+1,sNx+Olx-1
309              Kux(i,j,k,bi,bj) = MAX( Kux(i,j,k,bi,bj), GM_Kmin_horiz )
310             ENDDO
311            ENDDO
312    #endif /* GM_NON_UNITY_DIAGONAL */
313    
314    #ifdef GM_EXTRA_DIAGONAL
315    
316    #ifdef ALLOW_AUTODIFF_TAMC
317    CADJ STORE SlopeX(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
318    CADJ STORE taperFct(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte
319    #endif /* ALLOW_AUTODIFF_TAMC */
320          IF (GM_ExtraDiag) THEN
321            DO j=1-Oly+1,sNy+Oly-1
322             DO i=1-Olx+1,sNx+Olx-1
323              Kuz(i,j,k,bi,bj) =
324         &     ( GM_isopycK - GM_skewflx*GM_background_K
325    #ifdef GM_VISBECK_VARIABLE_K
326         &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))*GM_advect
327    #endif
328         &     )*SlopeX(i,j)*taperFct(i,j)
329             ENDDO
330            ENDDO
331          ENDIF
332    #endif /* GM_EXTRA_DIAGONAL */
333    
334  C     Gradient of Sigma at V points  C     Gradient of Sigma at V points
335        DO j=1-Oly+1,sNy+Oly-1        DO j=1-Oly+1,sNy+Oly-1
336         DO i=1-Olx+1,sNx+Olx-1         DO i=1-Olx+1,sNx+Olx-1
337          SlopeX(i,j)=0.25*( sigmaX(i, j ,k) +sigmaX(i+1, j ,k)          dSigmaDx(i,j)=op25*( sigmaX(i, j ,k) +sigmaX(i+1, j ,k)
338       &                    +sigmaX(i,j-1,k) +sigmaX(i+1,j-1,k) )       &                    +sigmaX(i,j-1,k) +sigmaX(i+1,j-1,k) )
339       &          *_maskS(i,j,k,bi,bj)       &          *_maskS(i,j,k,bi,bj)
340          SlopeY(i,j)=sigmaY(i,j,km1)          dSigmaDy(i,j)=sigmaY(i,j,k)
      &          *_maskS(i,j,k,bi,bj)  
         dSigmaDrReal(i,j)=0.25*( sigmaR(i,j-1, k ) +sigmaR(i,j, k )  
      &                          +sigmaR(i,j-1,kp1) +sigmaR(i,j,kp1) )  
341       &          *_maskS(i,j,k,bi,bj)       &          *_maskS(i,j,k,bi,bj)
342            dSigmaDrReal(i,j)=op25*( sigmaR(i,j-1, k ) +sigmaR(i,j, k )
343         &                  +maskp1*(sigmaR(i,j-1,kp1) +sigmaR(i,j,kp1)) )
344         &          *_maskS(i,j,k,bi,bj)*maskp1
345         ENDDO         ENDDO
346        ENDDO        ENDDO
347    
348    #ifdef ALLOW_AUTODIFF_TAMC
349    CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
350    CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
351    CADJ STORE dsigmadrreal(:,:)   = comlev1_bibj_k, key=kkey, byte=isbyte
352    #endif /* ALLOW_AUTODIFF_TAMC */
353    
354  C     Calculate slopes for use in tensor, taper and/or clip  C     Calculate slopes for use in tensor, taper and/or clip
355        CALL GMREDI_SLOPE_LIMIT(        CALL GMREDI_SLOPE_LIMIT(
356       I             dSigmadRReal,       U             dSigmadRReal,
357       I             rF(K),       I             rF(K),K,
358       U             SlopeX, SlopeY,       U             SlopeX, SlopeY,
359       O             dRdSigmaLtd,       U             dSigmaDx, dSigmaDy,
360         O             SlopeSqr, taperFct,
361       I             bi, bj, myThid )       I             bi, bj, myThid )
362    
363    #ifdef GM_NON_UNITY_DIAGONAL
364            DO j=1-Oly+1,sNy+Oly-1
365             DO i=1-Olx+1,sNx+Olx-1
366              Kvy(i,j,k,bi,bj) =
367         &     ( GM_isopycK
368    #ifdef GM_VISBECK_VARIABLE_K
369         &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))
370    #endif
371         &     )
372         &     *taperFct(i,j)
373             ENDDO
374            ENDDO
375    #ifdef ALLOW_AUTODIFF_TAMC
376    # ifndef GM_TAPER_ORIG_CLIPPING
377    CADJ STORE Kvy(:,:,k,bi,bj)  = comlev1_bibj_k, key=kkey, byte=isbyte
378    # endif
379    #endif
380            DO j=1-Oly+1,sNy+Oly-1
381             DO i=1-Olx+1,sNx+Olx-1
382              Kvy(i,j,k,bi,bj) = MAX( Kvy(i,j,k,bi,bj), GM_Kmin_horiz )
383             ENDDO
384            ENDDO
385    #endif /* GM_NON_UNITY_DIAGONAL */
386    
387    #ifdef GM_EXTRA_DIAGONAL
388    
389    #ifdef ALLOW_AUTODIFF_TAMC
390    CADJ STORE SlopeY(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
391    CADJ STORE taperFct(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte
392    #endif /* ALLOW_AUTODIFF_TAMC */
393          IF (GM_ExtraDiag) THEN
394            DO j=1-Oly+1,sNy+Oly-1
395             DO i=1-Olx+1,sNx+Olx-1
396              Kvz(i,j,k,bi,bj) =
397         &     ( GM_isopycK - GM_skewflx*GM_background_K
398    #ifdef GM_VISBECK_VARIABLE_K
399         &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))*GM_advect
400    #endif
401         &     )*SlopeY(i,j)*taperFct(i,j)
402             ENDDO
403            ENDDO
404          ENDIF
405    #endif /* GM_EXTRA_DIAGONAL */
406    
407    #endif /* GM_NON_UNITY_DIAGONAL || GM_EXTRA_DIAGONAL */
408    
409    #ifdef ALLOW_TIMEAVE
410    C--   Time-average
411        DO j=1-Oly+1,sNy+Oly-1        DO j=1-Oly+1,sNy+Oly-1
412         DO i=1-Olx+1,sNx+Olx-1         DO i=1-Olx+1,sNx+Olx-1
413          Kvy(i,j,k,myThid)=(dSigmaDrReal(i,j)*dRdSigmaLtd(i,j))**2          GM_Kwx_T(i,j,k,bi,bj)=GM_Kwx_T(i,j,k,bi,bj)
414         &                       +Kwx(i,j,k,bi,bj)*deltaTclock
415            GM_Kwy_T(i,j,k,bi,bj)=GM_Kwy_T(i,j,k,bi,bj)
416         &                       +Kwy(i,j,k,bi,bj)*deltaTclock
417            GM_Kwz_T(i,j,k,bi,bj)=GM_Kwz_T(i,j,k,bi,bj)
418         &                       +Kwz(i,j,k,bi,bj)*deltaTclock
419         ENDDO         ENDDO
420        ENDDO        ENDDO
421          GM_TimeAve(k,bi,bj)=GM_TimeAve(k,bi,bj)+deltaTclock
422    #endif /* ALLOW_TIMEAVE */
423    
424  #endif /* GM_NON_UNITY_DIAGONAL */  C-- end 2nd loop on vertical level index k
425          ENDDO
426    
427    
428    #ifdef GM_BOLUS_ADVEC
429          IF (GM_AdvForm) THEN
430            CALL GMREDI_CALC_PSI_B(
431         I             bi, bj, iMin, iMax, jMin, jMax,
432         I             sigmaX, sigmaY, sigmaR,
433         I             myThid )
434          ENDIF
435    #endif
436    
437  #endif /* ALLOW_GMREDI */  #endif /* ALLOW_GMREDI */
438    
# Line 219  C     Calculate slopes for use in tensor Line 441  C     Calculate slopes for use in tensor
441    
442    
443        SUBROUTINE GMREDI_CALC_TENSOR_DUMMY(        SUBROUTINE GMREDI_CALC_TENSOR_DUMMY(
444       I             bi, bj, iMin, iMax, jMin, jMax, K,       I             bi, bj, iMin, iMax, jMin, jMax,
445       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
446       I             myThid )       I             myThid )
447  C     /==========================================================\  C     /==========================================================\
# Line 242  C Line 464  C
464        _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
465        _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
466        _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
467        INTEGER bi,bj,iMin,iMax,jMin,jMax,K        INTEGER bi,bj,iMin,iMax,jMin,jMax
468        INTEGER myThid        INTEGER myThid
469  CEndOfInterface  CEndOfInterface
470    
471        INTEGER i, j        INTEGER i, j, k
472    
473  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
474    
475        DO j=1-Oly+1,sNy+Oly-1        DO k=1,Nr
476         DO i=1-Olx+1,sNx+Olx-1         DO j=1-Oly+1,sNy+Oly-1
477          Kwx(i,j,k,myThid) = 0.0          DO i=1-Olx+1,sNx+Olx-1
478          Kwy(i,j,k,myThid) = 0.0           Kwx(i,j,k,bi,bj) = 0.0
479          Kwz(i,j,k,myThid) = 0.0           Kwy(i,j,k,bi,bj) = 0.0
480             Kwz(i,j,k,bi,bj) = 0.0
481            ENDDO
482         ENDDO         ENDDO
483        ENDDO        ENDDO
484  #endif /* ALLOW_GMREDI */  #endif /* ALLOW_GMREDI */
485    
486        end        RETURN
487          END
Legend:
Removed from v.1.2  
changed lines
  Added in v.1.14
  ViewVC Help
Powered by ViewVC 1.1.22