/[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.7 by heimbach, Tue Aug 21 15:27:19 2001 UTC revision 1.18 by edhill, Mon Sep 29 19:24:31 2003 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,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, 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    
85  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
86  !HPF$ INDEPENDENT         kkey = (igmkey-1)*Nr + k
87           DO j=1-Oly,sNy+Oly
88            DO i=1-Olx,sNx+Olx
89             SlopeX(i,j)       = 0. _d 0
90             SlopeY(i,j)       = 0. _d 0
91             dSigmaDx(i,j)     = 0. _d 0
92             dSigmaDy(i,j)     = 0. _d 0
93             dSigmaDrReal(i,j) = 0. _d 0
94             SlopeSqr(i,j)     = 0. _d 0
95             taperFct(i,j)     = 0. _d 0
96             Kwx(i,j,k,bi,bj)  = 0. _d 0
97             Kwy(i,j,k,bi,bj)  = 0. _d 0
98             Kwz(i,j,k,bi,bj)  = 0. _d 0
99    # ifdef GM_NON_UNITY_DIAGONAL
100             Kux(i,j,k,bi,bj)  = 0. _d 0
101             Kvy(i,j,k,bi,bj)  = 0. _d 0
102    # endif
103    # ifdef GM_EXTRA_DIAGONAL
104             Kuz(i,j,k,bi,bj)  = 0. _d 0
105             Kvz(i,j,k,bi,bj)  = 0. _d 0
106    # endif
107    # ifdef GM_BOLUS_ADVEC
108             GM_PsiX(i,j,k,bi,bj)  = 0. _d 0
109             GM_PsiY(i,j,k,bi,bj)  = 0. _d 0
110    # endif
111            ENDDO
112           ENDDO
113  #endif  #endif
114    
115        DO j=1-Oly+1,sNy+Oly-1        DO j=1-Oly+1,sNy+Oly-1
 #ifdef ALLOW_AUTODIFF_TAMC  
 !HPF$ INDEPENDENT  
 #endif  
116         DO i=1-Olx+1,sNx+Olx-1         DO i=1-Olx+1,sNx+Olx-1
   
117  C      Gradient of Sigma at rVel points  C      Gradient of Sigma at rVel points
118          SlopeX(i,j)=0.25*( sigmaX(i+1, j ,km1) +sigmaX(i,j,km1)          dSigmaDx(i,j)=op25*( sigmaX(i+1, j ,k-1) +sigmaX(i,j,k-1)
119       &                    +sigmaX(i+1, j , k ) +sigmaX(i,j, k ) )       &                    +sigmaX(i+1, j , k ) +sigmaX(i,j, k ) )
120          SlopeY(i,j)=0.25*( sigmaY( i ,j+1,km1) +sigmaY(i,j,km1)       &                  *maskC(i,j,k,bi,bj)
121            dSigmaDy(i,j)=op25*( sigmaY( i ,j+1,k-1) +sigmaY(i,j,k-1)
122       &                    +sigmaY( i ,j+1, k ) +sigmaY(i,j, k ) )       &                    +sigmaY( i ,j+1, k ) +sigmaY(i,j, k ) )
123         &                  *maskC(i,j,k,bi,bj)
124          dSigmaDrReal(i,j)=sigmaR(i,j,k)          dSigmaDrReal(i,j)=sigmaR(i,j,k)
   
         if (hFacC(i,j,k,bi,bj).eq.0.) then  
          SlopeX(i,j)=0.  
          SlopeY(i,j)=0.  
         endif  
   
125         ENDDO         ENDDO
126        ENDDO        ENDDO
127    
128    #ifdef ALLOW_AUTODIFF_TAMC
129    CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
130    CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
131    CADJ STORE dsigmadrreal(:,:)   = comlev1_bibj_k, key=kkey, byte=isbyte
132    #endif /* ALLOW_AUTODIFF_TAMC */
133    
134  C     Calculate slopes for use in tensor, taper and/or clip  C     Calculate slopes for use in tensor, taper and/or clip
135        CALL GMREDI_SLOPE_LIMIT(        CALL GMREDI_SLOPE_LIMIT(
136       I             dSigmadRReal,       U             dSigmadRReal,
137       I             rF(K),       I             rF(K),K,
138       U             SlopeX, SlopeY,       U             SlopeX, SlopeY,
139       O             dRdSigmaLtd,       U             dSigmaDx, dSigmaDy,
140         O             SlopeSqr, taperFct,
141       I             bi, bj, myThid )       I             bi, bj, myThid )
142    
143        DO j=1-Oly+1,sNy+Oly-1        DO j=1-Oly+1,sNy+Oly-1
144         DO i=1-Olx+1,sNx+Olx-1         DO i=1-Olx+1,sNx+Olx-1
145    
146  C       Mask Iso-neutral slopes  C       Mask Iso-neutral slopes
147          if (hFacC(i,j,k,bi,bj).eq.0.) then          SlopeX(i,j)=SlopeX(i,j)*maskC(i,j,k,bi,bj)
148           SlopeX(i,j)=0.          SlopeY(i,j)=SlopeY(i,j)*maskC(i,j,k,bi,bj)
149           SlopeY(i,j)=0.          SlopeSqr(i,j)=SlopeSqr(i,j)*maskC(i,j,k,bi,bj)
150          endif  
151          Ssq=SlopeX(i,j)*SlopeX(i,j)+SlopeY(i,j)*SlopeY(i,j)         ENDDO
152          ENDDO
153  C       Components of Redi/GM tensor  
154          Kwx(i,j,k,bi,bj)=2.*SlopeX(i,j)  #ifdef ALLOW_AUTODIFF_TAMC
155          Kwy(i,j,k,bi,bj)=2.*SlopeY(i,j)  CADJ STORE SlopeX(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
156          Kwz(i,j,k,bi,bj)=Ssq  CADJ STORE SlopeY(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
157    CADJ STORE SlopeSqr(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte
158    CADJ STORE dsigmadrreal(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
159    CADJ STORE taperfct(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte
160    #endif /* ALLOW_AUTODIFF_TAMC */
161    
162          DO j=1-Oly+1,sNy+Oly-1
163           DO i=1-Olx+1,sNx+Olx-1
164    
165    C       Components of Redi/GM tensor
166            Kwx(i,j,k,bi,bj)= SlopeX(i,j)*taperFct(i,j)
167            Kwy(i,j,k,bi,bj)= SlopeY(i,j)*taperFct(i,j)
168            Kwz(i,j,k,bi,bj)= SlopeSqr(i,j)*taperFct(i,j)
169    
170  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
171    
172    C- note (jmc) : moved here since only used in VISBECK_VARIABLE_K
173    C           but do not know if *taperFct (or **2 ?) is necessary
174            Ssq(i,j)=SlopeSqr(i,j)*taperFct(i,j)
175    
176  C--     Depth average of M^2/N^2 * N  C--     Depth average of M^2/N^2 * N
177    
178  C       Calculate terms for mean Richardson number  C       Calculate terms for mean Richardson number
# Line 115  C       If negative then we are below th Line 186  C       If negative then we are below th
186  C       Now we convert deltaH to a non-dimensional fraction  C       Now we convert deltaH to a non-dimensional fraction
187          deltaH=deltaH/GM_Visbeck_depth          deltaH=deltaH/GM_Visbeck_depth
188    
189          if (K.eq.2) VisbeckK(i,j,bi,bj)=0.          IF (K.eq.2) VisbeckK(i,j,bi,bj)=0.
190  Calt?   if (dSigmaDrReal(i,j).NE.0.) then          IF ( Ssq(i,j).NE.0. .AND. dSigmaDrReal(i,j).NE.0. ) THEN
191  Calt?    N2=(-Gravity*recip_Rhonil)*dSigmaDrReal(i,j)           N2= -Gravity*recip_RhoConst*dSigmaDrReal(i,j)
192          if ( dRdSigmaLtd(i,j).NE.0. .AND. Ssq.NE.0. ) then           SN=sqrt(Ssq(i,j)*N2)
          N2=(-Gravity*recip_Rhonil)/dRdSigmaLtd(i,j)  
          SN=sqrt(Ssq*N2)  
193           VisbeckK(i,j,bi,bj)=VisbeckK(i,j,bi,bj)+deltaH           VisbeckK(i,j,bi,bj)=VisbeckK(i,j,bi,bj)+deltaH
194       &      *GM_Visbeck_alpha*GM_Visbeck_length*GM_Visbeck_length*SN       &      *GM_Visbeck_alpha*GM_Visbeck_length*GM_Visbeck_length*SN
195          endif          ENDIF
   
 C       Limit range that KapGM can take  
         VisbeckK(i,j,bi,bj)=  
      &     min(VisbeckK(i,j,bi,bj),GM_Visbeck_maxval_K)  
196    
197  #endif /* GM_VISBECK_VARIABLE_K */  #endif /* GM_VISBECK_VARIABLE_K */
198    
199           ENDDO
200          ENDDO
201    
202    C-- end 1rst loop on vertical level index k
203          ENDDO
204    
205    
 #ifdef ALLOW_TIMEAVE  
 C--     Time-average  
         GM_Kwx_T(i,j,k,bi,bj)=GM_Kwx_T(i,j,k,bi,bj)  
      &                       +Kwx(i,j,k,bi,bj)*deltaTclock  
         GM_Kwy_T(i,j,k,bi,bj)=GM_Kwy_T(i,j,k,bi,bj)  
      &                       +Kwy(i,j,k,bi,bj)*deltaTclock  
         GM_Kwz_T(i,j,k,bi,bj)=GM_Kwz_T(i,j,k,bi,bj)  
      &                       +Kwz(i,j,k,bi,bj)*deltaTclock  
206  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
207          IF (K.EQ.Nr)  #ifdef ALLOW_AUTODIFF_TAMC
208       &  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
      &                       +VisbeckK(i,j,bi,bj)*deltaTclock  
209  #endif  #endif
210  #endif /* ALLOW_TIMEAVE */        IF ( GM_Visbeck_alpha.NE.0. ) THEN
211    C-     Limit range that KapGM can take
212           DO j=1-Oly+1,sNy+Oly-1
213            DO i=1-Olx+1,sNx+Olx-1
214             VisbeckK(i,j,bi,bj)=
215         &       MIN(VisbeckK(i,j,bi,bj),GM_Visbeck_maxval_K)
216    #ifdef ALLOW_TIMEAVE
217             Visbeck_K_T(i,j,bi,bj)=Visbeck_K_T(i,j,bi,bj)
218         &                         +VisbeckK(i,j,bi,bj)*deltaTclock
219    #endif
220            ENDDO
221         ENDDO         ENDDO
222        ENDDO        ENDIF
223    cph( NEW
224    #ifdef ALLOW_AUTODIFF_TAMC
225    CADJ STORE VisbeckK(:,:,bi,bj) = comlev1_bibj, key=igmkey, byte=isbyte
226    #endif
227    cph)
228    #endif /* GM_VISBECK_VARIABLE_K */
229    
230  #ifdef ALLOW_TIMEAVE  
231        GM_TimeAve(k,bi,bj)=GM_TimeAve(k,bi,bj)+deltaTclock  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    #if (defined (GM_NON_UNITY_DIAGONAL) || \
242         defined (GM_VISBECK_VARIABLE_K))
243    CADJ STORE Kwx(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
244    CADJ STORE Kwy(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
245    CADJ STORE Kwz(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
246    #endif
247    #endif
248    
249    C-    express the Tensor in term of Diffusivity (= m**2 / s )
250          DO j=1-Oly+1,sNy+Oly-1
251           DO i=1-Olx+1,sNx+Olx-1
252            Kgm_tmp = GM_isopycK + GM_skewflx*GM_background_K
253    #ifdef GM_VISBECK_VARIABLE_K
254         &          + VisbeckK(i,j,bi,bj)*(1. _d 0 + GM_skewflx)    
255    #endif
256            Kwx(i,j,k,bi,bj)= Kgm_tmp*Kwx(i,j,k,bi,bj)
257            Kwy(i,j,k,bi,bj)= Kgm_tmp*Kwy(i,j,k,bi,bj)
258            Kwz(i,j,k,bi,bj)= ( GM_isopycK
259    #ifdef GM_VISBECK_VARIABLE_K
260         &                    + VisbeckK(i,j,bi,bj)
261  #endif  #endif
262         &                    )*Kwz(i,j,k,bi,bj)
263           ENDDO
264          ENDDO
265    
266    #if ( defined (GM_NON_UNITY_DIAGONAL) || defined (GM_EXTRA_DIAGONAL) )
267    
 #ifdef GM_NON_UNITY_DIAGONAL  
268  C     Gradient of Sigma at U points  C     Gradient of Sigma at U points
269        DO j=1-Oly+1,sNy+Oly-1        DO j=1-Oly+1,sNy+Oly-1
270         DO i=1-Olx+1,sNx+Olx-1         DO i=1-Olx+1,sNx+Olx-1
271          SlopeX(i,j)=sigmaX(i,j,km1)          dSigmaDx(i,j)=sigmaX(i,j,k)
272       &          *_maskW(i,j,k,bi,bj)       &          *_maskW(i,j,k,bi,bj)
273          SlopeY(i,j)=0.25*( sigmaY(i-1,j+1,k) +sigmaY(i,j+1,k)          dSigmaDy(i,j)=op25*( sigmaY(i-1,j+1,k) +sigmaY(i,j+1,k)
274       &                    +sigmaY(i-1, j ,k) +sigmaY(i, j ,k) )       &                      +sigmaY(i-1, j ,k) +sigmaY(i, j ,k) )
275       &          *_maskW(i,j,k,bi,bj)       &          *_maskW(i,j,k,bi,bj)
276          dSigmaDrReal(i,j)=0.25*( sigmaR(i-1,j, k ) +sigmaR(i,j, k )          dSigmaDrReal(i,j)=op25*( sigmaR(i-1,j, k ) +sigmaR(i,j, k )
277       &                          +sigmaR(i-1,j,kp1) +sigmaR(i,j,kp1) )       &                  +maskp1*(sigmaR(i-1,j,kp1) +sigmaR(i,j,kp1)) )
278       &          *_maskW(i,j,k,bi,bj)       &          *_maskW(i,j,k,bi,bj)
279         ENDDO         ENDDO
280        ENDDO        ENDDO
281    
282    #ifdef ALLOW_AUTODIFF_TAMC
283    CADJ STORE SlopeSqr(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
284    CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
285    CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
286    CADJ STORE dsigmadrreal(:,:)   = comlev1_bibj_k, key=kkey, byte=isbyte
287    #endif /* ALLOW_AUTODIFF_TAMC */
288    
289  C     Calculate slopes for use in tensor, taper and/or clip  C     Calculate slopes for use in tensor, taper and/or clip
290        CALL GMREDI_SLOPE_LIMIT(        CALL GMREDI_SLOPE_LIMIT(
291       I             dSigmadRReal,       U             dSigmadRReal,
292       I             rF(K),       I             rF(K),K,
293       U             SlopeX, SlopeY,       U             SlopeX, SlopeY,
294       O             dRdSigmaLtd,       U             dSigmaDx, dSigmaDy,
295         O             SlopeSqr, taperFct,
296       I             bi, bj, myThid )       I             bi, bj, myThid )
297    
298        DO j=1-Oly+1,sNy+Oly-1  cph( NEW
299         DO i=1-Olx+1,sNx+Olx-1  #ifdef ALLOW_AUTODIFF_TAMC
300          Kux(i,j,k,bi,bj)=(dSigmaDrReal(i,j)*dRdSigmaLtd(i,j))**2  cph(
301         ENDDO  CADJ STORE SlopeSqr(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
302        ENDDO  CADJ STORE taperfct(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
303    cph)
304    #endif /* ALLOW_AUTODIFF_TAMC */
305    cph)
306    
307    #ifdef GM_NON_UNITY_DIAGONAL
308            DO j=1-Oly+1,sNy+Oly-1
309             DO i=1-Olx+1,sNx+Olx-1
310              Kux(i,j,k,bi,bj) =
311         &     ( GM_isopycK
312    #ifdef GM_VISBECK_VARIABLE_K
313         &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))
314    #endif
315         &     )
316         &     *taperFct(i,j)
317             ENDDO
318            ENDDO
319    #ifdef ALLOW_AUTODIFF_TAMC
320    # ifdef GM_EXCLUDE_CLIPPING
321    CADJ STORE Kux(:,:,k,bi,bj)  = comlev1_bibj_k, key=kkey, byte=isbyte
322    # endif
323    #endif
324            DO j=1-Oly+1,sNy+Oly-1
325             DO i=1-Olx+1,sNx+Olx-1
326              Kux(i,j,k,bi,bj) = MAX( Kux(i,j,k,bi,bj), GM_Kmin_horiz )
327             ENDDO
328            ENDDO
329    #endif /* GM_NON_UNITY_DIAGONAL */
330    
331    #ifdef GM_EXTRA_DIAGONAL
332    
333    #ifdef ALLOW_AUTODIFF_TAMC
334    CADJ STORE SlopeX(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
335    CADJ STORE taperFct(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte
336    #endif /* ALLOW_AUTODIFF_TAMC */
337          IF (GM_ExtraDiag) THEN
338            DO j=1-Oly+1,sNy+Oly-1
339             DO i=1-Olx+1,sNx+Olx-1
340              Kuz(i,j,k,bi,bj) =
341         &     ( GM_isopycK - GM_skewflx*GM_background_K
342    #ifdef GM_VISBECK_VARIABLE_K
343         &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))*GM_advect
344    #endif
345         &     )*SlopeX(i,j)*taperFct(i,j)
346             ENDDO
347            ENDDO
348          ENDIF
349    #endif /* GM_EXTRA_DIAGONAL */
350    
351  C     Gradient of Sigma at V points  C     Gradient of Sigma at V points
352        DO j=1-Oly+1,sNy+Oly-1        DO j=1-Oly+1,sNy+Oly-1
353         DO i=1-Olx+1,sNx+Olx-1         DO i=1-Olx+1,sNx+Olx-1
354          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)
355       &                    +sigmaX(i,j-1,k) +sigmaX(i+1,j-1,k) )       &                    +sigmaX(i,j-1,k) +sigmaX(i+1,j-1,k) )
356       &          *_maskS(i,j,k,bi,bj)       &          *_maskS(i,j,k,bi,bj)
357          SlopeY(i,j)=sigmaY(i,j,km1)          dSigmaDy(i,j)=sigmaY(i,j,k)
358       &          *_maskS(i,j,k,bi,bj)       &          *_maskS(i,j,k,bi,bj)
359          dSigmaDrReal(i,j)=0.25*( sigmaR(i,j-1, k ) +sigmaR(i,j, k )          dSigmaDrReal(i,j)=op25*( sigmaR(i,j-1, k ) +sigmaR(i,j, k )
360       &                          +sigmaR(i,j-1,kp1) +sigmaR(i,j,kp1) )       &                  +maskp1*(sigmaR(i,j-1,kp1) +sigmaR(i,j,kp1)) )
361       &          *_maskS(i,j,k,bi,bj)       &          *_maskS(i,j,k,bi,bj)
362         ENDDO         ENDDO
363        ENDDO        ENDDO
364    
365    #ifdef ALLOW_AUTODIFF_TAMC
366    CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
367    CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
368    CADJ STORE dsigmadrreal(:,:)   = comlev1_bibj_k, key=kkey, byte=isbyte
369    #endif /* ALLOW_AUTODIFF_TAMC */
370    
371  C     Calculate slopes for use in tensor, taper and/or clip  C     Calculate slopes for use in tensor, taper and/or clip
372        CALL GMREDI_SLOPE_LIMIT(        CALL GMREDI_SLOPE_LIMIT(
373       I             dSigmadRReal,       U             dSigmadRReal,
374       I             rF(K),       I             rF(K),K,
375       U             SlopeX, SlopeY,       U             SlopeX, SlopeY,
376       O             dRdSigmaLtd,       U             dSigmaDx, dSigmaDy,
377         O             SlopeSqr, taperFct,
378       I             bi, bj, myThid )       I             bi, bj, myThid )
379    
380    cph(
381    #ifdef ALLOW_AUTODIFF_TAMC
382    cph(
383    CADJ STORE taperfct(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
384    cph)
385    #endif /* ALLOW_AUTODIFF_TAMC */
386    cph)
387    
388    #ifdef GM_NON_UNITY_DIAGONAL
389            DO j=1-Oly+1,sNy+Oly-1
390             DO i=1-Olx+1,sNx+Olx-1
391              Kvy(i,j,k,bi,bj) =
392         &     ( GM_isopycK
393    #ifdef GM_VISBECK_VARIABLE_K
394         &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))
395    #endif
396         &     )
397         &     *taperFct(i,j)
398             ENDDO
399            ENDDO
400    #ifdef ALLOW_AUTODIFF_TAMC
401    # ifdef GM_EXCLUDE_CLIPPING
402    CADJ STORE Kvy(:,:,k,bi,bj)  = comlev1_bibj_k, key=kkey, byte=isbyte
403    # endif
404    #endif
405            DO j=1-Oly+1,sNy+Oly-1
406             DO i=1-Olx+1,sNx+Olx-1
407              Kvy(i,j,k,bi,bj) = MAX( Kvy(i,j,k,bi,bj), GM_Kmin_horiz )
408             ENDDO
409            ENDDO
410    #endif /* GM_NON_UNITY_DIAGONAL */
411    
412    #ifdef GM_EXTRA_DIAGONAL
413    
414    #ifdef ALLOW_AUTODIFF_TAMC
415    CADJ STORE SlopeY(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
416    CADJ STORE taperFct(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte
417    #endif /* ALLOW_AUTODIFF_TAMC */
418          IF (GM_ExtraDiag) THEN
419            DO j=1-Oly+1,sNy+Oly-1
420             DO i=1-Olx+1,sNx+Olx-1
421              Kvz(i,j,k,bi,bj) =
422         &     ( GM_isopycK - GM_skewflx*GM_background_K
423    #ifdef GM_VISBECK_VARIABLE_K
424         &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))*GM_advect
425    #endif
426         &     )*SlopeY(i,j)*taperFct(i,j)
427             ENDDO
428            ENDDO
429          ENDIF
430    #endif /* GM_EXTRA_DIAGONAL */
431    
432    #endif /* GM_NON_UNITY_DIAGONAL || GM_EXTRA_DIAGONAL */
433    
434    #ifdef ALLOW_TIMEAVE
435    C--   Time-average
436        DO j=1-Oly+1,sNy+Oly-1        DO j=1-Oly+1,sNy+Oly-1
437         DO i=1-Olx+1,sNx+Olx-1         DO i=1-Olx+1,sNx+Olx-1
438          Kvy(i,j,k,bi,bj)=(dSigmaDrReal(i,j)*dRdSigmaLtd(i,j))**2          GM_Kwx_T(i,j,k,bi,bj)=GM_Kwx_T(i,j,k,bi,bj)
439         &                       +Kwx(i,j,k,bi,bj)*deltaTclock
440            GM_Kwy_T(i,j,k,bi,bj)=GM_Kwy_T(i,j,k,bi,bj)
441         &                       +Kwy(i,j,k,bi,bj)*deltaTclock
442            GM_Kwz_T(i,j,k,bi,bj)=GM_Kwz_T(i,j,k,bi,bj)
443         &                       +Kwz(i,j,k,bi,bj)*deltaTclock
444         ENDDO         ENDDO
445        ENDDO        ENDDO
446          GM_TimeAve(k,bi,bj)=GM_TimeAve(k,bi,bj)+deltaTclock
447    #endif /* ALLOW_TIMEAVE */
448    
449  #endif /* GM_NON_UNITY_DIAGONAL */  C-- end 2nd loop on vertical level index k
450          ENDDO
451    
452    
453    #ifdef GM_BOLUS_ADVEC
454          IF (GM_AdvForm) THEN
455            CALL GMREDI_CALC_PSI_B(
456         I             bi, bj, iMin, iMax, jMin, jMax,
457         I             sigmaX, sigmaY, sigmaR,
458         I             myThid )
459          ENDIF
460    #endif
461    
462  #endif /* ALLOW_GMREDI */  #endif /* ALLOW_GMREDI */
463    
# Line 222  C     Calculate slopes for use in tensor Line 466  C     Calculate slopes for use in tensor
466    
467    
468        SUBROUTINE GMREDI_CALC_TENSOR_DUMMY(        SUBROUTINE GMREDI_CALC_TENSOR_DUMMY(
469       I             bi, bj, iMin, iMax, jMin, jMax, K,       I             bi, bj, iMin, iMax, jMin, jMax,
470       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
471       I             myThid )       I             myThid )
472  C     /==========================================================\  C     /==========================================================\
# Line 245  C Line 489  C
489        _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
490        _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
491        _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
492        INTEGER bi,bj,iMin,iMax,jMin,jMax,K        INTEGER bi,bj,iMin,iMax,jMin,jMax
493        INTEGER myThid        INTEGER myThid
494  CEndOfInterface  CEndOfInterface
495    
496        INTEGER i, j        INTEGER i, j, k
497    
498  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
499    
500        DO j=1-Oly+1,sNy+Oly-1        DO k=1,Nr
501         DO i=1-Olx+1,sNx+Olx-1         DO j=1-Oly+1,sNy+Oly-1
502          Kwx(i,j,k,bi,bj) = 0.0          DO i=1-Olx+1,sNx+Olx-1
503          Kwy(i,j,k,bi,bj) = 0.0           Kwx(i,j,k,bi,bj) = 0.0
504          Kwz(i,j,k,bi,bj) = 0.0           Kwy(i,j,k,bi,bj) = 0.0
505             Kwz(i,j,k,bi,bj) = 0.0
506            ENDDO
507         ENDDO         ENDDO
508        ENDDO        ENDDO
509  #endif /* ALLOW_GMREDI */  #endif /* ALLOW_GMREDI */
510    
511        end        RETURN
512          END
Legend:
Removed from v.1.7  
changed lines
  Added in v.1.18
  ViewVC Help
Powered by ViewVC 1.1.22