/[MITgcm]/MITgcm/model/src/calc_viscosity.F
ViewVC logotype

Contents of /MITgcm/model/src/calc_viscosity.F

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


Revision 1.12 - (show annotations) (download)
Thu Nov 2 18:00:52 2017 UTC (6 years, 5 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint66o, checkpoint66n, checkpoint66m, HEAD
Changes since 1.11: +230 -1 lines
- new options 1) to account for true vertical distance (including hFac)
     in vertical viscous flux and diffusive flux ;
  2) to increase vertical mixing near surface and/or bottom where partial
     cell is too thin ;
- for now, both additions above are within: #ifndef EXCLUDE_PCELL_MIX_CODE ;

1 C $Header: /u/gcmpack/MITgcm/model/src/calc_viscosity.F,v 1.11 2014/12/24 19:09:33 jmc Exp $
2 C $Name: $
3
4 #include "PACKAGES_CONFIG.h"
5 #include "CPP_OPTIONS.h"
6
7 CBOP
8 C !ROUTINE: CALC_VISCOSITY
9 C !INTERFACE:
10 SUBROUTINE CALC_VISCOSITY(
11 I bi,bj, iMin,iMax,jMin,jMax,
12 O kappaRU, kappaRV,
13 I myThid )
14
15 C !DESCRIPTION: \bv
16 C *==========================================================*
17 C | SUBROUTINE CALC_VISCOSITY
18 C | o Calculate net vertical viscosity
19 C *==========================================================*
20 C \ev
21
22 C !USES:
23 IMPLICIT NONE
24 C == GLobal variables ==
25 #include "SIZE.h"
26 #include "EEPARAMS.h"
27 #include "PARAMS.h"
28 #include "DYNVARS.h"
29 #include "GRID.h"
30
31 C !INPUT/OUTPUT PARAMETERS:
32 C == Routine arguments ==
33 C iMin,iMax,jMin,jMax :: Range of points for which calculation
34 C bi,bj :: current tile indices
35 C kappaRU :: Total vertical viscosity for zonal flow.
36 C kappaRV :: Total vertical viscosity for meridional flow.
37 C myThid :: my Thread Id number
38 INTEGER iMin,iMax,jMin,jMax
39 INTEGER bi,bj
40 _RL kappaRU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1)
41 _RL kappaRV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1)
42 INTEGER myThid
43
44 C !LOCAL VARIABLES:
45 C == Local variables ==
46 C i, j, k :: Loop counters
47 INTEGER i,j,k
48 INTEGER ki
49 #ifndef EXCLUDE_PCELL_MIX_CODE
50 INTEGER km, mixSurf, mixBott
51 _RL pC_kFac
52 _RL tmpFac(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
53 #endif
54 CEOP
55
56 DO k = 1,Nr+1
57 ki = MIN(k,Nr)
58
59 DO j = 1-OLy, sNy+OLy
60 DO i = 1-OLx, sNx+OLx
61 kappaRU(i,j,k) = viscArNr(ki)
62 kappaRV(i,j,k) = viscArNr(ki)
63 ENDDO
64 ENDDO
65
66 #ifdef ALLOW_KPP
67 IF ( useKPP .AND. k.LE.Nr ) THEN
68 CALL KPP_CALC_VISC(
69 I bi,bj, iMin,iMax,jMin,jMax, k,
70 O kappaRU, kappaRV,
71 I myThid)
72 ENDIF
73 #endif
74
75 #ifdef ALLOW_PP81
76 IF ( usePP81 .AND. k.LE.Nr ) THEN
77 CALL PP81_CALC_VISC(
78 I bi,bj, iMin,iMax,jMin,jMax, k,
79 O kappaRU, kappaRV,
80 I myThid)
81 ENDIF
82 #endif
83
84 #ifdef ALLOW_KL10
85 IF ( useKL10 .AND. k.LE.Nr ) THEN
86 CALL KL10_CALC_VISC(
87 I bi,bj, iMin,iMax,jMin,jMax, k,
88 O kappaRU, kappaRV,
89 I myThid)
90 ENDIF
91 #endif
92
93 #ifdef ALLOW_MY82
94 IF ( useMY82 .AND. k.LE.Nr ) THEN
95 CALL MY82_CALC_VISC(
96 I bi,bj, iMin,iMax,jMin,jMax, k,
97 O kappaRU, kappaRV,
98 I myThid)
99 ENDIF
100 #endif
101
102 #ifdef ALLOW_GGL90
103 IF ( useGGL90 .AND. k.LE.Nr ) THEN
104 CALL GGL90_CALC_VISC(
105 I bi,bj, iMin,iMax,jMin,jMax, k,
106 O kappaRU, kappaRV,
107 I myThid)
108 ENDIF
109 #endif
110
111 IF ( k.EQ.Nr+1 .AND.
112 & ( usePP81 .OR. useKL10 .OR. useMY82 .OR. useGGL90 )
113 & ) THEN
114 DO j = 1-OLy, sNy+OLy
115 DO i = 1-OLx, sNx+OLx
116 kappaRU(i,j,k) = kappaRU(i,j,ki)
117 kappaRV(i,j,k) = kappaRV(i,j,ki)
118 ENDDO
119 ENDDO
120 ENDIF
121
122 C-- end of k loop
123 ENDDO
124
125 #ifndef EXCLUDE_PCELL_MIX_CODE
126 IF ( interViscAr_pCell ) THEN
127 C-- This is a hack: alter vertical viscosity (instead of changing many S/R)
128 C in order to account for missing hFac in viscous term
129 DO k = 2,Nr
130 km = k - 1
131 C- account for true distance (including hFac) in vertical gradient
132 DO j = 2-OLy, sNy+OLy
133 DO i = 2-OLx, sNx+OLx
134 IF ( k.GT.kSurfW(i,j,bi,bj) .AND.
135 & k.LE.MIN( kLowC(i,j,bi,bj), kLowC(i-1,j,bi,bj) )
136 & ) THEN
137 kappaRU(i,j,k) = kappaRU(i,j,k)
138 & *twoRL/(hFacW(i,j,km,bi,bj)+hFacW(i,j,k,bi,bj))
139 ENDIF
140 ENDDO
141 ENDDO
142 DO j = 2-OLy, sNy+OLy
143 DO i = 2-OLx, sNx+OLx
144 IF ( k.GT.kSurfS(i,j,bi,bj) .AND.
145 & k.LE.MIN( kLowC(i,j,bi,bj), kLowC(i,j-1,bi,bj) )
146 & ) THEN
147 kappaRV(i,j,k) = kappaRV(i,j,k)
148 & *twoRL/(hFacS(i,j,km,bi,bj)+hFacS(i,j,k,bi,bj))
149 ENDIF
150 ENDDO
151 ENDDO
152 ENDDO
153 ENDIF
154
155 IF ( pCellMix_select.GT.0 ) THEN
156 C-- This is a hack: alter vertical viscosity (instead of changing many S/R)
157 C in order to to increase mixing for too thin surface/bottom partial cell
158 mixSurf = pCellMix_select/10
159 mixBott = MOD(pCellMix_select,10)
160 DO k = 2,Nr
161 km = k - 1
162 pC_kFac = 1.
163 IF ( pCellMix_delR.LT.drF(k) )
164 & pC_kFac = pCellMix_delR*recip_drF(k)
165
166 C- Increase KappaRU above bottom level:
167 IF ( mixBott.GE.1 ) THEN
168 DO j = 2-OLy, sNy+OLy
169 DO i = 2-OLx, sNx+OLx
170 tmpFac(i,j) = 0. _d 0
171 IF ( k.EQ.MIN( kLowC(i,j,bi,bj), kLowC(i-1,j,bi,bj) )
172 & .AND. k.GT.kSurfW(i,j,bi,bj) ) THEN
173 tmpFac(i,j) = pC_kFac*_recip_hFacW(i,j,k,bi,bj)
174 ENDIF
175 ENDDO
176 ENDDO
177 ENDIF
178 IF ( mixBott.EQ.2 ) THEN
179 DO j = 2-OLy, sNy+OLy
180 DO i = 2-OLx, sNx+OLx
181 tmpFac(i,j) = tmpFac(i,j)*tmpFac(i,j)
182 ENDDO
183 ENDDO
184 ELSEIF ( mixBott.EQ.3 ) THEN
185 DO j = 2-OLy, sNy+OLy
186 DO i = 2-OLx, sNx+OLx
187 tmpFac(i,j) = tmpFac(i,j)*tmpFac(i,j)*tmpFac(i,j)
188 ENDDO
189 ENDDO
190 ELSEIF ( mixBott.EQ.4 ) THEN
191 DO j = 2-OLy, sNy+OLy
192 DO i = 2-OLx, sNx+OLx
193 tmpFac(i,j) = tmpFac(i,j)*tmpFac(i,j)
194 & *tmpFac(i,j)*tmpFac(i,j)
195 ENDDO
196 ENDDO
197 ENDIF
198 IF ( mixBott.GE.1 ) THEN
199 C- increase mixing above bottom (by ~(1/hFac)^mixBott) if too thin p-cell
200 DO j = 2-OLy, sNy+OLy
201 DO i = 2-OLx, sNx+OLx
202 tmpFac(i,j) = MIN( tmpFac(i,j), pCellMix_maxFac )
203 kappaRU(i,j,k) = MAX( kappaRU(i,j,k),
204 & pCellMix_viscAr(k)*tmpFac(i,j) )
205 ENDDO
206 ENDDO
207 ENDIF
208
209 C- Increase KappaRV above bottom level:
210 IF ( mixBott.GE.1 ) THEN
211 DO j = 2-OLy, sNy+OLy
212 DO i = 2-OLx, sNx+OLx
213 tmpFac(i,j) = 0. _d 0
214 IF ( k.EQ.MIN( kLowC(i,j,bi,bj), kLowC(i,j-1,bi,bj) )
215 & .AND. k.GT.kSurfS(i,j,bi,bj) ) THEN
216 tmpFac(i,j) = pC_kFac*_recip_hFacS(i,j,k,bi,bj)
217 ENDIF
218 ENDDO
219 ENDDO
220 ENDIF
221 IF ( mixBott.EQ.2 ) THEN
222 DO j = 2-OLy, sNy+OLy
223 DO i = 2-OLx, sNx+OLx
224 tmpFac(i,j) = tmpFac(i,j)*tmpFac(i,j)
225 ENDDO
226 ENDDO
227 ELSEIF ( mixBott.EQ.3 ) THEN
228 DO j = 2-OLy, sNy+OLy
229 DO i = 2-OLx, sNx+OLx
230 tmpFac(i,j) = tmpFac(i,j)*tmpFac(i,j)*tmpFac(i,j)
231 ENDDO
232 ENDDO
233 ELSEIF ( mixBott.EQ.4 ) THEN
234 DO j = 2-OLy, sNy+OLy
235 DO i = 2-OLx, sNx+OLx
236 tmpFac(i,j) = tmpFac(i,j)*tmpFac(i,j)
237 & *tmpFac(i,j)*tmpFac(i,j)
238 ENDDO
239 ENDDO
240 ENDIF
241 IF ( mixBott.GE.1 ) THEN
242 C- increase mixing above bottom (by ~(1/hFac)^mixBott) if too thin p-cell
243 DO j = 2-OLy, sNy+OLy
244 DO i = 2-OLx, sNx+OLx
245 tmpFac(i,j) = MIN( tmpFac(i,j), pCellMix_maxFac )
246 kappaRV(i,j,k) = MAX( kappaRV(i,j,k),
247 & pCellMix_viscAr(k)*tmpFac(i,j) )
248 ENDDO
249 ENDDO
250 ENDIF
251
252 pC_kFac = 1.
253 IF ( pCellMix_delR.LT.drF(km) )
254 & pC_kFac = pCellMix_delR*recip_drF(km)
255
256 C- Increase KappaRU below surface level:
257 IF ( mixSurf.GE.1 ) THEN
258 DO j = 2-OLy, sNy+OLy
259 DO i = 2-OLx, sNx+OLx
260 tmpFac(i,j) = 0. _d 0
261 IF ( km.EQ.kSurfW(i,j,bi,bj) .AND.
262 & km.LT.MIN( kLowC(i,j,bi,bj), kLowC(i-1,j,bi,bj) )
263 & ) THEN
264 tmpFac(i,j) = pC_kFac*_recip_hFacW(i,j,km,bi,bj)
265 ENDIF
266 ENDDO
267 ENDDO
268 ENDIF
269 IF ( mixSurf.EQ.2 ) THEN
270 DO j = 2-OLy, sNy+OLy
271 DO i = 2-OLx, sNx+OLx
272 tmpFac(i,j) = tmpFac(i,j)*tmpFac(i,j)
273 ENDDO
274 ENDDO
275 ELSEIF ( mixSurf.EQ.3 ) THEN
276 DO j = 2-OLy, sNy+OLy
277 DO i = 2-OLx, sNx+OLx
278 tmpFac(i,j) = tmpFac(i,j)*tmpFac(i,j)*tmpFac(i,j)
279 ENDDO
280 ENDDO
281 ELSEIF ( mixSurf.EQ.4 ) THEN
282 DO j = 2-OLy, sNy+OLy
283 DO i = 2-OLx, sNx+OLx
284 tmpFac(i,j) = tmpFac(i,j)*tmpFac(i,j)
285 & *tmpFac(i,j)*tmpFac(i,j)
286 ENDDO
287 ENDDO
288 ENDIF
289 IF ( mixSurf.GE.1 ) THEN
290 C- increase mixing below surface (by ~(1/hFac)^mixSurf) if too thin p-cell
291 DO j = 2-OLy, sNy+OLy
292 DO i = 2-OLx, sNx+OLx
293 tmpFac(i,j) = MIN( tmpFac(i,j), pCellMix_maxFac )
294 kappaRU(i,j,k) = MAX( kappaRU(i,j,k),
295 & pCellMix_viscAr(k)*tmpFac(i,j) )
296 ENDDO
297 ENDDO
298 ENDIF
299
300 C- Increase KappaRV below surface level:
301 IF ( mixSurf.GE.1 ) THEN
302 DO j = 2-OLy, sNy+OLy
303 DO i = 2-OLx, sNx+OLx
304 tmpFac(i,j) = 0. _d 0
305 IF ( km.EQ.kSurfS(i,j,bi,bj) .AND.
306 & km.LT.MIN( kLowC(i,j,bi,bj), kLowC(i,j-1,bi,bj) )
307 & ) THEN
308 tmpFac(i,j) = pC_kFac*_recip_hFacS(i,j,km,bi,bj)
309 ENDIF
310 ENDDO
311 ENDDO
312 ENDIF
313 IF ( mixSurf.EQ.2 ) THEN
314 DO j = 2-OLy, sNy+OLy
315 DO i = 2-OLx, sNx+OLx
316 tmpFac(i,j) = tmpFac(i,j)*tmpFac(i,j)
317 ENDDO
318 ENDDO
319 ELSEIF ( mixSurf.EQ.3 ) THEN
320 DO j = 2-OLy, sNy+OLy
321 DO i = 2-OLx, sNx+OLx
322 tmpFac(i,j) = tmpFac(i,j)*tmpFac(i,j)*tmpFac(i,j)
323 ENDDO
324 ENDDO
325 ELSEIF ( mixSurf.EQ.4 ) THEN
326 DO j = 2-OLy, sNy+OLy
327 DO i = 2-OLx, sNx+OLx
328 tmpFac(i,j) = tmpFac(i,j)*tmpFac(i,j)
329 & *tmpFac(i,j)*tmpFac(i,j)
330 ENDDO
331 ENDDO
332 ENDIF
333 IF ( mixSurf.GE.1 ) THEN
334 C- increase mixing below surface (by ~(1/hFac)^mixSurf) if too thin p-cell
335 DO j = 2-OLy, sNy+OLy
336 DO i = 2-OLx, sNx+OLx
337 tmpFac(i,j) = MIN( tmpFac(i,j), pCellMix_maxFac )
338 kappaRV(i,j,k) = MAX( kappaRV(i,j,k),
339 & pCellMix_viscAr(k)*tmpFac(i,j) )
340 ENDDO
341 ENDDO
342 ENDIF
343
344 C-- end of k loop
345 ENDDO
346 ENDIF
347 #endif /* ndef EXCLUDE_PCELL_MIX_CODE */
348
349 RETURN
350 END

  ViewVC Help
Powered by ViewVC 1.1.22