25 SUBROUTINE u3dbc (ng, tile, nout)
34 integer,
intent(in) :: ng, tile, nout
41 & lbi, ubi, lbj, ubj,
n(ng), &
42 & imins, imaxs, jmins, jmaxs, &
51 & LBi, UBi, LBj, UBj, UBk, &
52 & IminS, ImaxS, JminS, JmaxS, &
66 integer,
intent(in) :: ng, tile
67 integer,
intent(in) :: LBi, UBi, LBj, UBj, UBk
68 integer,
intent(in) :: IminS, ImaxS, JminS, JmaxS
69 integer,
intent(in) :: nstp, nout
72 real(r8),
intent(inout) :: u(LBi:,LBj:,:,:)
74 real(r8),
intent(inout) :: u(LBi:UBi,LBj:UBj,UBk,2)
82 real(r8),
parameter :: eps = 1.0e-20_r8
84 real(r8) :: Ce, Cx, cff, dUde, dUdt, dUdx
85 real(r8) :: obc_in, obc_out, tau
87 real(r8),
dimension(IminS:ImaxS,JminS:JmaxS) :: grad
89# include "set_bounds.h"
95 IF (
domain(ng)%Western_Edge(tile))
THEN
102 grad(istr ,j)=u(istr ,j ,k,nstp)- &
103 & u(istr ,j-1,k,nstp)
104 grad(istr+1,j)=u(istr+1,j ,k,nstp)- &
105 & u(istr+1,j-1,k,nstp)
109 dudt=u(istr+1,j,k,nstp)-u(istr+1,j,k,nout)
110 dudx=u(istr+1,j,k,nstp)-u(istr+2,j,k,nstp)
115 & (
clima(ng)%M3nudgcof(istr-1,j,k)+ &
116 &
clima(ng)%M3nudgcof(istr ,j,k))
117 obc_in =
obcfac(ng)*obc_out
122 IF ((dudt*dudx).lt.0.0_r8)
THEN
130 IF ((dudt*dudx).lt.0.0_r8) dudt=0.0_r8
131 IF ((dudt*(grad(istr+1,j )+ &
132 & grad(istr+1,j+1))).gt.0.0_r8)
THEN
135 dude=grad(istr+1,j+1)
137 cff=dudt/max(dudx*dudx+dude*dude,eps)
138 cx=min(1.0_r8,cff*dudx)
140 ce=min(1.0_r8,max(cff*dude,-1.0_r8))
144# if defined CELERITY_WRITE && defined FORWARD_WRITE
149 u(istr,j,k,nout)=(1.0_r8-cx)*u(istr,j,k,nstp)+ &
150 & cx*u(istr+1,j,k,nstp)- &
151 & max(ce,0.0_r8)*grad(istr,j )- &
152 & min(ce,0.0_r8)*grad(istr,j+1)
155 u(istr,j,k,nout)=u(istr,j,k,nout)+ &
160 u(istr,j,k,nout)=u(istr,j,k,nout)* &
161 &
grid(ng)%umask(istr,j)
164 u(istr,j,k,nout)=u(istr,j,k,nout)* &
165 &
grid(ng)%umask_wet(istr,j)
177 u(istr,j,k,nout)=
boundary(ng)%u_west(j,k)
179 u(istr,j,k,nout)=u(istr,j,k,nout)* &
180 &
grid(ng)%umask(istr,j)
183 u(istr,j,k,nout)=u(istr,j,k,nout)* &
184 &
grid(ng)%umask_wet(istr,j)
196 u(istr,j,k,nout)=u(istr+1,j,k,nout)
198 u(istr,j,k,nout)=u(istr,j,k,nout)* &
199 &
grid(ng)%umask(istr,j)
202 u(istr,j,k,nout)=u(istr,j,k,nout)* &
203 &
grid(ng)%umask_wet(istr,j)
215 u(istr,j,k,nout)=0.0_r8
226 IF (
domain(ng)%Eastern_Edge(tile))
THEN
233 grad(iend ,j)=u(iend ,j ,k,nstp)- &
234 & u(iend ,j-1,k,nstp)
235 grad(iend+1,j)=u(iend+1,j ,k,nstp)- &
236 & u(iend+1,j-1,k,nstp)
240 dudt=u(iend,j,k,nstp)-u(iend ,j,k,nout)
241 dudx=u(iend,j,k,nstp)-u(iend-1,j,k,nstp)
246 & (
clima(ng)%M3nudgcof(iend ,j,k)+ &
247 &
clima(ng)%M3nudgcof(iend+1,j,k))
248 obc_in =
obcfac(ng)*obc_out
253 IF ((dudt*dudx).lt.0.0_r8)
THEN
261 IF ((dudt*dudx).lt.0.0_r8) dudt=0.0_r8
262 IF ((dudt*(grad(iend,j )+ &
263 & grad(iend,j+1))).gt.0.0_r8)
THEN
268 cff=dudt/max(dudx*dudx+dude*dude,eps)
269 cx=min(1.0_r8,cff*dudx)
271 ce=min(1.0_r8,max(cff*dude,-1.0_r8))
275# if defined CELERITY_WRITE && defined FORWARD_WRITE
280 u(iend+1,j,k,nout)=(1.0_r8-cx)*u(iend+1,j,k,nstp)+ &
281 & cx*u(iend,j,k,nstp)- &
282 & max(ce,0.0_r8)*grad(iend+1,j )- &
283 & min(ce,0.0_r8)*grad(iend+1,j+1)
286 u(iend+1,j,k,nout)=u(iend+1,j,k,nout)+ &
288 & u(iend+1,j,k,nstp))
291 u(iend+1,j,k,nout)=u(iend+1,j,k,nout)* &
292 &
grid(ng)%umask(iend+1,j)
295 u(iend+1,j,k,nout)=u(iend+1,j,k,nout)* &
296 &
grid(ng)%umask_wet(iend+1,j)
308 u(iend+1,j,k,nout)=
boundary(ng)%u_east(j,k)
310 u(iend+1,j,k,nout)=u(iend+1,j,k,nout)* &
311 &
grid(ng)%umask(iend+1,j)
314 u(iend+1,j,k,nout)=u(iend+1,j,k,nout)* &
315 &
grid(ng)%umask_wet(iend+1,j)
327 u(iend+1,j,k,nout)=u(iend,j,k,nout)
329 u(iend+1,j,k,nout)=u(iend+1,j,k,nout)* &
330 &
grid(ng)%umask(iend+1,j)
333 u(iend+1,j,k,nout)=u(iend+1,j,k,nout)* &
334 &
grid(ng)%umask_wet(iend+1,j)
346 u(iend+1,j,k,nout)=0.0_r8
357 IF (
domain(ng)%Southern_Edge(tile))
THEN
364 grad(i,jstr-1)=u(i+1,jstr-1,k,nstp)- &
365 & u(i ,jstr-1,k,nstp)
366 grad(i,jstr )=u(i+1,jstr ,k,nstp)- &
371 dudt=u(i,jstr,k,nstp)-u(i,jstr ,k,nout)
372 dude=u(i,jstr,k,nstp)-u(i,jstr+1,k,nstp)
377 & (
clima(ng)%M3nudgcof(i-1,jstr-1,k)+ &
378 &
clima(ng)%M3nudgcof(i ,jstr-1,k))
379 obc_in =
obcfac(ng)*obc_out
384 IF ((dudt*dude).lt.0.0_r8)
THEN
392 IF ((dudt*dude).lt.0.0_r8) dudt=0.0_r8
393 IF ((dudt*(grad(i-1,jstr)+ &
394 & grad(i ,jstr))).gt.0.0_r8)
THEN
399 cff=dudt/max(dudx*dudx+dude*dude,eps)
401 cx=min(1.0_r8,max(cff*dudx,-1.0_r8))
405 ce=min(1.0_r8,cff*dude)
406# if defined CELERITY_WRITE && defined FORWARD_WRITE
411 u(i,jstr-1,k,nout)=(1.0_r8-ce)*u(i,jstr-1,k,nstp)+ &
412 & ce*u(i,jstr,k,nstp)- &
413 & max(cx,0.0_r8)*grad(i-1,jstr-1)- &
414 & min(cx,0.0_r8)*grad(i ,jstr-1)
417 u(i,jstr-1,k,nout)=u(i,jstr-1,k,nout)+ &
419 & u(i,jstr-1,k,nstp))
422 u(i,jstr-1,k,nout)=u(i,jstr-1,k,nout)* &
423 &
grid(ng)%umask(i,jstr-1)
426 u(i,jstr-1,k,nout)=u(i,jstr-1,k,nout)* &
427 &
grid(ng)%umask_wet(i,jstr-1)
439 u(i,jstr-1,k,nout)=
boundary(ng)%u_south(i,k)
441 u(i,jstr-1,k,nout)=u(i,jstr-1,k,nout)* &
442 &
grid(ng)%umask(i,jstr-1)
445 u(i,jstr-1,k,nout)=u(i,jstr-1,k,nout)* &
446 &
grid(ng)%umask_wet(i,jstr-1)
458 u(i,jstr-1,k,nout)=u(i,jstr,k,nout)
460 u(i,jstr-1,k,nout)=u(i,jstr-1,k,nout)* &
461 &
grid(ng)%umask(i,jstr-1)
464 u(i,jstr-1,k,nout)=u(i,jstr-1,k,nout)* &
465 &
grid(ng)%umask_wet(i,jstr-1)
485 u(i,jstr-1,k,nout)=
gamma2(ng)*u(i,jstr,k,nout)
487 u(i,jstr-1,k,nout)=u(i,jstr-1,k,nout)* &
488 &
grid(ng)%umask(i,jstr-1)
491 u(i,jstr-1,k,nout)=u(i,jstr-1,k,nout)* &
492 &
grid(ng)%umask_wet(i,jstr-1)
504 IF (
domain(ng)%Northern_Edge(tile))
THEN
511 grad(i,jend )=u(i+1,jend ,k,nstp)- &
513 grad(i,jend+1)=u(i+1,jend+1,k,nstp)- &
514 & u(i ,jend+1,k,nstp)
518 dudt=u(i,jend,k,nstp)-u(i,jend ,k,nout)
519 dude=u(i,jend,k,nstp)-u(i,jend-1,k,nstp)
524 & (
clima(ng)%M3nudgcof(i-1,jend+1,k)+ &
525 &
clima(ng)%M3nudgcof(i ,jend+1,k))
526 obc_in =
obcfac(ng)*obc_out
531 IF ((dudt*dude).lt.0.0_r8)
THEN
539 IF ((dudt*dude).lt.0.0_r8) dudt=0.0_r8
540 IF ((dudt*(grad(i-1,jend)+ &
541 & grad(i ,jend))).gt.0.0_r8)
THEN
546 cff=dudt/max(dudx*dudx+dude*dude,eps)
548 cx=min(1.0_r8,max(cff*dudx,-1.0_r8))
552 ce=min(1.0_r8,cff*dude)
553# if defined CELERITY_WRITE && defined FORWARD_WRITE
558 u(i,jend+1,k,nout)=(1.0_r8-ce)*u(i,jend+1,k,nstp)+ &
559 & ce*u(i,jend,k,nstp)- &
560 & max(cx,0.0_r8)*grad(i-1,jend+1)- &
561 & min(cx,0.0_r8)*grad(i ,jend+1)
564 u(i,jend+1,k,nout)=u(i,jend+1,k,nout)+ &
566 & u(i,jend+1,k,nstp))
569 u(i,jend+1,k,nout)=u(i,jend+1,k,nout)* &
570 &
grid(ng)%umask(i,jend+1)
573 u(i,jend+1,k,nout)=u(i,jend+1,k,nout)* &
574 &
grid(ng)%umask_wet(i,jend+1)
586 u(i,jend+1,k,nout)=
boundary(ng)%u_north(i,k)
588 u(i,jend+1,k,nout)=u(i,jend+1,k,nout)* &
589 &
grid(ng)%umask(i,jend+1)
592 u(i,jend+1,k,nout)=u(i,jend+1,k,nout)* &
593 &
grid(ng)%umask_wet(i,jend+1)
605 u(i,jend+1,k,nout)=u(i,jend,k,nout)
607 u(i,jend+1,k,nout)=u(i,jend+1,k,nout)* &
608 &
grid(ng)%umask(i,jend+1)
611 u(i,jend+1,k,nout)=u(i,jend+1,k,nout)* &
612 &
grid(ng)%umask_wet(i,jend+1)
632 u(i,jend+1,k,nout)=
gamma2(ng)*u(i,jend,k,nout)
634 u(i,jend+1,k,nout)=u(i,jend+1,k,nout)* &
635 &
grid(ng)%umask(i,jend+1)
638 u(i,jend+1,k,nout)=u(i,jend+1,k,nout)* &
639 &
grid(ng)%umask_wet(i,jend+1)
652 IF (
domain(ng)%SouthWest_Corner(tile))
THEN
656 u(istr,jstr-1,k,nout)=0.5_r8*(u(istr+1,jstr-1,k,nout)+ &
657 & u(istr ,jstr ,k,nout))
661 IF (
domain(ng)%SouthEast_Corner(tile))
THEN
665 u(iend+1,jstr-1,k,nout)=0.5_r8*(u(iend ,jstr-1,k,nout)+ &
666 & u(iend+1,jstr ,k,nout))
670 IF (
domain(ng)%NorthWest_Corner(tile))
THEN
674 u(istr,jend+1,k,nout)=0.5_r8*(u(istr ,jend ,k,nout)+ &
675 & u(istr+1,jend+1,k,nout))
679 IF (
domain(ng)%NorthEast_Corner(tile))
THEN
683 u(iend+1,jend+1,k,nout)=0.5_r8*(u(iend+1,jend ,k,nout)+ &
684 & u(iend ,jend+1,k,nout))
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_ocean), dimension(:), allocatable ocean
integer, dimension(:), allocatable n
type(t_lbc), dimension(:,:,:), allocatable lbc
type(t_domain), dimension(:), allocatable domain
real(dp), dimension(:,:), allocatable m3obc_out
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
integer, dimension(:), allocatable nstp
subroutine u3dbc(ng, tile, nout)
subroutine, public u3dbc_tile(ng, tile, lbi, ubi, lbj, ubj, ubk, imins, imaxs, jmins, jmaxs, nstp, nout, u)