ROMS
Loading...
Searching...
No Matches
ice_bc2d.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 lateral boundary conditions for any 2D ice field with array !
14! index, ifld.
15! !
16!=======================================================================
17!
18 USE mod_param
19 USE mod_grid
20 USE mod_ice
21 USE mod_scalars
22!
23 implicit none
24!
25 PRIVATE
26 PUBLIC ice_bc2d_tile
27!
28 CONTAINS
29!
30!***********************************************************************
31 SUBROUTINE ice_bc2d_tile (ng, tile, model, ifld, &
32 & LBi, UBi, LBj, UBj, &
33 & IminS, ImaxS, JminS, JmaxS, &
34 & liold, linew, &
35 & ui, vi, field, S)
36!***********************************************************************
37!
38! Imported variable declarations.
39!
40 integer, intent(in) :: ng, tile, model, ifld
41 integer, intent(in) :: LBi, UBi, LBj, UBj
42 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
43 integer, intent(in) :: liold, linew
44!
45 TYPE(T_LBC), intent(in) :: S(4)
46!
47# ifdef ASSUMED_SHAPE
48 real(r8), intent(in) :: ui(LBi:,LBj:,:)
49 real(r8), intent(in) :: vi(LBi:,LBj:,:)
50 real(r8), intent(inout) :: field(LBi:,LBj:,:)
51# else
52 real(r8), intent(in) :: ui(LBi:UBi,LBj:UBj,2)
53 real(r8), intent(in) :: vi(LBi:UBi,LBj:UBj,2)
54 real(r8), intent(inout) :: field(LBi:UBi,LBj:UBj,2)
55# endif
56!
57! Local variable declarations.
58!
59 integer :: i, j, know
60!
61 real(r8), parameter :: eps =1.0e-20_r8
62 real(r8) :: Ce, Cx, cff, dTde, dTdt, dTdx, tau
63
64 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: grad
65
66#include "set_bounds.h"
67!
68!-----------------------------------------------------------------------
69! Set time-indices.
70!-----------------------------------------------------------------------
71!
72 know=liold
73!
74!-----------------------------------------------------------------------
75! Lateral boundary conditions at the western edge.
76!-----------------------------------------------------------------------
77!
78 IF (domain(ng)%Western_Edge(tile)) THEN
79!
80! Western edge, implicit upstream radiation condition.
81!
82 IF (s(iwest)%radiation) THEN
83 DO j=jstr,jend+1
84 grad(istr-1,j)=field(istr-1,j,know)-field(istr-1,j-1,know)
85# ifdef MASKING
86 grad(istr-1,j)=grad(istr-1,j)*grid(ng)%vmask(istr-1,j)
87# endif
88 grad(istr,j)=field(istr,j,know)-field(istr,j-1,know)
89# ifdef MASKING
90 grad(istr,j)=grad(istr,j)*grid(ng)%vmask(istr,j)
91# endif
92 END DO
93 DO j=jstr,jend
94 dtdt=field(istr,j,know )-field(istr ,j,linew)
95 dtdx=field(istr,j,linew)-field(istr+1,j,linew)
96 IF (s(iwest)%nudging) THEN
97 tau=tobc_out(1,ng,iwest)
98 IF ((dtdt*dtdx).lt.0.0_r8) tau=tobc_in(1,ng,iwest)
99 tau=tau*dt(ng)
100 END IF
101 IF ((dtdt*dtdx).lt.0.0_r8) dtdt=0.0_r8
102 IF ((dtdt*(grad(istr,j)+grad(istr,j+1))).gt.0.0_r8) THEN
103 dtde=grad(istr,j )
104 ELSE
105 dtde=grad(istr,j+1)
106 END IF
107 cff=max(dtdx*dtdx+dtde*dtde,eps)
108 cx=dtdt*dtdx
109# ifdef RADIATION_2D
110 ce=min(cff,max(dtdt*dtde,-cff))
111# else
112 ce=0.0_r8
113# endif
114 field(istr-1,j,linew)=(cff*field(istr-1,j,know)+ &
115 & cx *field(istr ,j,linew)- &
116 & max(ce,0.0_r8)*grad(istr-1,j )- &
117 & min(ce,0.0_r8)*grad(istr-1,j+1))/ &
118 & (cff+cx)
119 IF (s(iwest)%nudging) THEN
120 field(istr-1,j,linew)=field(istr-1,j,linew)+ &
121 & tau*(ice_lobc(ifld,ng)%ice_west(j)- &
122 & field(istr-1,j,know))
123 END IF
124# ifdef MASKING
125 field(istr-1,j,linew)=field(istr-1,j,linew)* &
126 & grid(ng)%rmask(istr-1,j)
127# endif
128 END DO
129!
130! Western edge, clamped boundary condition.
131!
132 ELSE IF (s(iwest)%clamped) THEN
133 DO j=jstr,jend
134 field(istr-1,j,linew)=ice_lobc(ifld,ng)%ice_west(j)
135# ifdef MASKING
136 field(istr-1,j,linew)=field(istr-1,j,linew)* &
137 & grid(ng)%rmask(istr-1,j)
138# endif
139# ifdef WET_DRY
140 field(istr-1,j,linew)=field(istr-1,j,linew)* &
141 & grid(ng)%rmask_wet(istr-1,j)
142# endif
143 END DO
144!
145! Western edge, clamped on inflow, gradient on outflow.
146!
147 ELSE IF (s(iwest)%mixed) THEN
148 DO j=jstr,jend
149 IF (ui(1,j,linew).ge.0.0_r8) THEN
150 field(istr-1,j,linew)=ice_lobc(ifld,ng)%ice_west(j)
151# ifdef MASKING
152 field(istr-1,j,linew)=field(istr-1,j,linew)* &
153 & grid(ng)%rmask(istr-1,j)
154# endif
155# ifdef WET_DRY
156 field(istr-1,j,linew)=field(istr-1,j,linew)* &
157 & grid(ng)%rmask_wet(istr-1,j)
158# endif
159 ELSE
160 field(istr-1,j,linew)=field(istr,j,liold)
161# ifdef MASKING
162 field(istr-1,j,linew)=field(istr-1,j,linew)* &
163 & grid(ng)%rmask(istr-1,j)
164# endif
165# ifdef WET_DRY
166 field(istr-1,j,linew)=field(istr-1,j,linew)* &
167 & grid(ng)%rmask_wet(istr-1,j)
168# endif
169 END IF
170 END DO
171!
172! Western edge, closed boundary condition.
173!
174 ELSE IF (s(iwest)%closed) THEN
175 DO j=jstr,jend
176 field(istr-1,j,linew)=field(istr,j,linew)
177# ifdef MASKING
178 field(istr-1,j,linew)=field(istr-1,j,linew)* &
179 & grid(ng)%rmask(istr-1,j)
180# endif
181# ifdef WET_DRY
182 field(istr-1,j,linew)=field(istr-1,j,linew)* &
183 & grid(ng)%rmask_wet(istr-1,j)
184# endif
185 END DO
186 END IF
187 END IF
188!
189!-----------------------------------------------------------------------
190! Lateral boundary conditions at the eastern edge.
191!-----------------------------------------------------------------------
192!
193 IF (domain(ng)%Eastern_Edge(tile)) THEN
194!
195! Eastern edge, implicit upstream radiation condition.
196!
197 IF (s(ieast)%radiation) THEN
198 DO j=jstr,jend+1
199 grad(iend,j)=field(iend,j,know)-field(iend,j-1,know)
200# ifdef MASKING
201 grad(iend,j)=grad(iend,j)*grid(ng)%vmask(iend ,j)
202# endif
203 grad(iend+1,j)=field(iend+1,j,know)-field(iend+1,j-1,know)
204# ifdef MASKING
205 grad(iend+1,j)=grad(iend+1,j)*grid(ng)%vmask(iend+1,j)
206# endif
207 END DO
208 DO j=jstr,jend
209 dtdt=field(iend,j,know )-field(iend ,j,linew)
210 dtdx=field(iend,j,linew)-field(iend-1,j,linew)
211 IF (s(ieast)%nudging) THEN
212 tau=tobc_out(1,ng,ieast)
213 IF ((dtdt*dtdx).lt.0.0_r8) tau=tobc_in(1,ng,ieast)
214 tau=tau*dt(ng)
215 END IF
216 IF ((dtdt*dtdx).lt.0.0_r8) dtdt=0.0_r8
217 IF ((dtdt*(grad(iend,j)+grad(iend,j+1))).gt.0.0_r8) THEN
218 dtde=grad(iend,j )
219 ELSE
220 dtde=grad(iend,j+1)
221 END IF
222 cff=max(dtdx*dtdx+dtde*dtde,eps)
223 cx=dtdt*dtdx
224# ifdef RADIATION_2D
225 ce=min(cff,max(dtdt*dtde,-cff))
226# else
227 ce=0.0_r8
228# endif
229 field(iend+1,j,linew)=(cff*field(iend+1,j,know)+ &
230 & cx *field(iend ,j,linew)- &
231 & max(ce,0.0_r8)*grad(iend+1,j )- &
232 & min(ce,0.0_r8)*grad(iend+1,j+1))/ &
233 & (cff+cx)
234 IF (s(ieast)%nudging) THEN
235 field(iend+1,j,linew)=field(iend+1,j,linew)+ &
236 & tau*(ice_lobc(ifld,ng)%ice_east(j)- &
237 & field(iend+1,j,know))
238 END IF
239# ifdef MASKING
240 field(iend+1,j,linew)=field(iend+1,j,linew)* &
241 & grid(ng)%rmask(iend+1,j)
242# endif
243 END DO
244!
245! Eastern edge, clamped boundary condition.
246!
247 ELSE IF (s(ieast)%clamped) THEN
248 DO j=jstr,jend
249 field(iend+1,j,linew)=ice_lobc(ifld,ng)%ice_east(j)
250# ifdef MASKING
251 field(iend+1,j,linew)=field(iend+1,j,linew)* &
252 & grid(ng)%rmask(iend+1,j)
253# endif
254# ifdef WET_DRY
255 field(iend+1,j,linew)=field(iend+1,j,linew)* &
256 & grid(ng)%rmask_wet(iend+1,j)
257# endif
258 END DO
259!
260! Eastern edge, clamped on inflow, gradient on outflow.
261!
262 ELSE IF (s(ieast)%mixed) THEN
263 DO j=jstr,jend
264 IF (ui(iend+1,j,linew).le.0.0_r8) THEN
265 field(iend+1,j,linew)=ice_lobc(ifld,ng)%ice_east(j)
266# ifdef MASKING
267 field(iend+1,j,linew)=field(iend+1,j,linew)* &
268 & grid(ng)%rmask(iend+1,j)
269# endif
270# ifdef WET_DRY
271 field(iend+1,j,linew)=field(iend+1,j,linew)* &
272 & grid(ng)%rmask_wet(iend+1,j)
273# endif
274 ELSE
275 field(iend+1,j,linew)=field(iend,j,liold)
276# ifdef MASKING
277 field(iend+1,j,linew)=field(iend+1,j,linew)* &
278 & grid(ng)%rmask(iend+1,j)
279# endif
280# ifdef WET_DRY
281 field(iend+1,j,linew)=field(iend+1,j,linew)* &
282 & grid(ng)%rmask_wet(iend+1,j)
283# endif
284 END IF
285 END DO
286!
287! Eastern edge, closed boundary condition.
288!
289 ELSE IF (s(ieast)%closed) THEN
290 DO j=jstr,jend
291 field(iend+1,j,linew)=field(iend,j,linew)
292# ifdef MASKING
293 field(iend+1,j,linew)=field(iend+1,j,linew)* &
294 & grid(ng)%rmask(iend+1,j)
295# endif
296# ifdef WET_DRY
297 field(iend+1,j,linew)=field(iend+1,j,linew)* &
298 & grid(ng)%rmask_wet(iend+1,j)
299# endif
300 END DO
301 END IF
302 END IF
303!
304!-----------------------------------------------------------------------
305! Lateral boundary conditions at the southern edge.
306!-----------------------------------------------------------------------
307!
308 IF (domain(ng)%Southern_Edge(tile)) THEN
309!
310! Southern edge, implicit upstream radiation condition.
311!
312 IF (s(isouth)%radiation) THEN
313 DO i=istr,iend+1
314 grad(i,jstr)=field(i,jstr,know)-field(i-1,jstr,know)
315# ifdef MASKING
316 grad(i,jstr)=grad(i,jstr)*grid(ng)%umask(i,jstr)
317# endif
318 grad(i,jstr-1)=field(i,jstr-1,know)-field(i-1,jstr-1,know)
319# ifdef MASKING
320 grad(i,jstr-1)=grad(i,jstr-1)*grid(ng)%umask(i,jstr-1)
321# endif
322 END DO
323 DO i=istr,iend
324 dtdt=field(i,jstr,know )-field(i,jstr ,linew)
325 dtde=field(i,jstr,linew)-field(i,jstr+1,linew)
326 IF (s(isouth)%nudging) THEN
327 tau=tobc_out(1,ng,isouth)
328 IF ((dtdt*dtde).lt.0.0_r8) tau=tobc_in(1,ng,isouth)
329 tau=tau*dt(ng)
330 END IF
331 IF ((dtdt*dtde).lt.0.0_r8) dtdt=0.0_r8
332 IF ((dtdt*(grad(i,jstr)+grad(i+1,jstr))).gt.0.0_r8) THEN
333 dtdx=grad(i ,jstr)
334 ELSE
335 dtdx=grad(i+1,jstr)
336 END IF
337 cff=max(dtdx*dtdx+dtde*dtde,eps)
338# ifdef RADIATION_2D
339 cx=min(cff,max(dtdt*dtdx,-cff))
340# else
341 cx=0.0_r8
342# endif
343 ce=dtdt*dtde
344 field(i,jstr-1,linew)=(cff*field(i,jstr-1,know)+ &
345 & ce *field(i,jstr ,linew)- &
346 & max(cx,0.0_r8)*grad(i ,jstr-1)- &
347 & min(cx,0.0_r8)*grad(i+1,jstr-1))/ &
348 & (cff+ce)
349 IF (s(isouth)%nudging) THEN
350 field(i,jstr-1,linew)=field(i,jstr-1,linew)+ &
351 & tau*(ice_lobc(ifld,ng)%ice_south(i)-&
352 & field(i,jstr-1,know))
353 END IF
354# ifdef MASKING
355 field(i,jstr-1,linew)=field(i,jstr-1,linew)* &
356 & grid(ng)%rmask(i,jstr-1)
357# endif
358 END DO
359!
360! Southern edge, clamped boundary condition.
361!
362 ELSE IF (s(isouth)%clamped) THEN
363 DO i=istr,iend
364 field(i,jstr-1,linew)=ice_lobc(ifld,ng)%ice_south(i)
365# ifdef MASKING
366 field(i,jstr-1,linew)=field(i,jstr-1,linew)* &
367 & grid(ng)%rmask(i,jstr-1)
368# endif
369# ifdef WET_DRY
370 field(i,jstr-1,linew)=field(i,jstr-1,linew)* &
371 & grid(ng)%rmask_wet(i,jstr-1)
372# endif
373 END DO
374!
375! Southern edge, clamped on inflow, gradient on outflow.
376!
377 ELSE IF (s(isouth)%mixed) THEN
378 DO i=istr,iend
379 IF (vi(i,1,linew).ge.0._r8) THEN
380 field(i,jstr-1,linew)=ice_lobc(ifld,ng)%ice_south(i)
381# ifdef MASKING
382 field(i,jstr-1,linew)=field(i,jstr-1,linew)* &
383 & grid(ng)%rmask(i,jstr-1)
384# endif
385# ifdef WET_DRY
386 field(i,jstr-1,linew)=field(i,jstr-1,linew)* &
387 & grid(ng)%rmask_wet(i,jstr-1)
388# endif
389 ELSE
390 field(i,jstr-1,linew)=field(i,jstr,liold)
391# ifdef MASKING
392 field(i,jstr-1,linew)=field(i,jstr-1,linew)* &
393 & grid(ng)%rmask(i,jstr-1)
394# endif
395# ifdef WET_DRY
396 field(i,jstr-1,linew)=field(i,jstr-1,linew)* &
397 & grid(ng)%rmask_wet(i,jstr-1)
398# endif
399 END IF
400 END DO
401!
402! Southern edge, closed boundary condition.
403!
404 ELSE IF (s(isouth)%closed) THEN
405 DO i=istr,iend
406 field(i,jstr-1,linew)=field(i,jstr,linew)
407# ifdef MASKING
408 field(i,jstr-1,linew)=field(i,jstr-1,linew)* &
409 & grid(ng)%rmask(i,jstr-1)
410# endif
411# ifdef WET_DRY
412 field(i,jstr-1,linew)=field(i,jstr-1,linew)* &
413 & grid(ng)%rmask_wet(i,jstr-1)
414# endif
415 END DO
416 END IF
417 END IF
418!
419!-----------------------------------------------------------------------
420! Lateral boundary conditions at the northern edge.
421!-----------------------------------------------------------------------
422!
423 IF (domain(ng)%Northern_Edge(tile)) THEN
424!
425! Northern edge, implicit upstream radiation condition.
426!
427 IF (s(inorth)%radiation) THEN
428 DO i=istr,iend+1
429 grad(i,jend)=field(i,jend,know)-field(i-1,jend,know)
430# ifdef MASKING
431 grad(i,jend)=grad(i,jend)*grid(ng)%umask(i,jend)
432# endif
433 grad(i,jend+1)=field(i,jend+1,know)-field(i-1,jend+1,know)
434# ifdef MASKING
435 grad(i,jend+1)=grad(i,jend+1)*grid(ng)%umask(i,jend+1)
436# endif
437 END DO
438 DO i=istr,iend
439 dtdt=field(i,jend,know )-field(i,jend ,linew)
440 dtde=field(i,jend,linew)-field(i,jend-1,linew)
441 IF (s(inorth)%nudging) THEN
442 tau=tobc_out(1,ng,inorth)
443 IF ((dtdt*dtde).lt.0.0_r8) tau=tobc_in(1,ng,inorth)
444 tau=tau*dt(ng)
445 END IF
446 IF ((dtdt*dtde).lt.0.0_r8) dtdt=0.0_r8
447 IF ((dtdt*(grad(i,jend)+grad(i+1,jend))).gt.0.0_r8) THEN
448 dtdx=grad(i ,jend)
449 ELSE
450 dtdx=grad(i+1,jend)
451 END IF
452 cff=max(dtdx*dtdx+dtde*dtde,eps)
453# ifdef RADIATION_2D
454 cx=min(cff,max(dtdt*dtdx,-cff))
455# else
456 cx=0.0_r8
457# endif
458 ce=dtdt*dtde
459 field(i,jend+1,linew)=(cff*field(i,jend+1,know)+ &
460 & ce *field(i,jend ,linew)- &
461 & max(cx,0.0_r8)*grad(i ,jend+1)- &
462 & min(cx,0.0_r8)*grad(i+1,jend+1))/ &
463 & (cff+ce)
464 IF (s(inorth)%nudging) THEN
465 field(i,jend+1,linew)=field(i,jend+1,linew)+ &
466 & tau*(ice_lobc(ifld,ng)%ice_north(i)-&
467 & field(i,jend+1,know))
468 END IF
469# ifdef MASKING
470 field(i,jend+1,linew)=field(i,jend+1,linew)* &
471 & grid(ng)%rmask(i,jend+1)
472# endif
473 END DO
474!
475! Northern edge, clamped boundary condition.
476!
477 ELSE IF (s(inorth)%clamped) THEN
478 DO i=istr,iend
479 field(i,jend+1,linew)=ice_lobc(ifld,ng)%ice_north(i)
480# ifdef MASKING
481 field(i,jend+1,linew)=field(i,jend+1,linew)* &
482 & grid(ng)%rmask(i,jend+1)
483# endif
484# ifdef WET_DRY
485 field(i,jend+1,linew)=field(i,jend+1,linew)* &
486 & grid(ng)%rmask_wet(i,jend+1)
487# endif
488 END DO
489!
490! Northern edge, clamped on inflow, gradient on outflow.
491!
492 ELSE IF (s(inorth)%mixed) THEN
493 DO i=istr,iend
494 IF (vi(i,jend+1,linew).le.0.0_r8) THEN
495 field(i,jend+1,linew)=ice_lobc(ifld,ng)%ice_north(i)
496# ifdef MASKING
497 field(i,jend+1,linew)=field(i,mm(ng)+1,linew)* &
498 & grid(ng)%rmask(i,jend+1)
499# endif
500# ifdef WET_DRY
501 field(i,jend+1,linew)=field(i,jend+1,linew)* &
502 & grid(ng)%rmask_wet(i,jend+1)
503# endif
504 ELSE
505 field(i,jend+1,linew)=field(i,jend,liold)
506# ifdef MASKING
507 field(i,jend+1,linew)=field(i,jend+1,linew)* &
508 & grid(ng)%rmask(i,jend+1)
509# endif
510# ifdef WET_DRY
511 field(i,jend+1,linew)=field(i,jend+1,linew)* &
512 & grid(ng)%rmask_wet(i,jend+1)
513# endif
514 END IF
515 END DO
516!
517! Northern edge, closed boundary condition.
518!
519 ELSE IF (s(inorth)%closed) THEN
520 DO i=istr,iend
521 field(i,jend+1,linew)=field(i,jend,linew)
522# ifdef MASKING
523 field(i,jend+1,linew)=field(i,jend+1,linew)* &
524 & grid(ng)%rmask(i,jend+1)
525# endif
526# ifdef WET_DRY
527 field(i,jend+1,linew)=field(i,jend+1,linew)* &
528 & grid(ng)%rmask_wet(i,jend+1)
529# endif
530 END DO
531 END IF
532 END IF
533!
534!-----------------------------------------------------------------------
535! Boundary corners.
536!-----------------------------------------------------------------------
537!
538 IF (.not.ewperiodic(ng).and. .not.nsperiodic(ng)) THEN
539 IF (domain(ng)%SouthWest_Corner(tile)) THEN
540 field(istr-1,jstr-1,linew)=0.5_r8* &
541 & (field(istr ,jstr-1,linew)+ &
542 & field(istr-1,jstr ,linew))
543 END IF
544 IF (domain(ng)%SouthEast_Corner(tile)) THEN
545 field(iend+1,jstr-1,linew)=0.5_r8* &
546 & (field(iend+1,jstr ,linew)+ &
547 & field(iend ,jstr-1,linew))
548 END IF
549 IF (domain(ng)%NorthWest_Corner(tile)) THEN
550 field(istr-1,jend+1,linew)=0.5_r8* &
551 & (field(istr-1,jend ,linew)+ &
552 & field(istr ,jend+1,linew))
553 END IF
554 IF (domain(ng)%NorthEast_Corner(tile)) THEN
555 field(iend+1,jend+1,linew)=0.5_r8* &
556 & (field(iend+1,jend ,linew)+ &
557 & field(iend ,jend+1,linew))
558 END IF
559 END IF
560 RETURN
561 END SUBROUTINE ice_bc2d_tile
562#endif
563
564 END MODULE ice_bc2d_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_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable mm
Definition mod_param.F:456
real(dp), dimension(:), allocatable dt
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
real(dp), dimension(:,:,:), allocatable tobc_out
real(dp), dimension(:,:,:), allocatable tobc_in
integer, parameter isouth
integer, parameter ieast
integer, parameter inorth