/[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.22 - (hide annotations) (download)
Thu Dec 8 03:29:32 2005 UTC (18 years, 6 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint58b_post, checkpoint57y_post, checkpoint58, checkpoint58f_post, checkpoint58d_post, checkpoint58a_post, checkpoint57z_post, checkpoint57y_pre, checkpoint58e_post, checkpoint58g_post, checkpoint58c_post
Changes since 1.21: +93 -1 lines 
add diagnostics of Redi Off-diagonal fluxes (X & Y directions).
1 jmc 1.22 C $Header: /u/gcmpack/MITgcm/pkg/gmredi/gmredi_calc_tensor.F,v 1.21 2005/10/26 20:53:14 heimbach Exp $
2 heimbach 1.13 C $Name: $
3 adcroft 1.1
4     #include "GMREDI_OPTIONS.h"
5    
6     CStartOfInterface
7     SUBROUTINE GMREDI_CALC_TENSOR(
8 jmc 1.9 I bi, bj, iMin, iMax, jMin, jMax,
9 adcroft 1.1 I sigmaX, sigmaY, sigmaR,
10     I myThid )
11     C /==========================================================\
12     C | SUBROUTINE GMREDI_CALC_TENSOR |
13     C | o Calculate tensor elements for GM/Redi tensor. |
14     C |==========================================================|
15     C \==========================================================/
16     IMPLICIT NONE
17    
18     C == Global variables ==
19     #include "SIZE.h"
20     #include "GRID.h"
21     #include "DYNVARS.h"
22     #include "EEPARAMS.h"
23     #include "PARAMS.h"
24     #include "GMREDI.h"
25 jmc 1.20 #include "GMREDI_TAVE.h"
26 adcroft 1.1
27 heimbach 1.10 #ifdef ALLOW_AUTODIFF_TAMC
28     #include "tamc.h"
29     #include "tamc_keys.h"
30     #endif /* ALLOW_AUTODIFF_TAMC */
31    
32 adcroft 1.1 C == Routine arguments ==
33     C
34     _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
35     _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
36     _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
37 jmc 1.9 INTEGER bi,bj,iMin,iMax,jMin,jMax
38 adcroft 1.1 INTEGER myThid
39     CEndOfInterface
40    
41     #ifdef ALLOW_GMREDI
42    
43     C == Local variables ==
44 jmc 1.15 INTEGER i,j,k,kp1
45 adcroft 1.1 _RL SlopeX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
46     _RL SlopeY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
47 heimbach 1.12 _RL dSigmaDx(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
48     _RL dSigmaDy(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
49 jmc 1.19 _RL dSigmaDr(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
50 jmc 1.8 _RL SlopeSqr(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
51     _RL taperFct(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
52 jmc 1.15 _RL maskp1, Kgm_tmp
53 jmc 1.19 _RL ldd97_LrhoC(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
54     _RL ldd97_LrhoW(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
55     _RL ldd97_LrhoS(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
56     _RL Cspd, LrhoInf, LrhoSup, fCoriLoc
57 adcroft 1.1
58     #ifdef GM_VISBECK_VARIABLE_K
59 heimbach 1.14 _RL deltaH,zero_rs
60     PARAMETER(zero_rs=0.D0)
61 adcroft 1.1 _RL N2,SN
62 heimbach 1.10 _RL Ssq(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
63 adcroft 1.1 #endif
64    
65 jmc 1.22 #ifdef ALLOW_DIAGNOSTICS
66     LOGICAL doDiagRediFlx
67     LOGICAL DIAGNOSTICS_IS_ON
68     EXTERNAL DIAGNOSTICS_IS_ON
69     INTEGER km1
70     _RL dTdz
71     _RL tmp1k(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
72     #endif
73    
74 jmc 1.19 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
75    
76 heimbach 1.10 #ifdef ALLOW_AUTODIFF_TAMC
77     act1 = bi - myBxLo(myThid)
78     max1 = myBxHi(myThid) - myBxLo(myThid) + 1
79     act2 = bj - myByLo(myThid)
80     max2 = myByHi(myThid) - myByLo(myThid) + 1
81     act3 = myThid - 1
82     max3 = nTx*nTy
83     act4 = ikey_dynamics - 1
84 heimbach 1.12 igmkey = (act1 + 1) + act2*max1
85 heimbach 1.10 & + act3*max1*max2
86     & + act4*max1*max2*max3
87     #endif /* ALLOW_AUTODIFF_TAMC */
88    
89 jmc 1.22 #ifdef ALLOW_DIAGNOSTICS
90     doDiagRediFlx = .FALSE.
91     IF ( useDiagnostics ) THEN
92     doDiagRediFlx = DIAGNOSTICS_IS_ON('GM_KuzTz', myThid )
93     doDiagRediFlx = doDiagRediFlx .OR.
94     & DIAGNOSTICS_IS_ON('GM_KvzTz', myThid )
95     ENDIF
96     #endif
97    
98 heimbach 1.12 #ifdef GM_VISBECK_VARIABLE_K
99     DO j=1-Oly,sNy+Oly
100     DO i=1-Olx,sNx+Olx
101     VisbeckK(i,j,bi,bj) = 0. _d 0
102     ENDDO
103     ENDDO
104     #endif
105    
106 jmc 1.19 C-- set ldd97_Lrho (for tapering scheme ldd97):
107     IF (GM_taper_scheme.EQ.'ldd97') THEN
108     Cspd = 2. _d 0
109     LrhoInf = 15. _d 3
110     LrhoSup = 100. _d 3
111     C- Tracer point location (center):
112     DO j=1-Oly,sNy+Oly
113     DO i=1-Olx,sNx+Olx
114     IF (fCori(i,j,bi,bj).NE.0.) THEN
115     ldd97_LrhoC(i,j) = Cspd/ABS(fCori(i,j,bi,bj))
116     ELSE
117     ldd97_LrhoC(i,j) = LrhoSup
118     ENDIF
119     ldd97_LrhoC(i,j) = MAX(LrhoInf,MIN(ldd97_LrhoC(i,j),LrhoSup))
120     ENDDO
121     ENDDO
122     C- U point location (West):
123     DO j=1-Oly,sNy+Oly
124     ldd97_LrhoW(1-Olx,j) = LrhoSup
125     DO i=1-Olx+1,sNx+Olx
126     fCoriLoc = op5*(fCori(i-1,j,bi,bj)+fCori(i,j,bi,bj))
127     IF (fCoriLoc.NE.0.) THEN
128     ldd97_LrhoW(i,j) = Cspd/ABS(fCoriLoc)
129     ELSE
130     ldd97_LrhoW(i,j) = LrhoSup
131     ENDIF
132     ldd97_LrhoW(i,j) = MAX(LrhoInf,MIN(ldd97_LrhoW(i,j),LrhoSup))
133     ENDDO
134     ENDDO
135     C- V point location (South):
136     DO i=1-Olx+1,sNx+Olx
137     ldd97_LrhoS(i,1-Oly) = LrhoSup
138     ENDDO
139     DO j=1-Oly+1,sNy+Oly
140     DO i=1-Olx,sNx+Olx
141     fCoriLoc = op5*(fCori(i,j-1,bi,bj)+fCori(i,j,bi,bj))
142     IF (fCoriLoc.NE.0.) THEN
143     ldd97_LrhoS(i,j) = Cspd/ABS(fCoriLoc)
144     ELSE
145     ldd97_LrhoS(i,j) = LrhoSup
146     ENDIF
147     ldd97_LrhoS(i,j) = MAX(LrhoInf,MIN(ldd97_LrhoS(i,j),LrhoSup))
148     ENDDO
149     ENDDO
150     ELSE
151     C- Just initialize to zero (not use anyway)
152     DO j=1-Oly,sNy+Oly
153     DO i=1-Olx,sNx+Olx
154     ldd97_LrhoC(i,j) = 0. _d 0
155     ldd97_LrhoW(i,j) = 0. _d 0
156     ldd97_LrhoS(i,j) = 0. _d 0
157     ENDDO
158     ENDDO
159     ENDIF
160     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
161    
162 jmc 1.9 DO k=2,Nr
163     C-- 1rst loop on k : compute Tensor Coeff. at W points.
164 adcroft 1.1
165     #ifdef ALLOW_AUTODIFF_TAMC
166 heimbach 1.12 kkey = (igmkey-1)*Nr + k
167 heimbach 1.10 DO j=1-Oly,sNy+Oly
168     DO i=1-Olx,sNx+Olx
169     SlopeX(i,j) = 0. _d 0
170     SlopeY(i,j) = 0. _d 0
171 heimbach 1.12 dSigmaDx(i,j) = 0. _d 0
172     dSigmaDy(i,j) = 0. _d 0
173 jmc 1.19 dSigmaDr(i,j) = 0. _d 0
174 heimbach 1.10 SlopeSqr(i,j) = 0. _d 0
175     taperFct(i,j) = 0. _d 0
176     Kwx(i,j,k,bi,bj) = 0. _d 0
177     Kwy(i,j,k,bi,bj) = 0. _d 0
178     Kwz(i,j,k,bi,bj) = 0. _d 0
179 heimbach 1.12 # ifdef GM_NON_UNITY_DIAGONAL
180     Kux(i,j,k,bi,bj) = 0. _d 0
181     Kvy(i,j,k,bi,bj) = 0. _d 0
182     # endif
183     # ifdef GM_EXTRA_DIAGONAL
184     Kuz(i,j,k,bi,bj) = 0. _d 0
185     Kvz(i,j,k,bi,bj) = 0. _d 0
186     # endif
187     # ifdef GM_BOLUS_ADVEC
188     GM_PsiX(i,j,k,bi,bj) = 0. _d 0
189     GM_PsiY(i,j,k,bi,bj) = 0. _d 0
190     # endif
191 heimbach 1.10 ENDDO
192     ENDDO
193 adcroft 1.1 #endif
194 heimbach 1.10
195 adcroft 1.1 DO j=1-Oly+1,sNy+Oly-1
196     DO i=1-Olx+1,sNx+Olx-1
197     C Gradient of Sigma at rVel points
198 jmc 1.15 dSigmaDx(i,j)=op25*( sigmaX(i+1, j ,k-1) +sigmaX(i,j,k-1)
199 adcroft 1.1 & +sigmaX(i+1, j , k ) +sigmaX(i,j, k ) )
200 jmc 1.15 & *maskC(i,j,k,bi,bj)
201     dSigmaDy(i,j)=op25*( sigmaY( i ,j+1,k-1) +sigmaY(i,j,k-1)
202 adcroft 1.1 & +sigmaY( i ,j+1, k ) +sigmaY(i,j, k ) )
203 jmc 1.15 & *maskC(i,j,k,bi,bj)
204 jmc 1.19 dSigmaDr(i,j)=sigmaR(i,j,k)
205 adcroft 1.1 ENDDO
206     ENDDO
207    
208 heimbach 1.10 #ifdef ALLOW_AUTODIFF_TAMC
209 heimbach 1.12 CADJ STORE dSigmaDx(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
210     CADJ STORE dSigmaDy(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
211 jmc 1.19 CADJ STORE dSigmaDr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
212 heimbach 1.10 #endif /* ALLOW_AUTODIFF_TAMC */
213    
214 adcroft 1.1 C Calculate slopes for use in tensor, taper and/or clip
215     CALL GMREDI_SLOPE_LIMIT(
216 jmc 1.19 O SlopeX, SlopeY,
217 jmc 1.8 O SlopeSqr, taperFct,
218 jmc 1.19 U dSigmaDr,
219     I dSigmaDx, dSigmaDy,
220     I ldd97_LrhoC,rF(k),k,
221 adcroft 1.1 I bi, bj, myThid )
222    
223     DO j=1-Oly+1,sNy+Oly-1
224     DO i=1-Olx+1,sNx+Olx-1
225    
226     C Mask Iso-neutral slopes
227 jmc 1.15 SlopeX(i,j)=SlopeX(i,j)*maskC(i,j,k,bi,bj)
228     SlopeY(i,j)=SlopeY(i,j)*maskC(i,j,k,bi,bj)
229     SlopeSqr(i,j)=SlopeSqr(i,j)*maskC(i,j,k,bi,bj)
230 heimbach 1.10
231     ENDDO
232     ENDDO
233    
234     #ifdef ALLOW_AUTODIFF_TAMC
235 heimbach 1.16 CADJ STORE SlopeX(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
236     CADJ STORE SlopeY(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
237 heimbach 1.14 CADJ STORE SlopeSqr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
238 jmc 1.19 CADJ STORE dSigmaDr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
239     CADJ STORE taperFct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
240 heimbach 1.10 #endif /* ALLOW_AUTODIFF_TAMC */
241    
242     DO j=1-Oly+1,sNy+Oly-1
243     DO i=1-Olx+1,sNx+Olx-1
244 adcroft 1.1
245 jmc 1.9 C Components of Redi/GM tensor
246     Kwx(i,j,k,bi,bj)= SlopeX(i,j)*taperFct(i,j)
247     Kwy(i,j,k,bi,bj)= SlopeY(i,j)*taperFct(i,j)
248 jmc 1.8 Kwz(i,j,k,bi,bj)= SlopeSqr(i,j)*taperFct(i,j)
249 adcroft 1.1
250     #ifdef GM_VISBECK_VARIABLE_K
251 jmc 1.8
252     C- note (jmc) : moved here since only used in VISBECK_VARIABLE_K
253 edhill 1.18 C but do not know if *taperFct (or **2 ?) is necessary
254 heimbach 1.10 Ssq(i,j)=SlopeSqr(i,j)*taperFct(i,j)
255 jmc 1.8
256 adcroft 1.1 C-- Depth average of M^2/N^2 * N
257    
258     C Calculate terms for mean Richardson number
259     C which is used in the "variable K" parameterisaton.
260     C Distance between interface above layer and the integration depth
261     deltaH=abs(GM_Visbeck_depth)-abs(rF(k))
262     C If positive we limit this to the layer thickness
263     deltaH=min(deltaH,drF(k))
264     C If negative then we are below the integration level
265     deltaH=max(deltaH,zero_rs)
266     C Now we convert deltaH to a non-dimensional fraction
267     deltaH=deltaH/GM_Visbeck_depth
268    
269 jmc 1.8 IF (K.eq.2) VisbeckK(i,j,bi,bj)=0.
270 jmc 1.19 IF ( Ssq(i,j).NE.0. .AND. dSigmaDr(i,j).NE.0. ) THEN
271     N2= -Gravity*recip_RhoConst*dSigmaDr(i,j)
272 heimbach 1.10 SN=sqrt(Ssq(i,j)*N2)
273 heimbach 1.3 VisbeckK(i,j,bi,bj)=VisbeckK(i,j,bi,bj)+deltaH
274 adcroft 1.1 & *GM_Visbeck_alpha*GM_Visbeck_length*GM_Visbeck_length*SN
275 jmc 1.8 ENDIF
276 adcroft 1.1
277 jmc 1.9 #endif /* GM_VISBECK_VARIABLE_K */
278    
279     ENDDO
280     ENDDO
281    
282     C-- end 1rst loop on vertical level index k
283     ENDDO
284    
285 adcroft 1.1
286 jmc 1.9 #ifdef GM_VISBECK_VARIABLE_K
287 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
288     CADJ STORE VisbeckK(:,:,bi,bj) = comlev1_bibj, key=igmkey, byte=isbyte
289     #endif
290 jmc 1.9 IF ( GM_Visbeck_alpha.NE.0. ) THEN
291     C- Limit range that KapGM can take
292     DO j=1-Oly+1,sNy+Oly-1
293     DO i=1-Olx+1,sNx+Olx-1
294     VisbeckK(i,j,bi,bj)=
295     & MIN(VisbeckK(i,j,bi,bj),GM_Visbeck_maxval_K)
296     ENDDO
297     ENDDO
298     ENDIF
299 heimbach 1.16 cph( NEW
300     #ifdef ALLOW_AUTODIFF_TAMC
301     CADJ STORE VisbeckK(:,:,bi,bj) = comlev1_bibj, key=igmkey, byte=isbyte
302     #endif
303     cph)
304 adcroft 1.1 #endif /* GM_VISBECK_VARIABLE_K */
305    
306    
307 jmc 1.9 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
308 heimbach 1.10
309 jmc 1.9 C-- 2nd loop on k : compute Tensor Coeff. at U,V levels.
310     DO k=1,Nr
311     kp1 = MIN(Nr,k+1)
312     maskp1 = 1. _d 0
313     IF (k.GE.Nr) maskp1 = 0. _d 0
314    
315 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
316     kkey = (igmkey-1)*Nr + k
317 heimbach 1.16 #if (defined (GM_NON_UNITY_DIAGONAL) || \
318     defined (GM_VISBECK_VARIABLE_K))
319 heimbach 1.14 CADJ STORE Kwx(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
320     CADJ STORE Kwy(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
321     CADJ STORE Kwz(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
322 heimbach 1.12 #endif
323     #endif
324    
325 jmc 1.9 C- express the Tensor in term of Diffusivity (= m**2 / s )
326     DO j=1-Oly+1,sNy+Oly-1
327     DO i=1-Olx+1,sNx+Olx-1
328     Kgm_tmp = GM_isopycK + GM_skewflx*GM_background_K
329     #ifdef GM_VISBECK_VARIABLE_K
330 heimbach 1.16 & + VisbeckK(i,j,bi,bj)*(1. _d 0 + GM_skewflx)
331 jmc 1.9 #endif
332     Kwx(i,j,k,bi,bj)= Kgm_tmp*Kwx(i,j,k,bi,bj)
333     Kwy(i,j,k,bi,bj)= Kgm_tmp*Kwy(i,j,k,bi,bj)
334     Kwz(i,j,k,bi,bj)= ( GM_isopycK
335 adcroft 1.1 #ifdef GM_VISBECK_VARIABLE_K
336 jmc 1.9 & + VisbeckK(i,j,bi,bj)
337 adcroft 1.1 #endif
338 jmc 1.9 & )*Kwz(i,j,k,bi,bj)
339 adcroft 1.1 ENDDO
340     ENDDO
341 adcroft 1.4
342 jmc 1.9 #if ( defined (GM_NON_UNITY_DIAGONAL) || defined (GM_EXTRA_DIAGONAL) )
343 adcroft 1.1
344     C Gradient of Sigma at U points
345     DO j=1-Oly+1,sNy+Oly-1
346     DO i=1-Olx+1,sNx+Olx-1
347 heimbach 1.12 dSigmaDx(i,j)=sigmaX(i,j,k)
348 adcroft 1.1 & *_maskW(i,j,k,bi,bj)
349 heimbach 1.14 dSigmaDy(i,j)=op25*( sigmaY(i-1,j+1,k) +sigmaY(i,j+1,k)
350 heimbach 1.12 & +sigmaY(i-1, j ,k) +sigmaY(i, j ,k) )
351 adcroft 1.1 & *_maskW(i,j,k,bi,bj)
352 jmc 1.19 dSigmaDr(i,j)=op25*( sigmaR(i-1,j, k ) +sigmaR(i,j, k )
353 jmc 1.9 & +maskp1*(sigmaR(i-1,j,kp1) +sigmaR(i,j,kp1)) )
354 jmc 1.15 & *_maskW(i,j,k,bi,bj)
355 adcroft 1.1 ENDDO
356     ENDDO
357    
358 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
359 heimbach 1.17 CADJ STORE SlopeSqr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
360 heimbach 1.12 CADJ STORE dSigmaDx(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
361     CADJ STORE dSigmaDy(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
362 jmc 1.19 CADJ STORE dSigmaDr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
363 heimbach 1.12 #endif /* ALLOW_AUTODIFF_TAMC */
364    
365 adcroft 1.1 C Calculate slopes for use in tensor, taper and/or clip
366     CALL GMREDI_SLOPE_LIMIT(
367 jmc 1.19 O SlopeX, SlopeY,
368 jmc 1.8 O SlopeSqr, taperFct,
369 jmc 1.19 U dSigmaDr,
370     I dSigmaDx, dSigmaDy,
371     I ldd97_LrhoW,rC(k),k,
372 adcroft 1.1 I bi, bj, myThid )
373    
374 heimbach 1.16 cph( NEW
375     #ifdef ALLOW_AUTODIFF_TAMC
376     cph(
377 heimbach 1.17 CADJ STORE SlopeSqr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
378 jmc 1.19 CADJ STORE taperFct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
379 heimbach 1.16 cph)
380     #endif /* ALLOW_AUTODIFF_TAMC */
381     cph)
382    
383 jmc 1.9 #ifdef GM_NON_UNITY_DIAGONAL
384     DO j=1-Oly+1,sNy+Oly-1
385     DO i=1-Olx+1,sNx+Olx-1
386     Kux(i,j,k,bi,bj) =
387     & ( GM_isopycK
388     #ifdef GM_VISBECK_VARIABLE_K
389 heimbach 1.14 & +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))
390 jmc 1.9 #endif
391 heimbach 1.10 & )
392     & *taperFct(i,j)
393     ENDDO
394     ENDDO
395 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
396 heimbach 1.16 # ifdef GM_EXCLUDE_CLIPPING
397 heimbach 1.12 CADJ STORE Kux(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
398     # endif
399     #endif
400 heimbach 1.10 DO j=1-Oly+1,sNy+Oly-1
401     DO i=1-Olx+1,sNx+Olx-1
402 jmc 1.9 Kux(i,j,k,bi,bj) = MAX( Kux(i,j,k,bi,bj), GM_Kmin_horiz )
403     ENDDO
404     ENDDO
405     #endif /* GM_NON_UNITY_DIAGONAL */
406    
407     #ifdef GM_EXTRA_DIAGONAL
408 heimbach 1.12
409     #ifdef ALLOW_AUTODIFF_TAMC
410     CADJ STORE SlopeX(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
411     CADJ STORE taperFct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
412     #endif /* ALLOW_AUTODIFF_TAMC */
413 jmc 1.9 IF (GM_ExtraDiag) THEN
414     DO j=1-Oly+1,sNy+Oly-1
415     DO i=1-Olx+1,sNx+Olx-1
416     Kuz(i,j,k,bi,bj) =
417     & ( GM_isopycK - GM_skewflx*GM_background_K
418     #ifdef GM_VISBECK_VARIABLE_K
419 heimbach 1.14 & +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))*GM_advect
420 jmc 1.9 #endif
421     & )*SlopeX(i,j)*taperFct(i,j)
422     ENDDO
423     ENDDO
424     ENDIF
425     #endif /* GM_EXTRA_DIAGONAL */
426 adcroft 1.1
427 jmc 1.22 #ifdef ALLOW_DIAGNOSTICS
428     IF (doDiagRediFlx) THEN
429     km1 = MAX(k-1,1)
430     DO j=1,sNy
431     DO i=1,sNx+1
432     C store in tmp1k Kuz_Redi
433     tmp1k(i,j) = ( GM_isopycK
434     #ifdef GM_VISBECK_VARIABLE_K
435     & +(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))*0.5 _d 0
436     #endif
437     & )*SlopeX(i,j)*taperFct(i,j)
438     ENDDO
439     ENDDO
440     DO j=1,sNy
441     DO i=1,sNx+1
442     C- Vertical gradients interpolated to U points
443     dTdz = (
444     & +recip_drC(k)*
445     & ( maskC(i-1,j,k,bi,bj)*
446     & (theta(i-1,j,km1,bi,bj)-theta(i-1,j,k,bi,bj))
447     & +maskC( i ,j,k,bi,bj)*
448     & (theta( i ,j,km1,bi,bj)-theta( i ,j,k,bi,bj))
449     & )
450     & +recip_drC(kp1)*
451     & ( maskC(i-1,j,kp1,bi,bj)*
452     & (theta(i-1,j,k,bi,bj)-theta(i-1,j,kp1,bi,bj))
453     & +maskC( i ,j,kp1,bi,bj)*
454     & (theta( i ,j,k,bi,bj)-theta( i ,j,kp1,bi,bj))
455     & ) ) * 0.25 _d 0
456     tmp1k(i,j) = dyG(i,j,bi,bj)*drF(k)*hFacW(i,j,k,bi,bj)
457     & * tmp1k(i,j) * dTdz
458     ENDDO
459     ENDDO
460     CALL DIAGNOSTICS_FILL(tmp1k, 'GM_KuzTz', k,1,2,bi,bj,myThid)
461     ENDIF
462     #endif /* ALLOW_DIAGNOSTICS */
463    
464 adcroft 1.1 C Gradient of Sigma at V points
465     DO j=1-Oly+1,sNy+Oly-1
466     DO i=1-Olx+1,sNx+Olx-1
467 heimbach 1.14 dSigmaDx(i,j)=op25*( sigmaX(i, j ,k) +sigmaX(i+1, j ,k)
468 adcroft 1.1 & +sigmaX(i,j-1,k) +sigmaX(i+1,j-1,k) )
469     & *_maskS(i,j,k,bi,bj)
470 heimbach 1.12 dSigmaDy(i,j)=sigmaY(i,j,k)
471 adcroft 1.1 & *_maskS(i,j,k,bi,bj)
472 jmc 1.19 dSigmaDr(i,j)=op25*( sigmaR(i,j-1, k ) +sigmaR(i,j, k )
473 jmc 1.9 & +maskp1*(sigmaR(i,j-1,kp1) +sigmaR(i,j,kp1)) )
474 jmc 1.15 & *_maskS(i,j,k,bi,bj)
475 adcroft 1.1 ENDDO
476     ENDDO
477    
478 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
479     CADJ STORE dSigmaDx(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
480     CADJ STORE dSigmaDy(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
481 jmc 1.19 CADJ STORE dSigmaDr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
482 heimbach 1.12 #endif /* ALLOW_AUTODIFF_TAMC */
483    
484 adcroft 1.1 C Calculate slopes for use in tensor, taper and/or clip
485     CALL GMREDI_SLOPE_LIMIT(
486 jmc 1.19 O SlopeX, SlopeY,
487 jmc 1.8 O SlopeSqr, taperFct,
488 jmc 1.19 U dSigmaDr,
489     I dSigmaDx, dSigmaDy,
490     I ldd97_LrhoS,rC(k),k,
491 adcroft 1.1 I bi, bj, myThid )
492    
493 heimbach 1.16 cph(
494     #ifdef ALLOW_AUTODIFF_TAMC
495     cph(
496     CADJ STORE taperfct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
497     cph)
498     #endif /* ALLOW_AUTODIFF_TAMC */
499     cph)
500    
501 jmc 1.9 #ifdef GM_NON_UNITY_DIAGONAL
502     DO j=1-Oly+1,sNy+Oly-1
503     DO i=1-Olx+1,sNx+Olx-1
504     Kvy(i,j,k,bi,bj) =
505     & ( GM_isopycK
506     #ifdef GM_VISBECK_VARIABLE_K
507 heimbach 1.14 & +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))
508 jmc 1.9 #endif
509 heimbach 1.10 & )
510     & *taperFct(i,j)
511     ENDDO
512     ENDDO
513 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
514 heimbach 1.16 # ifdef GM_EXCLUDE_CLIPPING
515 heimbach 1.12 CADJ STORE Kvy(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
516     # endif
517     #endif
518 heimbach 1.10 DO j=1-Oly+1,sNy+Oly-1
519     DO i=1-Olx+1,sNx+Olx-1
520 jmc 1.9 Kvy(i,j,k,bi,bj) = MAX( Kvy(i,j,k,bi,bj), GM_Kmin_horiz )
521     ENDDO
522     ENDDO
523     #endif /* GM_NON_UNITY_DIAGONAL */
524    
525     #ifdef GM_EXTRA_DIAGONAL
526 heimbach 1.12
527     #ifdef ALLOW_AUTODIFF_TAMC
528     CADJ STORE SlopeY(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
529     CADJ STORE taperFct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
530     #endif /* ALLOW_AUTODIFF_TAMC */
531 jmc 1.9 IF (GM_ExtraDiag) THEN
532     DO j=1-Oly+1,sNy+Oly-1
533     DO i=1-Olx+1,sNx+Olx-1
534     Kvz(i,j,k,bi,bj) =
535     & ( GM_isopycK - GM_skewflx*GM_background_K
536     #ifdef GM_VISBECK_VARIABLE_K
537 heimbach 1.14 & +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))*GM_advect
538 jmc 1.9 #endif
539     & )*SlopeY(i,j)*taperFct(i,j)
540     ENDDO
541     ENDDO
542     ENDIF
543     #endif /* GM_EXTRA_DIAGONAL */
544    
545 jmc 1.22 #ifdef ALLOW_DIAGNOSTICS
546     IF (doDiagRediFlx) THEN
547     c km1 = MAX(k-1,1)
548     DO j=1,sNy+1
549     DO i=1,sNx
550     C store in tmp1k Kvz_Redi
551     tmp1k(i,j) = ( GM_isopycK
552     #ifdef GM_VISBECK_VARIABLE_K
553     & +(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))*0.5 _d 0
554     #endif
555     & )*SlopeY(i,j)*taperFct(i,j)
556     ENDDO
557     ENDDO
558     DO j=1,sNy+1
559     DO i=1,sNx
560     C- Vertical gradients interpolated to U points
561     dTdz = (
562     & +recip_drC(k)*
563     & ( maskC(i,j-1,k,bi,bj)*
564     & (theta(i,j-1,km1,bi,bj)-theta(i,j-1,k,bi,bj))
565     & +maskC(i, j ,k,bi,bj)*
566     & (theta(i, j ,km1,bi,bj)-theta(i, j ,k,bi,bj))
567     & )
568     & +recip_drC(kp1)*
569     & ( maskC(i,j-1,kp1,bi,bj)*
570     & (theta(i,j-1,k,bi,bj)-theta(i,j-1,kp1,bi,bj))
571     & +maskC(i, j ,kp1,bi,bj)*
572     & (theta(i, j ,k,bi,bj)-theta(i, j ,kp1,bi,bj))
573     & ) ) * 0.25 _d 0
574     tmp1k(i,j) = dxG(i,j,bi,bj)*drF(k)*hFacS(i,j,k,bi,bj)
575     & * tmp1k(i,j) * dTdz
576     ENDDO
577     ENDDO
578     CALL DIAGNOSTICS_FILL(tmp1k, 'GM_KvzTz', k,1,2,bi,bj,myThid)
579     ENDIF
580     #endif /* ALLOW_DIAGNOSTICS */
581    
582 jmc 1.9 #endif /* GM_NON_UNITY_DIAGONAL || GM_EXTRA_DIAGONAL */
583    
584     C-- end 2nd loop on vertical level index k
585     ENDDO
586 adcroft 1.1
587    
588 jmc 1.9 #ifdef GM_BOLUS_ADVEC
589     IF (GM_AdvForm) THEN
590     CALL GMREDI_CALC_PSI_B(
591     I bi, bj, iMin, iMax, jMin, jMax,
592     I sigmaX, sigmaY, sigmaR,
593 jmc 1.19 I ldd97_LrhoW, ldd97_LrhoS,
594 jmc 1.9 I myThid )
595     ENDIF
596     #endif
597 adcroft 1.1
598 jmc 1.19 #ifdef ALLOW_TIMEAVE
599     C-- Time-average
600     IF ( taveFreq.GT.0. ) THEN
601    
602     CALL TIMEAVE_CUMULATE( GM_Kwx_T, Kwx, Nr,
603     & deltaTclock, bi, bj, myThid )
604     CALL TIMEAVE_CUMULATE( GM_Kwy_T, Kwy, Nr,
605     & deltaTclock, bi, bj, myThid )
606     CALL TIMEAVE_CUMULATE( GM_Kwz_T, Kwz, Nr,
607     & deltaTclock, bi, bj, myThid )
608     #ifdef GM_VISBECK_VARIABLE_K
609     IF ( GM_Visbeck_alpha.NE.0. ) THEN
610     CALL TIMEAVE_CUMULATE( Visbeck_K_T, VisbeckK, 1,
611     & deltaTclock, bi, bj, myThid )
612     ENDIF
613     #endif
614     #ifdef GM_BOLUS_ADVEC
615     IF ( GM_AdvForm ) THEN
616     CALL TIMEAVE_CUMULATE( GM_PsiXtave, GM_PsiX, Nr,
617     & deltaTclock, bi, bj, myThid )
618     CALL TIMEAVE_CUMULATE( GM_PsiYtave, GM_PsiY, Nr,
619     & deltaTclock, bi, bj, myThid )
620     ENDIF
621     #endif
622     DO k=1,Nr
623     GM_TimeAve(k,bi,bj)=GM_TimeAve(k,bi,bj)+deltaTclock
624     ENDDO
625    
626     ENDIF
627     #endif /* ALLOW_TIMEAVE */
628    
629 jmc 1.20 #ifdef ALLOW_DIAGNOSTICS
630     IF ( useDiagnostics ) THEN
631 heimbach 1.21 CALL GMREDI_DIAGNOSTICS_DRIVER(bi,bj,myThid)
632 jmc 1.20 ENDIF
633     #endif /* ALLOW_DIAGNOSTICS */
634    
635 adcroft 1.1 #endif /* ALLOW_GMREDI */
636    
637     RETURN
638     END
639 heimbach 1.2
640    
641     SUBROUTINE GMREDI_CALC_TENSOR_DUMMY(
642 jmc 1.9 I bi, bj, iMin, iMax, jMin, jMax,
643 heimbach 1.2 I sigmaX, sigmaY, sigmaR,
644     I myThid )
645     C /==========================================================\
646     C | SUBROUTINE GMREDI_CALC_TENSOR |
647     C | o Calculate tensor elements for GM/Redi tensor. |
648     C |==========================================================|
649     C \==========================================================/
650     IMPLICIT NONE
651    
652     C == Global variables ==
653     #include "SIZE.h"
654     #include "EEPARAMS.h"
655     #include "GMREDI.h"
656    
657     C == Routine arguments ==
658     C
659     _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
660     _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
661     _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
662 jmc 1.9 INTEGER bi,bj,iMin,iMax,jMin,jMax
663 heimbach 1.2 INTEGER myThid
664     CEndOfInterface
665    
666 jmc 1.9 INTEGER i, j, k
667 heimbach 1.2
668     #ifdef ALLOW_GMREDI
669    
670 jmc 1.9 DO k=1,Nr
671     DO j=1-Oly+1,sNy+Oly-1
672     DO i=1-Olx+1,sNx+Olx-1
673     Kwx(i,j,k,bi,bj) = 0.0
674     Kwy(i,j,k,bi,bj) = 0.0
675     Kwz(i,j,k,bi,bj) = 0.0
676     ENDDO
677 heimbach 1.2 ENDDO
678     ENDDO
679     #endif /* ALLOW_GMREDI */
680    
681 jmc 1.9 RETURN
682     END
  ViewVC Help
Powered by ViewVC 1.1.22