30 SUBROUTINE ice_uibc (ng, tile, model)
37 integer,
intent(in) :: ng, tile, model
41 character (len=*),
parameter :: MyFile = &
47 CALL wclock_on (ng, model, 42, __line__, myfile)
49 CALL ice_uibc_tile (ng, tile, model, &
50 & lbi, ubi, lbj, ubj, &
51 & imins, imaxs, jmins, jmaxs, &
52 & liuol(ng), liunw(ng), &
55 CALL wclock_off (ng, model, 42, __line__, myfile)
59 END SUBROUTINE ice_uibc
62 SUBROUTINE ice_uibc_tile (ng, tile, model, &
63 & LBi, UBi, LBj, UBj, &
64 & IminS, ImaxS, JminS, JmaxS, &
71 integer,
intent(in) :: ng, tile, model
72 integer,
intent(in) :: LBi, UBi, LBj, UBj
73 integer,
intent(in) :: IminS, ImaxS, JminS, JmaxS
74 integer,
intent(in) :: liuol, liunw
77 real(r8),
intent(inout) :: ui(LBi:,LBj:,:)
79 real(r8),
intent(inout) :: ui(LBi:UBi,LBj:UBj,2)
84 integer :: i, Imax, Imin, j, know
86 real(r8),
parameter :: eps =1.0e-20_r8
87 real(r8) :: Ce, Cx, cff, dUde, dUdt, dUdx, tau
89 real(r8),
dimension(IminS:ImaxS,JminS:JmaxS) :: grad
91#include "set_bounds.h"
103 IF (
domain(ng)%Western_Edge(tile))
THEN
109 grad(istr ,j)=ui(istr ,j ,know)- &
111 grad(istr+1,j)=ui(istr+1,j ,know)- &
112 & ui(istr+1,j-1,know)
115 dudt=ui(istr+1,j,know )-ui(istr+1,j,liunw)
116 dudx=ui(istr+1,j,liunw)-ui(istr+2,j,liunw)
118 IF ((dudt*dudx).lt.0.0_r8)
THEN
125 IF ((dudt*dudx).lt.0.0_r8) dudt=0.0_r8
126 IF ((dudt*(grad(istr+1,j)+grad(istr+1,j+1))).gt.0.0_r8)
THEN
129 dude=grad(istr+1,j+1)
131 cff=max(dudx*dudx+dude*dude,eps)
134 ce=min(cff,max(dudt*dude,-cff))
138 ui(istr,j,liunw)=(cff*ui(istr ,j,know)+ &
139 & cx *ui(istr+1,j,liunw)- &
140 & max(ce,0.0_r8)*grad(istr,j )- &
141 & min(ce,0.0_r8)*grad(istr,j+1))/ &
144 ui(istr,j,liunw)=ui(istr,j,liunw)+ &
149 ui(istr,j,liunw)=ui(istr,j,liunw)* &
150 &
grid(ng)%umask(istr,j)
160 ui(istr,j,liunw)=ui(istr,j,liunw)* &
161 &
grid(ng)%umask(istr,j)
164 ui(istr,j,liunw)=ui(istr,j,liunw)* &
165 &
grid(ng)%umask_wet(istr,j)
173 ui(istr,j,liunw)=ui(istr+1,j,liunw)
175 ui(istr,j,liunw)=ui(istr,j,liunw)* &
176 &
grid(ng)%umask(istr,j)
179 ui(istr,j,liunw)=ui(istr,j,liunw)* &
180 &
grid(ng)%umask_wet(istr,j)
188 ui(istr,j,liunw)=0.0_r8
197 IF (
domain(ng)%Eastern_Edge(tile))
THEN
203 grad(iend ,j)=ui(iend ,j ,know)- &
205 grad(iend+1,j)=ui(iend+1,j ,know)- &
206 & ui(iend+1,j-1,know)
209 dudt=ui(iend,j,know )-ui(iend ,j,liunw)
210 dudx=ui(iend,j,liunw)-ui(iend-1,j,liunw)
212 IF ((dudt*dudx).lt.0.0_r8)
THEN
219 IF ((dudt*dudx).lt.0.0_r8) dudt=0.0_r8
220 IF ((dudt*(grad(iend,j)+grad(iend,j+1))).gt.0.0_r8)
THEN
225 cff=max(dudx*dudx+dude*dude,eps)
228 ce=min(cff,max(dudt*dude,-cff))
232 ui(iend+1,j,liunw)=(cff*ui(iend+1,j,know)+ &
233 & cx *ui(iend ,j,liunw)- &
234 & max(ce,0.0_r8)*grad(iend+1,j )- &
235 & min(ce,0.0_r8)*grad(iend+1,j+1))/ &
238 ui(iend+1,j,liunw)=ui(iend+1,j,liunw)+ &
243 ui(iend+1,j,liunw)=ui(iend+1,j,liunw)* &
244 &
grid(ng)%umask(iend+1,j)
254 ui(iend+1,j,liunw)=ui(iend+1,j,liunw)* &
255 &
grid(ng)%umask(iend+1,j)
258 ui(iend+1,j,liunw)=ui(iend+1,j,liunw)* &
259 &
grid(ng)%umask_wet(iend+1,j)
267 ui(iend+1,j,liunw)=ui(iend,j,liunw)
269 ui(iend+1,j,liunw)=ui(iend+1,j,liunw)* &
270 &
grid(ng)%umask(iend+1,j)
273 ui(iend+1,j,liunw)=ui(iend+1,j,liunw)* &
274 &
grid(ng)%umask_wet(iend+1,j)
282 ui(iend+1,j,liunw)=0.0_r8
291 IF (
domain(ng)%Southern_Edge(tile))
THEN
297 grad(i,jstr-1)=ui(i+1,jstr-1,know)- &
299 grad(i,jstr )=ui(i+1,jstr ,know)- &
303 dudt=ui(i,jstr,know )-ui(i,jstr ,liunw)
304 dude=ui(i,jstr,liunw)-ui(i,jstr+1,liunw)
306 IF ((dudt*dude).lt.0.0_r8)
THEN
313 IF ((dudt*dude).lt.0.0_r8) dudt=0.0_r8
314 IF ((dudt*(grad(i-1,jstr)+grad(i,jstr))).gt.0.0_r8)
THEN
319 cff=max(dudx*dudx+dude*dude,eps)
321 cx=min(cff,max(dudt*dudx,-cff))
326 ui(i,jstr-1,liunw)=(cff*ui(i,jstr-1,know)+ &
327 & ce *ui(i,jstr ,liunw)- &
328 & max(cx,0.0_r8)*grad(i-1,jstr-1)- &
329 & min(cx,0.0_r8)*grad(i ,jstr-1))/ &
332 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)+ &
337 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)* &
338 &
grid(ng)%umask(i,jstr-1)
348 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)* &
349 &
grid(ng)%umask(i,jstr-1)
352 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)* &
353 &
grid(ng)%umask_wet(i,jstr-1)
361 ui(i,jstr-1,liunw)=ui(i,1,liunw)
363 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)* &
364 &
grid(ng)%umask(i,jstr-1)
367 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)* &
368 &
grid(ng)%umask_wet(i,jstr-1)
384 ui(i,jstr-1,liunw)=
gamma2(ng)*ui(i,1,liunw)
386 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)* &
387 &
grid(ng)%umask(i,jstr-1)
390 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)* &
391 &
grid(ng)%umask_wet(i,jstr-1)
401 IF (
domain(ng)%Northern_Edge(tile))
THEN
407 grad(i,jend )=ui(i+1,jend ,know)- &
409 grad(i,jend+1)=ui(i+1,jend+1,know)- &
413 dudt=ui(i,jend,know )-ui(i,jend ,liunw)
414 dude=ui(i,jend,liunw)-ui(i,jend-1,liunw)
416 IF ((dudt*dude).lt.0.0_r8)
THEN
423 IF ((dudt*dude).lt.0.0_r8) dudt=0.0_r8
424 IF ((dudt*(grad(i-1,jend)+grad(i,jend))).gt.0.0_r8)
THEN
429 cff=max(dudx*dudx+dude*dude,eps)
431 cx=min(cff,max(dudt*dudx,-cff))
436 ui(i,jend+1,liunw)=(cff*ui(i,jend+1,know)+ &
437 & ce *ui(i,jend ,liunw)- &
438 & max(cx,0.0_r8)*grad(i-1,jend+1)- &
439 & min(cx,0.0_r8)*grad(i ,jend+1))/ &
441# ifdef NORTH_MINUDGING
443 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)+ &
448 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)* &
449 &
grid(ng)%umask(i,jend+1)
459 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)* &
460 &
grid(ng)%umask(i,jend+1)
463 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)* &
464 &
grid(ng)%umask_wet(i,jend+1)
472 ui(i,jend+1,liunw)=ui(i,jend,liunw)
474 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)* &
475 &
grid(ng)%umask(i,jend+1)
478 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)* &
479 &
grid(ng)%umask_wet(i,jend+1)
495 ui(i,jend+1,liunw)=
gamma2(ng)*ui(i,jend,liunw)
497 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)* &
498 &
grid(ng)%umask(i,jend+1)
501 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)* &
502 &
grid(ng)%umask_wet(i,jend+1)
513 IF (
domain(ng)%SouthWest_Corner(tile))
THEN
514 ui(istr,jstr-1,liunw)=0.5_r8*(ui(istr+1,jstr-1,liunw)+ &
515 & ui(istr ,jstr ,liunw))
517 IF (
domain(ng)%SouthEast_Corner(tile))
THEN
518 ui(iend+1,jstr-1,liunw)=0.5_r8*(ui(iend ,jstr-1,liunw)+ &
519 & ui(iend+1,jstr ,liunw))
521 IF (
domain(ng)%NorthWest_Corner(tile))
THEN
522 ui(istr,jend+1,liunw)=0.5_r8*(ui(istr+1,jend+1,liunw)+ &
523 & ui(istr ,jend ,liunw))
525 IF (
domain(ng)%NorthEast_Corner(tile))
THEN
526 ui(iend+1,jend+1,liunw)=0.5_r8*(ui(iend ,jend+1,liunw)+ &
527 & ui(iend+1,jend ,liunw))
532 END SUBROUTINE ice_uibc_tile
type(t_grid), dimension(:), allocatable grid
type(t_ice_lobc), dimension(:,:), allocatable ice_lobc
type(t_ice), dimension(:), allocatable ice
integer, dimension(nices) ibice
integer, parameter isuice
type(t_lbc), dimension(:,:,:), allocatable lbc
type(t_domain), dimension(:), allocatable domain
real(dp), dimension(:), allocatable dt
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
real(r8), dimension(:), allocatable gamma2
integer, parameter isouth
real(dp), dimension(:,:), allocatable m2obc_out
integer, parameter inorth
real(dp), dimension(:,:), allocatable m2obc_in
recursive subroutine wclock_off(ng, model, region, line, routine)
recursive subroutine wclock_on(ng, model, region, line, routine)