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

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

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


Revision 1.31 - (hide annotations) (download)
Sat Feb 2 02:35:53 2008 UTC (16 years, 3 months ago) by gforget
Branch: MAIN
CVS Tags: checkpoint59o, checkpoint59n
Changes since 1.30: +42 -7 lines
introduce isopycnal diffusion coefficient control.

1 gforget 1.31 C $Header: /u/gcmpack/MITgcm/pkg/gmredi/gmredi_calc_tensor.F,v 1.30 2007/09/07 14:27:19 dfer Exp $
2 heimbach 1.13 C $Name: $
3 adcroft 1.1
4     #include "GMREDI_OPTIONS.h"
5 jmc 1.27 #ifdef ALLOW_KPP
6     # include "KPP_OPTIONS.h"
7     #endif
8 jmc 1.26 #undef OLD_VISBECK_CALC
9 adcroft 1.1
10 jmc 1.27 CBOP
11     C !ROUTINE: GMREDI_CALC_TENSOR
12     C !INTERFACE:
13 adcroft 1.1 SUBROUTINE GMREDI_CALC_TENSOR(
14 jmc 1.27 I iMin, iMax, jMin, jMax,
15 adcroft 1.1 I sigmaX, sigmaY, sigmaR,
16 jmc 1.27 I bi, bj, myTime, myIter, myThid )
17    
18     C !DESCRIPTION: \bv
19     C *==========================================================*
20     C | SUBROUTINE GMREDI_CALC_TENSOR
21     C | o Calculate tensor elements for GM/Redi tensor.
22     C *==========================================================*
23     C *==========================================================*
24     C \ev
25    
26     C !USES:
27 adcroft 1.1 IMPLICIT NONE
28    
29     C == Global variables ==
30     #include "SIZE.h"
31     #include "GRID.h"
32     #include "DYNVARS.h"
33     #include "EEPARAMS.h"
34     #include "PARAMS.h"
35     #include "GMREDI.h"
36 jmc 1.20 #include "GMREDI_TAVE.h"
37 jmc 1.27 #ifdef ALLOW_KPP
38     # include "KPP.h"
39     #endif
40 adcroft 1.1
41 heimbach 1.10 #ifdef ALLOW_AUTODIFF_TAMC
42     #include "tamc.h"
43     #include "tamc_keys.h"
44     #endif /* ALLOW_AUTODIFF_TAMC */
45    
46 jmc 1.27 C !INPUT/OUTPUT PARAMETERS:
47 adcroft 1.1 C == Routine arguments ==
48 jmc 1.27 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 adcroft 1.1 C
53 jmc 1.27 INTEGER iMin,iMax,jMin,jMax
54 adcroft 1.1 _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
55     _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
56     _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
57 jmc 1.27 INTEGER bi, bj
58     _RL myTime
59     INTEGER myIter
60 adcroft 1.1 INTEGER myThid
61 jmc 1.27 CEOP
62 adcroft 1.1
63     #ifdef ALLOW_GMREDI
64    
65 jmc 1.27 C !LOCAL VARIABLES:
66 adcroft 1.1 C == Local variables ==
67 jmc 1.15 INTEGER i,j,k,kp1
68 adcroft 1.1 _RL SlopeX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
69     _RL SlopeY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
70 heimbach 1.12 _RL dSigmaDx(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
71     _RL dSigmaDy(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
72 jmc 1.19 _RL dSigmaDr(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
73 jmc 1.8 _RL SlopeSqr(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
74     _RL taperFct(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
75 jmc 1.15 _RL maskp1, Kgm_tmp
76 jmc 1.19 _RL ldd97_LrhoC(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
77     _RL ldd97_LrhoW(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
78     _RL ldd97_LrhoS(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
79     _RL Cspd, LrhoInf, LrhoSup, fCoriLoc
80 adcroft 1.1
81 jmc 1.27 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 adcroft 1.1 #ifdef GM_VISBECK_VARIABLE_K
89 jmc 1.26 #ifdef OLD_VISBECK_CALC
90 heimbach 1.14 _RL deltaH,zero_rs
91     PARAMETER(zero_rs=0.D0)
92 adcroft 1.1 _RL N2,SN
93 heimbach 1.10 _RL Ssq(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
94 jmc 1.26 #else
95     _RL dSigmaH
96     _RL deltaH, integrDepth
97     _RL Sloc, M2loc, SNloc
98     #endif
99 adcroft 1.1 #endif
100    
101 jmc 1.22 #ifdef ALLOW_DIAGNOSTICS
102     LOGICAL doDiagRediFlx
103     LOGICAL DIAGNOSTICS_IS_ON
104     EXTERNAL DIAGNOSTICS_IS_ON
105     INTEGER km1
106     _RL dTdz
107     _RL tmp1k(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
108     #endif
109    
110 jmc 1.19 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
111    
112 heimbach 1.10 #ifdef ALLOW_AUTODIFF_TAMC
113     act1 = bi - myBxLo(myThid)
114     max1 = myBxHi(myThid) - myBxLo(myThid) + 1
115     act2 = bj - myByLo(myThid)
116     max2 = myByHi(myThid) - myByLo(myThid) + 1
117     act3 = myThid - 1
118     max3 = nTx*nTy
119     act4 = ikey_dynamics - 1
120 heimbach 1.12 igmkey = (act1 + 1) + act2*max1
121 heimbach 1.10 & + act3*max1*max2
122     & + act4*max1*max2*max3
123     #endif /* ALLOW_AUTODIFF_TAMC */
124    
125 jmc 1.22 #ifdef ALLOW_DIAGNOSTICS
126     doDiagRediFlx = .FALSE.
127     IF ( useDiagnostics ) THEN
128     doDiagRediFlx = DIAGNOSTICS_IS_ON('GM_KuzTz', myThid )
129 jmc 1.26 doDiagRediFlx = doDiagRediFlx .OR.
130 jmc 1.22 & DIAGNOSTICS_IS_ON('GM_KvzTz', myThid )
131     ENDIF
132     #endif
133 jmc 1.26
134 heimbach 1.12 #ifdef GM_VISBECK_VARIABLE_K
135     DO j=1-Oly,sNy+Oly
136     DO i=1-Olx,sNx+Olx
137     VisbeckK(i,j,bi,bj) = 0. _d 0
138     ENDDO
139     ENDDO
140     #endif
141    
142 jmc 1.19 C-- set ldd97_Lrho (for tapering scheme ldd97):
143 jmc 1.27 IF ( GM_taper_scheme.EQ.'ldd97' .OR.
144     & GM_taper_scheme.EQ.'fm07' ) THEN
145 jmc 1.19 Cspd = 2. _d 0
146     LrhoInf = 15. _d 3
147     LrhoSup = 100. _d 3
148     C- Tracer point location (center):
149     DO j=1-Oly,sNy+Oly
150     DO i=1-Olx,sNx+Olx
151     IF (fCori(i,j,bi,bj).NE.0.) THEN
152     ldd97_LrhoC(i,j) = Cspd/ABS(fCori(i,j,bi,bj))
153     ELSE
154     ldd97_LrhoC(i,j) = LrhoSup
155     ENDIF
156     ldd97_LrhoC(i,j) = MAX(LrhoInf,MIN(ldd97_LrhoC(i,j),LrhoSup))
157     ENDDO
158     ENDDO
159     C- U point location (West):
160     DO j=1-Oly,sNy+Oly
161 jmc 1.27 kLow_W(1-Olx,j) = 0
162 jmc 1.19 ldd97_LrhoW(1-Olx,j) = LrhoSup
163     DO i=1-Olx+1,sNx+Olx
164 jmc 1.27 kLow_W(i,j) = MIN(kLowC(i-1,j,bi,bj),kLowC(i,j,bi,bj))
165 jmc 1.19 fCoriLoc = op5*(fCori(i-1,j,bi,bj)+fCori(i,j,bi,bj))
166     IF (fCoriLoc.NE.0.) THEN
167     ldd97_LrhoW(i,j) = Cspd/ABS(fCoriLoc)
168     ELSE
169     ldd97_LrhoW(i,j) = LrhoSup
170     ENDIF
171     ldd97_LrhoW(i,j) = MAX(LrhoInf,MIN(ldd97_LrhoW(i,j),LrhoSup))
172     ENDDO
173     ENDDO
174     C- V point location (South):
175     DO i=1-Olx+1,sNx+Olx
176 jmc 1.27 kLow_S(i,1-Oly) = 0
177 jmc 1.19 ldd97_LrhoS(i,1-Oly) = LrhoSup
178     ENDDO
179     DO j=1-Oly+1,sNy+Oly
180     DO i=1-Olx,sNx+Olx
181 jmc 1.27 kLow_S(i,j) = MIN(kLowC(i,j-1,bi,bj),kLowC(i,j,bi,bj))
182 jmc 1.19 fCoriLoc = op5*(fCori(i,j-1,bi,bj)+fCori(i,j,bi,bj))
183     IF (fCoriLoc.NE.0.) THEN
184     ldd97_LrhoS(i,j) = Cspd/ABS(fCoriLoc)
185     ELSE
186     ldd97_LrhoS(i,j) = LrhoSup
187     ENDIF
188     ldd97_LrhoS(i,j) = MAX(LrhoInf,MIN(ldd97_LrhoS(i,j),LrhoSup))
189     ENDDO
190     ENDDO
191     ELSE
192     C- Just initialize to zero (not use anyway)
193     DO j=1-Oly,sNy+Oly
194     DO i=1-Olx,sNx+Olx
195     ldd97_LrhoC(i,j) = 0. _d 0
196     ldd97_LrhoW(i,j) = 0. _d 0
197     ldd97_LrhoS(i,j) = 0. _d 0
198     ENDDO
199     ENDDO
200     ENDIF
201 jmc 1.27
202 jmc 1.19 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
203 jmc 1.27 C-- 1rst loop on k : compute Tensor Coeff. at W points.
204    
205 heimbach 1.28 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 jmc 1.27 #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 jmc 1.19
231 jmc 1.27 DO k=Nr,2,-1
232 adcroft 1.1
233     #ifdef ALLOW_AUTODIFF_TAMC
234 heimbach 1.12 kkey = (igmkey-1)*Nr + k
235 heimbach 1.10 DO j=1-Oly,sNy+Oly
236     DO i=1-Olx,sNx+Olx
237     SlopeX(i,j) = 0. _d 0
238     SlopeY(i,j) = 0. _d 0
239 heimbach 1.12 dSigmaDx(i,j) = 0. _d 0
240     dSigmaDy(i,j) = 0. _d 0
241 jmc 1.19 dSigmaDr(i,j) = 0. _d 0
242 heimbach 1.10 SlopeSqr(i,j) = 0. _d 0
243     taperFct(i,j) = 0. _d 0
244     Kwx(i,j,k,bi,bj) = 0. _d 0
245     Kwy(i,j,k,bi,bj) = 0. _d 0
246     Kwz(i,j,k,bi,bj) = 0. _d 0
247 heimbach 1.12 # ifdef GM_NON_UNITY_DIAGONAL
248     Kux(i,j,k,bi,bj) = 0. _d 0
249     Kvy(i,j,k,bi,bj) = 0. _d 0
250     # endif
251     # ifdef GM_EXTRA_DIAGONAL
252     Kuz(i,j,k,bi,bj) = 0. _d 0
253     Kvz(i,j,k,bi,bj) = 0. _d 0
254     # endif
255     # ifdef GM_BOLUS_ADVEC
256     GM_PsiX(i,j,k,bi,bj) = 0. _d 0
257     GM_PsiY(i,j,k,bi,bj) = 0. _d 0
258     # endif
259 heimbach 1.10 ENDDO
260     ENDDO
261 adcroft 1.1 #endif
262 heimbach 1.10
263 jmc 1.26 DO j=1-Oly+1,sNy+Oly-1
264     DO i=1-Olx+1,sNx+Olx-1
265 adcroft 1.1 C Gradient of Sigma at rVel points
266 jmc 1.26 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 )
268     & )*maskC(i,j,k,bi,bj)
269     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 )
271     & )*maskC(i,j,k,bi,bj)
272     dSigmaDr(i,j)=sigmaR(i,j,k)
273     ENDDO
274 adcroft 1.1 ENDDO
275    
276 heimbach 1.10 #ifdef ALLOW_AUTODIFF_TAMC
277 heimbach 1.12 CADJ STORE dSigmaDx(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
278     CADJ STORE dSigmaDy(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
279 jmc 1.19 CADJ STORE dSigmaDr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
280 heimbach 1.28 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 heimbach 1.10 #endif /* ALLOW_AUTODIFF_TAMC */
284    
285 jmc 1.26 #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 adcroft 1.1 C Calculate slopes for use in tensor, taper and/or clip
340 jmc 1.26 CALL GMREDI_SLOPE_LIMIT(
341 jmc 1.19 O SlopeX, SlopeY,
342 jmc 1.8 O SlopeSqr, taperFct,
343 jmc 1.27 U hTransLay, baseSlope, recipLambda,
344 jmc 1.19 U dSigmaDr,
345     I dSigmaDx, dSigmaDy,
346 jmc 1.27 I ldd97_LrhoC, locMixLayer, rF,
347     I kLowC(1-Olx,1-Oly,bi,bj),
348     I k, bi, bj, myTime, myIter, myThid )
349 adcroft 1.1
350 jmc 1.26 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 heimbach 1.10 ENDDO
358    
359     #ifdef ALLOW_AUTODIFF_TAMC
360 heimbach 1.16 CADJ STORE SlopeX(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
361     CADJ STORE SlopeY(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
362 heimbach 1.14 CADJ STORE SlopeSqr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
363 jmc 1.19 CADJ STORE dSigmaDr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
364     CADJ STORE taperFct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
365 heimbach 1.10 #endif /* ALLOW_AUTODIFF_TAMC */
366    
367 jmc 1.27 C Components of Redi/GM tensor
368 jmc 1.26 DO j=1-Oly+1,sNy+Oly-1
369     DO i=1-Olx+1,sNx+Olx-1
370 jmc 1.27 Kwx(i,j,k,bi,bj)= SlopeX(i,j)*taperFct(i,j)
371     Kwy(i,j,k,bi,bj)= SlopeY(i,j)*taperFct(i,j)
372     Kwz(i,j,k,bi,bj)= SlopeSqr(i,j)*taperFct(i,j)
373     ENDDO
374     ENDDO
375 adcroft 1.1
376     #ifdef GM_VISBECK_VARIABLE_K
377 jmc 1.26 #ifdef OLD_VISBECK_CALC
378 jmc 1.27 DO j=1-Oly+1,sNy+Oly-1
379     DO i=1-Olx+1,sNx+Olx-1
380 jmc 1.8
381 jmc 1.27 C- note (jmc) : moved here since only used in VISBECK_VARIABLE_K
382     C but do not know if *taperFct (or **2 ?) is necessary
383 heimbach 1.10 Ssq(i,j)=SlopeSqr(i,j)*taperFct(i,j)
384 jmc 1.8
385 adcroft 1.1 C-- Depth average of M^2/N^2 * N
386    
387     C Calculate terms for mean Richardson number
388     C which is used in the "variable K" parameterisaton.
389     C Distance between interface above layer and the integration depth
390     deltaH=abs(GM_Visbeck_depth)-abs(rF(k))
391     C If positive we limit this to the layer thickness
392     deltaH=min(deltaH,drF(k))
393     C If negative then we are below the integration level
394     deltaH=max(deltaH,zero_rs)
395     C Now we convert deltaH to a non-dimensional fraction
396     deltaH=deltaH/GM_Visbeck_depth
397    
398 jmc 1.8 IF (K.eq.2) VisbeckK(i,j,bi,bj)=0.
399 jmc 1.19 IF ( Ssq(i,j).NE.0. .AND. dSigmaDr(i,j).NE.0. ) THEN
400     N2= -Gravity*recip_RhoConst*dSigmaDr(i,j)
401 heimbach 1.10 SN=sqrt(Ssq(i,j)*N2)
402 heimbach 1.3 VisbeckK(i,j,bi,bj)=VisbeckK(i,j,bi,bj)+deltaH
403 adcroft 1.1 & *GM_Visbeck_alpha*GM_Visbeck_length*GM_Visbeck_length*SN
404 jmc 1.8 ENDIF
405 adcroft 1.1
406 jmc 1.27 ENDDO
407     ENDDO
408 jmc 1.26 #endif /* OLD_VISBECK_CALC */
409 jmc 1.9 #endif /* GM_VISBECK_VARIABLE_K */
410    
411     C-- end 1rst loop on vertical level index k
412     ENDDO
413    
414 adcroft 1.1
415 jmc 1.9 #ifdef GM_VISBECK_VARIABLE_K
416 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
417     CADJ STORE VisbeckK(:,:,bi,bj) = comlev1_bibj, key=igmkey, byte=isbyte
418     #endif
419 jmc 1.27 IF ( GM_Visbeck_alpha.GT.0. ) THEN
420 jmc 1.9 C- Limit range that KapGM can take
421     DO j=1-Oly+1,sNy+Oly-1
422     DO i=1-Olx+1,sNx+Olx-1
423     VisbeckK(i,j,bi,bj)=
424     & MIN(VisbeckK(i,j,bi,bj),GM_Visbeck_maxval_K)
425     ENDDO
426     ENDDO
427     ENDIF
428 heimbach 1.16 cph( NEW
429     #ifdef ALLOW_AUTODIFF_TAMC
430     CADJ STORE VisbeckK(:,:,bi,bj) = comlev1_bibj, key=igmkey, byte=isbyte
431     #endif
432     cph)
433 adcroft 1.1 #endif /* GM_VISBECK_VARIABLE_K */
434    
435 heimbach 1.29 C- express the Tensor in term of Diffusivity (= m**2 / s )
436     DO k=1,Nr
437 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
438     kkey = (igmkey-1)*Nr + k
439 heimbach 1.29 # if (defined (GM_NON_UNITY_DIAGONAL) || \
440 heimbach 1.16 defined (GM_VISBECK_VARIABLE_K))
441 heimbach 1.14 CADJ STORE Kwx(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
442     CADJ STORE Kwy(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
443     CADJ STORE Kwz(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
444 heimbach 1.29 # endif
445 heimbach 1.12 #endif
446 jmc 1.26 DO j=1-Oly+1,sNy+Oly-1
447     DO i=1-Olx+1,sNx+Olx-1
448 gforget 1.31 #ifdef ALLOW_KAPREDI_CONTROL
449     Kgm_tmp = kapredi(i,j,k,bi,bj)
450     #else
451     Kgm_tmp = GM_isopycK
452     #endif
453 heimbach 1.25 #ifdef ALLOW_KAPGM_CONTROL
454 gforget 1.31 & + GM_skewflx*kapgm(i,j,k,bi,bj)
455 heimbach 1.25 #else
456 gforget 1.31 & + GM_skewflx*GM_background_K
457 heimbach 1.25 #endif
458 jmc 1.9 #ifdef GM_VISBECK_VARIABLE_K
459 jmc 1.26 & + VisbeckK(i,j,bi,bj)*(1. _d 0 + GM_skewflx)
460 jmc 1.9 #endif
461 jmc 1.26 Kwx(i,j,k,bi,bj)= Kgm_tmp*Kwx(i,j,k,bi,bj)
462     Kwy(i,j,k,bi,bj)= Kgm_tmp*Kwy(i,j,k,bi,bj)
463 gforget 1.31 #ifdef ALLOW_KAPREDI_CONTROL
464     Kwz(i,j,k,bi,bj)= ( kapredi(i,j,k,bi,bj)
465     #else
466 jmc 1.26 Kwz(i,j,k,bi,bj)= ( GM_isopycK
467 gforget 1.31 #endif
468 adcroft 1.1 #ifdef GM_VISBECK_VARIABLE_K
469 jmc 1.26 & + VisbeckK(i,j,bi,bj)
470 adcroft 1.1 #endif
471 jmc 1.26 & )*Kwz(i,j,k,bi,bj)
472     ENDDO
473 adcroft 1.1 ENDDO
474 jmc 1.27 ENDDO
475    
476     #ifdef ALLOW_DIAGNOSTICS
477     IF ( useDiagnostics .AND. GM_taper_scheme.EQ.'fm07' ) THEN
478     CALL DIAGNOSTICS_FILL( hTransLay, 'GM_hTrsL', 0,1,2,bi,bj,myThid)
479     CALL DIAGNOSTICS_FILL( baseSlope, 'GM_baseS', 0,1,2,bi,bj,myThid)
480     CALL DIAGNOSTICS_FILL(recipLambda,'GM_rLamb', 0,1,2,bi,bj,myThid)
481     ENDIF
482     #endif /* ALLOW_DIAGNOSTICS */
483    
484 adcroft 1.4
485 jmc 1.9 #if ( defined (GM_NON_UNITY_DIAGONAL) || defined (GM_EXTRA_DIAGONAL) )
486 jmc 1.27 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
487     C-- 2nd k loop : compute Tensor Coeff. at U point
488    
489     #ifdef ALLOW_KPP
490     IF ( useKPP ) THEN
491     DO j=1-Oly,sNy+Oly
492     DO i=2-Olx,sNx+Olx
493     locMixLayer(i,j) = ( KPPhbl(i-1,j,bi,bj)
494     & + KPPhbl( i ,j,bi,bj) )*op5
495     ENDDO
496     ENDDO
497     ELSE
498     #else
499     IF ( .TRUE. ) THEN
500     #endif
501     DO j=1-Oly,sNy+Oly
502     DO i=2-Olx,sNx+Olx
503     locMixLayer(i,j) = ( hMixLayer(i-1,j,bi,bj)
504     & + hMixLayer( i ,j,bi,bj) )*op5
505     ENDDO
506     ENDDO
507     ENDIF
508     DO j=1-Oly,sNy+Oly
509     DO i=1-Olx,sNx+Olx
510     hTransLay(i,j) = 0.
511     baseSlope(i,j) = 0.
512     recipLambda(i,j)= 0.
513     ENDDO
514     DO i=2-Olx,sNx+Olx
515     hTransLay(i,j) = MAX( R_low(i-1,j,bi,bj), R_low(i,j,bi,bj) )
516     ENDDO
517     ENDDO
518    
519     DO k=Nr,1,-1
520     kp1 = MIN(Nr,k+1)
521     maskp1 = 1. _d 0
522     IF (k.GE.Nr) maskp1 = 0. _d 0
523 heimbach 1.29 #ifdef ALLOW_AUTODIFF_TAMC
524     kkey = (igmkey-1)*Nr + k
525     #endif
526 adcroft 1.1
527     C Gradient of Sigma at U points
528 jmc 1.26 DO j=1-Oly+1,sNy+Oly-1
529     DO i=1-Olx+1,sNx+Olx-1
530     dSigmaDx(i,j)=sigmaX(i,j,k)
531     & *_maskW(i,j,k,bi,bj)
532     dSigmaDy(i,j)=op25*( sigmaY(i-1,j+1,k)+sigmaY(i,j+1,k)
533     & +sigmaY(i-1, j ,k)+sigmaY(i, j ,k)
534     & )*_maskW(i,j,k,bi,bj)
535     dSigmaDr(i,j)=op25*( sigmaR(i-1,j, k )+sigmaR(i,j, k )
536     & +(sigmaR(i-1,j,kp1)+sigmaR(i,j,kp1))*maskp1
537     & )*_maskW(i,j,k,bi,bj)
538     ENDDO
539 adcroft 1.1 ENDDO
540    
541 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
542 heimbach 1.17 CADJ STORE SlopeSqr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
543 heimbach 1.12 CADJ STORE dSigmaDx(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
544     CADJ STORE dSigmaDy(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
545 heimbach 1.28 CADJ STORE dSigmaDr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
546     CADJ STORE locMixLayer(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
547     CADJ STORE baseSlope(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
548     CADJ STORE hTransLay(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
549     CADJ STORE recipLambda(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
550 heimbach 1.12 #endif /* ALLOW_AUTODIFF_TAMC */
551    
552 adcroft 1.1 C Calculate slopes for use in tensor, taper and/or clip
553 jmc 1.26 CALL GMREDI_SLOPE_LIMIT(
554 jmc 1.19 O SlopeX, SlopeY,
555 jmc 1.8 O SlopeSqr, taperFct,
556 jmc 1.27 U hTransLay, baseSlope, recipLambda,
557 jmc 1.19 U dSigmaDr,
558     I dSigmaDx, dSigmaDy,
559 jmc 1.27 I ldd97_LrhoW, locMixLayer, rC,
560     I kLow_W,
561     I k, bi, bj, myTime, myIter, myThid )
562 adcroft 1.1
563 heimbach 1.16 #ifdef ALLOW_AUTODIFF_TAMC
564 heimbach 1.17 CADJ STORE SlopeSqr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
565 jmc 1.19 CADJ STORE taperFct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
566 heimbach 1.16 #endif /* ALLOW_AUTODIFF_TAMC */
567    
568 jmc 1.9 #ifdef GM_NON_UNITY_DIAGONAL
569 jmc 1.26 c IF ( GM_nonUnitDiag ) THEN
570 jmc 1.9 DO j=1-Oly+1,sNy+Oly-1
571     DO i=1-Olx+1,sNx+Olx-1
572     Kux(i,j,k,bi,bj) =
573 gforget 1.31 #ifdef ALLOW_KAPREDI_CONTROL
574     & ( kapredi(i,j,k,bi,bj)
575     #else
576 jmc 1.9 & ( GM_isopycK
577 gforget 1.31 #endif
578 jmc 1.9 #ifdef GM_VISBECK_VARIABLE_K
579 heimbach 1.14 & +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))
580 jmc 1.9 #endif
581 jmc 1.27 & )*taperFct(i,j)
582 heimbach 1.10 ENDDO
583     ENDDO
584 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
585 heimbach 1.16 # ifdef GM_EXCLUDE_CLIPPING
586 heimbach 1.12 CADJ STORE Kux(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
587     # endif
588     #endif
589 heimbach 1.10 DO j=1-Oly+1,sNy+Oly-1
590     DO i=1-Olx+1,sNx+Olx-1
591 jmc 1.9 Kux(i,j,k,bi,bj) = MAX( Kux(i,j,k,bi,bj), GM_Kmin_horiz )
592     ENDDO
593     ENDDO
594 jmc 1.26 c ENDIF
595 jmc 1.9 #endif /* GM_NON_UNITY_DIAGONAL */
596    
597     #ifdef GM_EXTRA_DIAGONAL
598 heimbach 1.12
599     #ifdef ALLOW_AUTODIFF_TAMC
600     CADJ STORE SlopeX(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
601     CADJ STORE taperFct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
602     #endif /* ALLOW_AUTODIFF_TAMC */
603 jmc 1.27 IF ( GM_ExtraDiag ) THEN
604 jmc 1.9 DO j=1-Oly+1,sNy+Oly-1
605     DO i=1-Olx+1,sNx+Olx-1
606     Kuz(i,j,k,bi,bj) =
607 gforget 1.31 #ifdef ALLOW_KAPREDI_CONTROL
608     & ( kapredi(i,j,k,bi,bj)
609     #else
610     & ( GM_isopycK
611     #endif
612 heimbach 1.25 #ifdef ALLOW_KAPGM_CONTROL
613 gforget 1.31 & - GM_skewflx*kapgm(i,j,k,bi,bj)
614 heimbach 1.25 #else
615 gforget 1.31 & - GM_skewflx*GM_background_K
616 heimbach 1.25 #endif
617 jmc 1.9 #ifdef GM_VISBECK_VARIABLE_K
618 heimbach 1.14 & +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))*GM_advect
619 jmc 1.9 #endif
620     & )*SlopeX(i,j)*taperFct(i,j)
621     ENDDO
622     ENDDO
623 jmc 1.26 ENDIF
624 jmc 1.9 #endif /* GM_EXTRA_DIAGONAL */
625 adcroft 1.1
626 jmc 1.22 #ifdef ALLOW_DIAGNOSTICS
627 jmc 1.26 IF (doDiagRediFlx) THEN
628 jmc 1.22 km1 = MAX(k-1,1)
629     DO j=1,sNy
630     DO i=1,sNx+1
631     C store in tmp1k Kuz_Redi
632 gforget 1.31 #ifdef ALLOW_KAPREDI_CONTROL
633     tmp1k(i,j) = ( kapredi(i,j,k,bi,bj)
634     #else
635 jmc 1.22 tmp1k(i,j) = ( GM_isopycK
636 gforget 1.31 #endif
637 jmc 1.22 #ifdef GM_VISBECK_VARIABLE_K
638     & +(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))*0.5 _d 0
639     #endif
640     & )*SlopeX(i,j)*taperFct(i,j)
641     ENDDO
642     ENDDO
643     DO j=1,sNy
644     DO i=1,sNx+1
645     C- Vertical gradients interpolated to U points
646     dTdz = (
647     & +recip_drC(k)*
648     & ( maskC(i-1,j,k,bi,bj)*
649     & (theta(i-1,j,km1,bi,bj)-theta(i-1,j,k,bi,bj))
650     & +maskC( i ,j,k,bi,bj)*
651     & (theta( i ,j,km1,bi,bj)-theta( i ,j,k,bi,bj))
652     & )
653     & +recip_drC(kp1)*
654     & ( maskC(i-1,j,kp1,bi,bj)*
655     & (theta(i-1,j,k,bi,bj)-theta(i-1,j,kp1,bi,bj))
656     & +maskC( i ,j,kp1,bi,bj)*
657     & (theta( i ,j,k,bi,bj)-theta( i ,j,kp1,bi,bj))
658     & ) ) * 0.25 _d 0
659 heimbach 1.23 tmp1k(i,j) = dyG(i,j,bi,bj)*drF(k)
660     & * _hFacW(i,j,k,bi,bj)
661 jmc 1.22 & * tmp1k(i,j) * dTdz
662     ENDDO
663     ENDDO
664     CALL DIAGNOSTICS_FILL(tmp1k, 'GM_KuzTz', k,1,2,bi,bj,myThid)
665 jmc 1.26 ENDIF
666 jmc 1.22 #endif /* ALLOW_DIAGNOSTICS */
667    
668 jmc 1.27 C-- end 2nd loop on vertical level index k
669     ENDDO
670    
671     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
672     C-- 3rd k loop : compute Tensor Coeff. at V point
673    
674     #ifdef ALLOW_KPP
675     IF ( useKPP ) THEN
676     DO j=2-Oly,sNy+Oly
677     DO i=1-Olx,sNx+Olx
678     locMixLayer(i,j) = ( KPPhbl(i,j-1,bi,bj)
679     & + KPPhbl(i, j ,bi,bj) )*op5
680     ENDDO
681     ENDDO
682     ELSE
683     #else
684     IF ( .TRUE. ) THEN
685     #endif
686     DO j=2-Oly,sNy+Oly
687     DO i=1-Olx,sNx+Olx
688     locMixLayer(i,j) = ( hMixLayer(i,j-1,bi,bj)
689     & + hMixLayer(i, j ,bi,bj) )*op5
690     ENDDO
691     ENDDO
692     ENDIF
693     DO j=1-Oly,sNy+Oly
694     DO i=1-Olx,sNx+Olx
695     hTransLay(i,j) = 0.
696     baseSlope(i,j) = 0.
697     recipLambda(i,j)= 0.
698     ENDDO
699     ENDDO
700     DO j=2-Oly,sNy+Oly
701     DO i=1-Olx,sNx+Olx
702     hTransLay(i,j) = MAX( R_low(i,j-1,bi,bj), R_low(i,j,bi,bj) )
703     ENDDO
704     ENDDO
705    
706 adcroft 1.1 C Gradient of Sigma at V points
707 jmc 1.27 DO k=Nr,1,-1
708     kp1 = MIN(Nr,k+1)
709     maskp1 = 1. _d 0
710     IF (k.GE.Nr) maskp1 = 0. _d 0
711 heimbach 1.29 #ifdef ALLOW_AUTODIFF_TAMC
712     kkey = (igmkey-1)*Nr + k
713     #endif
714 jmc 1.27
715 jmc 1.26 DO j=1-Oly+1,sNy+Oly-1
716     DO i=1-Olx+1,sNx+Olx-1
717     dSigmaDx(i,j)=op25*( sigmaX(i, j ,k) +sigmaX(i+1, j ,k)
718     & +sigmaX(i,j-1,k) +sigmaX(i+1,j-1,k)
719     & )*_maskS(i,j,k,bi,bj)
720     dSigmaDy(i,j)=sigmaY(i,j,k)
721     & *_maskS(i,j,k,bi,bj)
722     dSigmaDr(i,j)=op25*( sigmaR(i,j-1, k )+sigmaR(i,j, k )
723     & +(sigmaR(i,j-1,kp1)+sigmaR(i,j,kp1))*maskp1
724     & )*_maskS(i,j,k,bi,bj)
725     ENDDO
726 adcroft 1.1 ENDDO
727    
728 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
729     CADJ STORE dSigmaDx(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
730     CADJ STORE dSigmaDy(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
731 heimbach 1.28 CADJ STORE dSigmaDr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
732     CADJ STORE baseSlope(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
733     CADJ STORE hTransLay(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
734     CADJ STORE recipLambda(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
735 heimbach 1.12 #endif /* ALLOW_AUTODIFF_TAMC */
736    
737 adcroft 1.1 C Calculate slopes for use in tensor, taper and/or clip
738 jmc 1.26 CALL GMREDI_SLOPE_LIMIT(
739 jmc 1.19 O SlopeX, SlopeY,
740 jmc 1.8 O SlopeSqr, taperFct,
741 jmc 1.27 U hTransLay, baseSlope, recipLambda,
742 jmc 1.19 U dSigmaDr,
743     I dSigmaDx, dSigmaDy,
744 jmc 1.27 I ldd97_LrhoS, locMixLayer, rC,
745     I kLow_S,
746     I k, bi, bj, myTime, myIter, myThid )
747 adcroft 1.1
748 heimbach 1.16 cph(
749     #ifdef ALLOW_AUTODIFF_TAMC
750     cph(
751     CADJ STORE taperfct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
752     cph)
753     #endif /* ALLOW_AUTODIFF_TAMC */
754     cph)
755    
756 jmc 1.9 #ifdef GM_NON_UNITY_DIAGONAL
757 jmc 1.26 c IF ( GM_nonUnitDiag ) THEN
758 jmc 1.9 DO j=1-Oly+1,sNy+Oly-1
759     DO i=1-Olx+1,sNx+Olx-1
760     Kvy(i,j,k,bi,bj) =
761 gforget 1.31 #ifdef ALLOW_KAPREDI_CONTROL
762     & ( kapredi(i,j,k,bi,bj)
763     #else
764 jmc 1.9 & ( GM_isopycK
765 gforget 1.31 #endif
766 jmc 1.9 #ifdef GM_VISBECK_VARIABLE_K
767 heimbach 1.14 & +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))
768 jmc 1.9 #endif
769 jmc 1.27 & )*taperFct(i,j)
770 heimbach 1.10 ENDDO
771     ENDDO
772 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
773 heimbach 1.16 # ifdef GM_EXCLUDE_CLIPPING
774 heimbach 1.12 CADJ STORE Kvy(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
775     # endif
776     #endif
777 heimbach 1.10 DO j=1-Oly+1,sNy+Oly-1
778     DO i=1-Olx+1,sNx+Olx-1
779 jmc 1.9 Kvy(i,j,k,bi,bj) = MAX( Kvy(i,j,k,bi,bj), GM_Kmin_horiz )
780     ENDDO
781     ENDDO
782 jmc 1.26 c ENDIF
783 jmc 1.9 #endif /* GM_NON_UNITY_DIAGONAL */
784    
785     #ifdef GM_EXTRA_DIAGONAL
786 heimbach 1.12
787     #ifdef ALLOW_AUTODIFF_TAMC
788     CADJ STORE SlopeY(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
789     CADJ STORE taperFct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
790     #endif /* ALLOW_AUTODIFF_TAMC */
791 jmc 1.27 IF ( GM_ExtraDiag ) THEN
792 jmc 1.9 DO j=1-Oly+1,sNy+Oly-1
793     DO i=1-Olx+1,sNx+Olx-1
794     Kvz(i,j,k,bi,bj) =
795 gforget 1.31 #ifdef ALLOW_KAPREDI_CONTROL
796     & ( kapredi(i,j,k,bi,bj)
797     #else
798     & ( GM_isopycK
799     #endif
800 heimbach 1.25 #ifdef ALLOW_KAPGM_CONTROL
801 gforget 1.31 & - GM_skewflx*kapgm(i,j,k,bi,bj)
802 heimbach 1.25 #else
803 gforget 1.31 & - GM_skewflx*GM_background_K
804 heimbach 1.25 #endif
805 jmc 1.9 #ifdef GM_VISBECK_VARIABLE_K
806 heimbach 1.14 & +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))*GM_advect
807 jmc 1.9 #endif
808     & )*SlopeY(i,j)*taperFct(i,j)
809     ENDDO
810     ENDDO
811 jmc 1.26 ENDIF
812 jmc 1.9 #endif /* GM_EXTRA_DIAGONAL */
813    
814 jmc 1.22 #ifdef ALLOW_DIAGNOSTICS
815 jmc 1.26 IF (doDiagRediFlx) THEN
816 dfer 1.30 km1 = MAX(k-1,1)
817 jmc 1.22 DO j=1,sNy+1
818     DO i=1,sNx
819     C store in tmp1k Kvz_Redi
820 gforget 1.31 #ifdef ALLOW_KAPREDI_CONTROL
821     tmp1k(i,j) = ( kapredi(i,j,k,bi,bj)
822     #else
823 jmc 1.22 tmp1k(i,j) = ( GM_isopycK
824 gforget 1.31 #endif
825 jmc 1.22 #ifdef GM_VISBECK_VARIABLE_K
826     & +(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))*0.5 _d 0
827     #endif
828     & )*SlopeY(i,j)*taperFct(i,j)
829     ENDDO
830     ENDDO
831     DO j=1,sNy+1
832     DO i=1,sNx
833     C- Vertical gradients interpolated to U points
834     dTdz = (
835     & +recip_drC(k)*
836     & ( maskC(i,j-1,k,bi,bj)*
837     & (theta(i,j-1,km1,bi,bj)-theta(i,j-1,k,bi,bj))
838     & +maskC(i, j ,k,bi,bj)*
839     & (theta(i, j ,km1,bi,bj)-theta(i, j ,k,bi,bj))
840     & )
841     & +recip_drC(kp1)*
842     & ( maskC(i,j-1,kp1,bi,bj)*
843     & (theta(i,j-1,k,bi,bj)-theta(i,j-1,kp1,bi,bj))
844     & +maskC(i, j ,kp1,bi,bj)*
845     & (theta(i, j ,k,bi,bj)-theta(i, j ,kp1,bi,bj))
846     & ) ) * 0.25 _d 0
847 heimbach 1.23 tmp1k(i,j) = dxG(i,j,bi,bj)*drF(k)
848     & * _hFacS(i,j,k,bi,bj)
849 jmc 1.22 & * tmp1k(i,j) * dTdz
850     ENDDO
851     ENDDO
852     CALL DIAGNOSTICS_FILL(tmp1k, 'GM_KvzTz', k,1,2,bi,bj,myThid)
853 jmc 1.26 ENDIF
854 jmc 1.22 #endif /* ALLOW_DIAGNOSTICS */
855    
856 jmc 1.27 C-- end 3rd loop on vertical level index k
857     ENDDO
858    
859 jmc 1.9 #endif /* GM_NON_UNITY_DIAGONAL || GM_EXTRA_DIAGONAL */
860    
861 adcroft 1.1
862 jmc 1.9 #ifdef GM_BOLUS_ADVEC
863     IF (GM_AdvForm) THEN
864 jmc 1.26 CALL GMREDI_CALC_PSI_B(
865 jmc 1.9 I bi, bj, iMin, iMax, jMin, jMax,
866     I sigmaX, sigmaY, sigmaR,
867 jmc 1.19 I ldd97_LrhoW, ldd97_LrhoS,
868 jmc 1.26 I myThid )
869 jmc 1.9 ENDIF
870     #endif
871 adcroft 1.1
872 jmc 1.19 #ifdef ALLOW_TIMEAVE
873     C-- Time-average
874     IF ( taveFreq.GT.0. ) THEN
875    
876     CALL TIMEAVE_CUMULATE( GM_Kwx_T, Kwx, Nr,
877     & deltaTclock, bi, bj, myThid )
878     CALL TIMEAVE_CUMULATE( GM_Kwy_T, Kwy, Nr,
879     & deltaTclock, bi, bj, myThid )
880     CALL TIMEAVE_CUMULATE( GM_Kwz_T, Kwz, Nr,
881     & deltaTclock, bi, bj, myThid )
882     #ifdef GM_VISBECK_VARIABLE_K
883     IF ( GM_Visbeck_alpha.NE.0. ) THEN
884     CALL TIMEAVE_CUMULATE( Visbeck_K_T, VisbeckK, 1,
885     & deltaTclock, bi, bj, myThid )
886     ENDIF
887     #endif
888     #ifdef GM_BOLUS_ADVEC
889     IF ( GM_AdvForm ) THEN
890     CALL TIMEAVE_CUMULATE( GM_PsiXtave, GM_PsiX, Nr,
891     & deltaTclock, bi, bj, myThid )
892     CALL TIMEAVE_CUMULATE( GM_PsiYtave, GM_PsiY, Nr,
893     & deltaTclock, bi, bj, myThid )
894     ENDIF
895     #endif
896     DO k=1,Nr
897     GM_TimeAve(k,bi,bj)=GM_TimeAve(k,bi,bj)+deltaTclock
898     ENDDO
899    
900     ENDIF
901     #endif /* ALLOW_TIMEAVE */
902    
903 jmc 1.20 #ifdef ALLOW_DIAGNOSTICS
904     IF ( useDiagnostics ) THEN
905 jmc 1.24 CALL GMREDI_DIAGNOSTICS_FILL(bi,bj,myThid)
906 jmc 1.20 ENDIF
907     #endif /* ALLOW_DIAGNOSTICS */
908    
909 adcroft 1.1 #endif /* ALLOW_GMREDI */
910    
911     RETURN
912     END
913 heimbach 1.2
914    
915     SUBROUTINE GMREDI_CALC_TENSOR_DUMMY(
916 jmc 1.27 I iMin, iMax, jMin, jMax,
917 heimbach 1.2 I sigmaX, sigmaY, sigmaR,
918 jmc 1.27 I bi, bj, myTime, myIter, myThid )
919 heimbach 1.2 C /==========================================================\
920     C | SUBROUTINE GMREDI_CALC_TENSOR |
921     C | o Calculate tensor elements for GM/Redi tensor. |
922     C |==========================================================|
923     C \==========================================================/
924     IMPLICIT NONE
925    
926     C == Global variables ==
927     #include "SIZE.h"
928     #include "EEPARAMS.h"
929     #include "GMREDI.h"
930    
931     C == Routine arguments ==
932     C
933     _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
934     _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
935     _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
936 jmc 1.27 INTEGER iMin,iMax,jMin,jMax
937     INTEGER bi, bj
938     _RL myTime
939     INTEGER myIter
940 heimbach 1.2 INTEGER myThid
941     CEndOfInterface
942    
943 jmc 1.27 #ifdef ALLOW_GMREDI
944    
945 jmc 1.9 INTEGER i, j, k
946 heimbach 1.2
947 jmc 1.9 DO k=1,Nr
948     DO j=1-Oly+1,sNy+Oly-1
949     DO i=1-Olx+1,sNx+Olx-1
950     Kwx(i,j,k,bi,bj) = 0.0
951     Kwy(i,j,k,bi,bj) = 0.0
952     Kwz(i,j,k,bi,bj) = 0.0
953     ENDDO
954 heimbach 1.2 ENDDO
955     ENDDO
956     #endif /* ALLOW_GMREDI */
957    
958 jmc 1.9 RETURN
959     END

  ViewVC Help
Powered by ViewVC 1.1.22