95 |
_RL a(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) |
_RL a(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) |
96 |
_RL b(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) |
_RL b(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) |
97 |
_RL c(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) |
_RL c(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) |
98 |
|
#ifdef ALLOW_GGL90_HORIZDIFF |
99 |
|
C xA, yA - area of lateral faces |
100 |
|
C dfx, dfy - diffusive flux across lateral faces |
101 |
|
_RL xA (1-OLx:sNx+OLx,1-OLy:sNy+OLy) |
102 |
|
_RL yA (1-OLx:sNx+OLx,1-OLy:sNy+OLy) |
103 |
|
_RL dfx(1-OLx:sNx+OLx,1-OLy:sNy+OLy) |
104 |
|
_RL dfy(1-OLx:sNx+OLx,1-OLy:sNy+OLy) |
105 |
|
#endif /* ALLOW_GGL90_HORIZDIFF */ |
106 |
CEOP |
CEOP |
107 |
iMin = 2-OLx |
iMin = 2-OLx |
108 |
iMax = sNx+OLx-1 |
iMax = sNx+OLx-1 |
200 |
& ) |
& ) |
201 |
ENDDO |
ENDDO |
202 |
ENDDO |
ENDDO |
203 |
ENDDO |
ENDDO |
204 |
C |
C horizontal diffusion of TKE (requires an exchange in |
205 |
C Implicit time step to update TKE for k=1,Nr; TKE(Nr+1)=0 by default |
C do_fields_blocking_exchanges) |
206 |
C |
#ifdef ALLOW_GGL90_HORIZDIFF |
207 |
|
IF ( GGL90diffTKEh .GT. 0. _d 0 ) THEN |
208 |
|
DO K=2,Nr |
209 |
|
C common factors |
210 |
|
DO j=1-Oly,sNy+Oly |
211 |
|
DO i=1-Olx,sNx+Olx |
212 |
|
xA(i,j) = _dyG(i,j,bi,bj) |
213 |
|
& *drF(k)*_hFacW(i,j,k,bi,bj) |
214 |
|
yA(i,j) = _dxG(i,j,bi,bj) |
215 |
|
& *drF(k)*_hFacS(i,j,k,bi,bj) |
216 |
|
ENDDO |
217 |
|
ENDDO |
218 |
|
C Compute diffusive fluxes |
219 |
|
C ... across x-faces |
220 |
|
DO j=1-Oly,sNy+Oly |
221 |
|
dfx(1-Olx,j)=0. |
222 |
|
DO i=1-Olx+1,sNx+Olx |
223 |
|
dfx(i,j) = -GGL90diffTKEh*xA(i,j) |
224 |
|
& *_recip_dxC(i,j,bi,bj) |
225 |
|
& *(GGL90TKE(i,j,k,bi,bj)-GGL90TKE(i-1,j,k,bi,bj)) |
226 |
|
& *CosFacU(j,bi,bj) |
227 |
|
ENDDO |
228 |
|
ENDDO |
229 |
|
C ... across y-faces |
230 |
|
DO i=1-Olx,sNx+Olx |
231 |
|
dfy(i,1-Oly)=0. |
232 |
|
ENDDO |
233 |
|
DO j=1-Oly+1,sNy+Oly |
234 |
|
DO i=1-Olx,sNx+Olx |
235 |
|
dfy(i,j) = -GGL90diffTKEh*yA(i,j) |
236 |
|
& *_recip_dyC(i,j,bi,bj) |
237 |
|
& *(GGL90TKE(i,j,k,bi,bj)-GGL90TKE(i,j-1,k,bi,bj)) |
238 |
|
#ifdef ISOTROPIC_COS_SCALING |
239 |
|
& *CosFacV(j,bi,bj) |
240 |
|
#endif /* ISOTROPIC_COS_SCALING */ |
241 |
|
ENDDO |
242 |
|
ENDDO |
243 |
|
C Compute divergence of fluxes |
244 |
|
DO j=1-Oly,sNy+Oly-1 |
245 |
|
DO i=1-Olx,sNx+Olx-1 |
246 |
|
gTKE(i,j,k)=gTKE(i,j,k) |
247 |
|
& -_recip_hFacC(i,j,k,bi,bj)*recip_drF(k)*recip_rA(i,j,bi,bj) |
248 |
|
& *( (dfx(i+1,j)-dfx(i,j)) |
249 |
|
& +(dfy(i,j+1)-dfy(i,j)) |
250 |
|
& ) |
251 |
|
ENDDO |
252 |
|
ENDDO |
253 |
|
C end of k-loop |
254 |
|
ENDDO |
255 |
|
C end if GGL90diffTKEh .eq. 0. |
256 |
|
ENDIF |
257 |
|
#endif /* ALLOW_GGL90_HORIZDIFF */ |
258 |
|
|
259 |
|
C ============================================ |
260 |
|
C Implicit time step to update TKE for k=1,Nr; |
261 |
|
C TKE(Nr+1)=0 by default |
262 |
|
C ============================================ |
263 |
C set up matrix |
C set up matrix |
264 |
C-- Old aLower |
C-- Lower diagonal |
265 |
DO j=jMin,jMax |
DO j=jMin,jMax |
266 |
DO i=iMin,iMax |
DO i=iMin,iMax |
267 |
a(i,j,1) = 0. _d 0 |
a(i,j,1) = 0. _d 0 |
279 |
ENDDO |
ENDDO |
280 |
ENDDO |
ENDDO |
281 |
ENDDO |
ENDDO |
282 |
|
C-- Upper diagonal |
|
C-- Old aUpper |
|
283 |
DO j=jMin,jMax |
DO j=jMin,jMax |
284 |
DO i=iMin,iMax |
DO i=iMin,iMax |
285 |
c(i,j,1) = 0. _d 0 |
c(i,j,1) = 0. _d 0 |
287 |
ENDDO |
ENDDO |
288 |
ENDDO |
ENDDO |
289 |
CML DO k=1,Nr-1 |
CML DO k=1,Nr-1 |
290 |
DO k=2,Nr |
DO k=2,Nr-1 |
291 |
kp1=min(Nr,k+1) |
kp1=min(Nr,k+1) |
292 |
DO j=jMin,jMax |
DO j=jMin,jMax |
293 |
DO i=iMin,iMax |
DO i=iMin,iMax |
295 |
& *recip_drF( k )*recip_hFacI(i,j,k,bi,bj) |
& *recip_drF( k )*recip_hFacI(i,j,k,bi,bj) |
296 |
& *.5*(KappaE(i,j,k)+KappaE(i,j,kp1)) |
& *.5*(KappaE(i,j,k)+KappaE(i,j,kp1)) |
297 |
& *recip_drC(k) |
& *recip_drC(k) |
298 |
C IF (recip_hFacI(i,j,kp1,bi,bj).EQ.0.) c(i,j,k)=0. |
IF (recip_hFacI(i,j,kp1,bi,bj).EQ.0.) c(i,j,k)=0. |
299 |
ENDDO |
ENDDO |
300 |
ENDDO |
ENDDO |
301 |
ENDDO |
ENDDO |
302 |
|
C-- Center diagonal |
|
C-- Old aCenter |
|
303 |
DO k=1,Nr |
DO k=1,Nr |
304 |
DO j=jMin,jMax |
DO j=jMin,jMax |
305 |
DO i=iMin,iMax |
DO i=iMin,iMax |
348 |
C impose minimum TKE to avoid numerical undershoots below zero |
C impose minimum TKE to avoid numerical undershoots below zero |
349 |
GGL90TKE(I,J,K,bi,bj) = MAX( gTKE(I,J,K), GGL90TKEmin ) |
GGL90TKE(I,J,K,bi,bj) = MAX( gTKE(I,J,K), GGL90TKEmin ) |
350 |
& * maskC(I,J,K,bi,bj) |
& * maskC(I,J,K,bi,bj) |
|
C |
|
|
C end of time step |
|
|
C |
|
351 |
ENDDO |
ENDDO |
352 |
ENDDO |
ENDDO |
353 |
ENDDO |
ENDDO |
354 |
C |
C |
355 |
|
C end of time step |
356 |
|
C =============================== |
357 |
C compute viscosity coefficients |
C compute viscosity coefficients |
358 |
C |
C |
359 |
DO K=2,Nr |
DO K=2,Nr |