58
59
66
67
68
69 integer, intent(in) :: ng, tile
70 integer, intent(in) :: LBi, UBi, LBj, UBj, UBk
71 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
72 integer, intent(in) :: nstp, nout
73
74# ifdef ASSUMED_SHAPE
75 real(r8), intent(inout) :: tl_u(LBi:,LBj:,:,:)
76# else
77 real(r8), intent(inout) :: tl_u(LBi:UBi,LBj:UBj,UBk,2)
78# endif
79
80
81
82 integer :: Imin, Imax
83 integer :: i, j, k
84
85 real(r8) :: Ce, Cx, cff
86 real(r8) :: obc_in, obc_out, tau
87
88 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_grad
89
90# include "set_bounds.h"
91
92
93
94
95
96 IF (
domain(ng)%Western_Edge(tile))
THEN
97
98
99
101 IF (
iic(ng).ne.0)
THEN
103 DO j=jstr,jend+1
104
105
106
107 tl_grad(istr,j)=0.0_r8
108 END DO
109 DO j=jstr,jend
111# if defined CELERITY_READ && defined FORWARD_READ
114 obc_out=0.5_r8* &
115 & (
clima(ng)%M3nudgcof(istr-1,j,k)+ &
116 &
clima(ng)%M3nudgcof(istr ,j,k))
117 obc_in =
obcfac(ng)*obc_out
118 ELSE
121 END IF
122 IF (
boundary(ng)%u_west_Cx(j,k).lt.0.0_r8)
THEN
123 tau=obc_in
124 ELSE
125 tau=obc_out
126 END IF
128 END IF
130# ifdef RADIATION_2D
132# else
133 ce=0.0_r8
134# endif
136# endif
137
138
139
140
141
142
143 tl_u(istr,j,k,nout)=(cff*tl_u(istr ,j,k,nstp)+ &
144 & cx *tl_u(istr+1,j,k,nout)- &
145 & max(ce,0.0_r8)* &
146 & tl_grad(istr,j )- &
147 & min(ce,0.0_r8)* &
148 & tl_grad(istr,j+1))/ &
149 & (cff+cx)
150
152
153
154
155
156 tl_u(istr,j,k,nout)=tl_u(istr,j,k,nout)- &
157 & tau*tl_u(istr,j,k,nstp)
158 END IF
159# ifdef MASKING
160
161
162
163 tl_u(istr,j,k,nout)=tl_u(istr,j,k,nout)* &
164 &
grid(ng)%umask(istr,j)
165# endif
166 END IF
167 END DO
168 END DO
169 END IF
170
171
172
175 DO j=jstr,jend
177
178
179 tl_u(istr,j,k,nout)=
boundary(ng)%tl_u_west(j,k)
180# ifdef MASKING
181
182
183
184 tl_u(istr,j,k,nout)=tl_u(istr,j,k,nout)* &
185 &
grid(ng)%umask(istr,j)
186# endif
187 END IF
188 END DO
189 END DO
190
191
192
195 DO j=jstr,jend
197
198
199 tl_u(istr,j,k,nout)=tl_u(istr+1,j,k,nout)
200# ifdef MASKING
201
202
203
204 tl_u(istr,j,k,nout)=tl_u(istr,j,k,nout)* &
205 &
grid(ng)%umask(istr,j)
206# endif
207 END IF
208 END DO
209 END DO
210
211
212
215 DO j=jstr,jend
217
218
219 tl_u(istr,j,k,nout)=0.0_r8
220 END IF
221 END DO
222 END DO
223 END IF
224 END IF
225
226
227
228
229
230 IF (
domain(ng)%Eastern_Edge(tile))
THEN
231
232
233
235 IF (
iic(ng).ne.0)
THEN
237 DO j=jstr,jend+1
238
239
240
241 tl_grad(iend+1,j)=0.0_r8
242 END DO
243 DO j=jstr,jend
245# if defined CELERITY_READ && defined FORWARD_READ
248 obc_out=0.5_r8* &
249 & (
clima(ng)%M3nudgcof(iend ,j,k)+ &
250 &
clima(ng)%M3nudgcof(iend+1,j,k))
251 obc_in =
obcfac(ng)*obc_out
252 ELSE
255 END IF
256 IF (
boundary(ng)%u_east_Cx(j,k).lt.0.0_r8)
THEN
257 tau=obc_in
258 ELSE
259 tau=obc_out
260 END IF
262 END IF
264# ifdef RADIATION_2D
266# else
267 ce=0.0_r8
268# endif
270# endif
271
272
273
274
275
276
277 tl_u(iend+1,j,k,nout)=(cff*tl_u(iend+1,j,k,nstp)+ &
278 & cx *tl_u(iend ,j,k,nout)- &
279 & max(ce,0.0_r8)* &
280 & tl_grad(iend+1,j )- &
281 & min(ce,0.0_r8)* &
282 & tl_grad(iend+1,j+1))/ &
283 & (cff+cx)
284
286
287
288
289
290 tl_u(iend+1,j,k,nout)=tl_u(iend+1,j,k,nout)- &
291 & tau*tl_u(iend+1,j,k,nstp)
292 END IF
293# ifdef MASKING
294
295
296
297 tl_u(iend+1,j,k,nout)=tl_u(iend+1,j,k,nout)* &
298 &
grid(ng)%umask(iend+1,j)
299# endif
300 END IF
301 END DO
302 END DO
303 END IF
304
305
306
309 DO j=jstr,jend
311
312
313 tl_u(iend+1,j,k,nout)=
boundary(ng)%tl_u_east(j,k)
314# ifdef MASKING
315
316
317
318 tl_u(iend+1,j,k,nout)=tl_u(iend+1,j,k,nout)* &
319 &
grid(ng)%umask(iend+1,j)
320# endif
321 END IF
322 END DO
323 END DO
324
325
326
329 DO j=jstr,jend
331
332
333 tl_u(iend+1,j,k,nout)=tl_u(iend,j,k,nout)
334# ifdef MASKING
335
336
337
338 tl_u(iend+1,j,k,nout)=tl_u(iend+1,j,k,nout)* &
339 &
grid(ng)%umask(iend+1,j)
340# endif
341 END IF
342 END DO
343 END DO
344
345
346
349 DO j=jstr,jend
351
352
353 tl_u(iend+1,j,k,nout)=0.0_r8
354 END IF
355 END DO
356 END DO
357 END IF
358 END IF
359
360
361
362
363
364 IF (
domain(ng)%Southern_Edge(tile))
THEN
365
366
367
369 IF (
iic(ng).ne.0)
THEN
371 DO i=istru-1,iend
372
373
374
375 tl_grad(i,jstr-1)=0.0_r8
376 END DO
377 DO i=istru,iend
379# if defined CELERITY_READ && defined FORWARD_READ
382 obc_out=0.5_r8* &
383 & (
clima(ng)%M3nudgcof(i-1,jstr-1,k)+ &
384 &
clima(ng)%M3nudgcof(i ,jstr-1,k))
385 obc_in =
obcfac(ng)*obc_out
386 ELSE
389 END IF
390 IF (
boundary(ng)%u_south_Ce(i,k).lt.0.0_r8)
THEN
391 tau=obc_in
392 ELSE
393 tau=obc_out
394 END IF
396 END IF
397# ifdef RADIATION_2D
399# else
400 cx=0.0_r8
401# endif
404# endif
405
406
407
408
409
410
411 tl_u(i,jstr-1,k,nout)=(cff*tl_u(i,jstr-1,k,nstp)+ &
412 & ce *tl_u(i,jstr ,k,nout)- &
413 & max(cx,0.0_r8)* &
414 & tl_grad(i-1,jstr-1)- &
415 & min(cx,0.0_r8)* &
416 & tl_grad(i ,jstr-1))/ &
417 & (cff+ce)
418
420
421
422
423
424 tl_u(i,jstr-1,k,nout)=tl_u(i,jstr-1,k,nout)- &
425 & tau*tl_u(i,jstr-1,k,nstp)
426 END IF
427# ifdef MASKING
428
429
430
431 tl_u(i,jstr-1,k,nout)=tl_u(i,jstr-1,k,nout)* &
432 &
grid(ng)%umask(i,jstr-1)
433# endif
434 END IF
435 END DO
436 END DO
437 END IF
438
439
440
443 DO i=istru,iend
445
446
447 tl_u(i,jstr-1,k,nout)=
boundary(ng)%tl_u_south(i,k)
448# ifdef MASKING
449
450
451
452 tl_u(i,jstr-1,k,nout)=tl_u(i,jstr-1,k,nout)* &
453 &
grid(ng)%umask(i,jstr-1)
454# endif
455 END IF
456 END DO
457 END DO
458
459
460
463 DO i=istru,iend
465
466
467 tl_u(i,jstr-1,k,nout)=tl_u(i,jstr,k,nout)
468# ifdef MASKING
469
470
471
472 tl_u(i,jstr-1,k,nout)=tl_u(i,jstr-1,k,nout)* &
473 &
grid(ng)%umask(i,jstr-1)
474# endif
475 END IF
476 END DO
477 END DO
478
479
480
481
484 imin=istru
485 imax=iend
486 ELSE
487 imin=istr
488 imax=iendr
489 END IF
491 DO i=imin,imax
493
494
495 tl_u(i,jstr-1,k,nout)=
gamma2(ng)*tl_u(i,jstr,k,nout)
496# ifdef MASKING
497
498
499
500 tl_u(i,jstr-1,k,nout)=tl_u(i,jstr-1,k,nout)* &
501 &
grid(ng)%umask(i,jstr-1)
502# endif
503 END IF
504 END DO
505 END DO
506 END IF
507 END IF
508
509
510
511
512
513 IF (
domain(ng)%Northern_Edge(tile))
THEN
514
515
516
518 IF (
iic(ng).ne.0)
THEN
520 DO i=istru-1,iend
521
522
523
524 tl_grad(i,jend+1)=0.0_r8
525 END DO
526 DO i=istru,iend
528# if defined CELERITY_READ && defined FORWARD_READ
531 obc_out=0.5_r8* &
532 & (
clima(ng)%M3nudgcof(i-1,jend+1,k)+ &
533 &
clima(ng)%M3nudgcof(i ,jend+1,k))
534 obc_in =
obcfac(ng)*obc_out
535 ELSE
538 END IF
539 IF (
boundary(ng)%u_north_Ce(i,k).lt.0.0_r8)
THEN
540 tau=obc_in
541 ELSE
542 tau=obc_out
543 END IF
545 END IF
546# ifdef RADIATION_2D
548# else
549 cx=0.0_r8
550# endif
553# endif
554
555
556
557
558
559
560 tl_u(i,jend+1,k,nout)=(cff*tl_u(i,jend+1,k,nstp)+ &
561 & ce *tl_u(i,jend ,k,nout)- &
562 & max(cx,0.0_r8)* &
563 & tl_grad(i-1,jend+1)- &
564 & min(cx,0.0_r8)* &
565 & tl_grad(i ,jend+1))/ &
566 & (cff+ce)
567
569
570
571
572
573 tl_u(i,jend+1,k,nout)=tl_u(i,jend+1,k,nout)- &
574 & tau*tl_u(i,jend+1,k,nstp)
575 END IF
576# ifdef MASKING
577
578
579
580 tl_u(i,jend+1,k,nout)=tl_u(i,jend+1,k,nout)* &
581 &
grid(ng)%umask(i,jend+1)
582# endif
583 END IF
584 END DO
585 END DO
586 END IF
587
588
589
592 DO i=istru,iend
594
595
596 tl_u(i,jend+1,k,nout)=
boundary(ng)%tl_u_north(i,k)
597# ifdef MASKING
598
599
600
601 tl_u(i,jend+1,k,nout)=tl_u(i,jend+1,k,nout)* &
602 &
grid(ng)%umask(i,jend+1)
603# endif
604 END IF
605 END DO
606 END DO
607
608
609
612 DO i=istru,iend
614
615
616 tl_u(i,jend+1,k,nout)=tl_u(i,jend,k,nout)
617# ifdef MASKING
618
619
620
621 tl_u(i,jend+1,k,nout)=tl_u(i,jend+1,k,nout)* &
622 &
grid(ng)%umask(i,jend+1)
623# endif
624 END IF
625 END DO
626 END DO
627
628
629
630
633 imin=istru
634 imax=iend
635 ELSE
636 imin=istr
637 imax=iendr
638 END IF
640 DO i=imin,imax
642
643
644 tl_u(i,jend+1,k,nout)=
gamma2(ng)*tl_u(i,jend,k,nout)
645# ifdef MASKING
646
647
648
649 tl_u(i,jend+1,k,nout)=tl_u(i,jend+1,k,nout)* &
650 &
grid(ng)%umask(i,jend+1)
651# endif
652 END IF
653 END DO
654 END DO
655 END IF
656 END IF
657
658
659
660
661
663 IF (
domain(ng)%SouthWest_Corner(tile))
THEN
667
668
669
670 tl_u(istr,jstr-1,k,nout)=0.5_r8* &
671 & (tl_u(istr+1,jstr-1,k,nout)+ &
672 & tl_u(istr ,jstr ,k,nout))
673 END DO
674 END IF
675 END IF
676 IF (
domain(ng)%SouthEast_Corner(tile))
THEN
680
681
682
683 tl_u(iend+1,jstr-1,k,nout)=0.5_r8* &
684 & (tl_u(iend ,jstr-1,k,nout)+ &
685 & tl_u(iend+1,jstr ,k,nout))
686 END DO
687 END IF
688 END IF
689 IF (
domain(ng)%NorthWest_Corner(tile))
THEN
693
694
695
696 tl_u(istr,jend+1,k,nout)=0.5_r8* &
697 & (tl_u(istr ,jend ,k,nout)+ &
698 & tl_u(istr+1,jend+1,k,nout))
699 END DO
700 END IF
701 END IF
702 IF (
domain(ng)%NorthEast_Corner(tile))
THEN
706
707
708
709 tl_u(iend+1,jend+1,k,nout)=0.5_r8* &
710 & (tl_u(iend+1,jend ,k,nout)+ &
711 & tl_u(iend ,jend+1,k,nout))
712 END DO
713 END IF
714 END IF
715 END IF
716
717 RETURN
type(t_boundary), dimension(:), allocatable boundary
type(t_apply), dimension(:), allocatable lbc_apply
type(t_clima), dimension(:), allocatable clima
type(t_grid), dimension(:), allocatable grid
type(t_lbc), dimension(:,:,:), allocatable tl_lbc
type(t_domain), dimension(:), allocatable domain
real(dp), dimension(:,:), allocatable m3obc_out
integer, dimension(:), allocatable iic
real(dp), dimension(:), allocatable dt
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
real(dp), dimension(:), allocatable obcfac
real(r8), dimension(:), allocatable gamma2
logical, dimension(:), allocatable lnudgem3clm
integer, parameter isouth
integer, parameter inorth
real(dp), dimension(:,:), allocatable m3obc_in