ROMS
Loading...
Searching...
No Matches
ice_uibc.F
Go to the documentation of this file.
1#include "cppdefs.h"
3
4#ifdef ICE_MODEL
5!
6!git $Id$
7!=======================================================================
8! Copyright (c) 2002-2025 The ROMS Group !
9! Licensed under a MIT/X style license W. Paul Budgell !
10! See License_ROMS.md Katherine Hedstrom !
11!================================================== Hernan G. Arango ===
12! !
13! Sets the lateral boundary conditions on the the ice U-velocity. !
14! !
15!=======================================================================
16!
17 USE mod_param
18 USE mod_grid
19 USE mod_ice
20 USE mod_scalars
21!
22 implicit none
23!
24 PRIVATE
25 PUBLIC ice_uibc_tile
26!
27 CONTAINS
28!
29!***********************************************************************
30 SUBROUTINE ice_uibc (ng, tile, model)
31!***********************************************************************
32!
33 USE mod_stepping
34!
35! Imported variable declarations.
36!
37 integer, intent(in) :: ng, tile, model
38!
39! Local variable declarations.
40!
41 character (len=*), parameter :: MyFile = &
42 & __FILE__
43!
44#include "tile.h"
45!
46# ifdef PROFILE
47 CALL wclock_on (ng, model, 42, __line__, myfile)
48# endif
49 CALL ice_uibc_tile (ng, tile, model, &
50 & lbi, ubi, lbj, ubj, &
51 & imins, imaxs, jmins, jmaxs, &
52 & liuol(ng), liunw(ng), &
53 & ice(ng) % Si(:,:,:,isuice))
54# ifdef PROFILE
55 CALL wclock_off (ng, model, 42, __line__, myfile)
56# endif
57!
58 RETURN
59 END SUBROUTINE ice_uibc
60!
61!***********************************************************************
62 SUBROUTINE ice_uibc_tile (ng, tile, model, &
63 & LBi, UBi, LBj, UBj, &
64 & IminS, ImaxS, JminS, JmaxS, &
65 & liuol, liunw, &
66 & ui)
67!***********************************************************************
68!
69! Imported variable declarations.
70!
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
75!
76# ifdef ASSUMED_SHAPE
77 real(r8), intent(inout) :: ui(LBi:,LBj:,:)
78# else
79 real(r8), intent(inout) :: ui(LBi:UBi,LBj:UBj,2)
80# endif
81!
82! Local variable declarations.
83!
84 integer :: i, Imax, Imin, j, know
85!
86 real(r8), parameter :: eps =1.0e-20_r8
87 real(r8) :: Ce, Cx, cff, dUde, dUdt, dUdx, tau
88
89 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: grad
90
91#include "set_bounds.h"
92!
93!-----------------------------------------------------------------------
94! Set time-indices
95!-----------------------------------------------------------------------
96!
97 know=liuol
98!
99!-----------------------------------------------------------------------
100! Lateral boundary conditions at the western edge.
101!-----------------------------------------------------------------------
102!
103 IF (domain(ng)%Western_Edge(tile)) THEN
104 IF (lbc(iwest,ibice(isuice),ng)%radiation) THEN
105!
106! Western edge, implicit upstream radiation condition.
107!
108 DO j=jstr,jend+1
109 grad(istr ,j)=ui(istr ,j ,know)- &
110 & ui(istr ,j-1,know)
111 grad(istr+1,j)=ui(istr+1,j ,know)- &
112 & ui(istr+1,j-1,know)
113 END DO
114 DO j=jstr,jend
115 dudt=ui(istr+1,j,know )-ui(istr+1,j,liunw)
116 dudx=ui(istr+1,j,liunw)-ui(istr+2,j,liunw)
117 IF (lbc(iwest,ibice(isuice),ng)%nudging) THEN
118 IF ((dudt*dudx).lt.0.0_r8) THEN
119 tau=m2obc_in(ng,iwest)
120 ELSE
121 tau=m2obc_out(ng,iwest)
122 END IF
123 tau=tau*dt(ng)
124 END IF
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
127 dude=grad(istr+1,j )
128 ELSE
129 dude=grad(istr+1,j+1)
130 END IF
131 cff=max(dudx*dudx+dude*dude,eps)
132 cx=dudt*dudx
133# ifdef RADIATION_2D
134 ce=min(cff,max(dudt*dude,-cff))
135# else
136 ce=0.0_r8
137# endif
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))/ &
142 & (cff+cx)
143 IF (lbc(iwest,ibice(isuice),ng)%nudging) THEN
144 ui(istr,j,liunw)=ui(istr,j,liunw)+ &
145 & tau*(ice_lobc(isuice,ng)%ice_west(j)- &
146 & ui(istr,j,know))
147 END IF
148# ifdef MASKING
149 ui(istr,j,liunw)=ui(istr,j,liunw)* &
150 & grid(ng)%umask(istr,j)
151# endif
152 END DO
153!
154! Western edge, clamped boundary condition.
155!
156 ELSE IF (lbc(iwest,ibice(isuice),ng)%clamped) THEN
157 DO j=jstr,jend
158 ui(istr,j,liunw)=ice_lobc(isuice,ng)%ice_west(j)
159# ifdef MASKING
160 ui(istr,j,liunw)=ui(istr,j,liunw)* &
161 & grid(ng)%umask(istr,j)
162# endif
163# ifdef WET_DRY
164 ui(istr,j,liunw)=ui(istr,j,liunw)* &
165 & grid(ng)%umask_wet(istr,j)
166# endif
167 END DO
168!
169! Western edge, gradient boundary condition.
170!
171 ELSE IF (lbc(iwest,ibice(isuice),ng)%gradient) THEN
172 DO j=jstr,jend
173 ui(istr,j,liunw)=ui(istr+1,j,liunw)
174# ifdef MASKING
175 ui(istr,j,liunw)=ui(istr,j,liunw)* &
176 & grid(ng)%umask(istr,j)
177# endif
178# ifdef WET_DRY
179 ui(istr,j,liunw)=ui(istr,j,liunw)* &
180 & grid(ng)%umask_wet(istr,j)
181# endif
182 END DO
183!
184! Western edge, closed boundary condition.
185!
186 ELSE IF (lbc(iwest,ibice(isuice),ng)%closed) THEN
187 DO j=jstr,jend
188 ui(istr,j,liunw)=0.0_r8
189 END DO
190 END IF
191 END IF
192!
193!-----------------------------------------------------------------------
194! Lateral boundary conditions at the eastern edge.
195!-----------------------------------------------------------------------
196!
197 IF (domain(ng)%Eastern_Edge(tile)) THEN
198 IF (lbc(ieast,ibice(isuice),ng)%radiation) THEN
199!
200! Eastern edge, implicit upstream radiation condition.
201!
202 DO j=jstr,jend+1
203 grad(iend ,j)=ui(iend ,j ,know)- &
204 & ui(iend ,j-1,know)
205 grad(iend+1,j)=ui(iend+1,j ,know)- &
206 & ui(iend+1,j-1,know)
207 END DO
208 DO j=jstr,jend
209 dudt=ui(iend,j,know )-ui(iend ,j,liunw)
210 dudx=ui(iend,j,liunw)-ui(iend-1,j,liunw)
211 IF (lbc(ieast,ibice(isuice),ng)%nudging) THEN
212 IF ((dudt*dudx).lt.0.0_r8) THEN
213 tau=m2obc_in(ng,ieast)
214 ELSE
215 tau=m2obc_out(ng,ieast)
216 END IF
217 tau=tau*dt(ng)
218 END IF
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
221 dude=grad(iend,j)
222 ELSE
223 dude=grad(iend,j+1)
224 END IF
225 cff=max(dudx*dudx+dude*dude,eps)
226 cx=dudt*dudx
227# ifdef RADIATION_2D
228 ce=min(cff,max(dudt*dude,-cff))
229# else
230 ce=0.0_r8
231# endif
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))/ &
236 & (cff+cx)
237 IF (lbc(ieast,ibice(isuice),ng)%nudging) THEN
238 ui(iend+1,j,liunw)=ui(iend+1,j,liunw)+ &
239 & tau*(ice_lobc(isuice,ng)%ice_east(j)- &
240 & ui(iend+1,j,know))
241 END IF
242# ifdef MASKING
243 ui(iend+1,j,liunw)=ui(iend+1,j,liunw)* &
244 & grid(ng)%umask(iend+1,j)
245# endif
246 END DO
247!
248! Eastern edge, clamped boundary condition.
249!
250 ELSE IF (lbc(ieast,ibice(isuice),ng)%clamped) THEN
251 DO j=jstr,jend
252 ui(iend+1,j,liunw)=ice_lobc(isuice,ng)%ice_east(j)
253# ifdef MASKING
254 ui(iend+1,j,liunw)=ui(iend+1,j,liunw)* &
255 & grid(ng)%umask(iend+1,j)
256# endif
257# ifdef WET_DRY
258 ui(iend+1,j,liunw)=ui(iend+1,j,liunw)* &
259 & grid(ng)%umask_wet(iend+1,j)
260# endif
261 END DO
262!
263! Eastern edge, gradient boundary condition.
264!
265 ELSE IF (lbc(ieast,ibice(isuice),ng)%gradient) THEN
266 DO j=jstr,jend
267 ui(iend+1,j,liunw)=ui(iend,j,liunw)
268# ifdef MASKING
269 ui(iend+1,j,liunw)=ui(iend+1,j,liunw)* &
270 & grid(ng)%umask(iend+1,j)
271# endif
272# ifdef WET_DRY
273 ui(iend+1,j,liunw)=ui(iend+1,j,liunw)* &
274 & grid(ng)%umask_wet(iend+1,j)
275# endif
276 END DO
277!
278! Eastern edge, closed boundary condition.
279!
280 ELSE IF (lbc(ieast,ibice(isuice),ng)%closed) THEN
281 DO j=jstr,jend
282 ui(iend+1,j,liunw)=0.0_r8
283 END DO
284 END IF
285 END IF
286!
287!-----------------------------------------------------------------------
288! Lateral boundary conditions at the southern edge.
289!-----------------------------------------------------------------------
290!
291 IF (domain(ng)%Southern_Edge(tile)) THEN
292 IF (lbc(isouth,ibice(isuice),ng)%radiation) THEN
293!
294! Southern edge, implicit upstream radiation condition.
295!
296 DO i=istrp-1,iend
297 grad(i,jstr-1)=ui(i+1,jstr-1,know)- &
298 & ui(i ,jstr-1,know)
299 grad(i,jstr )=ui(i+1,jstr ,know)- &
300 & ui(i ,jstr ,know)
301 END DO
302 DO i=istrp,iend
303 dudt=ui(i,jstr,know )-ui(i,jstr ,liunw)
304 dude=ui(i,jstr,liunw)-ui(i,jstr+1,liunw)
305 IF (lbc(isouth,ibice(isuice),ng)%nudging) THEN
306 IF ((dudt*dude).lt.0.0_r8) THEN
307 tau=m2obc_in(ng,isouth)
308 ELSE
309 tau=m2obc_out(ng,isouth)
310 END IF
311 tau=tau*dt(ng)
312 END IF
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
315 dudx=grad(i-1,jstr)
316 ELSE
317 dudx=grad(i ,jstr)
318 END IF
319 cff=max(dudx*dudx+dude*dude,eps)
320# ifdef RADIATION_2D
321 cx=min(cff,max(dudt*dudx,-cff))
322# else
323 cx=0.0_r8
324# endif
325 ce=dudt*dude
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))/ &
330 & (cff+ce)
331 IF (lbc(isouth,ibice(isuice),ng)%nudging) THEN
332 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)+ &
333 & tau*(ice_lobc(isuice,ng)%ice_south(i)- &
334 & ui(i,jstr-1,know))
335 END IF
336# ifdef MASKING
337 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)* &
338 & grid(ng)%umask(i,jstr-1)
339# endif
340 END DO
341!
342! Southern edge, clamped boundary condition.
343!
344 ELSE IF (lbc(isouth,ibice(isuice),ng)%clamped) THEN
345 DO i=istrp,iend
346 ui(i,jstr-1,liunw)=ice_lobc(isuice,ng)%ice_south(i)
347# ifdef MASKING
348 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)* &
349 & grid(ng)%umask(i,jstr-1)
350# endif
351# ifdef WET_DRY
352 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)* &
353 & grid(ng)%umask_wet(i,jstr-1)
354# endif
355 END DO
356!
357! Southern edge, gradient boundary condition.
358!
359 ELSE IF (lbc(isouth,ibice(isuice),ng)%gradient) THEN
360 DO i=istrp,iend
361 ui(i,jstr-1,liunw)=ui(i,1,liunw)
362# ifdef MASKING
363 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)* &
364 & grid(ng)%umask(i,jstr-1)
365# endif
366# ifdef WET_DRY
367 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)* &
368 & grid(ng)%umask_wet(i,jstr-1)
369# endif
370 END DO
371!
372! Southern edge, closed boundary condition: free slip (gamma2=1) or
373! no slip (gamma2=-1).
374!
375 ELSE IF (lbc(isouth,ibice(isuice),ng)%closed) THEN
376 IF (ewperiodic(ng)) THEN
377 imin=istrp
378 imax=iend
379 ELSE
380 imin=istr
381 imax=iendt
382 END IF
383 DO i=imin,imax
384 ui(i,jstr-1,liunw)=gamma2(ng)*ui(i,1,liunw)
385# ifdef MASKING
386 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)* &
387 & grid(ng)%umask(i,jstr-1)
388# endif
389# ifdef WET_DRY
390 ui(i,jstr-1,liunw)=ui(i,jstr-1,liunw)* &
391 & grid(ng)%umask_wet(i,jstr-1)
392# endif
393 END DO
394 END IF
395 END IF
396!
397!-----------------------------------------------------------------------
398! Lateral boundary conditions at the northern edge.
399!-----------------------------------------------------------------------
400!
401 IF (domain(ng)%Northern_Edge(tile)) THEN
402 IF (lbc(inorth,ibice(isuice),ng)%radiation) THEN
403!
404! Northern edge, implicit upstream radiation condition.
405!
406 DO i=istrp-1,iend
407 grad(i,jend )=ui(i+1,jend ,know)- &
408 & ui(i ,jend ,know)
409 grad(i,jend+1)=ui(i+1,jend+1,know)- &
410 & ui(i ,jend+1,know)
411 END DO
412 DO i=istrp,iend
413 dudt=ui(i,jend,know )-ui(i,jend ,liunw)
414 dude=ui(i,jend,liunw)-ui(i,jend-1,liunw)
415 IF (lbc(inorth,ibice(isuice),ng)%nudging) THEN
416 IF ((dudt*dude).lt.0.0_r8) THEN
417 tau=m2obc_in(ng,inorth)
418 ELSE
419 tau=m2obc_out(ng,inorth)
420 END IF
421 tau=tau*dt(ng)
422 END IF
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
425 dudx=grad(i-1,jend)
426 ELSE
427 dudx=grad(i ,jend)
428 END IF
429 cff=max(dudx*dudx+dude*dude,eps)
430# ifdef RADIATION_2D
431 cx=min(cff,max(dudt*dudx,-cff))
432# else
433 cx=0.0_r8
434# endif
435 ce=dudt*dude
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))/ &
440 & (cff+ce)
441# ifdef NORTH_MINUDGING
442 IF (lbc(inorth,ibice(isuice),ng)%nudging) THEN
443 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)+ &
444 & tau*(ice_lobc(isuice.ng)%ice_north(i)- &
445 & ui(i,jend+1,know))
446# endif
447# ifdef MASKING
448 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)* &
449 & grid(ng)%umask(i,jend+1)
450# endif
451 END DO
452!
453! Northern edge, clamped boundary condition.
454!
455 ELSE IF (lbc(inorth,ibice(isuice),ng)%clamped) THEN
456 DO i=istrp,iend
457 ui(i,jend+1,liunw)=ice_lobc(isuice,ng)%ice_north(i)
458# ifdef MASKING
459 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)* &
460 & grid(ng)%umask(i,jend+1)
461# endif
462# ifdef WET_DRY
463 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)* &
464 & grid(ng)%umask_wet(i,jend+1)
465# endif
466 END DO
467!
468! Northern edge, gradient boundary condition.
469!
470 ELSE IF (lbc(inorth,ibice(isuice),ng)%gradient) THEN
471 DO i=istrp,iend
472 ui(i,jend+1,liunw)=ui(i,jend,liunw)
473# ifdef MASKING
474 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)* &
475 & grid(ng)%umask(i,jend+1)
476# endif
477# ifdef WET_DRY
478 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)* &
479 & grid(ng)%umask_wet(i,jend+1)
480# endif
481 END DO
482!
483! Northern edge, closed boundary condition: free slip (gamma2=1) or
484! no slip (gamma2=-1).
485!
486 ELSE IF (lbc(inorth,ibice(isuice),ng)%closed) THEN
487 IF (ewperiodic(ng)) THEN
488 imin=istrp
489 imax=iend
490 ELSE
491 imin=istr
492 imax=iendt
493 END IF
494 DO i=imin,imax
495 ui(i,jend+1,liunw)=gamma2(ng)*ui(i,jend,liunw)
496# ifdef MASKING
497 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)* &
498 & grid(ng)%umask(i,jend+1)
499# endif
500# ifdef WET_DRY
501 ui(i,jend+1,liunw)=ui(i,jend+1,liunw)* &
502 & grid(ng)%umask_wet(i,jend+1)
503# endif
504 END DO
505 END IF
506 END IF
507!
508!-----------------------------------------------------------------------
509! Boundary corners.
510!-----------------------------------------------------------------------
511!
512 IF (.not.(ewperiodic(ng).or.nsperiodic(ng))) THEN
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))
516 END IF
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))
520 END IF
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))
524 END IF
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))
528 END IF
529 END IF
530!
531 RETURN
532 END SUBROUTINE ice_uibc_tile
533#endif
534 END MODULE ice_uibc_mod
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_ice_lobc), dimension(:,:), allocatable ice_lobc
Definition ice_mod.h:303
type(t_ice), dimension(:), allocatable ice
Definition ice_mod.h:283
integer, dimension(nices) ibice
Definition ice_mod.h:162
integer, parameter isuice
Definition ice_mod.h:146
type(t_lbc), dimension(:,:,:), allocatable lbc
Definition mod_param.F:375
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
real(dp), dimension(:), allocatable dt
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
real(r8), dimension(:), allocatable gamma2
integer, parameter isouth
real(dp), dimension(:,:), allocatable m2obc_out
integer, parameter ieast
integer, parameter inorth
real(dp), dimension(:,:), allocatable m2obc_in
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3