/[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.23 by heimbach, Wed Jun 7 01:55:14 2006 UTC revision 1.28 by heimbach, Tue Jun 26 15:34:15 2007 UTC
# Line 2  C $Header$ Line 2  C $Header$
2  C $Name$  C $Name$
3    
4  #include "GMREDI_OPTIONS.h"  #include "GMREDI_OPTIONS.h"
5    #ifdef ALLOW_KPP
6    # include "KPP_OPTIONS.h"
7    #endif
8    #undef OLD_VISBECK_CALC
9    
10  CStartOfInterface  CBOP
11    C     !ROUTINE: GMREDI_CALC_TENSOR
12    C     !INTERFACE:
13        SUBROUTINE GMREDI_CALC_TENSOR(        SUBROUTINE GMREDI_CALC_TENSOR(
14       I             bi, bj, iMin, iMax, jMin, jMax,       I             iMin, iMax, jMin, jMax,
15       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
16       I             myThid )       I             bi, bj, myTime, myIter, myThid )
17  C     /==========================================================\  
18  C     | SUBROUTINE GMREDI_CALC_TENSOR                            |  C     !DESCRIPTION: \bv
19  C     | o Calculate tensor elements for GM/Redi tensor.          |  C     *==========================================================*
20  C     |==========================================================|  C     | SUBROUTINE GMREDI_CALC_TENSOR
21  C     \==========================================================/  C     | o Calculate tensor elements for GM/Redi tensor.
22    C     *==========================================================*
23    C     *==========================================================*
24    C     \ev
25    
26    C     !USES:
27        IMPLICIT NONE        IMPLICIT NONE
28    
29  C     == Global variables ==  C     == Global variables ==
# Line 23  C     == Global variables == Line 34  C     == Global variables ==
34  #include "PARAMS.h"  #include "PARAMS.h"
35  #include "GMREDI.h"  #include "GMREDI.h"
36  #include "GMREDI_TAVE.h"  #include "GMREDI_TAVE.h"
37    #ifdef ALLOW_KPP
38    # include "KPP.h"
39    #endif
40    
41  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
42  #include "tamc.h"  #include "tamc.h"
43  #include "tamc_keys.h"  #include "tamc_keys.h"
44  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
45    
46    C     !INPUT/OUTPUT PARAMETERS:
47  C     == Routine arguments ==  C     == Routine arguments ==
48    C     bi, bj    :: tile indices
49    C     myTime    :: Current time in simulation
50    C     myIter    :: Current iteration number in simulation
51    C     myThid    :: My Thread Id. number
52  C  C
53          INTEGER iMin,iMax,jMin,jMax
54        _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
55        _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
56        _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
57        INTEGER bi,bj,iMin,iMax,jMin,jMax        INTEGER bi, bj
58          _RL     myTime
59          INTEGER myIter
60        INTEGER myThid        INTEGER myThid
61  CEndOfInterface  CEOP
62    
63  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
64    
65    C     !LOCAL VARIABLES:
66  C     == Local variables ==  C     == Local variables ==
67        INTEGER i,j,k,kp1        INTEGER i,j,k,kp1
68        _RL SlopeX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL SlopeX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
# Line 55  C     == Local variables == Line 78  C     == Local variables ==
78        _RL ldd97_LrhoS(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL ldd97_LrhoS(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
79        _RL Cspd, LrhoInf, LrhoSup, fCoriLoc        _RL Cspd, LrhoInf, LrhoSup, fCoriLoc
80    
81          INTEGER kLow_W (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
82          INTEGER kLow_S (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
83          _RL locMixLayer(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
84          _RL baseSlope  (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
85          _RL hTransLay  (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
86          _RL recipLambda(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
87    
88  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
89    #ifdef OLD_VISBECK_CALC
90        _RL deltaH,zero_rs        _RL deltaH,zero_rs
91        PARAMETER(zero_rs=0.D0)        PARAMETER(zero_rs=0.D0)
92        _RL N2,SN        _RL N2,SN
93        _RL Ssq(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL Ssq(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
94    #else
95          _RL dSigmaH
96          _RL deltaH, integrDepth
97          _RL Sloc, M2loc, SNloc
98    #endif
99  #endif  #endif
100    
101  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
# Line 90  C---+----1----+----2----+----3----+----4 Line 126  C---+----1----+----2----+----3----+----4
126        doDiagRediFlx = .FALSE.        doDiagRediFlx = .FALSE.
127        IF ( useDiagnostics ) THEN        IF ( useDiagnostics ) THEN
128          doDiagRediFlx = DIAGNOSTICS_IS_ON('GM_KuzTz', myThid )          doDiagRediFlx = DIAGNOSTICS_IS_ON('GM_KuzTz', myThid )
129          doDiagRediFlx = doDiagRediFlx .OR.          doDiagRediFlx = doDiagRediFlx .OR.
130       &                  DIAGNOSTICS_IS_ON('GM_KvzTz', myThid )       &                  DIAGNOSTICS_IS_ON('GM_KvzTz', myThid )
131        ENDIF        ENDIF
132  #endif  #endif
133        
134  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
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
# Line 104  C---+----1----+----2----+----3----+----4 Line 140  C---+----1----+----2----+----3----+----4
140  #endif  #endif
141    
142  C--   set ldd97_Lrho (for tapering scheme ldd97):  C--   set ldd97_Lrho (for tapering scheme ldd97):
143        IF (GM_taper_scheme.EQ.'ldd97') THEN        IF ( GM_taper_scheme.EQ.'ldd97' .OR.
144         &     GM_taper_scheme.EQ.'fm07' ) THEN
145         Cspd = 2. _d 0         Cspd = 2. _d 0
146         LrhoInf = 15. _d 3         LrhoInf = 15. _d 3
147         LrhoSup = 100. _d 3         LrhoSup = 100. _d 3
# Line 121  C-     Tracer point location (center): Line 158  C-     Tracer point location (center):
158         ENDDO         ENDDO
159  C-     U point location (West):  C-     U point location (West):
160         DO j=1-Oly,sNy+Oly         DO j=1-Oly,sNy+Oly
161            kLow_W(1-Olx,j) = 0
162          ldd97_LrhoW(1-Olx,j) = LrhoSup          ldd97_LrhoW(1-Olx,j) = LrhoSup
163          DO i=1-Olx+1,sNx+Olx          DO i=1-Olx+1,sNx+Olx
164             kLow_W(i,j) = MIN(kLowC(i-1,j,bi,bj),kLowC(i,j,bi,bj))
165           fCoriLoc = op5*(fCori(i-1,j,bi,bj)+fCori(i,j,bi,bj))           fCoriLoc = op5*(fCori(i-1,j,bi,bj)+fCori(i,j,bi,bj))
166           IF (fCoriLoc.NE.0.) THEN           IF (fCoriLoc.NE.0.) THEN
167             ldd97_LrhoW(i,j) = Cspd/ABS(fCoriLoc)             ldd97_LrhoW(i,j) = Cspd/ABS(fCoriLoc)
# Line 134  C-     U point location (West): Line 173  C-     U point location (West):
173         ENDDO         ENDDO
174  C-     V point location (South):  C-     V point location (South):
175         DO i=1-Olx+1,sNx+Olx         DO i=1-Olx+1,sNx+Olx
176             kLow_S(i,1-Oly) = 0
177           ldd97_LrhoS(i,1-Oly) = LrhoSup           ldd97_LrhoS(i,1-Oly) = LrhoSup
178         ENDDO         ENDDO
179         DO j=1-Oly+1,sNy+Oly         DO j=1-Oly+1,sNy+Oly
180          DO i=1-Olx,sNx+Olx          DO i=1-Olx,sNx+Olx
181             kLow_S(i,j) = MIN(kLowC(i,j-1,bi,bj),kLowC(i,j,bi,bj))
182           fCoriLoc = op5*(fCori(i,j-1,bi,bj)+fCori(i,j,bi,bj))           fCoriLoc = op5*(fCori(i,j-1,bi,bj)+fCori(i,j,bi,bj))
183           IF (fCoriLoc.NE.0.) THEN           IF (fCoriLoc.NE.0.) THEN
184             ldd97_LrhoS(i,j) = Cspd/ABS(fCoriLoc)             ldd97_LrhoS(i,j) = Cspd/ABS(fCoriLoc)
# Line 157  C-    Just initialize to zero (not use a Line 198  C-    Just initialize to zero (not use a
198          ENDDO          ENDDO
199         ENDDO         ENDDO
200        ENDIF        ENDIF
 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|  
201    
202        DO k=2,Nr  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
203  C-- 1rst loop on k : compute Tensor Coeff. at W points.  C-- 1rst loop on k : compute Tensor Coeff. at W points.
204    
205          DO j=1-Oly,sNy+Oly
206           DO i=1-Olx,sNx+Olx
207             hTransLay(i,j) = R_low(i,j,bi,bj)
208             baseSlope(i,j) =  0. _d 0
209             recipLambda(i,j) = 0. _d 0
210             locMixLayer(i,j) = 0. _d 0
211           ENDDO
212          ENDDO
213    #ifdef ALLOW_KPP
214          IF ( useKPP ) THEN
215           DO j=1-Oly,sNy+Oly
216            DO i=1-Olx,sNx+Olx
217             locMixLayer(i,j) = KPPhbl(i,j,bi,bj)
218            ENDDO
219           ENDDO
220          ELSE
221    #else
222          IF ( .TRUE. ) THEN
223    #endif
224           DO j=1-Oly,sNy+Oly
225            DO i=1-Olx,sNx+Olx
226             locMixLayer(i,j) = hMixLayer(i,j,bi,bj)
227            ENDDO
228           ENDDO
229          ENDIF
230    
231          DO k=Nr,2,-1
232    
233  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
234         kkey = (igmkey-1)*Nr + k         kkey = (igmkey-1)*Nr + k
235         DO j=1-Oly,sNy+Oly         DO j=1-Oly,sNy+Oly
# Line 192  C-- 1rst loop on k : compute Tensor Coef Line 260  C-- 1rst loop on k : compute Tensor Coef
260         ENDDO         ENDDO
261  #endif  #endif
262    
263        DO j=1-Oly+1,sNy+Oly-1         DO j=1-Oly+1,sNy+Oly-1
264         DO i=1-Olx+1,sNx+Olx-1          DO i=1-Olx+1,sNx+Olx-1
265  C      Gradient of Sigma at rVel points  C      Gradient of Sigma at rVel points
266          dSigmaDx(i,j)=op25*( sigmaX(i+1, j ,k-1) +sigmaX(i,j,k-1)           dSigmaDx(i,j)=op25*( sigmaX(i+1,j,k-1)+sigmaX(i,j,k-1)
267       &                    +sigmaX(i+1, j , k ) +sigmaX(i,j, k ) )       &                       +sigmaX(i+1,j, k )+sigmaX(i,j, k )
268       &                  *maskC(i,j,k,bi,bj)       &                      )*maskC(i,j,k,bi,bj)
269          dSigmaDy(i,j)=op25*( sigmaY( i ,j+1,k-1) +sigmaY(i,j,k-1)           dSigmaDy(i,j)=op25*( sigmaY(i,j+1,k-1)+sigmaY(i,j,k-1)
270       &                    +sigmaY( i ,j+1, k ) +sigmaY(i,j, k ) )       &                       +sigmaY(i,j+1, k )+sigmaY(i,j, k )
271       &                  *maskC(i,j,k,bi,bj)       &                      )*maskC(i,j,k,bi,bj)
272          dSigmaDr(i,j)=sigmaR(i,j,k)           dSigmaDr(i,j)=sigmaR(i,j,k)
273            ENDDO
274         ENDDO         ENDDO
       ENDDO  
275    
276  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
277  CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
278  CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
279  CADJ STORE dSigmaDr(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE dSigmaDr(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
280    CADJ STORE baseSlope(:,:)      = comlev1_bibj_k, key=kkey, byte=isbyte
281    CADJ STORE hTransLay(:,:)      = comlev1_bibj_k, key=kkey, byte=isbyte
282    CADJ STORE recipLambda(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
283  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
284    
285    #ifdef GM_VISBECK_VARIABLE_K
286    #ifndef OLD_VISBECK_CALC
287           IF ( GM_Visbeck_alpha.GT.0. .AND.
288         &      -rC(k-1).LT.GM_Visbeck_depth ) THEN
289    
290    C--     Depth average of f/sqrt(Ri) = M^2/N^2 * N
291    C       M^2 and N^2 are horizontal & vertical gradient of buoyancy.
292    
293    C       Calculate terms for mean Richardson number which is used
294    C       in the "variable K" parameterisaton:
295    C       compute depth average from surface down to the bottom or
296    C       GM_Visbeck_depth, whatever is the shallower.
297    
298            DO j=1-Oly+1,sNy+Oly-1
299             DO i=1-Olx+1,sNx+Olx-1
300              IF ( maskC(i,j,k,bi,bj).NE.0. ) THEN
301               integrDepth = -rC( kLowC(i,j,bi,bj) )
302    C-      in 2 steps to avoid mix of RS & RL type in min fct. arguments
303               integrDepth = MIN( integrDepth, GM_Visbeck_depth )
304    C       Distance between level center above and the integration depth
305               deltaH = integrDepth + rC(k-1)
306    C       If negative then we are below the integration level
307    C       (cannot be the case with 2 conditions on maskC & -rC(k-1))
308    C       If positive we limit this to the distance from center above
309               deltaH = MIN( deltaH, drC(k) )
310    C       Now we convert deltaH to a non-dimensional fraction
311               deltaH = deltaH/( integrDepth+rC(1) )
312    
313    C--      compute: ( M^2 * S )^1/2   (= M^2 / N since S=M^2/N^2 )
314               dSigmaH = dSigmaDx(i,j)*dSigmaDx(i,j)
315         &             + dSigmaDy(i,j)*dSigmaDy(i,j)
316               IF ( dSigmaH .GT. 0. _d 0 ) THEN
317                 dSigmaH = SQRT( dSigmaH )
318    C-       compute slope, limited by GM_maxSlope:
319                 IF ( -dSigmaDr(i,j).GT.dSigmaH*GM_rMaxSlope ) THEN
320                  Sloc = dSigmaH / ( -dSigmaDr(i,j) )
321                 ELSE
322                  Sloc = GM_maxSlope
323                 ENDIF
324                 M2loc = Gravity*recip_RhoConst*dSigmaH
325                 SNloc = SQRT( Sloc*M2loc )
326               ELSE
327                 SNloc = 0. _d 0
328               ENDIF
329               VisbeckK(i,j,bi,bj) = VisbeckK(i,j,bi,bj)
330         &       +deltaH*GM_Visbeck_alpha
331         &              *GM_Visbeck_length*GM_Visbeck_length*SNloc
332              ENDIF
333             ENDDO
334            ENDDO
335           ENDIF
336    #endif /* ndef OLD_VISBECK_CALC */
337    #endif /* GM_VISBECK_VARIABLE_K */
338    
339  C     Calculate slopes for use in tensor, taper and/or clip  C     Calculate slopes for use in tensor, taper and/or clip
340        CALL GMREDI_SLOPE_LIMIT(         CALL GMREDI_SLOPE_LIMIT(
341       O             SlopeX, SlopeY,       O             SlopeX, SlopeY,
342       O             SlopeSqr, taperFct,       O             SlopeSqr, taperFct,
343         U             hTransLay, baseSlope, recipLambda,
344       U             dSigmaDr,       U             dSigmaDr,
345       I             dSigmaDx, dSigmaDy,       I             dSigmaDx, dSigmaDy,
346       I             ldd97_LrhoC,rF(k),k,       I             ldd97_LrhoC, locMixLayer, rF,
347       I             bi, bj, myThid )       I             kLowC(1-Olx,1-Oly,bi,bj),
348         I             k, bi, bj, myTime, myIter, myThid )
       DO j=1-Oly+1,sNy+Oly-1  
        DO i=1-Olx+1,sNx+Olx-1  
   
 C       Mask Iso-neutral slopes  
         SlopeX(i,j)=SlopeX(i,j)*maskC(i,j,k,bi,bj)  
         SlopeY(i,j)=SlopeY(i,j)*maskC(i,j,k,bi,bj)  
         SlopeSqr(i,j)=SlopeSqr(i,j)*maskC(i,j,k,bi,bj)  
349    
350           DO j=1-Oly+1,sNy+Oly-1
351            DO i=1-Olx+1,sNx+Olx-1
352    C      Mask Iso-neutral slopes
353             SlopeX(i,j)=SlopeX(i,j)*maskC(i,j,k,bi,bj)
354             SlopeY(i,j)=SlopeY(i,j)*maskC(i,j,k,bi,bj)
355             SlopeSqr(i,j)=SlopeSqr(i,j)*maskC(i,j,k,bi,bj)
356            ENDDO
357         ENDDO         ENDDO
       ENDDO  
358    
359  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
360  CADJ STORE SlopeX(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE SlopeX(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
# Line 239  CADJ STORE dSigmaDr(:,:)     = comlev1_b Line 364  CADJ STORE dSigmaDr(:,:)     = comlev1_b
364  CADJ STORE taperFct(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE taperFct(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte
365  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
366    
367        DO j=1-Oly+1,sNy+Oly-1  C      Components of Redi/GM tensor
368         DO i=1-Olx+1,sNx+Olx-1         DO j=1-Oly+1,sNy+Oly-1
369            DO i=1-Olx+1,sNx+Olx-1
370  C       Components of Redi/GM tensor            Kwx(i,j,k,bi,bj)= SlopeX(i,j)*taperFct(i,j)
371          Kwx(i,j,k,bi,bj)= SlopeX(i,j)*taperFct(i,j)            Kwy(i,j,k,bi,bj)= SlopeY(i,j)*taperFct(i,j)
372          Kwy(i,j,k,bi,bj)= SlopeY(i,j)*taperFct(i,j)            Kwz(i,j,k,bi,bj)= SlopeSqr(i,j)*taperFct(i,j)
373          Kwz(i,j,k,bi,bj)= SlopeSqr(i,j)*taperFct(i,j)          ENDDO
374           ENDDO
375    
376  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
377    #ifdef OLD_VISBECK_CALC
378           DO j=1-Oly+1,sNy+Oly-1
379            DO i=1-Olx+1,sNx+Olx-1
380    
381  C- note (jmc) : moved here since only used in VISBECK_VARIABLE_K  C- note (jmc) : moved here since only used in VISBECK_VARIABLE_K
382  C           but do not know if *taperFct (or **2 ?) is necessary  C           but do not know if *taperFct (or **2 ?) is necessary
383          Ssq(i,j)=SlopeSqr(i,j)*taperFct(i,j)          Ssq(i,j)=SlopeSqr(i,j)*taperFct(i,j)
384    
385  C--     Depth average of M^2/N^2 * N  C--     Depth average of M^2/N^2 * N
# Line 274  C       Now we convert deltaH to a non-d Line 403  C       Now we convert deltaH to a non-d
403       &      *GM_Visbeck_alpha*GM_Visbeck_length*GM_Visbeck_length*SN       &      *GM_Visbeck_alpha*GM_Visbeck_length*GM_Visbeck_length*SN
404          ENDIF          ENDIF
405    
406  #endif /* GM_VISBECK_VARIABLE_K */          ENDDO
   
407         ENDDO         ENDDO
408        ENDDO  #endif /* OLD_VISBECK_CALC */
409    #endif /* GM_VISBECK_VARIABLE_K */
410    
411  C-- end 1rst loop on vertical level index k  C-- end 1rst loop on vertical level index k
412        ENDDO        ENDDO
# Line 287  C-- end 1rst loop on vertical level inde Line 416  C-- end 1rst loop on vertical level inde
416  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
417  CADJ STORE VisbeckK(:,:,bi,bj) = comlev1_bibj, key=igmkey, byte=isbyte  CADJ STORE VisbeckK(:,:,bi,bj) = comlev1_bibj, key=igmkey, byte=isbyte
418  #endif  #endif
419        IF ( GM_Visbeck_alpha.NE.0. ) THEN        IF ( GM_Visbeck_alpha.GT.0. ) THEN
420  C-     Limit range that KapGM can take  C-     Limit range that KapGM can take
421         DO j=1-Oly+1,sNy+Oly-1         DO j=1-Oly+1,sNy+Oly-1
422          DO i=1-Olx+1,sNx+Olx-1          DO i=1-Olx+1,sNx+Olx-1
# Line 303  CADJ STORE VisbeckK(:,:,bi,bj) = comlev1 Line 432  CADJ STORE VisbeckK(:,:,bi,bj) = comlev1
432  cph)  cph)
433  #endif /* GM_VISBECK_VARIABLE_K */  #endif /* GM_VISBECK_VARIABLE_K */
434    
   
 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|  
   
 C-- 2nd loop on k : compute Tensor Coeff. at U,V levels.  
       DO k=1,Nr  
        kp1 = MIN(Nr,k+1)  
        maskp1 = 1. _d 0  
        IF (k.GE.Nr) maskp1 = 0. _d 0  
   
435  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
436         kkey = (igmkey-1)*Nr + k         kkey = (igmkey-1)*Nr + k
437  #if (defined (GM_NON_UNITY_DIAGONAL) || \  #if (defined (GM_NON_UNITY_DIAGONAL) || \
# Line 323  CADJ STORE Kwz(:,:,k,bi,bj) = comlev1_bi Line 443  CADJ STORE Kwz(:,:,k,bi,bj) = comlev1_bi
443  #endif  #endif
444    
445  C-    express the Tensor in term of Diffusivity (= m**2 / s )  C-    express the Tensor in term of Diffusivity (= m**2 / s )
446        DO j=1-Oly+1,sNy+Oly-1        DO k=1,Nr
447         DO i=1-Olx+1,sNx+Olx-1         DO j=1-Oly+1,sNy+Oly-1
448          Kgm_tmp = GM_isopycK + GM_skewflx*GM_background_K          DO i=1-Olx+1,sNx+Olx-1
449    #ifdef ALLOW_KAPGM_CONTROL
450             Kgm_tmp = GM_isopycK + GM_skewflx*kapgm(i,j,k,bi,bj)
451    #else
452             Kgm_tmp = GM_isopycK + GM_skewflx*GM_background_K
453    #endif
454  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
455       &          + VisbeckK(i,j,bi,bj)*(1. _d 0 + GM_skewflx)           &           + VisbeckK(i,j,bi,bj)*(1. _d 0 + GM_skewflx)
456  #endif  #endif
457          Kwx(i,j,k,bi,bj)= Kgm_tmp*Kwx(i,j,k,bi,bj)           Kwx(i,j,k,bi,bj)= Kgm_tmp*Kwx(i,j,k,bi,bj)
458          Kwy(i,j,k,bi,bj)= Kgm_tmp*Kwy(i,j,k,bi,bj)           Kwy(i,j,k,bi,bj)= Kgm_tmp*Kwy(i,j,k,bi,bj)
459          Kwz(i,j,k,bi,bj)= ( GM_isopycK           Kwz(i,j,k,bi,bj)= ( GM_isopycK
460  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
461       &                    + VisbeckK(i,j,bi,bj)       &                     + VisbeckK(i,j,bi,bj)
462  #endif  #endif
463       &                    )*Kwz(i,j,k,bi,bj)       &                     )*Kwz(i,j,k,bi,bj)
464            ENDDO
465         ENDDO         ENDDO
466        ENDDO        ENDDO
467    
468    #ifdef ALLOW_DIAGNOSTICS
469          IF ( useDiagnostics .AND. GM_taper_scheme.EQ.'fm07' ) THEN
470           CALL DIAGNOSTICS_FILL( hTransLay, 'GM_hTrsL', 0,1,2,bi,bj,myThid)
471           CALL DIAGNOSTICS_FILL( baseSlope, 'GM_baseS', 0,1,2,bi,bj,myThid)
472           CALL DIAGNOSTICS_FILL(recipLambda,'GM_rLamb', 0,1,2,bi,bj,myThid)
473          ENDIF
474    #endif /* ALLOW_DIAGNOSTICS */
475    
476    
477  #if ( defined (GM_NON_UNITY_DIAGONAL) || defined (GM_EXTRA_DIAGONAL) )  #if ( defined (GM_NON_UNITY_DIAGONAL) || defined (GM_EXTRA_DIAGONAL) )
478    C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
479    C-- 2nd  k loop : compute Tensor Coeff. at U point
480    
481  C     Gradient of Sigma at U points  #ifdef ALLOW_KPP
482        DO j=1-Oly+1,sNy+Oly-1        IF ( useKPP ) THEN
483         DO i=1-Olx+1,sNx+Olx-1         DO j=1-Oly,sNy+Oly
484          dSigmaDx(i,j)=sigmaX(i,j,k)          DO i=2-Olx,sNx+Olx
485       &          *_maskW(i,j,k,bi,bj)           locMixLayer(i,j) = ( KPPhbl(i-1,j,bi,bj)
486          dSigmaDy(i,j)=op25*( sigmaY(i-1,j+1,k) +sigmaY(i,j+1,k)       &                      + KPPhbl( i ,j,bi,bj) )*op5
487       &                      +sigmaY(i-1, j ,k) +sigmaY(i, j ,k) )          ENDDO
488       &          *_maskW(i,j,k,bi,bj)         ENDDO
489          dSigmaDr(i,j)=op25*( sigmaR(i-1,j, k ) +sigmaR(i,j, k )        ELSE
490       &                  +maskp1*(sigmaR(i-1,j,kp1) +sigmaR(i,j,kp1)) )  #else
491       &          *_maskW(i,j,k,bi,bj)        IF ( .TRUE. ) THEN
492    #endif
493           DO j=1-Oly,sNy+Oly
494            DO i=2-Olx,sNx+Olx
495             locMixLayer(i,j) = ( hMixLayer(i-1,j,bi,bj)
496         &                      + hMixLayer( i ,j,bi,bj) )*op5
497            ENDDO
498           ENDDO
499          ENDIF
500          DO j=1-Oly,sNy+Oly
501           DO i=1-Olx,sNx+Olx
502             hTransLay(i,j) =  0.
503             baseSlope(i,j) =  0.
504             recipLambda(i,j)= 0.
505           ENDDO
506           DO i=2-Olx,sNx+Olx
507             hTransLay(i,j) = MAX( R_low(i-1,j,bi,bj), R_low(i,j,bi,bj) )
508         ENDDO         ENDDO
509        ENDDO        ENDDO
510    
511          DO k=Nr,1,-1
512           kp1 = MIN(Nr,k+1)
513           maskp1 = 1. _d 0
514           IF (k.GE.Nr) maskp1 = 0. _d 0
515    
516    C     Gradient of Sigma at U points
517           DO j=1-Oly+1,sNy+Oly-1
518            DO i=1-Olx+1,sNx+Olx-1
519             dSigmaDx(i,j)=sigmaX(i,j,k)
520         &                       *_maskW(i,j,k,bi,bj)
521             dSigmaDy(i,j)=op25*( sigmaY(i-1,j+1,k)+sigmaY(i,j+1,k)
522         &                       +sigmaY(i-1, j ,k)+sigmaY(i, j ,k)
523         &                      )*_maskW(i,j,k,bi,bj)
524             dSigmaDr(i,j)=op25*( sigmaR(i-1,j, k )+sigmaR(i,j, k )
525         &                      +(sigmaR(i-1,j,kp1)+sigmaR(i,j,kp1))*maskp1
526         &                      )*_maskW(i,j,k,bi,bj)
527            ENDDO
528           ENDDO
529    
530  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
531  CADJ STORE SlopeSqr(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE SlopeSqr(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
532  CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
533  CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
534  CADJ STORE dSigmaDr(:,:)   = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE dSigmaDr(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
535    CADJ STORE locMixLayer(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
536    CADJ STORE baseSlope(:,:)      = comlev1_bibj_k, key=kkey, byte=isbyte
537    CADJ STORE hTransLay(:,:)      = comlev1_bibj_k, key=kkey, byte=isbyte
538    CADJ STORE recipLambda(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
539  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
540    
541  C     Calculate slopes for use in tensor, taper and/or clip  C     Calculate slopes for use in tensor, taper and/or clip
542        CALL GMREDI_SLOPE_LIMIT(         CALL GMREDI_SLOPE_LIMIT(
543       O             SlopeX, SlopeY,       O             SlopeX, SlopeY,
544       O             SlopeSqr, taperFct,       O             SlopeSqr, taperFct,
545         U             hTransLay, baseSlope, recipLambda,
546       U             dSigmaDr,       U             dSigmaDr,
547       I             dSigmaDx, dSigmaDy,       I             dSigmaDx, dSigmaDy,
548       I             ldd97_LrhoW,rC(k),k,       I             ldd97_LrhoW, locMixLayer, rC,
549       I             bi, bj, myThid )       I             kLow_W,
550         I             k, bi, bj, myTime, myIter, myThid )
551    
552  cph( NEW  cph( NEW
553  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 381  cph) Line 559  cph)
559  cph)  cph)
560    
561  #ifdef GM_NON_UNITY_DIAGONAL  #ifdef GM_NON_UNITY_DIAGONAL
562    c      IF ( GM_nonUnitDiag ) THEN
563          DO j=1-Oly+1,sNy+Oly-1          DO j=1-Oly+1,sNy+Oly-1
564           DO i=1-Olx+1,sNx+Olx-1           DO i=1-Olx+1,sNx+Olx-1
565            Kux(i,j,k,bi,bj) =            Kux(i,j,k,bi,bj) =
# Line 388  cph) Line 567  cph)
567  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
568       &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))       &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))
569  #endif  #endif
570       &     )       &     )*taperFct(i,j)
      &     *taperFct(i,j)  
571           ENDDO           ENDDO
572          ENDDO          ENDDO
573  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 402  CADJ STORE Kux(:,:,k,bi,bj)  = comlev1_b Line 580  CADJ STORE Kux(:,:,k,bi,bj)  = comlev1_b
580            Kux(i,j,k,bi,bj) = MAX( Kux(i,j,k,bi,bj), GM_Kmin_horiz )            Kux(i,j,k,bi,bj) = MAX( Kux(i,j,k,bi,bj), GM_Kmin_horiz )
581           ENDDO           ENDDO
582          ENDDO          ENDDO
583    c      ENDIF
584  #endif /* GM_NON_UNITY_DIAGONAL */  #endif /* GM_NON_UNITY_DIAGONAL */
585    
586  #ifdef GM_EXTRA_DIAGONAL  #ifdef GM_EXTRA_DIAGONAL
# Line 410  CADJ STORE Kux(:,:,k,bi,bj)  = comlev1_b Line 589  CADJ STORE Kux(:,:,k,bi,bj)  = comlev1_b
589  CADJ STORE SlopeX(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE SlopeX(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
590  CADJ STORE taperFct(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE taperFct(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte
591  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
592        IF (GM_ExtraDiag) THEN         IF ( GM_ExtraDiag ) THEN
593          DO j=1-Oly+1,sNy+Oly-1          DO j=1-Oly+1,sNy+Oly-1
594           DO i=1-Olx+1,sNx+Olx-1           DO i=1-Olx+1,sNx+Olx-1
595            Kuz(i,j,k,bi,bj) =            Kuz(i,j,k,bi,bj) =
596    #ifdef ALLOW_KAPGM_CONTROL
597         &     ( GM_isopycK - GM_skewflx*kapgm(i,j,k,bi,bj)
598    #else
599       &     ( GM_isopycK - GM_skewflx*GM_background_K       &     ( GM_isopycK - GM_skewflx*GM_background_K
600    #endif
601  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
602       &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))*GM_advect       &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))*GM_advect
603  #endif  #endif
604       &     )*SlopeX(i,j)*taperFct(i,j)       &     )*SlopeX(i,j)*taperFct(i,j)
605           ENDDO           ENDDO
606          ENDDO          ENDDO
607        ENDIF         ENDIF
608  #endif /* GM_EXTRA_DIAGONAL */  #endif /* GM_EXTRA_DIAGONAL */
609    
610  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
611        IF (doDiagRediFlx) THEN         IF (doDiagRediFlx) THEN
612          km1 = MAX(k-1,1)          km1 = MAX(k-1,1)
613          DO j=1,sNy          DO j=1,sNy
614           DO i=1,sNx+1           DO i=1,sNx+1
# Line 459  C-        Vertical gradients interpolate Line 642  C-        Vertical gradients interpolate
642           ENDDO           ENDDO
643          ENDDO          ENDDO
644          CALL DIAGNOSTICS_FILL(tmp1k, 'GM_KuzTz', k,1,2,bi,bj,myThid)          CALL DIAGNOSTICS_FILL(tmp1k, 'GM_KuzTz', k,1,2,bi,bj,myThid)
645        ENDIF         ENDIF
646  #endif /* ALLOW_DIAGNOSTICS */  #endif /* ALLOW_DIAGNOSTICS */
647    
648  C     Gradient of Sigma at V points  C-- end 2nd  loop on vertical level index k
649        DO j=1-Oly+1,sNy+Oly-1        ENDDO
650         DO i=1-Olx+1,sNx+Olx-1  
651          dSigmaDx(i,j)=op25*( sigmaX(i, j ,k) +sigmaX(i+1, j ,k)  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
652       &                    +sigmaX(i,j-1,k) +sigmaX(i+1,j-1,k) )  C-- 3rd  k loop : compute Tensor Coeff. at V point
653       &          *_maskS(i,j,k,bi,bj)  
654          dSigmaDy(i,j)=sigmaY(i,j,k)  #ifdef ALLOW_KPP
655       &          *_maskS(i,j,k,bi,bj)        IF ( useKPP ) THEN
656          dSigmaDr(i,j)=op25*( sigmaR(i,j-1, k ) +sigmaR(i,j, k )         DO j=2-Oly,sNy+Oly
657       &                  +maskp1*(sigmaR(i,j-1,kp1) +sigmaR(i,j,kp1)) )          DO i=1-Olx,sNx+Olx
658       &          *_maskS(i,j,k,bi,bj)           locMixLayer(i,j) = ( KPPhbl(i,j-1,bi,bj)
659         &                      + KPPhbl(i, j ,bi,bj) )*op5
660            ENDDO
661           ENDDO
662          ELSE
663    #else
664          IF ( .TRUE. ) THEN
665    #endif
666           DO j=2-Oly,sNy+Oly
667            DO i=1-Olx,sNx+Olx
668             locMixLayer(i,j) = ( hMixLayer(i,j-1,bi,bj)
669         &                      + hMixLayer(i, j ,bi,bj) )*op5
670            ENDDO
671           ENDDO
672          ENDIF
673          DO j=1-Oly,sNy+Oly
674           DO i=1-Olx,sNx+Olx
675             hTransLay(i,j) =  0.
676             baseSlope(i,j) =  0.
677             recipLambda(i,j)= 0.
678         ENDDO         ENDDO
679        ENDDO        ENDDO
680          DO j=2-Oly,sNy+Oly
681           DO i=1-Olx,sNx+Olx
682             hTransLay(i,j) = MAX( R_low(i,j-1,bi,bj), R_low(i,j,bi,bj) )
683           ENDDO
684          ENDDO
685    
686    C     Gradient of Sigma at V points
687          DO k=Nr,1,-1
688           kp1 = MIN(Nr,k+1)
689           maskp1 = 1. _d 0
690           IF (k.GE.Nr) maskp1 = 0. _d 0
691    
692           DO j=1-Oly+1,sNy+Oly-1
693            DO i=1-Olx+1,sNx+Olx-1
694             dSigmaDx(i,j)=op25*( sigmaX(i, j ,k) +sigmaX(i+1, j ,k)
695         &                       +sigmaX(i,j-1,k) +sigmaX(i+1,j-1,k)
696         &                      )*_maskS(i,j,k,bi,bj)
697             dSigmaDy(i,j)=sigmaY(i,j,k)
698         &                       *_maskS(i,j,k,bi,bj)
699             dSigmaDr(i,j)=op25*( sigmaR(i,j-1, k )+sigmaR(i,j, k )
700         &                      +(sigmaR(i,j-1,kp1)+sigmaR(i,j,kp1))*maskp1
701         &                      )*_maskS(i,j,k,bi,bj)
702            ENDDO
703           ENDDO
704    
705  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
706  CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
707  CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
708  CADJ STORE dSigmaDr(:,:)   = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE dSigmaDr(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
709    CADJ STORE baseSlope(:,:)      = comlev1_bibj_k, key=kkey, byte=isbyte
710    CADJ STORE hTransLay(:,:)      = comlev1_bibj_k, key=kkey, byte=isbyte
711    CADJ STORE recipLambda(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
712  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
713    
714  C     Calculate slopes for use in tensor, taper and/or clip  C     Calculate slopes for use in tensor, taper and/or clip
715        CALL GMREDI_SLOPE_LIMIT(         CALL GMREDI_SLOPE_LIMIT(
716       O             SlopeX, SlopeY,       O             SlopeX, SlopeY,
717       O             SlopeSqr, taperFct,       O             SlopeSqr, taperFct,
718         U             hTransLay, baseSlope, recipLambda,
719       U             dSigmaDr,       U             dSigmaDr,
720       I             dSigmaDx, dSigmaDy,       I             dSigmaDx, dSigmaDy,
721       I             ldd97_LrhoS,rC(k),k,       I             ldd97_LrhoS, locMixLayer, rC,
722       I             bi, bj, myThid )       I             kLow_S,
723         I             k, bi, bj, myTime, myIter, myThid )
724    
725  cph(  cph(
726  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 500  cph) Line 731  cph)
731  cph)  cph)
732    
733  #ifdef GM_NON_UNITY_DIAGONAL  #ifdef GM_NON_UNITY_DIAGONAL
734    c      IF ( GM_nonUnitDiag ) THEN
735          DO j=1-Oly+1,sNy+Oly-1          DO j=1-Oly+1,sNy+Oly-1
736           DO i=1-Olx+1,sNx+Olx-1           DO i=1-Olx+1,sNx+Olx-1
737            Kvy(i,j,k,bi,bj) =            Kvy(i,j,k,bi,bj) =
# Line 507  cph) Line 739  cph)
739  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
740       &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))       &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))
741  #endif  #endif
742       &     )       &     )*taperFct(i,j)
      &     *taperFct(i,j)  
743           ENDDO           ENDDO
744          ENDDO          ENDDO
745  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 521  CADJ STORE Kvy(:,:,k,bi,bj)  = comlev1_b Line 752  CADJ STORE Kvy(:,:,k,bi,bj)  = comlev1_b
752            Kvy(i,j,k,bi,bj) = MAX( Kvy(i,j,k,bi,bj), GM_Kmin_horiz )            Kvy(i,j,k,bi,bj) = MAX( Kvy(i,j,k,bi,bj), GM_Kmin_horiz )
753           ENDDO           ENDDO
754          ENDDO          ENDDO
755    c      ENDIF
756  #endif /* GM_NON_UNITY_DIAGONAL */  #endif /* GM_NON_UNITY_DIAGONAL */
757    
758  #ifdef GM_EXTRA_DIAGONAL  #ifdef GM_EXTRA_DIAGONAL
# Line 529  CADJ STORE Kvy(:,:,k,bi,bj)  = comlev1_b Line 761  CADJ STORE Kvy(:,:,k,bi,bj)  = comlev1_b
761  CADJ STORE SlopeY(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE SlopeY(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
762  CADJ STORE taperFct(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE taperFct(:,:)     = comlev1_bibj_k, key=kkey, byte=isbyte
763  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
764        IF (GM_ExtraDiag) THEN         IF ( GM_ExtraDiag ) THEN
765          DO j=1-Oly+1,sNy+Oly-1          DO j=1-Oly+1,sNy+Oly-1
766           DO i=1-Olx+1,sNx+Olx-1           DO i=1-Olx+1,sNx+Olx-1
767            Kvz(i,j,k,bi,bj) =            Kvz(i,j,k,bi,bj) =
768    #ifdef ALLOW_KAPGM_CONTROL
769         &     ( GM_isopycK - GM_skewflx*kapgm(i,j,k,bi,bj)
770    #else
771       &     ( GM_isopycK - GM_skewflx*GM_background_K       &     ( GM_isopycK - GM_skewflx*GM_background_K
772    #endif
773  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
774       &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))*GM_advect       &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))*GM_advect
775  #endif  #endif
776       &     )*SlopeY(i,j)*taperFct(i,j)       &     )*SlopeY(i,j)*taperFct(i,j)
777           ENDDO           ENDDO
778          ENDDO          ENDDO
779        ENDIF         ENDIF
780  #endif /* GM_EXTRA_DIAGONAL */  #endif /* GM_EXTRA_DIAGONAL */
781    
782  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
783        IF (doDiagRediFlx) THEN         IF (doDiagRediFlx) THEN
784  c       km1 = MAX(k-1,1)  c       km1 = MAX(k-1,1)
785          DO j=1,sNy+1          DO j=1,sNy+1
786           DO i=1,sNx           DO i=1,sNx
# Line 578  C-        Vertical gradients interpolate Line 814  C-        Vertical gradients interpolate
814           ENDDO           ENDDO
815          ENDDO          ENDDO
816          CALL DIAGNOSTICS_FILL(tmp1k, 'GM_KvzTz', k,1,2,bi,bj,myThid)          CALL DIAGNOSTICS_FILL(tmp1k, 'GM_KvzTz', k,1,2,bi,bj,myThid)
817        ENDIF         ENDIF
818  #endif /* ALLOW_DIAGNOSTICS */  #endif /* ALLOW_DIAGNOSTICS */
819    
820  #endif /* GM_NON_UNITY_DIAGONAL || GM_EXTRA_DIAGONAL */  C-- end 3rd  loop on vertical level index k
   
 C-- end 2nd loop on vertical level index k  
821        ENDDO        ENDDO
822    
823    #endif /* GM_NON_UNITY_DIAGONAL || GM_EXTRA_DIAGONAL */
824    
825    
826  #ifdef GM_BOLUS_ADVEC  #ifdef GM_BOLUS_ADVEC
827        IF (GM_AdvForm) THEN        IF (GM_AdvForm) THEN
828          CALL GMREDI_CALC_PSI_B(         CALL GMREDI_CALC_PSI_B(
829       I             bi, bj, iMin, iMax, jMin, jMax,       I             bi, bj, iMin, iMax, jMin, jMax,
830       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
831       I             ldd97_LrhoW, ldd97_LrhoS,       I             ldd97_LrhoW, ldd97_LrhoS,
832       I             myThid )       I             myThid )
833        ENDIF        ENDIF
834  #endif  #endif
835    
# Line 630  C--   Time-average Line 866  C--   Time-average
866    
867  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
868        IF ( useDiagnostics ) THEN        IF ( useDiagnostics ) THEN
869         CALL GMREDI_DIAGNOSTICS_DRIVER(bi,bj,myThid)          CALL GMREDI_DIAGNOSTICS_FILL(bi,bj,myThid)
870        ENDIF        ENDIF
871  #endif /* ALLOW_DIAGNOSTICS */  #endif /* ALLOW_DIAGNOSTICS */
872    
# Line 641  C--   Time-average Line 877  C--   Time-average
877    
878    
879        SUBROUTINE GMREDI_CALC_TENSOR_DUMMY(        SUBROUTINE GMREDI_CALC_TENSOR_DUMMY(
880       I             bi, bj, iMin, iMax, jMin, jMax,       I             iMin, iMax, jMin, jMax,
881       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
882       I             myThid )       I             bi, bj, myTime, myIter, myThid )
883  C     /==========================================================\  C     /==========================================================\
884  C     | SUBROUTINE GMREDI_CALC_TENSOR                            |  C     | SUBROUTINE GMREDI_CALC_TENSOR                            |
885  C     | o Calculate tensor elements for GM/Redi tensor.          |  C     | o Calculate tensor elements for GM/Redi tensor.          |
# Line 661  C Line 897  C
897        _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
898        _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
899        _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
900        INTEGER bi,bj,iMin,iMax,jMin,jMax        INTEGER iMin,iMax,jMin,jMax
901          INTEGER bi, bj
902          _RL     myTime
903          INTEGER myIter
904        INTEGER myThid        INTEGER myThid
905  CEndOfInterface  CEndOfInterface
906    
       INTEGER i, j, k  
   
907  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
908    
909          INTEGER i, j, k
910    
911        DO k=1,Nr        DO k=1,Nr
912         DO j=1-Oly+1,sNy+Oly-1         DO j=1-Oly+1,sNy+Oly-1
913          DO i=1-Olx+1,sNx+Olx-1          DO i=1-Olx+1,sNx+Olx-1
Legend:
Removed from v.1.23  
changed lines
  Added in v.1.28
  ViewVC Help
Powered by ViewVC 1.1.22