/[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.21 - (hide annotations) (download)
Wed Oct 26 20:53:14 2005 UTC (18 years, 7 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint57w_post, checkpoint57x_post
Changes since 1.20: +2 -26 lines
Separate diagnostics_fill from active routine gmredi_calc_tensor.

1 heimbach 1.21 C $Header: /u/gcmpack/MITgcm/pkg/gmredi/gmredi_calc_tensor.F,v 1.20 2005/01/04 00:19:38 jmc 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.19 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
66    
67 heimbach 1.10 #ifdef ALLOW_AUTODIFF_TAMC
68     act1 = bi - myBxLo(myThid)
69     max1 = myBxHi(myThid) - myBxLo(myThid) + 1
70     act2 = bj - myByLo(myThid)
71     max2 = myByHi(myThid) - myByLo(myThid) + 1
72     act3 = myThid - 1
73     max3 = nTx*nTy
74     act4 = ikey_dynamics - 1
75 heimbach 1.12 igmkey = (act1 + 1) + act2*max1
76 heimbach 1.10 & + act3*max1*max2
77     & + act4*max1*max2*max3
78     #endif /* ALLOW_AUTODIFF_TAMC */
79    
80 heimbach 1.12 #ifdef GM_VISBECK_VARIABLE_K
81     DO j=1-Oly,sNy+Oly
82     DO i=1-Olx,sNx+Olx
83     VisbeckK(i,j,bi,bj) = 0. _d 0
84     ENDDO
85     ENDDO
86     #endif
87    
88 jmc 1.19 C-- set ldd97_Lrho (for tapering scheme ldd97):
89     IF (GM_taper_scheme.EQ.'ldd97') THEN
90     Cspd = 2. _d 0
91     LrhoInf = 15. _d 3
92     LrhoSup = 100. _d 3
93     C- Tracer point location (center):
94     DO j=1-Oly,sNy+Oly
95     DO i=1-Olx,sNx+Olx
96     IF (fCori(i,j,bi,bj).NE.0.) THEN
97     ldd97_LrhoC(i,j) = Cspd/ABS(fCori(i,j,bi,bj))
98     ELSE
99     ldd97_LrhoC(i,j) = LrhoSup
100     ENDIF
101     ldd97_LrhoC(i,j) = MAX(LrhoInf,MIN(ldd97_LrhoC(i,j),LrhoSup))
102     ENDDO
103     ENDDO
104     C- U point location (West):
105     DO j=1-Oly,sNy+Oly
106     ldd97_LrhoW(1-Olx,j) = LrhoSup
107     DO i=1-Olx+1,sNx+Olx
108     fCoriLoc = op5*(fCori(i-1,j,bi,bj)+fCori(i,j,bi,bj))
109     IF (fCoriLoc.NE.0.) THEN
110     ldd97_LrhoW(i,j) = Cspd/ABS(fCoriLoc)
111     ELSE
112     ldd97_LrhoW(i,j) = LrhoSup
113     ENDIF
114     ldd97_LrhoW(i,j) = MAX(LrhoInf,MIN(ldd97_LrhoW(i,j),LrhoSup))
115     ENDDO
116     ENDDO
117     C- V point location (South):
118     DO i=1-Olx+1,sNx+Olx
119     ldd97_LrhoS(i,1-Oly) = LrhoSup
120     ENDDO
121     DO j=1-Oly+1,sNy+Oly
122     DO i=1-Olx,sNx+Olx
123     fCoriLoc = op5*(fCori(i,j-1,bi,bj)+fCori(i,j,bi,bj))
124     IF (fCoriLoc.NE.0.) THEN
125     ldd97_LrhoS(i,j) = Cspd/ABS(fCoriLoc)
126     ELSE
127     ldd97_LrhoS(i,j) = LrhoSup
128     ENDIF
129     ldd97_LrhoS(i,j) = MAX(LrhoInf,MIN(ldd97_LrhoS(i,j),LrhoSup))
130     ENDDO
131     ENDDO
132     ELSE
133     C- Just initialize to zero (not use anyway)
134     DO j=1-Oly,sNy+Oly
135     DO i=1-Olx,sNx+Olx
136     ldd97_LrhoC(i,j) = 0. _d 0
137     ldd97_LrhoW(i,j) = 0. _d 0
138     ldd97_LrhoS(i,j) = 0. _d 0
139     ENDDO
140     ENDDO
141     ENDIF
142     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
143    
144 jmc 1.9 DO k=2,Nr
145     C-- 1rst loop on k : compute Tensor Coeff. at W points.
146 adcroft 1.1
147     #ifdef ALLOW_AUTODIFF_TAMC
148 heimbach 1.12 kkey = (igmkey-1)*Nr + k
149 heimbach 1.10 DO j=1-Oly,sNy+Oly
150     DO i=1-Olx,sNx+Olx
151     SlopeX(i,j) = 0. _d 0
152     SlopeY(i,j) = 0. _d 0
153 heimbach 1.12 dSigmaDx(i,j) = 0. _d 0
154     dSigmaDy(i,j) = 0. _d 0
155 jmc 1.19 dSigmaDr(i,j) = 0. _d 0
156 heimbach 1.10 SlopeSqr(i,j) = 0. _d 0
157     taperFct(i,j) = 0. _d 0
158     Kwx(i,j,k,bi,bj) = 0. _d 0
159     Kwy(i,j,k,bi,bj) = 0. _d 0
160     Kwz(i,j,k,bi,bj) = 0. _d 0
161 heimbach 1.12 # ifdef GM_NON_UNITY_DIAGONAL
162     Kux(i,j,k,bi,bj) = 0. _d 0
163     Kvy(i,j,k,bi,bj) = 0. _d 0
164     # endif
165     # ifdef GM_EXTRA_DIAGONAL
166     Kuz(i,j,k,bi,bj) = 0. _d 0
167     Kvz(i,j,k,bi,bj) = 0. _d 0
168     # endif
169     # ifdef GM_BOLUS_ADVEC
170     GM_PsiX(i,j,k,bi,bj) = 0. _d 0
171     GM_PsiY(i,j,k,bi,bj) = 0. _d 0
172     # endif
173 heimbach 1.10 ENDDO
174     ENDDO
175 adcroft 1.1 #endif
176 heimbach 1.10
177 adcroft 1.1 DO j=1-Oly+1,sNy+Oly-1
178     DO i=1-Olx+1,sNx+Olx-1
179     C Gradient of Sigma at rVel points
180 jmc 1.15 dSigmaDx(i,j)=op25*( sigmaX(i+1, j ,k-1) +sigmaX(i,j,k-1)
181 adcroft 1.1 & +sigmaX(i+1, j , k ) +sigmaX(i,j, k ) )
182 jmc 1.15 & *maskC(i,j,k,bi,bj)
183     dSigmaDy(i,j)=op25*( sigmaY( i ,j+1,k-1) +sigmaY(i,j,k-1)
184 adcroft 1.1 & +sigmaY( i ,j+1, k ) +sigmaY(i,j, k ) )
185 jmc 1.15 & *maskC(i,j,k,bi,bj)
186 jmc 1.19 dSigmaDr(i,j)=sigmaR(i,j,k)
187 adcroft 1.1 ENDDO
188     ENDDO
189    
190 heimbach 1.10 #ifdef ALLOW_AUTODIFF_TAMC
191 heimbach 1.12 CADJ STORE dSigmaDx(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
192     CADJ STORE dSigmaDy(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
193 jmc 1.19 CADJ STORE dSigmaDr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
194 heimbach 1.10 #endif /* ALLOW_AUTODIFF_TAMC */
195    
196 adcroft 1.1 C Calculate slopes for use in tensor, taper and/or clip
197     CALL GMREDI_SLOPE_LIMIT(
198 jmc 1.19 O SlopeX, SlopeY,
199 jmc 1.8 O SlopeSqr, taperFct,
200 jmc 1.19 U dSigmaDr,
201     I dSigmaDx, dSigmaDy,
202     I ldd97_LrhoC,rF(k),k,
203 adcroft 1.1 I bi, bj, myThid )
204    
205     DO j=1-Oly+1,sNy+Oly-1
206     DO i=1-Olx+1,sNx+Olx-1
207    
208     C Mask Iso-neutral slopes
209 jmc 1.15 SlopeX(i,j)=SlopeX(i,j)*maskC(i,j,k,bi,bj)
210     SlopeY(i,j)=SlopeY(i,j)*maskC(i,j,k,bi,bj)
211     SlopeSqr(i,j)=SlopeSqr(i,j)*maskC(i,j,k,bi,bj)
212 heimbach 1.10
213     ENDDO
214     ENDDO
215    
216     #ifdef ALLOW_AUTODIFF_TAMC
217 heimbach 1.16 CADJ STORE SlopeX(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
218     CADJ STORE SlopeY(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
219 heimbach 1.14 CADJ STORE SlopeSqr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
220 jmc 1.19 CADJ STORE dSigmaDr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
221     CADJ STORE taperFct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
222 heimbach 1.10 #endif /* ALLOW_AUTODIFF_TAMC */
223    
224     DO j=1-Oly+1,sNy+Oly-1
225     DO i=1-Olx+1,sNx+Olx-1
226 adcroft 1.1
227 jmc 1.9 C Components of Redi/GM tensor
228     Kwx(i,j,k,bi,bj)= SlopeX(i,j)*taperFct(i,j)
229     Kwy(i,j,k,bi,bj)= SlopeY(i,j)*taperFct(i,j)
230 jmc 1.8 Kwz(i,j,k,bi,bj)= SlopeSqr(i,j)*taperFct(i,j)
231 adcroft 1.1
232     #ifdef GM_VISBECK_VARIABLE_K
233 jmc 1.8
234     C- note (jmc) : moved here since only used in VISBECK_VARIABLE_K
235 edhill 1.18 C but do not know if *taperFct (or **2 ?) is necessary
236 heimbach 1.10 Ssq(i,j)=SlopeSqr(i,j)*taperFct(i,j)
237 jmc 1.8
238 adcroft 1.1 C-- Depth average of M^2/N^2 * N
239    
240     C Calculate terms for mean Richardson number
241     C which is used in the "variable K" parameterisaton.
242     C Distance between interface above layer and the integration depth
243     deltaH=abs(GM_Visbeck_depth)-abs(rF(k))
244     C If positive we limit this to the layer thickness
245     deltaH=min(deltaH,drF(k))
246     C If negative then we are below the integration level
247     deltaH=max(deltaH,zero_rs)
248     C Now we convert deltaH to a non-dimensional fraction
249     deltaH=deltaH/GM_Visbeck_depth
250    
251 jmc 1.8 IF (K.eq.2) VisbeckK(i,j,bi,bj)=0.
252 jmc 1.19 IF ( Ssq(i,j).NE.0. .AND. dSigmaDr(i,j).NE.0. ) THEN
253     N2= -Gravity*recip_RhoConst*dSigmaDr(i,j)
254 heimbach 1.10 SN=sqrt(Ssq(i,j)*N2)
255 heimbach 1.3 VisbeckK(i,j,bi,bj)=VisbeckK(i,j,bi,bj)+deltaH
256 adcroft 1.1 & *GM_Visbeck_alpha*GM_Visbeck_length*GM_Visbeck_length*SN
257 jmc 1.8 ENDIF
258 adcroft 1.1
259 jmc 1.9 #endif /* GM_VISBECK_VARIABLE_K */
260    
261     ENDDO
262     ENDDO
263    
264     C-- end 1rst loop on vertical level index k
265     ENDDO
266    
267 adcroft 1.1
268 jmc 1.9 #ifdef GM_VISBECK_VARIABLE_K
269 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
270     CADJ STORE VisbeckK(:,:,bi,bj) = comlev1_bibj, key=igmkey, byte=isbyte
271     #endif
272 jmc 1.9 IF ( GM_Visbeck_alpha.NE.0. ) THEN
273     C- Limit range that KapGM can take
274     DO j=1-Oly+1,sNy+Oly-1
275     DO i=1-Olx+1,sNx+Olx-1
276     VisbeckK(i,j,bi,bj)=
277     & MIN(VisbeckK(i,j,bi,bj),GM_Visbeck_maxval_K)
278     ENDDO
279     ENDDO
280     ENDIF
281 heimbach 1.16 cph( NEW
282     #ifdef ALLOW_AUTODIFF_TAMC
283     CADJ STORE VisbeckK(:,:,bi,bj) = comlev1_bibj, key=igmkey, byte=isbyte
284     #endif
285     cph)
286 adcroft 1.1 #endif /* GM_VISBECK_VARIABLE_K */
287    
288    
289 jmc 1.9 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
290 heimbach 1.10
291 jmc 1.9 C-- 2nd loop on k : compute Tensor Coeff. at U,V levels.
292     DO k=1,Nr
293     kp1 = MIN(Nr,k+1)
294     maskp1 = 1. _d 0
295     IF (k.GE.Nr) maskp1 = 0. _d 0
296    
297 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
298     kkey = (igmkey-1)*Nr + k
299 heimbach 1.16 #if (defined (GM_NON_UNITY_DIAGONAL) || \
300     defined (GM_VISBECK_VARIABLE_K))
301 heimbach 1.14 CADJ STORE Kwx(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
302     CADJ STORE Kwy(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
303     CADJ STORE Kwz(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
304 heimbach 1.12 #endif
305     #endif
306    
307 jmc 1.9 C- express the Tensor in term of Diffusivity (= m**2 / s )
308     DO j=1-Oly+1,sNy+Oly-1
309     DO i=1-Olx+1,sNx+Olx-1
310     Kgm_tmp = GM_isopycK + GM_skewflx*GM_background_K
311     #ifdef GM_VISBECK_VARIABLE_K
312 heimbach 1.16 & + VisbeckK(i,j,bi,bj)*(1. _d 0 + GM_skewflx)
313 jmc 1.9 #endif
314     Kwx(i,j,k,bi,bj)= Kgm_tmp*Kwx(i,j,k,bi,bj)
315     Kwy(i,j,k,bi,bj)= Kgm_tmp*Kwy(i,j,k,bi,bj)
316     Kwz(i,j,k,bi,bj)= ( GM_isopycK
317 adcroft 1.1 #ifdef GM_VISBECK_VARIABLE_K
318 jmc 1.9 & + VisbeckK(i,j,bi,bj)
319 adcroft 1.1 #endif
320 jmc 1.9 & )*Kwz(i,j,k,bi,bj)
321 adcroft 1.1 ENDDO
322     ENDDO
323 adcroft 1.4
324 jmc 1.9 #if ( defined (GM_NON_UNITY_DIAGONAL) || defined (GM_EXTRA_DIAGONAL) )
325 adcroft 1.1
326     C Gradient of Sigma at U points
327     DO j=1-Oly+1,sNy+Oly-1
328     DO i=1-Olx+1,sNx+Olx-1
329 heimbach 1.12 dSigmaDx(i,j)=sigmaX(i,j,k)
330 adcroft 1.1 & *_maskW(i,j,k,bi,bj)
331 heimbach 1.14 dSigmaDy(i,j)=op25*( sigmaY(i-1,j+1,k) +sigmaY(i,j+1,k)
332 heimbach 1.12 & +sigmaY(i-1, j ,k) +sigmaY(i, j ,k) )
333 adcroft 1.1 & *_maskW(i,j,k,bi,bj)
334 jmc 1.19 dSigmaDr(i,j)=op25*( sigmaR(i-1,j, k ) +sigmaR(i,j, k )
335 jmc 1.9 & +maskp1*(sigmaR(i-1,j,kp1) +sigmaR(i,j,kp1)) )
336 jmc 1.15 & *_maskW(i,j,k,bi,bj)
337 adcroft 1.1 ENDDO
338     ENDDO
339    
340 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
341 heimbach 1.17 CADJ STORE SlopeSqr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
342 heimbach 1.12 CADJ STORE dSigmaDx(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
343     CADJ STORE dSigmaDy(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
344 jmc 1.19 CADJ STORE dSigmaDr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
345 heimbach 1.12 #endif /* ALLOW_AUTODIFF_TAMC */
346    
347 adcroft 1.1 C Calculate slopes for use in tensor, taper and/or clip
348     CALL GMREDI_SLOPE_LIMIT(
349 jmc 1.19 O SlopeX, SlopeY,
350 jmc 1.8 O SlopeSqr, taperFct,
351 jmc 1.19 U dSigmaDr,
352     I dSigmaDx, dSigmaDy,
353     I ldd97_LrhoW,rC(k),k,
354 adcroft 1.1 I bi, bj, myThid )
355    
356 heimbach 1.16 cph( NEW
357     #ifdef ALLOW_AUTODIFF_TAMC
358     cph(
359 heimbach 1.17 CADJ STORE SlopeSqr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
360 jmc 1.19 CADJ STORE taperFct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
361 heimbach 1.16 cph)
362     #endif /* ALLOW_AUTODIFF_TAMC */
363     cph)
364    
365 jmc 1.9 #ifdef GM_NON_UNITY_DIAGONAL
366     DO j=1-Oly+1,sNy+Oly-1
367     DO i=1-Olx+1,sNx+Olx-1
368     Kux(i,j,k,bi,bj) =
369     & ( GM_isopycK
370     #ifdef GM_VISBECK_VARIABLE_K
371 heimbach 1.14 & +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))
372 jmc 1.9 #endif
373 heimbach 1.10 & )
374     & *taperFct(i,j)
375     ENDDO
376     ENDDO
377 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
378 heimbach 1.16 # ifdef GM_EXCLUDE_CLIPPING
379 heimbach 1.12 CADJ STORE Kux(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
380     # endif
381     #endif
382 heimbach 1.10 DO j=1-Oly+1,sNy+Oly-1
383     DO i=1-Olx+1,sNx+Olx-1
384 jmc 1.9 Kux(i,j,k,bi,bj) = MAX( Kux(i,j,k,bi,bj), GM_Kmin_horiz )
385     ENDDO
386     ENDDO
387     #endif /* GM_NON_UNITY_DIAGONAL */
388    
389     #ifdef GM_EXTRA_DIAGONAL
390 heimbach 1.12
391     #ifdef ALLOW_AUTODIFF_TAMC
392     CADJ STORE SlopeX(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
393     CADJ STORE taperFct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
394     #endif /* ALLOW_AUTODIFF_TAMC */
395 jmc 1.9 IF (GM_ExtraDiag) THEN
396     DO j=1-Oly+1,sNy+Oly-1
397     DO i=1-Olx+1,sNx+Olx-1
398     Kuz(i,j,k,bi,bj) =
399     & ( GM_isopycK - GM_skewflx*GM_background_K
400     #ifdef GM_VISBECK_VARIABLE_K
401 heimbach 1.14 & +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i-1,j,bi,bj))*GM_advect
402 jmc 1.9 #endif
403     & )*SlopeX(i,j)*taperFct(i,j)
404     ENDDO
405     ENDDO
406     ENDIF
407     #endif /* GM_EXTRA_DIAGONAL */
408 adcroft 1.1
409     C Gradient of Sigma at V points
410     DO j=1-Oly+1,sNy+Oly-1
411     DO i=1-Olx+1,sNx+Olx-1
412 heimbach 1.14 dSigmaDx(i,j)=op25*( sigmaX(i, j ,k) +sigmaX(i+1, j ,k)
413 adcroft 1.1 & +sigmaX(i,j-1,k) +sigmaX(i+1,j-1,k) )
414     & *_maskS(i,j,k,bi,bj)
415 heimbach 1.12 dSigmaDy(i,j)=sigmaY(i,j,k)
416 adcroft 1.1 & *_maskS(i,j,k,bi,bj)
417 jmc 1.19 dSigmaDr(i,j)=op25*( sigmaR(i,j-1, k ) +sigmaR(i,j, k )
418 jmc 1.9 & +maskp1*(sigmaR(i,j-1,kp1) +sigmaR(i,j,kp1)) )
419 jmc 1.15 & *_maskS(i,j,k,bi,bj)
420 adcroft 1.1 ENDDO
421     ENDDO
422    
423 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
424     CADJ STORE dSigmaDx(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
425     CADJ STORE dSigmaDy(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
426 jmc 1.19 CADJ STORE dSigmaDr(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
427 heimbach 1.12 #endif /* ALLOW_AUTODIFF_TAMC */
428    
429 adcroft 1.1 C Calculate slopes for use in tensor, taper and/or clip
430     CALL GMREDI_SLOPE_LIMIT(
431 jmc 1.19 O SlopeX, SlopeY,
432 jmc 1.8 O SlopeSqr, taperFct,
433 jmc 1.19 U dSigmaDr,
434     I dSigmaDx, dSigmaDy,
435     I ldd97_LrhoS,rC(k),k,
436 adcroft 1.1 I bi, bj, myThid )
437    
438 heimbach 1.16 cph(
439     #ifdef ALLOW_AUTODIFF_TAMC
440     cph(
441     CADJ STORE taperfct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
442     cph)
443     #endif /* ALLOW_AUTODIFF_TAMC */
444     cph)
445    
446 jmc 1.9 #ifdef GM_NON_UNITY_DIAGONAL
447     DO j=1-Oly+1,sNy+Oly-1
448     DO i=1-Olx+1,sNx+Olx-1
449     Kvy(i,j,k,bi,bj) =
450     & ( GM_isopycK
451     #ifdef GM_VISBECK_VARIABLE_K
452 heimbach 1.14 & +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))
453 jmc 1.9 #endif
454 heimbach 1.10 & )
455     & *taperFct(i,j)
456     ENDDO
457     ENDDO
458 heimbach 1.12 #ifdef ALLOW_AUTODIFF_TAMC
459 heimbach 1.16 # ifdef GM_EXCLUDE_CLIPPING
460 heimbach 1.12 CADJ STORE Kvy(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
461     # endif
462     #endif
463 heimbach 1.10 DO j=1-Oly+1,sNy+Oly-1
464     DO i=1-Olx+1,sNx+Olx-1
465 jmc 1.9 Kvy(i,j,k,bi,bj) = MAX( Kvy(i,j,k,bi,bj), GM_Kmin_horiz )
466     ENDDO
467     ENDDO
468     #endif /* GM_NON_UNITY_DIAGONAL */
469    
470     #ifdef GM_EXTRA_DIAGONAL
471 heimbach 1.12
472     #ifdef ALLOW_AUTODIFF_TAMC
473     CADJ STORE SlopeY(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
474     CADJ STORE taperFct(:,:) = comlev1_bibj_k, key=kkey, byte=isbyte
475     #endif /* ALLOW_AUTODIFF_TAMC */
476 jmc 1.9 IF (GM_ExtraDiag) THEN
477     DO j=1-Oly+1,sNy+Oly-1
478     DO i=1-Olx+1,sNx+Olx-1
479     Kvz(i,j,k,bi,bj) =
480     & ( GM_isopycK - GM_skewflx*GM_background_K
481     #ifdef GM_VISBECK_VARIABLE_K
482 heimbach 1.14 & +op5*(VisbeckK(i,j,bi,bj)+VisbeckK(i,j-1,bi,bj))*GM_advect
483 jmc 1.9 #endif
484     & )*SlopeY(i,j)*taperFct(i,j)
485     ENDDO
486     ENDDO
487     ENDIF
488     #endif /* GM_EXTRA_DIAGONAL */
489    
490     #endif /* GM_NON_UNITY_DIAGONAL || GM_EXTRA_DIAGONAL */
491    
492     C-- end 2nd loop on vertical level index k
493     ENDDO
494 adcroft 1.1
495    
496 jmc 1.9 #ifdef GM_BOLUS_ADVEC
497     IF (GM_AdvForm) THEN
498     CALL GMREDI_CALC_PSI_B(
499     I bi, bj, iMin, iMax, jMin, jMax,
500     I sigmaX, sigmaY, sigmaR,
501 jmc 1.19 I ldd97_LrhoW, ldd97_LrhoS,
502 jmc 1.9 I myThid )
503     ENDIF
504     #endif
505 adcroft 1.1
506 jmc 1.19 #ifdef ALLOW_TIMEAVE
507     C-- Time-average
508     IF ( taveFreq.GT.0. ) THEN
509    
510     CALL TIMEAVE_CUMULATE( GM_Kwx_T, Kwx, Nr,
511     & deltaTclock, bi, bj, myThid )
512     CALL TIMEAVE_CUMULATE( GM_Kwy_T, Kwy, Nr,
513     & deltaTclock, bi, bj, myThid )
514     CALL TIMEAVE_CUMULATE( GM_Kwz_T, Kwz, Nr,
515     & deltaTclock, bi, bj, myThid )
516     #ifdef GM_VISBECK_VARIABLE_K
517     IF ( GM_Visbeck_alpha.NE.0. ) THEN
518     CALL TIMEAVE_CUMULATE( Visbeck_K_T, VisbeckK, 1,
519     & deltaTclock, bi, bj, myThid )
520     ENDIF
521     #endif
522     #ifdef GM_BOLUS_ADVEC
523     IF ( GM_AdvForm ) THEN
524     CALL TIMEAVE_CUMULATE( GM_PsiXtave, GM_PsiX, Nr,
525     & deltaTclock, bi, bj, myThid )
526     CALL TIMEAVE_CUMULATE( GM_PsiYtave, GM_PsiY, Nr,
527     & deltaTclock, bi, bj, myThid )
528     ENDIF
529     #endif
530     DO k=1,Nr
531     GM_TimeAve(k,bi,bj)=GM_TimeAve(k,bi,bj)+deltaTclock
532     ENDDO
533    
534     ENDIF
535     #endif /* ALLOW_TIMEAVE */
536    
537 jmc 1.20 #ifdef ALLOW_DIAGNOSTICS
538     IF ( useDiagnostics ) THEN
539 heimbach 1.21 CALL GMREDI_DIAGNOSTICS_DRIVER(bi,bj,myThid)
540 jmc 1.20 ENDIF
541     #endif /* ALLOW_DIAGNOSTICS */
542    
543 adcroft 1.1 #endif /* ALLOW_GMREDI */
544    
545     RETURN
546     END
547 heimbach 1.2
548    
549     SUBROUTINE GMREDI_CALC_TENSOR_DUMMY(
550 jmc 1.9 I bi, bj, iMin, iMax, jMin, jMax,
551 heimbach 1.2 I sigmaX, sigmaY, sigmaR,
552     I myThid )
553     C /==========================================================\
554     C | SUBROUTINE GMREDI_CALC_TENSOR |
555     C | o Calculate tensor elements for GM/Redi tensor. |
556     C |==========================================================|
557     C \==========================================================/
558     IMPLICIT NONE
559    
560     C == Global variables ==
561     #include "SIZE.h"
562     #include "EEPARAMS.h"
563     #include "GMREDI.h"
564    
565     C == Routine arguments ==
566     C
567     _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
568     _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
569     _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
570 jmc 1.9 INTEGER bi,bj,iMin,iMax,jMin,jMax
571 heimbach 1.2 INTEGER myThid
572     CEndOfInterface
573    
574 jmc 1.9 INTEGER i, j, k
575 heimbach 1.2
576     #ifdef ALLOW_GMREDI
577    
578 jmc 1.9 DO k=1,Nr
579     DO j=1-Oly+1,sNy+Oly-1
580     DO i=1-Olx+1,sNx+Olx-1
581     Kwx(i,j,k,bi,bj) = 0.0
582     Kwy(i,j,k,bi,bj) = 0.0
583     Kwz(i,j,k,bi,bj) = 0.0
584     ENDDO
585 heimbach 1.2 ENDDO
586     ENDDO
587     #endif /* ALLOW_GMREDI */
588    
589 jmc 1.9 RETURN
590     END

  ViewVC Help
Powered by ViewVC 1.1.22