ROMS
Loading...
Searching...
No Matches
tkebc_ex.F
Go to the documentation of this file.
1#include "cppdefs.h"
2 MODULE tkebc_mod
3#if defined SOLVE3D && (defined MY25_MIXING || defined GLS_MIXING)
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This subroutine sets lateral boundary conditions for turbulent !
13! kinetic energy and turbulent length scale variables associated !
14! with the Mellor and Yamada or GOTM closures. !
15! !
16!=======================================================================
17!
18 implicit none
19
20 PRIVATE
21 PUBLIC :: tkebc_tile
22
23 CONTAINS
24!
25!***********************************************************************
26 SUBROUTINE tkebc (ng, tile, nout)
27!***********************************************************************
28!
29 USE mod_param
30 USE mod_mixing
31 USE mod_stepping
32!
33! Imported variable declarations.
34!
35 integer, intent(in) :: ng, tile, nout
36!
37! Local variable declarations.
38!
39# include "tile.h"
40!
41 CALL tkebc_tile (ng, tile, &
42 & lbi, ubi, lbj, ubj, n(ng), &
43 & imins, imaxs, jmins, jmaxs, &
44 & nout, nstp(ng), &
45 & mixing(ng)% gls, &
46 & mixing(ng)% tke)
47 RETURN
48 END SUBROUTINE tkebc
49!
50!***********************************************************************
51 SUBROUTINE tkebc_tile (ng, tile, &
52 & LBi, UBi, LBj, UBj, UBk, &
53 & IminS, ImaxS, JminS, JmaxS, &
54 & nout, nstp, &
55 & gls, tke)
56!***********************************************************************
57!
58 USE mod_param
59 USE mod_grid
60 USE mod_ncparam
61 USE mod_scalars
62!
63! Imported variable declarations.
64!
65 integer, intent(in) :: ng, tile
66 integer, intent(in) :: LBi, UBi, LBj, UBj, UBk
67 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
68 integer, intent(in) :: nout, nstp
69!
70# ifdef ASSUMED_SHAPE
71 real(r8), intent(inout) :: gls(LBi:,LBj:,0:,:)
72 real(r8), intent(inout) :: tke(LBi:,LBj:,0:,:)
73# else
74 real(r8), intent(inout) :: gls(LBi:UBi,LBj:UBj,0:UBk,3)
75 real(r8), intent(inout) :: tke(LBi:UBi,LBj:UBj,0:UBk,3)
76# endif
77!
78! Local variable declarations.
79!
80 integer :: i, j, k
81
82 real(r8), parameter :: eps = 1.0e-20_r8
83
84 real(r8) :: Ce, Cx, cff, dKde, dKdt, dKdx
85
86 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: grad
87 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: gradL
88
89# include "set_bounds.h"
90!
91!-----------------------------------------------------------------------
92! Lateral boundary conditions at the western edge.
93!-----------------------------------------------------------------------
94!
95 IF (domain(ng)%Western_Edge(tile)) THEN
96!
97! Western edge, implicit upstream radiation condition.
98!
99 IF (lbc(iwest,ismtke,ng)%radiation) THEN
100 DO k=0,n(ng)
101 DO j=jstr,jend+1
102 grad(istr-1,j)=tke(istr-1,j ,k,nstp)- &
103 & tke(istr-1,j-1,k,nstp)
104# ifdef MASKING
105 grad(istr-1,j)=grad(istr-1,j)* &
106 & grid(ng)%vmask(istr-1,j)
107# endif
108 grad(istr ,j)=tke(istr ,j ,k,nstp)- &
109 & tke(istr ,j-1,k,nstp)
110# ifdef MASKING
111 grad(istr ,j)=grad(istr ,j)* &
112 & grid(ng)%vmask(istr ,j)
113# endif
114 gradl(istr-1,j)=gls(istr-1,j ,k,nstp)- &
115 & gls(istr-1,j-1,k,nstp)
116# ifdef MASKING
117 gradl(istr-1,j)=gradl(istr-1,j)* &
118 & grid(ng)%vmask(istr-1,j)
119# endif
120 gradl(istr ,j)=gls(istr ,j ,k,nstp)- &
121 & gls(istr ,j-1,k,nstp)
122# ifdef MASKING
123 gradl(istr ,j)=gradl(istr ,j)* &
124 & grid(ng)%vmask(istr ,j)
125# endif
126 END DO
127 DO j=jstr,jend
128 IF (lbc_apply(ng)%west(j)) THEN
129 dkdt=tke(istr,j,k,nstp)-tke(istr ,j,k,nout)
130 dkdx=tke(istr,j,k,nstp)-tke(istr+1,j,k,nstp)
131 IF ((dkdt*dkdx).lt.0.0_r8) dkdt=0.0_r8
132 IF ((dkdt*(grad(istr,j )+ &
133 & grad(istr,j+1))).gt.0.0_r8) THEN
134 dkde=grad(istr,j )
135 ELSE
136 dkde=grad(istr,j+1)
137 END IF
138 cff=dkdt/max(dkdx*dkdx+dkde*dkde,eps)
139 cx=min(1.0_r8,cff*dkdx)
140# ifdef RADIATION_2D
141 ce=min(1.0_r8,max(cff*dkde,-1.0_r8))
142# else
143 ce=0.0_r8
144# endif
145 tke(istr-1,j,k,nout)=(1.0_r8-cx)*tke(istr-1,j,k,nstp)+ &
146 & cx*tke(istr,j,k,nstp)- &
147 & max(ce,0.0_r8)*grad(istr-1,j )- &
148 & min(ce,0.0_r8)*grad(istr-1,j+1)
149# ifdef MASKING
150 tke(istr-1,j,k,nout)=tke(istr-1,j,k,nout)* &
151 & grid(ng)%rmask(istr-1,j)
152# endif
153 dkdt=gls(istr,j,k,nstp)-gls(istr ,j,k,nout)
154 dkdx=gls(istr,j,k,nstp)-gls(istr+1,j,k,nstp)
155 IF ((dkdt*dkdx).lt.0.0_r8) dkdt=0.0_r8
156 IF ((dkdt*(gradl(istr,j )+ &
157 & gradl(istr,j+1))).gt.0.0_r8) THEN
158 dkde=gradl(istr,j )
159 ELSE
160 dkde=gradl(istr,j+1)
161 END IF
162 cff=dkdt/max(dkdx*dkdx+dkde*dkde,eps)
163 cx=min(1.0_r8,cff*dkdx)
164# ifdef RADIATION_2D
165 ce=min(1.0_r8,max(cff*dkde,-1.0_r8))
166# else
167 ce=0.0_r8
168# endif
169 gls(istr-1,j,k,nout)=(1.0_r8-cx)*gls(istr-1,j,k,nstp)+ &
170 & cx*gls(istr,j,k,nstp)- &
171 & max(ce,0.0_r8)*gradl(istr-1,j )- &
172 & min(ce,0.0_r8)*gradl(istr-1,j+1)
173# ifdef MASKING
174 gls(istr-1,j,k,nout)=gls(istr-1,j,k,nout)* &
175 & grid(ng)%rmask(istr-1,j)
176# endif
177 END IF
178 END DO
179 END DO
180!
181! Western edge, gradient boundary condition.
182!
183 ELSE IF (lbc(iwest,ismtke,ng)%gradient) THEN
184 DO k=0,n(ng)
185 DO j=jstr,jend
186 IF (lbc_apply(ng)%west(j)) THEN
187 tke(istr-1,j,k,nout)=tke(istr,j,k,nout)
188# ifdef MASKING
189 tke(istr-1,j,k,nout)=tke(istr-1,j,k,nout)* &
190 & grid(ng)%rmask(istr-1,j)
191# endif
192 gls(istr-1,j,k,nout)=gls(istr,j,k,nout)
193# ifdef MASKING
194 gls(istr-1,j,k,nout)=gls(istr-1,j,k,nout)* &
195 & grid(ng)%rmask(istr-1,j)
196# endif
197 END IF
198 END DO
199 END DO
200!
201! Western edge, closed boundary condition.
202!
203 ELSE IF (lbc(iwest,ismtke,ng)%closed) THEN
204 DO k=0,n(ng)
205 DO j=jstr,jend
206 IF (lbc_apply(ng)%west(j)) THEN
207 tke(istr-1,j,k,nout)=tke(istr,j,k,nout)
208# ifdef MASKING
209 tke(istr-1,j,k,nout)=tke(istr-1,j,k,nout)* &
210 & grid(ng)%rmask(istr-1,j)
211# endif
212 gls(istr-1,j,k,nout)=gls(istr,j,k,nout)
213# ifdef MASKING
214 gls(istr-1,j,k,nout)=gls(istr-1,j,k,nout)* &
215 & grid(ng)%rmask(istr-1,j)
216# endif
217 END IF
218 END DO
219 END DO
220 END IF
221 END IF
222!
223!-----------------------------------------------------------------------
224! Lateral boundary conditions at the eastern edge.
225!-----------------------------------------------------------------------
226!
227 IF (domain(ng)%Eastern_Edge(tile)) THEN
228!
229! Eastern edge, implicit upstream radiation condition.
230!
231 IF (lbc(ieast,ismtke,ng)%radiation) THEN
232 DO k=0,n(ng)
233 DO j=jstr,jend+1
234 grad(iend ,j)=tke(iend ,j ,k,nstp)- &
235 & tke(iend ,j-1,k,nstp)
236# ifdef MASKING
237 grad(iend ,j)=grad(iend ,j)* &
238 & grid(ng)%vmask(iend ,j)
239# endif
240 grad(iend+1,j)=tke(iend+1,j ,k,nstp)- &
241 & tke(iend+1,j-1,k,nstp)
242# ifdef MASKING
243 grad(iend+1,j)=grad(iend+1,j)* &
244 & grid(ng)%vmask(iend+1,j)
245# endif
246 gradl(iend ,j)=gls(iend ,j ,k,nstp)- &
247 & gls(iend ,j-1,k,nstp)
248# ifdef MASKING
249 gradl(iend ,j)=gradl(iend ,j)* &
250 & grid(ng)%vmask(iend ,j)
251# endif
252 gradl(iend+1,j)=gls(iend+1,j ,k,nstp)- &
253 & gls(iend+1,j-1,k,nstp)
254# ifdef MASKING
255 gradl(iend+1,j)=gradl(iend+1,j)* &
256 & grid(ng)%vmask(iend+1,j)
257# endif
258 END DO
259 DO j=jstr,jend
260 IF (lbc_apply(ng)%east(j)) THEN
261 dkdt=tke(iend,j,k,nstp)-tke(iend ,j,k,nout)
262 dkdx=tke(iend,j,k,nstp)-tke(iend-1,j,k,nstp)
263 IF ((dkdt*dkdx).lt.0.0_r8) dkdt=0.0_r8
264 IF ((dkdt*(grad(iend,j )+ &
265 & grad(iend,j+1))).gt.0.0_r8) THEN
266 dkde=grad(iend,j )
267 ELSE
268 dkde=grad(iend,j+1)
269 END IF
270 cff=dkdt/max(dkdx*dkdx+dkde*dkde,eps)
271 cx=min(1.0_r8,cff*dkdx)
272# ifdef RADIATION_2D
273 ce=min(1.0_r8,max(cff*dkde,-1.0_r8))
274# else
275 ce=0.0_r8
276# endif
277 tke(iend+1,j,k,nout)=(1.0_r8-cx)*tke(iend+1,j,k,nstp)+ &
278 & cx*tke(iend,j,k,nstp)- &
279 & max(ce,0.0_r8)*grad(iend+1,j )- &
280 & min(ce,0.0_r8)*grad(iend+1,j+1)
281# ifdef MASKING
282 tke(iend+1,j,k,nout)=tke(iend+1,j,k,nout)* &
283 & grid(ng)%rmask(iend+1,j)
284# endif
285 dkdt=gls(iend,j,k,nstp)-gls(iend ,j,k,nout)
286 dkdx=gls(iend,j,k,nstp)-gls(iend-1,j,k,nstp)
287 IF ((dkdt*dkdx).lt.0.0_r8) dkdt=0.0_r8
288 IF ((dkdt*(gradl(iend,j )+ &
289 & gradl(iend,j+1))).gt.0.0_r8) THEN
290 dkde=gradl(iend,j )
291 ELSE
292 dkde=gradl(iend,j+1)
293 END IF
294 cff=dkdt/max(dkdx*dkdx+dkde*dkde,eps)
295 cx=min(1.0_r8,cff*dkdx)
296# ifdef RADIATION_2D
297 ce=min(1.0_r8,max(cff*dkde,-1.0_r8))
298# else
299 ce=0.0_r8
300# endif
301 gls(iend+1,j,k,nout)=(1.0_r8-cx)*gls(iend+1,j,k,nstp)+ &
302 & cx*gls(iend,j,k,nstp)- &
303 & max(ce,0.0_r8)*gradl(iend+1,j )- &
304 & min(ce,0.0_r8)*gradl(iend+1,j+1)
305# ifdef MASKING
306 gls(iend+1,j,k,nout)=gls(iend+1,j,k,nout)* &
307 & grid(ng)%rmask(iend+1,j)
308# endif
309 END IF
310 END DO
311 END DO
312!
313! Eastern edge, gradient boundary condition.
314!
315 ELSE IF (lbc(ieast,ismtke,ng)%gradient) THEN
316 DO k=0,n(ng)
317 DO j=jstr,jend
318 IF (lbc_apply(ng)%east(j)) THEN
319 tke(iend+1,j,k,nout)=tke(iend,j,k,nout)
320# ifdef MASKING
321 tke(iend+1,j,k,nout)=tke(iend+1,j,k,nout)* &
322 & grid(ng)%rmask(iend+1,j)
323# endif
324 gls(iend+1,j,k,nout)=gls(iend,j,k,nout)
325# ifdef MASKING
326 gls(iend+1,j,k,nout)=gls(iend+1,j,k,nout)* &
327 & grid(ng)%rmask(iend+1,j)
328# endif
329 END IF
330 END DO
331 END DO
332!
333! Eastern edge, closed boundary condition.
334!
335 ELSE IF (lbc(ieast,ismtke,ng)%closed) THEN
336 DO k=0,n(ng)
337 DO j=jstr,jend
338 IF (lbc_apply(ng)%east(j)) THEN
339 tke(iend+1,j,k,nout)=tke(iend,j,k,nout)
340# ifdef MASKING
341 tke(iend+1,j,k,nout)=tke(iend+1,j,k,nout)* &
342 & grid(ng)%rmask(iend+1,j)
343# endif
344 gls(iend+1,j,k,nout)=gls(iend,j,k,nout)
345# ifdef MASKING
346 gls(iend+1,j,k,nout)=gls(iend+1,j,k,nout)* &
347 & grid(ng)%rmask(iend+1,j)
348# endif
349 END IF
350 END DO
351 END DO
352 END IF
353 END IF
354!
355!-----------------------------------------------------------------------
356! Lateral boundary conditions at the southern edge.
357!-----------------------------------------------------------------------
358!
359 IF (domain(ng)%Southern_Edge(tile)) THEN
360!
361! Southern edge, implicit upstream radiation condition.
362!
363 IF (lbc(isouth,ismtke,ng)%radiation) THEN
364 DO k=0,n(ng)
365 DO i=istr,iend+1
366 grad(i,jstr )=tke(i ,jstr ,k,nstp)- &
367 & tke(i-1,jstr ,k,nstp)
368# ifdef MASKING
369 grad(i,jstr )=grad(i,jstr )*grid(ng)%umask(i,jstr )
370# endif
371 grad(i,jstr-1)=tke(i ,jstr-1,k,nstp)- &
372 & tke(i-1,jstr-1,k,nstp)
373# ifdef MASKING
374 grad(i,jstr-1)=grad(i,jstr-1)*grid(ng)%umask(i,jstr-1)
375# endif
376 gradl(i,jstr )=gls(i ,jstr ,k,nstp)- &
377 & gls(i-1,jstr ,k,nstp)
378# ifdef MASKING
379 gradl(i,jstr )=gradl(i,jstr )*grid(ng)%umask(i,jstr )
380# endif
381 gradl(i,jstr-1)=gls(i ,jstr-1,k,nstp)- &
382 & gls(i-1,jstr-1,k,nstp)
383# ifdef MASKING
384 gradl(i,jstr-1)=gradl(i,jstr-1)*grid(ng)%umask(i,jstr-1)
385# endif
386 END DO
387 DO i=istr,iend
388 IF (lbc_apply(ng)%south(i)) THEN
389 dkdt=tke(i,jstr,k,nstp)-tke(i,jstr ,k,nout)
390 dkde=tke(i,jstr,k,nstp)-tke(i,jstr+1,k,nstp)
391 IF ((dkdt*dkde).lt.0.0_r8) dkdt=0.0_r8
392 IF ((dkdt*(grad(i ,jstr)+ &
393 & grad(i+1,jstr))).gt.0.0_r8) THEN
394 dkdx=grad(i ,jstr)
395 ELSE
396 dkdx=grad(i+1,jstr)
397 END IF
398 cff=dkdt/max(dkdx*dkdx+dkde*dkde,eps)
399# ifdef RADIATION_2D
400 cx=min(1.0_r8,max(cff*dkdx,-1.0_r8))
401# else
402 cx=0.0_r8
403# endif
404 ce=min(1.0_r8,cff*dkde)
405 tke(i,jstr-1,k,nout)=(1.0_r8-ce)*tke(i,jstr-1,k,nstp)+ &
406 & ce*tke(i,jstr,k,nstp)- &
407 & max(cx,0.0_r8)*grad(i ,jstr-1)- &
408 & min(cx,0.0_r8)*grad(i+1,jstr-1)
409# ifdef MASKING
410 tke(i,jstr-1,k,nout)=tke(i,jstr-1,k,nout)* &
411 & grid(ng)%rmask(i,jstr-1)
412# endif
413 dkdt=gls(i,jstr,k,nstp)-gls(i,jstr ,k,nout)
414 dkde=gls(i,jstr,k,nstp)-gls(i,jstr+1,k,nstp)
415 IF ((dkdt*dkde).lt.0.0_r8) dkdt=0.0_r8
416 IF ((dkdt*(gradl(i ,jstr)+ &
417 & gradl(i+1,jstr))).gt.0.0_r8) THEN
418 dkdx=gradl(i ,jstr)
419 ELSE
420 dkdx=gradl(i+1,jstr)
421 END IF
422 cff=dkdt/max(dkdx*dkdx+dkde*dkde,eps)
423# ifdef RADIATION_2D
424 cx=min(1.0_r8,max(cff*dkdx,-1.0_r8))
425# else
426 cx=0.0_r8
427# endif
428 ce=min(1.0_r8,cff*dkde)
429 gls(i,jstr-1,k,nout)=(1.0_r8-ce)*gls(i,jstr-1,k,nstp)+ &
430 & ce*gls(i,jstr,k,nstp)- &
431 & max(cx,0.0_r8)*gradl(i ,jstr-1)- &
432 & min(cx,0.0_r8)*gradl(i+1,jstr-1)
433# ifdef MASKING
434 gls(i,jstr-1,k,nout)=gls(i,jstr-1,k,nout)* &
435 & grid(ng)%rmask(i,jstr-1)
436# endif
437 END IF
438 END DO
439 END DO
440!
441! Southern edge, gradient boundary condition.
442!
443 ELSE IF (lbc(isouth,ismtke,ng)%gradient) THEN
444 DO k=0,n(ng)
445 DO i=istr,iend
446 IF (lbc_apply(ng)%south(i)) THEN
447 tke(i,jstr-1,k,nout)=tke(i,jstr,k,nout)
448# ifdef MASKING
449 tke(i,jstr-1,k,nout)=tke(i,jstr-1,k,nout)* &
450 & grid(ng)%rmask(i,jstr-1)
451# endif
452 gls(i,jstr-1,k,nout)=gls(i,jstr,k,nout)
453# ifdef MASKING
454 gls(i,jstr-1,k,nout)=gls(i,jstr-1,k,nout)* &
455 & grid(ng)%rmask(i,jstr-1)
456# endif
457 END IF
458 END DO
459 END DO
460!
461! Southern edge, closed boundary condition.
462!
463 ELSE IF (lbc(isouth,ismtke,ng)%closed) THEN
464 DO k=0,n(ng)
465 DO i=istr,iend
466 IF (lbc_apply(ng)%south(i)) THEN
467 tke(i,jstr-1,k,nout)=tke(i,jstr,k,nout)
468# ifdef MASKING
469 tke(i,jstr-1,k,nout)=tke(i,jstr-1,k,nout)* &
470 & grid(ng)%rmask(i,jstr-1)
471# endif
472 gls(i,jstr-1,k,nout)=gls(i,jstr,k,nout)
473# ifdef MASKING
474 gls(i,jstr-1,k,nout)=gls(i,jstr-1,k,nout)* &
475 & grid(ng)%rmask(i,jstr-1)
476# endif
477 END IF
478 END DO
479 END DO
480 END IF
481 END IF
482!
483!-----------------------------------------------------------------------
484! Lateral boundary conditions at the northern edge.
485!-----------------------------------------------------------------------
486!
487 IF (domain(ng)%Northern_Edge(tile)) THEN
488!
489! Northern edge, implicit upstream radiation condition.
490!
491 IF (lbc(inorth,ismtke,ng)%radiation) THEN
492 DO k=0,n(ng)
493 DO i=istr,iend+1
494 grad(i,jend )=tke(i ,jend ,k,nstp)- &
495 & tke(i-1,jend ,k,nstp)
496# ifdef MASKING
497 grad(i,jend )=grad(i,jend )* &
498 & grid(ng)%umask(i,jend )
499# endif
500 grad(i,jend+1)=tke(i ,jend+1,k,nstp)- &
501 & tke(i-1,jend+1,k,nstp)
502# ifdef MASKING
503 grad(i,jend+1)=grad(i,jend+1)* &
504 & grid(ng)%umask(i,jend+1)
505# endif
506 gradl(i,jend )=gls(i ,jend ,k,nstp)- &
507 & gls(i-1,jend ,k,nstp)
508# ifdef MASKING
509 gradl(i,jend )=gradl(i,jend )* &
510 & grid(ng)%umask(i,jend )
511# endif
512 gradl(i,jend+1)=gls(i ,jend+1,k,nstp)- &
513 & gls(i-1,jend+1,k,nstp)
514# ifdef MASKING
515 gradl(i,jend+1)=gradl(i,jend+1)* &
516 & grid(ng)%umask(i,jend+1)
517# endif
518 END DO
519 DO i=istr,iend
520 IF (lbc_apply(ng)%north(i)) THEN
521 dkdt=tke(i,jend,k,nstp)-tke(i,jend ,k,nout)
522 dkde=tke(i,jend,k,nstp)-tke(i,jend-1,k,nstp)
523 IF ((dkdt*dkde).lt.0.0_r8) dkdt=0.0_r8
524 IF ((dkdt*(grad(i ,jend)+ &
525 & grad(i+1,jend))).gt.0.0_r8) THEN
526 dkdx=grad(i ,jend)
527 ELSE
528 dkdx=grad(i+1,jend)
529 END IF
530 cff=dkdt/max(dkdx*dkdx+dkde*dkde,eps)
531# ifdef RADIATION_2D
532 cx=min(1.0_r8,max(cff*dkdx,-1.0_r8))
533# else
534 cx=0.0_r8
535# endif
536 ce=min(1.0_r8,cff*dkde)
537 tke(i,jend+1,k,nout)=(1.0_r8-ce)*tke(i,jend+1,k,nstp)+ &
538 & ce*tke(i,jend,k,nstp)- &
539 & max(cx,0.0_r8)*grad(i ,jend+1)- &
540 & min(cx,0.0_r8)*grad(i+1,jend+1)
541# ifdef MASKING
542 tke(i,jend+1,k,nout)=tke(i,jend+1,k,nout)* &
543 & grid(ng)%rmask(i,jend+1)
544# endif
545 dkdt=gls(i,jend,k,nstp)-gls(i,jend ,k,nout)
546 dkde=gls(i,jend,k,nstp)-gls(i,jend-1,k,nstp)
547 IF ((dkdt*dkde).lt.0.0_r8) dkdt=0.0_r8
548 IF ((dkdt*(gradl(i ,jend)+ &
549 & gradl(i+1,jend))).gt.0.0_r8) THEN
550 dkdx=gradl(i ,jend)
551 ELSE
552 dkdx=gradl(i+1,jend)
553 END IF
554 cff=dkdt/max(dkdx*dkdx+dkde*dkde,eps)
555# ifdef RADIATION_2D
556 cx=min(1.0_r8,max(cff*dkdx,-1.0_r8))
557# else
558 cx=0.0_r8
559# endif
560 ce=min(1.0_r8,cff*dkde)
561 gls(i,jend+1,k,nout)=(1.0_r8-ce)*gls(i,jend+1,k,nstp)+ &
562 & ce*gls(i,jend,k,nstp)- &
563 & max(cx,0.0_r8)*gradl(i ,jend+1)- &
564 & min(cx,0.0_r8)*gradl(i+1,jend+1)
565# ifdef MASKING
566 gls(i,jend+1,k,nout)=gls(i,jend+1,k,nout)* &
567 & grid(ng)%rmask(i,jend+1)
568# endif
569 END IF
570 END DO
571 END DO
572!
573! Northern edge, gradient boundary condition.
574!
575 ELSE IF (lbc(inorth,ismtke,ng)%gradient) THEN
576 DO k=0,n(ng)
577 DO i=istr,iend
578 IF (lbc_apply(ng)%north(i)) THEN
579 tke(i,jend+1,k,nout)=tke(i,jend,k,nout)
580# ifdef MASKING
581 tke(i,jend+1,k,nout)=tke(i,jend+1,k,nout)* &
582 & grid(ng)%rmask(i,jend+1)
583# endif
584 gls(i,jend+1,k,nout)=gls(i,jend,k,nout)
585# ifdef MASKING
586 gls(i,jend+1,k,nout)=gls(i,jend+1,k,nout)* &
587 & grid(ng)%rmask(i,jend+1)
588# endif
589 END IF
590 END DO
591 END DO
592!
593! Northern edge, closed boundary condition.
594!
595 ELSE IF (lbc(inorth,ismtke,ng)%closed) THEN
596 DO k=0,n(ng)
597 DO i=istr,iend
598 IF (lbc_apply(ng)%north(i)) THEN
599 tke(i,jend+1,k,nout)=tke(i,jend,k,nout)
600# ifdef MASKING
601 tke(i,jend+1,k,nout)=tke(i,jend+1,k,nout)* &
602 & grid(ng)%rmask(i,jend+1)
603# endif
604 gls(i,jend+1,k,nout)=gls(i,jend,k,nout)
605# ifdef MASKING
606 gls(i,jend+1,k,nout)=gls(i,jend+1,k,nout)* &
607 & grid(ng)%rmask(i,jend+1)
608# endif
609 END IF
610 END DO
611 END DO
612 END IF
613 END IF
614!
615!-----------------------------------------------------------------------
616! Boundary corners.
617!-----------------------------------------------------------------------
618!
619 IF (.not.(ewperiodic(ng).or.nsperiodic(ng))) THEN
620 IF (domain(ng)%SouthWest_Corner(tile)) THEN
621 IF (lbc_apply(ng)%south(istr-1).and. &
622 & lbc_apply(ng)%west (jstr-1)) THEN
623 DO k=0,n(ng)
624 tke(istr-1,jstr-1,k,nout)=0.5_r8* &
625 & (tke(istr ,jstr-1,k,nout)+ &
626 & tke(istr-1,jstr ,k,nout))
627 gls(istr-1,jstr-1,k,nout)=0.5_r8* &
628 & (gls(istr ,jstr-1,k,nout)+ &
629 & gls(istr-1,jstr ,k,nout))
630 END DO
631 END IF
632 END IF
633 IF (domain(ng)%SouthEast_Corner(tile)) THEN
634 IF (lbc_apply(ng)%south(iend+1).and. &
635 & lbc_apply(ng)%east (jstr-1)) THEN
636 DO k=0,n(ng)
637 tke(iend+1,jstr-1,k,nout)=0.5_r8* &
638 & (tke(iend ,jstr-1,k,nout)+ &
639 & tke(iend+1,jstr ,k,nout))
640 gls(iend+1,jstr-1,k,nout)=0.5_r8* &
641 & (gls(iend ,jstr-1,k,nout)+ &
642 & gls(iend+1,jstr ,k,nout))
643 END DO
644 END IF
645 END IF
646 IF (domain(ng)%NorthWest_Corner(tile)) THEN
647 IF (lbc_apply(ng)%north(istr-1).and. &
648 & lbc_apply(ng)%west (jend+1)) THEN
649 DO k=0,n(ng)
650 tke(istr-1,jend+1,k,nout)=0.5_r8* &
651 & (tke(istr ,jend+1,k,nout)+ &
652 & tke(istr-1,jend ,k,nout))
653 gls(istr-1,jend+1,k,nout)=0.5_r8* &
654 & (gls(istr ,jend+1,k,nout)+ &
655 & gls(istr-1,jend ,k,nout))
656 END DO
657 END IF
658 END IF
659 IF (domain(ng)%NorthEast_Corner(tile)) THEN
660 IF (lbc_apply(ng)%north(iend+1).and. &
661 & lbc_apply(ng)%east (jend+1)) THEN
662 DO k=0,n(ng)
663 tke(iend+1,jend+1,k,nout)=0.5_r8* &
664 & (tke(iend ,jend+1,k,nout)+ &
665 & tke(iend+1,jend ,k,nout))
666 gls(iend+1,jend+1,k,nout)=0.5_r8* &
667 & (gls(iend ,jend+1,k,nout)+ &
668 & gls(iend+1,jend ,k,nout))
669 END DO
670 END IF
671 END IF
672 END IF
673
674 RETURN
675 END SUBROUTINE tkebc_tile
676#endif
677 END MODULE tkebc_mod
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
integer ismtke
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_lbc), dimension(:,:,:), allocatable lbc
Definition mod_param.F:375
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
integer, parameter isouth
integer, parameter ieast
integer, parameter inorth
integer, dimension(:), allocatable nstp
subroutine tkebc(ng, tile, nout)
Definition tkebc_im.F:27
subroutine, public tkebc_tile(ng, tile, lbi, ubi, lbj, ubj, ubk, imins, imaxs, jmins, jmaxs, nout, nstp, gls, tke)
Definition tkebc_im.F:56