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

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

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


Revision 1.8 - (show annotations) (download)
Tue Dec 4 15:09:34 2001 UTC (22 years, 5 months ago) by jmc
Branch: MAIN
Changes since 1.7: +28 -28 lines
Modifications related to GM-Redi tapering options:
1) avoid changing the slope sign in negatively stratified region
2) tapered GM-Redi tensor was not truly adiabatic: fixed.
3) add one (simple) tapering function.

1 C $Header: /u/gcmpack/models/MITgcmUV/pkg/gmredi/gmredi_calc_tensor.F,v 1.7 2001/08/21 15:27:19 heimbach Exp $
2 C $Name: $
3
4 #include "GMREDI_OPTIONS.h"
5
6 CStartOfInterface
7 SUBROUTINE GMREDI_CALC_TENSOR(
8 I bi, bj, iMin, iMax, jMin, jMax, K,
9 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 #include "GMREDI_DIAGS.h"
26
27 C == Routine arguments ==
28 C
29 _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
30 _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
31 _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
32 INTEGER bi,bj,iMin,iMax,jMin,jMax,K
33 INTEGER myThid
34 CEndOfInterface
35
36 #ifdef ALLOW_GMREDI
37
38 C == Local variables ==
39 INTEGER i,j,km1,kp1
40 _RL SlopeX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
41 _RL SlopeY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
42 _RL dSigmaDrReal(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
43 _RL SlopeSqr(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
44 _RL taperFct(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
45
46 #ifdef GM_VISBECK_VARIABLE_K
47 _RS deltaH,zero_rs
48 PARAMETER(zero_rs=0.)
49 _RL N2,SN
50 _RL Ssq
51 #endif
52
53
54 km1=max(1,K-1)
55 kp1=min(Nr,K)
56
57
58 #ifdef ALLOW_AUTODIFF_TAMC
59 !HPF$ INDEPENDENT
60 #endif
61 DO j=1-Oly+1,sNy+Oly-1
62 #ifdef ALLOW_AUTODIFF_TAMC
63 !HPF$ INDEPENDENT
64 #endif
65 DO i=1-Olx+1,sNx+Olx-1
66
67 C Gradient of Sigma at rVel points
68 SlopeX(i,j)=0.25*( sigmaX(i+1, j ,km1) +sigmaX(i,j,km1)
69 & +sigmaX(i+1, j , k ) +sigmaX(i,j, k ) )
70 & *maskC(i,j,k,bi,bj)
71 SlopeY(i,j)=0.25*( sigmaY( i ,j+1,km1) +sigmaY(i,j,km1)
72 & +sigmaY( i ,j+1, k ) +sigmaY(i,j, k ) )
73 & *maskC(i,j,k,bi,bj)
74 dSigmaDrReal(i,j)=sigmaR(i,j,k)
75
76 ENDDO
77 ENDDO
78
79 C Calculate slopes for use in tensor, taper and/or clip
80 CALL GMREDI_SLOPE_LIMIT(
81 I dSigmadRReal,
82 I rF(K),
83 U SlopeX, SlopeY,
84 O SlopeSqr, taperFct,
85 I bi, bj, myThid )
86
87 DO j=1-Oly+1,sNy+Oly-1
88 DO i=1-Olx+1,sNx+Olx-1
89
90 C Mask Iso-neutral slopes
91 SlopeX(i,j)=SlopeX(i,j)*maskC(i,j,k,bi,bj)
92 SlopeY(i,j)=SlopeY(i,j)*maskC(i,j,k,bi,bj)
93 SlopeSqr(i,j)=SlopeSqr(i,j)*maskC(i,j,k,bi,bj)
94 c Ssq=SlopeX(i,j)*SlopeX(i,j)+SlopeY(i,j)*SlopeY(i,j)
95
96 C Components of Redi/GM tensor
97 Kwx(i,j,k,bi,bj)=2.*SlopeX(i,j)*taperFct(i,j)
98 Kwy(i,j,k,bi,bj)=2.*SlopeY(i,j)*taperFct(i,j)
99 Kwz(i,j,k,bi,bj)= SlopeSqr(i,j)*taperFct(i,j)
100
101 #ifdef GM_VISBECK_VARIABLE_K
102
103 C- note (jmc) : moved here since only used in VISBECK_VARIABLE_K
104 C but don't know if *taperFct (or **2 ?) is necessary
105 Ssq=SlopeSqr(i,j)*taperFct(i,j)
106
107 C-- Depth average of M^2/N^2 * N
108
109 C Calculate terms for mean Richardson number
110 C which is used in the "variable K" parameterisaton.
111 C Distance between interface above layer and the integration depth
112 deltaH=abs(GM_Visbeck_depth)-abs(rF(k))
113 C If positive we limit this to the layer thickness
114 deltaH=min(deltaH,drF(k))
115 C If negative then we are below the integration level
116 deltaH=max(deltaH,zero_rs)
117 C Now we convert deltaH to a non-dimensional fraction
118 deltaH=deltaH/GM_Visbeck_depth
119
120 IF (K.eq.2) VisbeckK(i,j,bi,bj)=0.
121 IF (Ssq.NE.0.) THEN
122 N2= -Gravity*recip_Rhonil*dSigmaDrReal(i,j)
123 SN=sqrt(Ssq*N2)
124 VisbeckK(i,j,bi,bj)=VisbeckK(i,j,bi,bj)+deltaH
125 & *GM_Visbeck_alpha*GM_Visbeck_length*GM_Visbeck_length*SN
126 ENDIF
127
128 C Limit range that KapGM can take
129 VisbeckK(i,j,bi,bj)=
130 & min(VisbeckK(i,j,bi,bj),GM_Visbeck_maxval_K)
131
132 #endif /* GM_VISBECK_VARIABLE_K */
133
134
135 #ifdef ALLOW_TIMEAVE
136 C-- Time-average
137 GM_Kwx_T(i,j,k,bi,bj)=GM_Kwx_T(i,j,k,bi,bj)
138 & +Kwx(i,j,k,bi,bj)*deltaTclock
139 GM_Kwy_T(i,j,k,bi,bj)=GM_Kwy_T(i,j,k,bi,bj)
140 & +Kwy(i,j,k,bi,bj)*deltaTclock
141 GM_Kwz_T(i,j,k,bi,bj)=GM_Kwz_T(i,j,k,bi,bj)
142 & +Kwz(i,j,k,bi,bj)*deltaTclock
143 #ifdef GM_VISBECK_VARIABLE_K
144 IF (K.EQ.Nr)
145 & Visbeck_K_T(i,j,bi,bj)=Visbeck_K_T(i,j,bi,bj)
146 & +VisbeckK(i,j,bi,bj)*deltaTclock
147 #endif
148 #endif /* ALLOW_TIMEAVE */
149 ENDDO
150 ENDDO
151
152 #ifdef ALLOW_TIMEAVE
153 GM_TimeAve(k,bi,bj)=GM_TimeAve(k,bi,bj)+deltaTclock
154 #endif
155
156
157 #ifdef GM_NON_UNITY_DIAGONAL
158 C Gradient of Sigma at U points
159 DO j=1-Oly+1,sNy+Oly-1
160 DO i=1-Olx+1,sNx+Olx-1
161 SlopeX(i,j)=sigmaX(i,j,km1)
162 & *_maskW(i,j,k,bi,bj)
163 SlopeY(i,j)=0.25*( sigmaY(i-1,j+1,k) +sigmaY(i,j+1,k)
164 & +sigmaY(i-1, j ,k) +sigmaY(i, j ,k) )
165 & *_maskW(i,j,k,bi,bj)
166 dSigmaDrReal(i,j)=0.25*( sigmaR(i-1,j, k ) +sigmaR(i,j, k )
167 & +sigmaR(i-1,j,kp1) +sigmaR(i,j,kp1) )
168 & *_maskW(i,j,k,bi,bj)
169 ENDDO
170 ENDDO
171
172 C Calculate slopes for use in tensor, taper and/or clip
173 CALL GMREDI_SLOPE_LIMIT(
174 I dSigmadRReal,
175 I rF(K),
176 U SlopeX, SlopeY,
177 O SlopeSqr, taperFct,
178 I bi, bj, myThid )
179
180 DO j=1-Oly+1,sNy+Oly-1
181 DO i=1-Olx+1,sNx+Olx-1
182 Kux(i,j,k,bi,bj)=taperFct(i,j)
183 ENDDO
184 ENDDO
185
186 C Gradient of Sigma at V points
187 DO j=1-Oly+1,sNy+Oly-1
188 DO i=1-Olx+1,sNx+Olx-1
189 SlopeX(i,j)=0.25*( sigmaX(i, j ,k) +sigmaX(i+1, j ,k)
190 & +sigmaX(i,j-1,k) +sigmaX(i+1,j-1,k) )
191 & *_maskS(i,j,k,bi,bj)
192 SlopeY(i,j)=sigmaY(i,j,km1)
193 & *_maskS(i,j,k,bi,bj)
194 dSigmaDrReal(i,j)=0.25*( sigmaR(i,j-1, k ) +sigmaR(i,j, k )
195 & +sigmaR(i,j-1,kp1) +sigmaR(i,j,kp1) )
196 & *_maskS(i,j,k,bi,bj)
197 ENDDO
198 ENDDO
199
200 C Calculate slopes for use in tensor, taper and/or clip
201 CALL GMREDI_SLOPE_LIMIT(
202 I dSigmadRReal,
203 I rF(K),
204 U SlopeX, SlopeY,
205 O SlopeSqr, taperFct,
206 I bi, bj, myThid )
207
208 DO j=1-Oly+1,sNy+Oly-1
209 DO i=1-Olx+1,sNx+Olx-1
210 Kvy(i,j,k,bi,bj)=taperFct(i,j)
211 ENDDO
212 ENDDO
213
214 #endif /* GM_NON_UNITY_DIAGONAL */
215
216
217
218 #endif /* ALLOW_GMREDI */
219
220 RETURN
221 END
222
223
224 SUBROUTINE GMREDI_CALC_TENSOR_DUMMY(
225 I bi, bj, iMin, iMax, jMin, jMax, K,
226 I sigmaX, sigmaY, sigmaR,
227 I myThid )
228 C /==========================================================\
229 C | SUBROUTINE GMREDI_CALC_TENSOR |
230 C | o Calculate tensor elements for GM/Redi tensor. |
231 C |==========================================================|
232 C \==========================================================/
233 IMPLICIT NONE
234
235 C == Global variables ==
236 #include "SIZE.h"
237 #include "GRID.h"
238 #include "DYNVARS.h"
239 #include "EEPARAMS.h"
240 #include "PARAMS.h"
241 #include "GMREDI.h"
242
243 C == Routine arguments ==
244 C
245 _RL sigmaX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
246 _RL sigmaY(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
247 _RL sigmaR(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
248 INTEGER bi,bj,iMin,iMax,jMin,jMax,K
249 INTEGER myThid
250 CEndOfInterface
251
252 INTEGER i, j
253
254 #ifdef ALLOW_GMREDI
255
256 DO j=1-Oly+1,sNy+Oly-1
257 DO i=1-Olx+1,sNx+Olx-1
258 Kwx(i,j,k,bi,bj) = 0.0
259 Kwy(i,j,k,bi,bj) = 0.0
260 Kwz(i,j,k,bi,bj) = 0.0
261 ENDDO
262 ENDDO
263 #endif /* ALLOW_GMREDI */
264
265 end

  ViewVC Help
Powered by ViewVC 1.1.22