/[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.31 by gforget, Sat Feb 2 02:35:53 2008 UTC revision 1.39 by jmc, Thu Feb 10 21:24:19 2011 UTC
# Line 5  C $Name$ Line 5  C $Name$
5  #ifdef ALLOW_KPP  #ifdef ALLOW_KPP
6  # include "KPP_OPTIONS.h"  # include "KPP_OPTIONS.h"
7  #endif  #endif
 #undef OLD_VISBECK_CALC  
8    
9  CBOP  CBOP
10  C     !ROUTINE: GMREDI_CALC_TENSOR  C     !ROUTINE: GMREDI_CALC_TENSOR
# Line 64  CEOP Line 63  CEOP
63    
64  C     !LOCAL VARIABLES:  C     !LOCAL VARIABLES:
65  C     == Local variables ==  C     == Local variables ==
66        INTEGER i,j,k,kp1        INTEGER i,j,k
67        _RL SlopeX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL SlopeX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
68        _RL SlopeY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL SlopeY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
69        _RL dSigmaDx(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL dSigmaDx(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
# Line 72  C     == Local variables == Line 71  C     == Local variables ==
71        _RL dSigmaDr(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL dSigmaDr(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
72        _RL SlopeSqr(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL SlopeSqr(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
73        _RL taperFct(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL taperFct(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
       _RL maskp1, Kgm_tmp  
74        _RL ldd97_LrhoC(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL ldd97_LrhoC(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
75        _RL ldd97_LrhoW(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL ldd97_LrhoW(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
76        _RL ldd97_LrhoS(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL ldd97_LrhoS(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
77        _RL Cspd, LrhoInf, LrhoSup, fCoriLoc        _RL Cspd, LrhoInf, LrhoSup, fCoriLoc
78          _RL Kgm_tmp, isopycK, bolus_K
79    
80        INTEGER kLow_W (1-Olx:sNx+Olx,1-Oly:sNy+Oly)        INTEGER kLow_W (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
81        INTEGER kLow_S (1-Olx:sNx+Olx,1-Oly:sNy+Oly)        INTEGER kLow_S (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
# Line 84  C     == Local variables == Line 83  C     == Local variables ==
83        _RL baseSlope  (1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL baseSlope  (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
84        _RL hTransLay  (1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL hTransLay  (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
85        _RL recipLambda(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL recipLambda(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
86          INTEGER  km1
87    #if ( defined (GM_NON_UNITY_DIAGONAL) || defined (GM_EXTRA_DIAGONAL) )
88          INTEGER kp1
89          _RL maskp1
90    #endif
91    
92  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
93  #ifdef OLD_VISBECK_CALC  #ifdef OLD_VISBECK_CALC
       _RL deltaH,zero_rs  
       PARAMETER(zero_rs=0.D0)  
       _RL N2,SN  
94        _RL Ssq(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL Ssq(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
95  #else  #else
96        _RL dSigmaH        _RL dSigmaH, dSigmaR
97        _RL deltaH, integrDepth        _RL Sloc, M2loc
       _RL Sloc, M2loc, SNloc  
 #endif  
98  #endif  #endif
99          _RL recipMaxSlope
100          _RL deltaH, integrDepth
101          _RL N2loc, SNloc
102    #endif /* GM_VISBECK_VARIABLE_K */
103    
104  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
105        LOGICAL  doDiagRediFlx        LOGICAL  doDiagRediFlx
106        LOGICAL  DIAGNOSTICS_IS_ON        LOGICAL  DIAGNOSTICS_IS_ON
107        EXTERNAL DIAGNOSTICS_IS_ON        EXTERNAL DIAGNOSTICS_IS_ON
108        INTEGER  km1  #if ( defined (GM_NON_UNITY_DIAGONAL) || defined (GM_EXTRA_DIAGONAL) )
109        _RL dTdz        _RL dTdz
110        _RL tmp1k(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL tmp1k(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
111  #endif  #endif
112    #endif /* ALLOW_DIAGNOSTICS */
113    
114  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
115    
# Line 132  C---+----1----+----2----+----3----+----4 Line 136  C---+----1----+----2----+----3----+----4
136  #endif  #endif
137    
138  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
139          recipMaxSlope = 0. _d 0
140          IF ( GM_Visbeck_maxSlope.GT.0. _d 0 ) THEN
141            recipMaxSlope = 1. _d 0 / GM_Visbeck_maxSlope
142          ENDIF
143        DO j=1-Oly,sNy+Oly        DO j=1-Oly,sNy+Oly
144         DO i=1-Olx,sNx+Olx         DO i=1-Olx,sNx+Olx
145          VisbeckK(i,j,bi,bj) = 0. _d 0          VisbeckK(i,j,bi,bj) = 0. _d 0
# Line 258  C-- 1rst loop on k : compute Tensor Coef Line 266  C-- 1rst loop on k : compute Tensor Coef
266  # endif  # endif
267          ENDDO          ENDDO
268         ENDDO         ENDDO
269  #endif  #endif /* ALLOW_AUTODIFF_TAMC */
270    
271         DO j=1-Oly+1,sNy+Oly-1         DO j=1-Oly+1,sNy+Oly-1
272          DO i=1-Olx+1,sNx+Olx-1          DO i=1-Olx+1,sNx+Olx-1
# Line 269  C      Gradient of Sigma at rVel points Line 277  C      Gradient of Sigma at rVel points
277           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)
278       &                       +sigmaY(i,j+1, k )+sigmaY(i,j, k )       &                       +sigmaY(i,j+1, k )+sigmaY(i,j, k )
279       &                      )*maskC(i,j,k,bi,bj)       &                      )*maskC(i,j,k,bi,bj)
280           dSigmaDr(i,j)=sigmaR(i,j,k)  c        dSigmaDr(i,j)=sigmaR(i,j,k)
281          ENDDO          ENDDO
282         ENDDO         ENDDO
283    
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte  
 CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte  
 CADJ STORE dSigmaDr(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte  
 CADJ STORE baseSlope(:,:)      = comlev1_bibj_k, key=kkey, byte=isbyte  
 CADJ STORE hTransLay(:,:)      = comlev1_bibj_k, key=kkey, byte=isbyte  
 CADJ STORE recipLambda(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte  
 #endif /* ALLOW_AUTODIFF_TAMC */  
   
284  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
285  #ifndef OLD_VISBECK_CALC  #ifndef OLD_VISBECK_CALC
286         IF ( GM_Visbeck_alpha.GT.0. .AND.         IF ( GM_Visbeck_alpha.GT.0. .AND.
287       &      -rC(k-1).LT.GM_Visbeck_depth ) THEN       &      -rC(k-1).LT.GM_Visbeck_depth ) THEN
288    
289            DO j=1-Oly,sNy+Oly
290             DO i=1-Olx,sNx+Olx
291              dSigmaDr(i,j) = MIN( sigmaR(i,j,k), 0. _d 0 )
292             ENDDO
293            ENDDO
294    
295  C--     Depth average of f/sqrt(Ri) = M^2/N^2 * N  C--     Depth average of f/sqrt(Ri) = M^2/N^2 * N
296  C       M^2 and N^2 are horizontal & vertical gradient of buoyancy.  C       M^2 and N^2 are horizontal & vertical gradient of buoyancy.
297    
# Line 301  C       GM_Visbeck_depth, whatever is th Line 306  C       GM_Visbeck_depth, whatever is th
306             integrDepth = -rC( kLowC(i,j,bi,bj) )             integrDepth = -rC( kLowC(i,j,bi,bj) )
307  C-      in 2 steps to avoid mix of RS & RL type in min fct. arguments  C-      in 2 steps to avoid mix of RS & RL type in min fct. arguments
308             integrDepth = MIN( integrDepth, GM_Visbeck_depth )             integrDepth = MIN( integrDepth, GM_Visbeck_depth )
309    C-      to recover "old-visbeck" form with Visbeck_minDepth = Visbeck_depth
310               integrDepth = MAX( integrDepth, GM_Visbeck_minDepth )
311  C       Distance between level center above and the integration depth  C       Distance between level center above and the integration depth
312             deltaH = integrDepth + rC(k-1)             deltaH = integrDepth + rC(k-1)
313  C       If negative then we are below the integration level  C       If negative then we are below the integration level
# Line 310  C       If positive we limit this to the Line 317  C       If positive we limit this to the
317  C       Now we convert deltaH to a non-dimensional fraction  C       Now we convert deltaH to a non-dimensional fraction
318             deltaH = deltaH/( integrDepth+rC(1) )             deltaH = deltaH/( integrDepth+rC(1) )
319    
320  C--      compute: ( M^2 * S )^1/2   (= M^2 / N since S=M^2/N^2 )  C--      compute: ( M^2 * S )^1/2   (= S*N since S=M^2/N^2 )
321    C        a 5 points average gives a more "homogeneous" formulation
322    C        (same stencil and same weights as for dSigmaH calculation)
323               dSigmaR = ( dSigmaDr(i,j)*4. _d 0
324         &               + dSigmaDr(i-1,j)
325         &               + dSigmaDr(i+1,j)
326         &               + dSigmaDr(i,j-1)
327         &               + dSigmaDr(i,j+1)
328         &               )/( 4. _d 0
329         &                 + maskC(i-1,j,k,bi,bj)
330         &                 + maskC(i+1,j,k,bi,bj)
331         &                 + maskC(i,j-1,k,bi,bj)
332         &                 + maskC(i,j+1,k,bi,bj)
333         &                 )
334             dSigmaH = dSigmaDx(i,j)*dSigmaDx(i,j)             dSigmaH = dSigmaDx(i,j)*dSigmaDx(i,j)
335       &             + dSigmaDy(i,j)*dSigmaDy(i,j)       &             + dSigmaDy(i,j)*dSigmaDy(i,j)
336             IF ( dSigmaH .GT. 0. _d 0 ) THEN             IF ( dSigmaH .GT. 0. _d 0 ) THEN
337               dSigmaH = SQRT( dSigmaH )               dSigmaH = SQRT( dSigmaH )
338  C-       compute slope, limited by GM_maxSlope:  C-       compute slope, limited by GM_Visbeck_maxSlope:
339               IF ( -dSigmaDr(i,j).GT.dSigmaH*GM_rMaxSlope ) THEN               IF ( -dSigmaR.GT.dSigmaH*recipMaxSlope ) THEN
340                Sloc = dSigmaH / ( -dSigmaDr(i,j) )                Sloc = dSigmaH / ( -dSigmaR )
341                 ELSE
342                  Sloc = GM_Visbeck_maxSlope
343                 ENDIF
344                 M2loc = gravity*recip_rhoConst*dSigmaH
345    c            SNloc = SQRT( Sloc*M2loc )
346                 N2loc = -gravity*recip_rhoConst*dSigmaR
347    c            N2loc = -gravity*recip_rhoConst*dSigmaDr(i,j)
348                 IF ( N2loc.GT.0. _d 0 ) THEN
349                   SNloc = Sloc*SQRT(N2loc)
350               ELSE               ELSE
351                Sloc = GM_maxSlope                 SNloc = 0. _d 0
352               ENDIF               ENDIF
              M2loc = Gravity*recip_RhoConst*dSigmaH  
              SNloc = SQRT( Sloc*M2loc )  
353             ELSE             ELSE
354               SNloc = 0. _d 0               SNloc = 0. _d 0
355             ENDIF             ENDIF
# Line 335  C-       compute slope, limited by GM_ma Line 362  C-       compute slope, limited by GM_ma
362         ENDIF         ENDIF
363  #endif /* ndef OLD_VISBECK_CALC */  #endif /* ndef OLD_VISBECK_CALC */
364  #endif /* GM_VISBECK_VARIABLE_K */  #endif /* GM_VISBECK_VARIABLE_K */
365           DO j=1-Oly,sNy+Oly
366            DO i=1-Olx,sNx+Olx
367             dSigmaDr(i,j)=sigmaR(i,j,k)
368            ENDDO
369           ENDDO
370    
371    #ifdef ALLOW_AUTODIFF_TAMC
372    CADJ STORE dSigmaDx(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
373    CADJ STORE dSigmaDy(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
374    CADJ STORE dSigmaDr(:,:)       = comlev1_bibj_k, key=kkey, byte=isbyte
375    CADJ STORE baseSlope(:,:)      = comlev1_bibj_k, key=kkey, byte=isbyte
376    CADJ STORE hTransLay(:,:)      = comlev1_bibj_k, key=kkey, byte=isbyte
377    CADJ STORE recipLambda(:,:)    = comlev1_bibj_k, key=kkey, byte=isbyte
378    #endif /* ALLOW_AUTODIFF_TAMC */
379    
380  C     Calculate slopes for use in tensor, taper and/or clip  C     Calculate slopes for use in tensor, taper and/or clip
381         CALL GMREDI_SLOPE_LIMIT(         CALL GMREDI_SLOPE_LIMIT(
# Line 389  C       which is used in the "variable K Line 430  C       which is used in the "variable K
430  C       Distance between interface above layer and the integration depth  C       Distance between interface above layer and the integration depth
431          deltaH=abs(GM_Visbeck_depth)-abs(rF(k))          deltaH=abs(GM_Visbeck_depth)-abs(rF(k))
432  C       If positive we limit this to the layer thickness  C       If positive we limit this to the layer thickness
433          deltaH=min(deltaH,drF(k))          integrDepth = drF(k)
434            deltaH=min(deltaH,integrDepth)
435  C       If negative then we are below the integration level  C       If negative then we are below the integration level
436          deltaH=max(deltaH,zero_rs)          deltaH=max(deltaH, 0. _d 0)
437  C       Now we convert deltaH to a non-dimensional fraction  C       Now we convert deltaH to a non-dimensional fraction
438          deltaH=deltaH/GM_Visbeck_depth          deltaH=deltaH/GM_Visbeck_depth
439    
         IF (K.eq.2) VisbeckK(i,j,bi,bj)=0.  
440          IF ( Ssq(i,j).NE.0. .AND. dSigmaDr(i,j).NE.0. ) THEN          IF ( Ssq(i,j).NE.0. .AND. dSigmaDr(i,j).NE.0. ) THEN
441           N2= -Gravity*recip_RhoConst*dSigmaDr(i,j)           N2loc = -gravity*recip_rhoConst*dSigmaDr(i,j)
442           SN=sqrt(Ssq(i,j)*N2)           SNloc = SQRT(Ssq(i,j)*N2loc )
443           VisbeckK(i,j,bi,bj)=VisbeckK(i,j,bi,bj)+deltaH           VisbeckK(i,j,bi,bj) = VisbeckK(i,j,bi,bj)
444       &      *GM_Visbeck_alpha*GM_Visbeck_length*GM_Visbeck_length*SN       &       +deltaH*GM_Visbeck_alpha
445         &              *GM_Visbeck_length*GM_Visbeck_length*SNloc
446          ENDIF          ENDIF
447    
448          ENDDO          ENDDO
# Line 421  C-     Limit range that KapGM can take Line 463  C-     Limit range that KapGM can take
463         DO j=1-Oly+1,sNy+Oly-1         DO j=1-Oly+1,sNy+Oly-1
464          DO i=1-Olx+1,sNx+Olx-1          DO i=1-Olx+1,sNx+Olx-1
465           VisbeckK(i,j,bi,bj)=           VisbeckK(i,j,bi,bj)=
466       &       MIN(VisbeckK(i,j,bi,bj),GM_Visbeck_maxval_K)       &       MIN( MAX( VisbeckK(i,j,bi,bj), GM_Visbeck_minVal_K ),
467         &            GM_Visbeck_maxVal_K )
468          ENDDO          ENDDO
469         ENDDO         ENDDO
470        ENDIF        ENDIF
# Line 437  C-    express the Tensor in term of Diff Line 480  C-    express the Tensor in term of Diff
480  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
481         kkey = (igmkey-1)*Nr + k         kkey = (igmkey-1)*Nr + k
482  # if (defined (GM_NON_UNITY_DIAGONAL) || \  # if (defined (GM_NON_UNITY_DIAGONAL) || \
483       defined (GM_VISBECK_VARIABLE_K))        defined (GM_VISBECK_VARIABLE_K))
484  CADJ STORE Kwx(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE Kwx(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
485  CADJ STORE Kwy(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE Kwy(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
486  CADJ STORE Kwz(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte  CADJ STORE Kwz(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
487  # endif  # endif
488  #endif  #endif
489           km1 = MAX(k-1,1)
490           isopycK = GM_isopycK
491         &         *(GM_isoFac1d(km1)+GM_isoFac1d(k))*op5
492           bolus_K = GM_background_K
493         &         *(GM_bolFac1d(km1)+GM_bolFac1d(k))*op5
494         DO j=1-Oly+1,sNy+Oly-1         DO j=1-Oly+1,sNy+Oly-1
495          DO i=1-Olx+1,sNx+Olx-1          DO i=1-Olx+1,sNx+Olx-1
496  #ifdef ALLOW_KAPREDI_CONTROL  #ifdef ALLOW_KAPREDI_CONTROL
497           Kgm_tmp = kapredi(i,j,k,bi,bj)           Kgm_tmp = kapredi(i,j,k,bi,bj)
498  #else  #else
499           Kgm_tmp = GM_isopycK           Kgm_tmp = isopycK*GM_isoFac2d(i,j,bi,bj)
500  #endif  #endif
501  #ifdef ALLOW_KAPGM_CONTROL  #ifdef ALLOW_KAPGM_CONTROL
502       &           + GM_skewflx*kapgm(i,j,k,bi,bj)       &           + GM_skewflx*kapgm(i,j,k,bi,bj)
503  #else  #else
504       &           + GM_skewflx*GM_background_K       &           + GM_skewflx*bolus_K*GM_bolFac2d(i,j,bi,bj)
505  #endif  #endif
506  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
507       &           + VisbeckK(i,j,bi,bj)*(1. _d 0 + GM_skewflx)       &           + VisbeckK(i,j,bi,bj)*(1. _d 0 + GM_skewflx)
# Line 461  CADJ STORE Kwz(:,:,k,bi,bj) = comlev1_bi Line 509  CADJ STORE Kwz(:,:,k,bi,bj) = comlev1_bi
509           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)
510           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)
511  #ifdef ALLOW_KAPREDI_CONTROL  #ifdef ALLOW_KAPREDI_CONTROL
512           Kwz(i,j,k,bi,bj)= ( kapredi(i,j,k,bi,bj)           Kwz(i,j,k,bi,bj)= ( kapredi(i,j,k,bi,bj)
513  #else  #else
514           Kwz(i,j,k,bi,bj)= ( GM_isopycK           Kwz(i,j,k,bi,bj)= ( isopycK*GM_isoFac2d(i,j,bi,bj)
515  #endif  #endif
516  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
517       &                     + VisbeckK(i,j,bi,bj)       &                     + VisbeckK(i,j,bi,bj)
# Line 573  c      IF ( GM_nonUnitDiag ) THEN Line 621  c      IF ( GM_nonUnitDiag ) THEN
621  #ifdef ALLOW_KAPREDI_CONTROL  #ifdef ALLOW_KAPREDI_CONTROL
622       &     ( kapredi(i,j,k,bi,bj)       &     ( kapredi(i,j,k,bi,bj)
623  #else  #else
624       &     ( GM_isopycK       &     ( GM_isopycK*GM_isoFac1d(k)
625         &        *op5*(GM_isoFac2d(i-1,j,bi,bj)+GM_isoFac2d(i,j,bi,bj))
626  #endif  #endif
627  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
628       &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))       &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))
# Line 607  CADJ STORE taperFct(:,:)     = comlev1_b Line 656  CADJ STORE taperFct(:,:)     = comlev1_b
656  #ifdef ALLOW_KAPREDI_CONTROL  #ifdef ALLOW_KAPREDI_CONTROL
657       &     ( kapredi(i,j,k,bi,bj)       &     ( kapredi(i,j,k,bi,bj)
658  #else  #else
659       &     ( GM_isopycK       &     ( GM_isopycK*GM_isoFac1d(k)
660         &        *op5*(GM_isoFac2d(i-1,j,bi,bj)+GM_isoFac2d(i,j,bi,bj))
661  #endif  #endif
662  #ifdef ALLOW_KAPGM_CONTROL  #ifdef ALLOW_KAPGM_CONTROL
663       &     - GM_skewflx*kapgm(i,j,k,bi,bj)       &     - GM_skewflx*kapgm(i,j,k,bi,bj)
664  #else  #else
665       &     - GM_skewflx*GM_background_K       &     - GM_skewflx*GM_background_K*GM_bolFac1d(k)
666         &        *op5*(GM_bolFac2d(i-1,j,bi,bj)+GM_bolFac2d(i,j,bi,bj))
667  #endif  #endif
668  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
669       &     +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
# Line 632  C         store in tmp1k Kuz_Redi Line 683  C         store in tmp1k Kuz_Redi
683  #ifdef ALLOW_KAPREDI_CONTROL  #ifdef ALLOW_KAPREDI_CONTROL
684            tmp1k(i,j) = ( kapredi(i,j,k,bi,bj)            tmp1k(i,j) = ( kapredi(i,j,k,bi,bj)
685  #else  #else
686            tmp1k(i,j) = ( GM_isopycK            tmp1k(i,j) = ( GM_isopycK*GM_isoFac1d(k)
687         &        *op5*(GM_isoFac2d(i-1,j,bi,bj)+GM_isoFac2d(i,j,bi,bj))
688  #endif  #endif
689  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
690       &     +(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))*0.5 _d 0       &     +(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))*0.5 _d 0
# Line 761  c      IF ( GM_nonUnitDiag ) THEN Line 813  c      IF ( GM_nonUnitDiag ) THEN
813  #ifdef ALLOW_KAPREDI_CONTROL  #ifdef ALLOW_KAPREDI_CONTROL
814       &     ( kapredi(i,j,k,bi,bj)       &     ( kapredi(i,j,k,bi,bj)
815  #else  #else
816       &     ( GM_isopycK       &     ( GM_isopycK*GM_isoFac1d(k)
817         &        *op5*(GM_isoFac2d(i,j-1,bi,bj)+GM_isoFac2d(i,j,bi,bj))
818  #endif  #endif
819  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
820       &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))       &     +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))
# Line 795  CADJ STORE taperFct(:,:)     = comlev1_b Line 848  CADJ STORE taperFct(:,:)     = comlev1_b
848  #ifdef ALLOW_KAPREDI_CONTROL  #ifdef ALLOW_KAPREDI_CONTROL
849       &     ( kapredi(i,j,k,bi,bj)       &     ( kapredi(i,j,k,bi,bj)
850  #else  #else
851       &     ( GM_isopycK       &     ( GM_isopycK*GM_isoFac1d(k)
852         &        *op5*(GM_isoFac2d(i,j-1,bi,bj)+GM_isoFac2d(i,j,bi,bj))
853  #endif  #endif
854  #ifdef ALLOW_KAPGM_CONTROL  #ifdef ALLOW_KAPGM_CONTROL
855       &     - GM_skewflx*kapgm(i,j,k,bi,bj)       &     - GM_skewflx*kapgm(i,j,k,bi,bj)
856  #else  #else
857       &     - GM_skewflx*GM_background_K       &     - GM_skewflx*GM_background_K*GM_bolFac1d(k)
858         &        *op5*(GM_bolFac2d(i,j-1,bi,bj)+GM_bolFac2d(i,j,bi,bj))
859  #endif  #endif
860  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
861       &     +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
# Line 820  C         store in tmp1k Kvz_Redi Line 875  C         store in tmp1k Kvz_Redi
875  #ifdef ALLOW_KAPREDI_CONTROL  #ifdef ALLOW_KAPREDI_CONTROL
876            tmp1k(i,j) = ( kapredi(i,j,k,bi,bj)            tmp1k(i,j) = ( kapredi(i,j,k,bi,bj)
877  #else  #else
878            tmp1k(i,j) = ( GM_isopycK            tmp1k(i,j) = ( GM_isopycK*GM_isoFac1d(k)
879         &        *op5*(GM_isoFac2d(i,j-1,bi,bj)+GM_isoFac2d(i,j,bi,bj))
880  #endif  #endif
881  #ifdef GM_VISBECK_VARIABLE_K  #ifdef GM_VISBECK_VARIABLE_K
882       &     +(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))*0.5 _d 0       &     +(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))*0.5 _d 0
# Line 861  C-- end 3rd  loop on vertical level inde Line 917  C-- end 3rd  loop on vertical level inde
917    
918  #ifdef GM_BOLUS_ADVEC  #ifdef GM_BOLUS_ADVEC
919        IF (GM_AdvForm) THEN        IF (GM_AdvForm) THEN
920         CALL GMREDI_CALC_PSI_B(  #ifdef GM_BOLUS_BVP
921           IF (GM_UseBVP) THEN
922            CALL GMREDI_CALC_PSI_BVP(
923       I             bi, bj, iMin, iMax, jMin, jMax,       I             bi, bj, iMin, iMax, jMin, jMax,
924       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
      I             ldd97_LrhoW, ldd97_LrhoS,  
925       I             myThid )       I             myThid )
926           ELSE
927    #endif
928            CALL GMREDI_CALC_PSI_B(
929         I              bi, bj, iMin, iMax, jMin, jMax,
930         I              sigmaX, sigmaY, sigmaR,
931         I              ldd97_LrhoW, ldd97_LrhoS,
932         I              myThid )
933    #ifdef GM_BOLUS_BVP
934           ENDIF
935    #endif
936        ENDIF        ENDIF
937  #endif  #endif
938    
# Line 893  C--   Time-average Line 960  C--   Time-average
960       &                          deltaTclock, bi, bj, myThid )       &                          deltaTclock, bi, bj, myThid )
961         ENDIF         ENDIF
962  #endif  #endif
963         DO k=1,Nr         GM_timeAve(bi,bj) = GM_timeAve(bi,bj)+deltaTclock
          GM_TimeAve(k,bi,bj)=GM_TimeAve(k,bi,bj)+deltaTclock  
        ENDDO  
964    
965        ENDIF        ENDIF
966  #endif /* ALLOW_TIMEAVE */  #endif /* ALLOW_TIMEAVE */
# Line 911  C--   Time-average Line 976  C--   Time-average
976        RETURN        RETURN
977        END        END
978    
979    C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
980    
981    CBOP
982    C     !ROUTINE: GMREDI_CALC_TENSOR_DUMMY
983    C     !INTERFACE:
984        SUBROUTINE GMREDI_CALC_TENSOR_DUMMY(        SUBROUTINE GMREDI_CALC_TENSOR_DUMMY(
985       I             iMin, iMax, jMin, jMax,       I             iMin, iMax, jMin, jMax,
986       I             sigmaX, sigmaY, sigmaR,       I             sigmaX, sigmaY, sigmaR,
987       I             bi, bj, myTime, myIter, myThid )       I             bi, bj, myTime, myIter, myThid )
988  C     /==========================================================\  
989  C     | SUBROUTINE GMREDI_CALC_TENSOR                            |  C     !DESCRIPTION: \bv
990  C     | o Calculate tensor elements for GM/Redi tensor.          |  C     *==========================================================*
991  C     |==========================================================|  C     | SUBROUTINE GMREDI_CALC_TENSOR_DUMMY
992  C     \==========================================================/  C     | o Calculate tensor elements for GM/Redi tensor.
993    C     *==========================================================*
994    C     \ev
995    
996    C     !USES:
997        IMPLICIT NONE        IMPLICIT NONE
998    
999  C     == Global variables ==  C     == Global variables ==
# Line 928  C     == Global variables == Line 1001  C     == Global variables ==
1001  #include "EEPARAMS.h"  #include "EEPARAMS.h"
1002  #include "GMREDI.h"  #include "GMREDI.h"
1003    
1004  C     == Routine arguments ==  C     !INPUT/OUTPUT PARAMETERS:
 C  
1005        _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
1006        _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
1007        _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
# Line 938  C Line 1010  C
1010        _RL     myTime        _RL     myTime
1011        INTEGER myIter        INTEGER myIter
1012        INTEGER myThid        INTEGER myThid
1013  CEndOfInterface  CEOP
1014    
1015  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
1016    C     !LOCAL VARIABLES:
1017        INTEGER i, j, k        INTEGER i, j, k
1018    
1019        DO k=1,Nr        DO k=1,Nr
Legend:
Removed from v.1.31  
changed lines
  Added in v.1.39
  ViewVC Help
Powered by ViewVC 1.1.22