ROMS
Loading...
Searching...
No Matches
ice_vibc.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 V-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_vibc_tile
26!
27 CONTAINS
28!
29!***********************************************************************
30 SUBROUTINE vibc (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_vibc_tile (ng, tile, model, &
50 & lbi, ubi, lbj, ubj, &
51 & imins, imaxs, jmins, jmaxs, &
52 & liuol(ng), liunw(ng), &
53 & ice(ng) % Si(:,:,:,isvice))
54# ifdef PROFILE
55 CALL wclock_off (ng, model, 42, __line__, myfile)
56# endif
57!
58 RETURN
59 END SUBROUTINE vibc
60!
61!***********************************************************************
62 SUBROUTINE ice_vibc_tile (ng, tile, model, &
63 & LBi, UBi, LBj, UBj, &
64 & IminS, ImaxS, JminS, JmaxS, &
65 & liuol, liunw, &
66 & vi)
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) :: vi(LBi:,LBj:,:)
78# else
79 real(r8), intent(inout) :: vi(LBi:UBi,LBj:UBj,2)
80# endif
81!
82! Local variable declarations.
83!
84 integer :: i, Jmax, Jmin, j, know
85!
86 real(r8), parameter :: eps =1.0e-20_r8
87 real(r8) :: Ce, Cx, cff, dVde, dVdt, dVdx, 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(isvice),ng)%radiation) THEN
105!
106! Western edge, implicit upstream radiation condition.
107!
108 DO j=jstrp-1,jend
109 grad(istr-1,j)=vi(istr-1,j+1,know)- &
110 & vi(istr-1,j ,know)
111 grad(istr ,j)=vi(istr ,j+1,know)- &
112 & vi(istr ,j ,know)
113 END DO
114 DO j=jstrp,jend
115 dvdt=vi(istr,j,know )-vi(istr ,j,liunw)
116 dvdx=vi(istr,j,liunw)-vi(istr+1,j,liunw)
117 IF (lbc(iwest,ibice(isvice),ng)%nudging) THEN
118 IF ((dvdt*dvdx).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 ((dvdt*dvdx).lt.0.0_r8) dvdt=0.0_r8
126 IF ((dvdt*(grad(istr,j-1)+grad(istr,j))).gt.0.0_r8) THEN
127 dvde=grad(istr,j-1)
128 ELSE
129 dvde=grad(istr,j )
130 END IF
131 cff=max(dvdx*dvdx+dvde*dvde,eps)
132 cx=dvdt*dvdx
133# ifdef RADIATION_2D
134 ce=min(cff,max(dvdt*dvde,-cff))
135# else
136 ce=0.0_r8
137# endif
138 vi(istr-1,j,liunw)=(cff*vi(istr-1,j,know)+ &
139 & cx *vi(istr ,j,liunw)- &
140 & max(ce,0.0_r8)*grad(istr-1,j-1)- &
141 & min(ce,0.0_r8)*grad(istr-1,j ))/ &
142 & (cff+cx)
143 IF (lbc(iwest,ibice(isvice),ng)%nudging) THEN
144 vi(istr-1,j,liunw)=vi(istr-1,j,liunw)+ &
145 & tau*(ice_lobc(isvice,ng)%ice_west(j)- &
146 & vi(istr-1,j,know))
147 END IF
148# ifdef MASKING
149 vi(istr-1,j,liunw)=vi(istr-1,j,liunw)* &
150 & grid(ng)%vmask(istr-1,j)
151# endif
152 END DO
153!
154! Western edge, clamped boundary condition.
155!
156 ELSE IF (lbc(iwest,ibice(isvice),ng)%clamped) THEN
157 DO j=jstrp,jend
158 vi(istr-1,j,liunw)=ice_lobc(isvice,ng)%ice_west(j)
159# ifdef MASKING
160 vi(istr-1,j,liunw)=vi(istr-1,j,liunw)* &
161 & grid(ng)%vmask(istr-1,j)
162# endif
163# ifdef WET_DRY
164 vi(istr-1,j,liunw)=vi(istr-1,j,liunw)* &
165 & grid(ng)%vmask_wet(istr-1,j)
166# endif
167 END DO
168!
169! Western edge, gradient boundary condition.
170!
171 ELSE IF (lbc(iwest,ibice(isvice),ng)%gradient) THEN
172 DO j=jstrp,jend
173 vi(istr-1,j,liunw)=vi(istr,j,liunw)
174# ifdef MASKING
175 vi(istr-1,j,liunw)=vi(istr-1,j,liunw)* &
176 & grid(ng)%vmask(istr-1,j)
177# endif
178# ifdef WET_DRY
179 vi(istr-1,j,liunw)=vi(istr-1,j,liunw)* &
180 & grid(ng)%vmask_wet(istr-1,j)
181# endif
182 END DO
183!
184! Western edge, closed boundary condition: free slip (gamma2=1) or
185! no slip (gamma2=-1).
186!
187 ELSE IF (lbc(iwest,ibice(isvice),ng)%closed) THEN
188 IF (nsperiodic(ng)) THEN
189 jmin=jstrp
190 jmax=jend
191 ELSE
192 jmin=jstr
193 jmax=jendt
194 END IF
195 DO j=jmin,jmax
196 vi(istr-1,j,liunw)=gamma2(ng)*vi(istr-1,j,liunw)
197# ifdef MASKING
198 vi(istr-1,j,liunw)=vi(istr-1,j,liunw)* &
199 & grid(ng)%vmask(istr-1,j)
200# endif
201# ifdef WET_DRY
202 vi(istr-1,j,liunw)=vi(istr-1,j,liunw)* &
203 & grid(ng)%vmask_wet(istr-1,j)
204# endif
205 END DO
206 END IF
207 END IF
208!
209!-----------------------------------------------------------------------
210! Lateral boundary conditions at the eastern edge.
211!-----------------------------------------------------------------------
212!
213 IF (domain(ng)%Eastern_Edge(tile)) THEN
214 IF (lbc(ieast,ibice(isvice),ng)%radiation) THEN
215!
216! Eastern edge, implicit upstream radiation condition.
217!
218 DO j=jstrp-1,jend
219 grad(iend ,j)=vi(iend ,j+1,know)- &
220 & vi(iend ,j ,know)
221 grad(iend+1,j)=vi(iend+1,j+1,know)- &
222 & vi(iend+1,j ,know)
223 END DO
224 DO j=jstrp,jend
225 dvdt=vi(iend,j,know )-vi(iend ,j,liunw)
226 dvdx=vi(iend,j,liunw)-vi(iend-1,j,liunw)
227 IF (lbc(ieast,ibice(isvice),ng)%nudging) THEN
228 IF ((dvdt*dvdx).lt.0.0_r8) THEN
229 tau=m2obc_in(ng,ieast)
230 ELSE
231 tau=m2obc_out(ng,ieast)
232 END IF
233 tau=tau*dt(ng)
234 END IF
235 IF ((dvdt*dvdx).lt.0.0_r8) dvdt=0.0_r8
236 IF ((dvdt*(grad(iend,j-1)+grad(iend,j))).gt.0.0_r8) THEN
237 dvde=grad(iend,j-1)
238 ELSE
239 dvde=grad(iend,j )
240 END IF
241 cff=max(dvdx*dvdx+dvde*dvde,eps)
242 cx=dvdt*dvdx
243# ifdef RADIATION_2D
244 ce=min(cff,max(dvdt*dvde,-cff))
245# else
246 ce=0.0_r8
247# endif
248 vi(iend+1,j,liunw)=(cff*vi(iend+1,j,know)+ &
249 & cx *vi(iend ,j,liunw)- &
250 & max(ce,0.0_r8)*grad(iend+1,j-1)- &
251 & min(ce,0.0_r8)*grad(iend+1,j ))/ &
252 & (cff+cx)
253 IF (lbc(ieast,ibice(isvice),ng)%nudging) THEN
254 vi(iend+1,j,liunw)=vi(iend+1,j,liunw)+ &
255 & tau*(ice_lobc(isvice,ng)%ice_east(j)- &
256 & vi(iend+1,j,know))
257 END IF
258# ifdef MASKING
259 vi(iend+1,j,liunw)=vi(iend+1,j,liunw)* &
260 & grid(ng)%vmask(iend+1,j)
261# endif
262 END DO
263!
264! Eastern edge, clamped boundary condition.
265!
266 ELSE IF (lbc(ieast,ibice(isvice),ng)%clamped) THEN
267 DO j=jstrp,jend
268 vi(iend+1,j,liunw)=ice_lobc(isvice,ng)%ice_east(j)
269# ifdef MASKING
270 vi(iend+1,j,liunw)=vi(iend+1,j,liunw)* &
271 & grid(ng)%vmask(iend+1,j)
272# endif
273# ifdef WET_DRY
274 vi(iend+1,j,liunw)=vi(iend+1,j,liunw)* &
275 & grid(ng)%vmask_wet(iend+1,j)
276# endif
277 END DO
278!
279! Eastern edge, gradient boundary condition.
280!
281 ELSE IF (lbc(ieast,ibice(isvice),ng)%gradient) THEN
282 DO j=jstrp,jend
283 vi(iend+1,j,liunw)=vi(iend,j,liunw)
284# ifdef MASKING
285 vi(iend+1,j,liunw)=vi(iend+1,j,liunw)* &
286 & grid(ng)%vmask(iend+1,j)
287# endif
288# ifdef WET_DRY
289 vi(iend+1,j,liunw)=vi(iend+1,j,liunw)* &
290 & grid(ng)%vmask_wet(iend+1,j)
291# endif
292 END DO
293!
294! Eastern edge, closed boundary condition: free slip (gamma2=1) or
295! no slip (gamma2=-1).
296!
297 ELSE IF (lbc(ieast,ibice(isvice),ng)%closed) THEN
298 IF (nsperiodic(ng)) THEN
299 jmin=jstrp
300 jmax=jend
301 ELSE
302 jmin=jstr
303 jmax=jendt
304 END IF
305 DO j=jmin,jmax
306 vi(iend+1,j,liunw)=gamma2(ng)*vi(iend,j,liunw)
307# ifdef MASKING
308 vi(iend+1,j,liunw)=vi(iend+1,j,liunw)* &
309 & grid(ng)%vmask(iend+1,j)
310# endif
311# ifdef WET_DRY
312 vi(iend+1,j,liunw)=vi(iend+1,j,liunw)* &
313 & grid(ng)%vmask_wet(iend+1,j)
314# endif
315 END DO
316 END IF
317 END IF
318!
319!-----------------------------------------------------------------------
320! Lateral boundary conditions at the southern edge.
321!-----------------------------------------------------------------------
322!
323 IF (domain(ng)%Southern_Edge(tile)) THEN
324 IF (lbc(isouth,isvice,ng)%radiation) THEN
325!
326! Southern edge, implicit upstream radiation condition.
327!
328 DO i=istr,iend+1
329 grad(i,jstr )=vi(i ,jstr ,know)- &
330 & vi(i-1,jstr ,know)
331 grad(i,jstr+1)=vi(i ,jstr+1,know)- &
332 & vi(i-1,jstr+1,know)
333 END DO
334 DO i=istr,iend
335 dvdt=vi(i,jstr+1,know )-vi(i,jstr+1,liunw)
336 dvde=vi(i,jstr+1,liunw)-vi(i,jstr+2,liunw)
337 IF (lbc(isouth,ibice(isvice),ng)%nudging) THEN
338 IF ((dvdt*dvde).lt.0.0_r8) THEN
339 tau=m2obc_in(ng,isouth)
340 ELSE
341 tau=m2obc_out(ng,isouth)
342 END IF
343 tau=tau*dt(ng)
344 END IF
345 IF ((dvdt*dvde).lt.0.0_r8) dvdt=0.0_r8
346 IF ((dvdt*(grad(i,jstr+1)+grad(i+1,jstr+1))).gt.0.0_r8) THEN
347 dvdx=grad(i ,jstr+1)
348 ELSE
349 dvdx=grad(i+1,jstr+1)
350 END IF
351 cff=max(dvdx*dvdx+dvde*dvde,eps)
352# ifdef RADIATION_2D
353 cx=min(cff,max(dvdt*dvdx,-cff))
354# else
355 cx=0.0_r8
356# endif
357 ce=dvdt*dvde
358 vi(i,jstr,liunw)=(cff*vi(i,jstr ,know)+ &
359 & ce *vi(i,jstr+1,liunw)- &
360 & max(cx,0.0_r8)*grad(i ,jstr)- &
361 & min(cx,0.0_r8)*grad(i+1,jstr))/ &
362 & (cff+ce)
363 IF (lbc(isouth,ibice(isvice),ng)%nudging) THEN
364 vi(i,jstr,liunw)=vi(i,jstr,liunw)+ &
365 & tau*(ice_lobc(isvice,ng)%ice_south(i)- &
366 & vi(i,jstr,know))
367 END IF
368# ifdef MASKING
369 vi(i,jstr,liunw)=vi(i,jstr,liunw)* &
370 & grid(ng)%vmask(i,jstr)
371# endif
372 END DO
373!
374! Southern edge, clamped boundary condition.
375!
376 ELSE IF (lbc(isouth,ibice(isvice),ng)%clamped) THEN
377 DO i=istr,iend
378 vi(i,jstr,liunw)=ice_lobc(isvice,ng)%ice_south(i)
379# ifdef MASKING
380 vi(i,jstr,liunw)=vi(i,jstr,liunw)* &
381 & grid(ng)%vmask(i,jstr)
382# endif
383# ifdef WET_DRY
384 vi(i,jstr,liunw)=vi(i,jstr,liunw)* &
385 & grid(ng)%vmask_wet(i,jstr)
386# endif
387 END DO
388!
389! Southern edge, gradient boundary condition.
390!
391 ELSE IF (lbc(isouth,ibice(isvice),ng)%gradient) THEN
392 DO i=istr,iend
393 vi(i,jstr,liunw)=vi(i,jstr+1,liunw)
394# ifdef MASKING
395 vi(i,jstr,liunw)=vi(i,jstr,liunw)* &
396 & grid(ng)%vmask(i,jstr)
397# endif
398# ifdef WET_DRY
399 vi(i,jstr,liunw)=vi(i,jstr,liunw)* &
400 & grid(ng)%vmask_wet(i,jstr)
401# endif
402 END DO
403!
404! Southern edge, closed boundary condition.
405!
406 ELSE IF (lbc(isouth,ibice(isvice),ng)%closed) THEN
407 DO i=istr,iend
408 vi(i,jstr,liunw)=0.0_r8
409 END DO
410 END IF
411 END IF
412!
413!-----------------------------------------------------------------------
414! Lateral boundary conditions at the northern edge.
415!-----------------------------------------------------------------------
416!
417 IF (domain(ng)%Northern_Edge(tile)) THEN
418 IF (lbc(inorth,ibice(isvice),ng)%radiation) THEN
419!
420! Northern edge, implicit upstream radiation condition.
421!
422 DO i=istr,iend+1
423 grad(i,jend )=vi(i ,jend ,know)- &
424 & vi(i-1,jend ,know)
425 grad(i,jend+1)=vi(i ,jend+1,know)- &
426 & vi(i-1,jend+1,know)
427 END DO
428 DO i=istr,iend
429 dvdt=vi(i,jend,know )-vi(i,jend ,liunw)
430 dvde=vi(i,jend,liunw)-vi(i,jend-1,liunw)
431 IF (lbc(inorth,ibice(isvice),ng)%nudging) THEN
432 IF ((dvdt*dvde).lt.0.0_r8) THEN
433 tau=m2obc_in(ng,inorth)
434 ELSE
435 tau=m2obc_out(ng,inorth)
436 END IF
437 tau=tau*dt(ng)
438 END IF
439 IF ((dvdt*dvde).lt.0.0_r8) dvdt=0.0_r8
440 IF ((dvdt*(grad(i,jend)+grad(i+1,jend))).gt.0.0_r8) THEN
441 dvdx=grad(i ,jend)
442 ELSE
443 dvdx=grad(i+1,jend)
444 END IF
445 cff=max(dvdx*dvdx+dvde*dvde,eps)
446# ifdef RADIATION_2D
447 cx=min(cff,max(dvdt*dvdx,-cff))
448# else
449 cx=0.0_r8
450# endif
451 ce=dvdt*dvde
452 vi(i,jend+1,liunw)=(cff*vi(i,jend+1,know)+ &
453 & ce *vi(i,jend ,liunw)- &
454 & max(cx,0.0_r8)*grad(i ,jend+1)- &
455 & min(cx,0.0_r8)*grad(i+1,jend+1))/ &
456 & (cff+ce)
457 IF (lbc(inorth,ibice(isvice),ng)%nudging) THEN
458 vi(i,jend+1,liunw)=vi(i,jend+1,liunw)+ &
459 & tau*(ice_lobc(isvice,ng)%ice_north(i)- &
460 & vi(i,jend+1,know))
461 END IF
462# ifdef MASKING
463 vi(i,jend+1,liunw)=vi(i,jend+1,liunw)* &
464 & grid(ng)%vmask(i,jend+1)
465# endif
466 END DO
467!
468! Northern edge, clamped boundary condition.
469!
470 ELSE IF (lbc(inorth,ibice(isvice),ng)%clamped) THEN
471 DO i=istr,iend
472 vi(i,jend+1,liunw)=ice_lobc(isvice,ng)%ice_north(i)
473# ifdef MASKING
474 vi(i,jend+1,liunw)=vi(i,jend+1,liunw)* &
475 & grid(ng)%vmask(i,jend)
476# endif
477# ifdef WET_DRY
478 vi(i,jend+1,liunw)=vi(i,jend+1,liunw)* &
479 & grid(ng)%vmask_wet(i,jend)
480# endif
481 END DO
482!
483! Northern edge, gradient boundary condition.
484!
485 ELSE IF (lbc(inorth,ibice(isvice),ng)%gradient) THEN
486 DO i=istr,iend
487 vi(i,jend+1,liunw)=vi(i,jend,liunw)
488# ifdef MASKING
489 vi(i,jend+1,liunw)=vi(i,jend+1,liunw)* &
490 & grid(ng)%vmask(i,jend)
491# endif
492# ifdef WET_DRY
493 vi(i,jend+1,liunw)=vi(i,jend+1,liunw)* &
494 & grid(ng)%vmask_wet(i,jend)
495# endif
496 END DO
497!
498! Northern edge, closed boundary condition.
499!
500 ELSE IF (lbc(inorth,ibice(isvice),ng)%closed) THEN
501 DO i=istr,iend
502 vi(i,jend+1,liunw)=0.0_r8
503 END DO
504 END IF
505 END IF
506!
507!-----------------------------------------------------------------------
508! Boundary corners.
509!-----------------------------------------------------------------------
510!
511 IF (.not.(ewperiodic(ng).or.nsperiodic(ng))) THEN
512 IF (domain(ng)%SouthWest_Corner(tile)) THEN
513 vi(istr-1,jstr,liunw)=0.5_r8*(vi(istr-1,jstr+1,liunw)+ &
514 & vi(istr ,jstr ,liunw))
515 END IF
516 IF (domain(ng)%SouthEast_Corner(tile)) THEN
517 vi(iend+1,jstr,liunw)=0.5_r8*(vi(iend ,jstr ,liunw)+ &
518 & vi(iend+1,jstr+1,liunw))
519 END IF
520 IF (domain(ng)%NorthWest_Corner(tile)) THEN
521 vi(istr-1,jend+1,liunw)=0.5_r8*(vi(istr-1,jend ,liunw)+ &
522 & vi(istr ,jend+1,liunw))
523 END IF
524 IF (domain(ng)%NorthEast_Corner(tile)) THEN
525 vi(iend+1,jend+1,liunw)=0.5_r8*(vi(iend+1,jend ,liunw)+ &
526 & vi(iend ,jend+1,liunw))
527 END IF
528 END IF
529!
530 RETURN
531 END SUBROUTINE ice_vibc_tile
532#endif
533 END MODULE ice_vibc_mod
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer, parameter isvice
Definition ice_mod.h:147
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
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