ROMS
Loading...
Searching...
No Matches
ice_tibc.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 internal ice !
14! temperature, ti. !
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_tibc_tile
27!
28 CONTAINS
29!
30!***********************************************************************
31 SUBROUTINE ice_tibc (ng, tile, model)
32!***********************************************************************
33!
34 USE mod_stepping
35!
36! Imported variable declarations.
37!
38 integer, intent(in) :: ng, tile, model
39!
40! Local variable declarations.
41!
42 character (len=*), parameter :: MyFile = &
43 & __FILE__
44!
45# include "tile.h"
46!
47# ifdef PROFILE
48 CALL wclock_on (ng, model, 42, __line__, myfile)
49# endif
50 CALL ice_tibc_tile (ng, tile, model, &
51 & lbi, ubi, lbj, ubj, &
52 & liold(ng), linew(ng), &
53 & ice(ng) % Si(:,:,:,isuice), &
54 & ice(ng) % Si(:,:,:,isvice), &
55 & ice(ng) % Si(:,:,:,ishice), &
56 & ice(ng) % Si(:,:,:,istice), &
57 & ice(ng) % Si(:,:,:,isenth))
58# ifdef PROFILE
59 CALL wclock_off (ng, model, 42, __line__, myfile)
60# endif
61!
62 RETURN
63 END SUBROUTINE ice_tibc
64!
65!***********************************************************************
66 SUBROUTINE ice_tibc_tile (ng, tile, model, &
67 & LBi, UBi, LBj, UBj, &
68 & liold, linew, &
69 & ui, vi, hi, ti, enthalpy)
70!***********************************************************************
71!
72! Imported variable declarations.
73!
74 integer, intent(in) :: ng, tile, model
75 integer, intent(in) :: LBi, UBi, LBj, UBj
76 integer, intent(in) :: liold, linew
77!
78# ifdef ASSUMED_SHAPE
79 real(r8), intent(in) :: ui(LBi:,LBj:,:)
80 real(r8), intent(in) :: vi(LBi:,LBj:,:)
81 real(r8), intent(in) :: hi(LBi:,LBj:,:)
82 real(r8), intent(inout) :: ti(LBi:,LBj:,:)
83 real(r8), intent(inout) :: enthalpy(LBi:,LBj:,:)
84# else
85 real(r8), intent(in) :: ui(LBi:UBi,LBj:UBj,2)
86 real(r8), intent(in) :: vi(LBi:UBi,LBj:UBj,2)
87 real(r8), intent(in) :: hi(LBi:UBi,LBj:UBj,2)
88 real(r8), intent(inout) :: ti(LBi:UBi,LBj:UBj,2)
89 real(r8), intent(inout) :: enthalpy(LBi:UBi,LBj:UBj,2)
90# endif
91!
92! Local variable declarations.
93!
94 integer :: i, j, know
95!
96 real(r8), parameter :: eps=1.0e-6_r8
97
98#include "set_bounds.h"
99!
100!-----------------------------------------------------------------------
101! Set time-indices
102!-----------------------------------------------------------------------
103!
104 know=liold
105!
106!-----------------------------------------------------------------------
107! Lateral boundary conditions at the western edge.
108!-----------------------------------------------------------------------
109!
110 IF (domain(ng)%Western_Edge(tile)) THEN
111!
112! Western edge, clamped boundary condition.
113!
114 IF (lbc(iwest,ibice(istice),ng)%clamped) THEN
115 DO j=jstr,jend
116 enthalpy(istr-1,j,linew)=ice_lobc(ishice,ng)%ice_west(j)* &
117 & ice_lobc(istice,ng)%ice_west(j)
118# ifdef MASKING
119 enthalpy(istr-1,j,linew)=enthalpy(istr-1,j,linew)* &
120 & grid(ng)%rmask(istr-1,j)
121# endif
122# ifdef WET_DRY
123 enthalpy(istr-1,j,linew)=enthalpy(istr-1,j,linew)* &
124 & grid(ng)%rmask_wet(istr-1,j)
125# endif
126 ti(istr-1,j,linew)=enthalpy(istr-1,j,linew)/ &
127 & max(hi(istr-1,j,linew),eps)
128 IF (hi(istr-1,j,linew).le.min_hi(ng)) THEN
129 enthalpy(istr-1,j,linew)=0.0_r8
130 ti(istr-1,j,linew)=0.0_r8
131 END IF
132 END DO
133!
134! Western edge, clamped on inflow, gradient on outflow.
135!
136 ELSE IF (lbc(iwest,ibice(istice),ng)%mixed) THEN
137 DO j=jstr,jend
138 IF (ui(istr,j,linew).ge.0.0_r8) THEN
139 enthalpy(istr-1,j,linew)=ice_lobc(ishice,ng)%ice_west(j)* &
140 & ice_lobc(istice,ng)%ice_west(j)
141# ifdef MASKING
142 enthalpy(istr-1,j,linew)=enthalpy(istr-1,j,linew)* &
143 & grid(ng)%rmask(istr-1,j)
144# endif
145# ifdef WET_DRY
146 enthalpy(istr-1,j,linew)=enthalpy(istr-1,j,linew)* &
147 & grid(ng)%rmask_wet(istr-1,j)
148# endif
149 ELSE
150 enthalpy(istr-1,j,linew)=enthalpy(istr,j,liold)
151# ifdef MASKING
152 enthalpy(istr-1,j,linew)=enthalpy(istr-1,j,linew)* &
153 & grid(ng)%rmask(istr-1,j)
154# endif
155# ifdef WET_DRY
156 enthalpy(istr-1,j,linew)=enthalpy(istr-1,j,linew)* &
157 & grid(ng)%rmask_wet(istr-1,j)
158# endif
159 END IF
160 ti(istr-1,j,linew)=enthalpy(istr-1,j,linew)/ &
161 & max(hi(istr-1,j,linew),eps)
162 IF (hi(istr-1,j,linew).le.min_hi(ng)) THEN
163 enthalpy(istr-1,j,linew)=0.0_r8
164 ti(istr-1,j,linew)=0.0_r8
165 END IF
166 END DO
167!
168! Western edge, gradient boundary condition.
169!
170 ELSE IF (lbc(iwest,ibice(istice),ng)%gradient) THEN
171 DO j=jstr,jend
172 enthalpy(istr-1,j,linew)=hi(istr,j,linew)* &
173 & ti(istr,j,linew)
174# ifdef MASKING
175 enthalpy(istr-1,j,linew)=enthalpy(istr-1,j,linew)* &
176 & grid(ng)%rmask(istr-1,j)
177# endif
178# ifdef WET_DRY
179 enthalpy(istr-1,j,linew)=enthalpy(istr-1,j,linew)* &
180 & grid(ng)%rmask_wet(istr-1,j)
181# endif
182 ti(istr-1,j,linew)=enthalpy(istr-1,j,linew)/ &
183 & max(hi(istr-1,j,linew),eps)
184 IF (hi(istr-1,j,linew).le.min_hi(ng)) THEN
185 enthalpy(istr-1,j,linew)=0.0_r8
186 ti(istr-1,j,linew)=0.0_r8
187 END IF
188 END DO
189!
190! Western edge, closed boundary condition.
191!
192 ELSE IF (lbc(iwest,ibice(istice),ng)%closed) THEN
193 DO j=jstr,jend
194 enthalpy(istr-1,j,linew)=hi(istr,j,linew)* &
195 & ti(istr,j,linew)
196# ifdef MASKING
197 enthalpy(istr-1,j,linew)=enthalpy(istr-1,j,linew)* &
198 & grid(ng)%rmask(istr-1,j)
199# endif
200# ifdef WET_DRY
201 enthalpy(istr-1,j,linew)=enthalpy(istr-1,j,linew)* &
202 & grid(ng)%rmask_wet(istr-1,j)
203# endif
204 ti(istr-1,j,linew)=enthalpy(istr-1,j,linew)/ &
205 & max(hi(istr-1,j,linew),eps)
206 IF (hi(istr-1,j,linew).le.min_hi(ng)) THEN
207 enthalpy(istr-1,j,linew)=0.0_r8
208 ti(istr-1,j,linew)=0.0_r8
209 END IF
210 END DO
211 END IF
212 END IF
213!
214!-----------------------------------------------------------------------
215! Lateral boundary conditions at the eastern edge.
216!-----------------------------------------------------------------------
217!
218 IF (domain(ng)%Eastern_Edge(tile)) THEN
219 IF (lbc(ieast,ibice(istice),ng)%clamped) THEN
220!
221! Eastern edge, clamped boundary condition.
222!
223 DO j=jstr,jend
224 enthalpy(iend+1,j,linew)=ice_lobc(ishice,ng)%ice_east(j)* &
225 & ice_lobc(istice,ng)%ice_east(j)
226# ifdef MASKING
227 enthalpy(iend+1,j,linew)=enthalpy(iend+1,j,linew)* &
228 & grid(ng)%rmask(iend+1,j)
229# endif
230# ifdef WET_DRY
231 enthalpy(iend+1,j,linew)=enthalpy(iend+1,j,linew)* &
232 & grid(ng)%rmask_wet(iend+1,j)
233# endif
234 ti(iend+1,j,linew)=enthalpy(iend+1,j,linew)/ &
235 & max(hi(iend+1,j,linew),eps)
236 IF (hi(iend+1,j,linew).le.min_hi(ng)) THEN
237 enthalpy(iend+1,j,linew)=0.0_r8
238 ti(iend+1,j,linew)=0.0_r8
239 END IF
240 END DO
241!
242! Eastern edge, clamped on inflow, gradient on outflow.
243!
244 ELSE IF (lbc(ieast,ibice(istice),ng)%mixed) THEN
245 DO j=jstr,jend
246 IF (ui(iend+1,j,linew).le.0.0_r8) THEN
247 enthalpy(iend+1,j,linew)=ice_lobc(ishice,ng)%ice_east(j)* &
248 & ice_lobc(istice,ng)%ice_east(j)
249# ifdef MASKING
250 enthalpy(iend+1,j,linew)=enthalpy(iend+1,j,linew)* &
251 & grid(ng)%rmask(iend+1,j)
252# endif
253# ifdef WET_DRY
254 enthalpy(iend+1,j,linew)=enthalpy(iend+1,j,linew)* &
255 & grid(ng)%rmask_wet(iend+1,j)
256# endif
257 ELSE
258 enthalpy(iend+1,j,linew)=hi(iend,j,liold)* &
259 & ti(iend,j,liold)
260# ifdef MASKING
261 enthalpy(iend+1,j,linew)=enthalpy(iend+1,j,linew)* &
262 & grid(ng)%rmask(iend+1,j)
263# endif
264# ifdef WET_DRY
265 enthalpy(iend+1,j,linew)=enthalpy(iend+1,j,linew)* &
266 & grid(ng)%rmask_wet(iend+1,j)
267# endif
268 ti(iend+1,j,linew)=enthalpy(iend+1,j,linew)/ &
269 & max(hi(iend+1,j,linew),eps)
270 IF (hi(iend+1,j,linew).le.min_hi(ng)) THEN
271 enthalpy(iend+1,j,linew)=0.0_r8
272 ti(iend+1,j,linew)=0.0_r8
273 END IF
274 END IF
275 END DO
276!
277! Eastern edge, gradient boundary condition.
278!
279 ELSE IF (lbc(ieast,ibice(istice),ng)%gradient) THEN
280 DO j=jstr,jend
281 enthalpy(iend+1,j,linew)=hi(iend,j,linew)* &
282 & ti(iend,j,linew)
283# ifdef MASKING
284 enthalpy(iend+1,j,linew)=enthalpy(iend+1,j,linew)* &
285 & grid(ng)%rmask(iend+1,j)
286# endif
287# ifdef WET_DRY
288 enthalpy(iend+1,j,linew)=enthalpy(iend+1,j,linew)* &
289 & grid(ng)%rmask_wet(iend+1,j)
290# endif
291 ti(iend+1,j,linew)=enthalpy(iend+1,j,linew)/ &
292 & max(hi(iend+1,j,linew),eps)
293 IF (hi(iend+1,j,linew).le.min_hi(ng)) THEN
294 enthalpy(iend+1,j,linew)=0.0_r8
295 ti(iend+1,j,linew)=0.0_r8
296 END IF
297 END DO
298!
299! Eastern edge, closed boundary condition.
300!
301 ELSE IF (lbc(ieast,ibice(istice),ng)%closed) THEN
302 DO j=jstr,jend
303 enthalpy(iend+1,j,linew)=hi(iend,j,linew)* &
304 & ti(iend,j,linew)
305# ifdef MASKING
306 enthalpy(iend+1,j,linew)=enthalpy(iend+1,j,linew)* &
307 & grid(ng)%rmask(iend+1,j)
308# endif
309# ifdef WET_DRY
310 enthalpy(iend+1,j,linew)=enthalpy(iend+1,j,linew)* &
311 & grid(ng)%rmask_wet(iend+1,j)
312# endif
313 ti(iend+1,j,linew)=enthalpy(iend+1,j,linew)/ &
314 & max(hi(iend+1,j,linew),eps)
315 IF (hi(iend+1,j,linew).le.min_hi(ng)) THEN
316 enthalpy(iend+1,j,linew)=0.0_r8
317 ti(iend+1,j,linew)=0.0_r8
318 END IF
319 END DO
320 END IF
321 END IF
322!
323!-----------------------------------------------------------------------
324! Lateral boundary conditions at the southern edge.
325!-----------------------------------------------------------------------
326!
327 IF (domain(ng)%Southern_Edge(tile)) THEN
328 IF (lbc(isouth,ibice(istice),ng)%clamped) THEN
329!
330! Southern edge, clamped boundary condition.
331!
332 DO i=istr,iend
333 enthalpy(i,jstr-1,linew)=ice_lobc(ishice,ng)%ice_south(i)* &
334 & ice_lobc(istice,ng)%ice_south(i)
335# ifdef MASKING
336 enthalpy(i,jstr-1,linew)=enthalpy(i,jstr-1,linew)* &
337 & grid(ng)%rmask(i,jstr-1)
338# endif
339# ifdef WET_DRY
340 enthalpy(i,jstr-1,linew)=enthalpy(i,jstr-1,linew)* &
341 & grid(ng)%rmask_wet(i,jstr-1)
342# endif
343 ti(i,jstr-1,linew)=enthalpy(i,jstr-1,linew)/ &
344 & max(hi(i,jstr-1,linew),eps)
345 IF (hi(i,jstr-1,linew).le.min_hi(ng)) THEN
346 enthalpy(i,jstr-1,linew)=0.0_r8
347 ti(i,jstr-1,linew)=0.0_r8
348 END IF
349 END DO
350!
351! Southern edge, clamped boundary condition.
352!
353 ELSE IF (lbc(isouth,ibice(istice),ng)%mixed) THEN
354 DO i=istr,iend
355 IF (vi(i,1,linew).ge.0.0_r8) THEN
356 enthalpy(i,jstr-1,linew)=ice_lobc(ishice,ng)%ice_south(i)*&
357 & ice_lobc(istice,ng)%ice_south(i)
358# ifdef MASKING
359 enthalpy(i,jstr-1,linew)=enthalpy(i,jstr-1,linew)* &
360 & grid(ng)%rmask(i,jstr-1)
361# endif
362# ifdef WET_DRY
363 enthalpy(i,jstr-1,linew)=enthalpy(i,jstr-1,linew)* &
364 & grid(ng)%rmask_wet(i,jstr-1)
365# endif
366 ELSE
367 enthalpy(i,jstr-1,linew)=enthalpy(i,jstr,liold)
368# ifdef MASKING
369 enthalpy(i,jstr-1,linew)=enthalpy(i,jstr-1,linew)* &
370 & grid(ng)%rmask(i,jstr-1)
371# endif
372# ifdef WET_DRY
373 enthalpy(i,jstr-1,linew)=enthalpy(i,jstr-1,linew)* &
374 & grid(ng)%rmask_wet(i,jstr-1)
375# endif
376 ti(i,jstr-1,linew)=enthalpy(i,jstr-1,linew)/ &
377 & max(hi(i,jstr-1,linew),eps)
378 IF (hi(i,jstr-1,linew).le.min_hi(ng)) THEN
379 enthalpy(i,jstr-1,linew)=0.0_r8
380 ti(i,jstr-1,linew)=0.0_r8
381 END IF
382 ENDIF
383 END DO
384!
385! Southern edge, gradient boundary condition.
386!
387 ELSE IF (lbc(isouth,ibice(istice),ng)%gradient) THEN
388 DO i=istr,iend
389 enthalpy(i,jstr-1,linew)=hi(i,jstr,linew)* &
390 & ti(i,jstr,linew)
391# ifdef MASKING
392 enthalpy(i,jstr-1,linew)=enthalpy(i,jstr-1,linew)* &
393 & grid(ng)%rmask(i,jstr-1)
394# endif
395# ifdef WET_DRY
396 enthalpy(i,jstr-1,linew)=enthalpy(i,jstr-1,linew)* &
397 & grid(ng)%rmask_wet(i,jstr-1)
398# endif
399 ti(i,jstr-1,linew)=enthalpy(i,jstr-1,linew)/ &
400 & max(hi(i,jstr-1,linew),eps)
401 IF (hi(i,jstr-1,linew).le.min_hi(ng)) THEN
402 enthalpy(i,jstr-1,linew)=0.0_r8
403 ti(i,jstr-1,linew)=0.0_r8
404 END IF
405 END DO
406!
407! Southern edge, closed boundary condition.
408!
409 ELSE IF (lbc(isouth,ibice(istice),ng)%closed) THEN
410 DO i=istr,iend
411 enthalpy(i,jstr-1,linew)=enthalpy(i,jstr,linew)
412# ifdef MASKING
413 enthalpy(i,jstr-1,linew)=enthalpy(i,jstr-1,linew)* &
414 & grid(ng)%rmask(i,jstr-1)
415# endif
416# ifdef WET_DRY
417 enthalpy(i,jstr-1,linew)=enthalpy(i,jstr-1,linew)* &
418 & grid(ng)%rmask_wet(i,jstr-1)
419# endif
420 ti(i,jstr-1,linew)=enthalpy(i,jstr-1,linew)/ &
421 & max(hi(i,jstr-1,linew),eps)
422 IF (hi(i,jstr-1,linew).le.min_hi(ng)) THEN
423 enthalpy(i,jstr-1,linew)=0.0_r8
424 ti(i,jstr-1,linew)=0.0_r8
425 END IF
426 END DO
427 END IF
428 END IF
429!
430!-----------------------------------------------------------------------
431! Lateral boundary conditions at the northern edge.
432!-----------------------------------------------------------------------
433!
434 IF (domain(ng)%Northern_Edge(tile)) THEN
435 IF (lbc(inorth,ibice(istice),ng)%clamped) THEN
436!
437! Northern edge, clamped boundary condition.
438!
439 DO i=istr,iend
440 enthalpy(i,jend+1,linew)=ice_lobc(ishice,ng)%ice_north(i)* &
441 & ice_lobc(istice,ng)%ice_north(i)
442# ifdef MASKING
443 enthalpy(i,jend+1,linew)=enthalpy(i,jend+1,linew)* &
444 & grid(ng)%rmask(i,jend+1)
445# endif
446# ifdef WET_DRY
447 enthalpy(i,jend+1,linew)=enthalpy(i,jend+1,linew)* &
448 & grid(ng)%rmask_wet(i,jend+1)
449# endif
450 ti(i,jend+1,linew)=enthalpy(i,jend+1,linew)/ &
451 & max(hi(i,jend+1,linew),eps)
452 IF (hi(i,jend+1,linew).le.min_hi(ng)) THEN
453 enthalpy(i,jend+1,linew)=0.0_r8
454 ti(i,jend+1,linew)=0.0_r8
455 END IF
456 END DO
457!
458! Northern edge, clamped on inflow, gradient on outflow.
459!
460 ELSE IF (lbc(inorth,ibice(istice),ng)%mixed) THEN
461 DO i=istr,iend
462 IF (vi(i,jend+1,linew).le.0.0_r8) THEN
463 enthalpy(i,jend+1,linew)=ice_lobc(ishice,ng)%ice_north(i)*&
464 & ice_lobc(istice,ng)%ice_north(i)
465# ifdef MASKING
466 enthalpy(i,jend+1,linew)=enthalpy(i,jend+1,linew)* &
467 & grid(ng)%rmask(i,jend+1)
468# endif
469# ifdef WET_DRY
470 enthalpy(i,jend+1,linew)=enthalpy(i,jend+1,linew)* &
471 & grid(ng)%rmask_wet(i,jend+1)
472# endif
473 ELSE
474 enthalpy(i,jend+1,linew)=enthalpy(i,jend,liold)
475# ifdef MASKING
476 enthalpy(i,jend+1,linew)=enthalpy(i,jend+1,linew)* &
477 & grid(ng)%rmask(i,jend+1)
478# endif
479# ifdef WET_DRY
480 enthalpy(i,jend+1,linew)=enthalpy(i,jend+1,linew)* &
481 & grid(ng)%rmask_wet(i,jend+1)
482# endif
483 ENDIF
484 ti(i,jend+1,linew)=enthalpy(i,jend+1,linew)/ &
485 & max(hi(i,jend+1,linew),eps)
486 IF (hi(i,jend+1,linew).le.min_hi(ng)) THEN
487 enthalpy(i,jend+1,linew)=0.0_r8
488 ti(i,jend+1,linew)=0.0_r8
489 END IF
490 END DO
491!
492! Northern edge, gradient boundary condition.
493!
494 ELSE IF (lbc(inorth,ibice(istice),ng)%gradient) THEN
495 DO i=istr,iend
496 enthalpy(i,jend+1,linew)=hi(i,jend,linew)* &
497 & ti(i,jend,linew)
498# ifdef MASKING
499 enthalpy(i,jend+1,linew)=enthalpy(i,jend+1,linew)* &
500 & grid(ng)%rmask(i,jend+1)
501# endif
502# ifdef WET_DRY
503 enthalpy(i,jend+1,linew)=enthalpy(i,jend+1,linew)* &
504 & grid(ng)%rmask_wet(i,jend+1)
505# endif
506# ifdef SOGLOBEC
507 ti(i,jend+1,linew)=enthalpy(i,jend+1,linew)/ &
508 & max(hi(i,jend,linew),eps)
509# else
510 ti(i,jend+1,linew)=enthalpy(i,jend+1,linew)/ &
511 & max(hi(i,jend+1,linew),eps)
512# endif
513 IF (hi(i,jend+1,linew).le.min_hi(ng)) THEN
514 enthalpy(i,jend+1,linew)=0.0_r8
515 ti(i,jend+1,linew)=0.0_r8
516 END IF
517 END DO
518!
519! Northern edge, closed boundary condition.
520!
521 ELSE IF (lbc(inorth,ibice(istice),ng)%closed) THEN
522 DO i=istr,iend
523 enthalpy(i,jend+1,linew)=hi(i,jend,linew)* &
524 & ti(i,jend,linew)
525# ifdef MASKING
526 enthalpy(i,jend+1,linew)=enthalpy(i,jend+1,linew)* &
527 & grid(ng)%rmask(i,jend+1)
528# endif
529# ifdef WET_DRY
530 enthalpy(i,jend+1,linew)=enthalpy(i,jend+1,linew)* &
531 & grid(ng)%rmask_wet(i,jend+1)
532# endif
533 ti(i,jend+1,linew)=enthalpy(i,jend+1,linew)/ &
534 & max(hi(i,jend+1,linew),eps)
535 IF (hi(i,jend+1,linew).le.min_hi(ng)) THEN
536 enthalpy(i,jend+1,linew)=0.0_r8
537 ti(i,jend+1,linew)=0.0_r8
538 END IF
539 END DO
540 END IF
541 END IF
542!
543!-----------------------------------------------------------------------
544! Boundary corners.
545!-----------------------------------------------------------------------
546!
547 IF (.not.(ewperiodic(ng).or.nsperiodic(ng))) THEN
548 IF (domain(ng)%SouthWest_Corner(tile)) THEN
549 enthalpy(istr-1,jstr-1,linew)=0.5_r8* &
550 & (enthalpy(istr ,jstr-1,linew)+ &
551 & enthalpy(istr-1,jstr ,linew))
552# ifdef MASKING
553 enthalpy(istr-1,jstr-1,linew)=enthalpy(istr-1,jstr-1,linew)* &
554 & grid(ng)%rmask(istr-1,jstr-1)
555# endif
556# ifdef WET_DRY
557 enthalpy(istr-1,jstr-1,linew)=enthalpy(istr-1,jstr-1,linew)* &
558 & grid(ng)%rmask_wet(istr-1,jstr-1)
559# endif
560 ti(istr-1,jstr-1,linew)=enthalpy(istr-1,jstr-1,linew)/ &
561 & max(hi(istr-1,jstr-1,linew),eps)
562 IF (hi(istr-1,jstr-1,linew).le.min_hi(ng)) THEN
563 enthalpy(istr-1,jstr-1,linew)=0.0_r8
564 ti(istr-1,jstr-1,linew)=0.0_r8
565 END IF
566 END IF
567 IF (domain(ng)%SouthEast_Corner(tile)) THEN
568 enthalpy(iend+1,jstr-1,linew)=0.5_r8* &
569 & (enthalpy(iend+1,jstr ,linew)+ &
570 & enthalpy(iend ,jstr-1,linew))
571# ifdef MASKING
572 enthalpy(iend+1,jstr-1,linew)=enthalpy(iend+1,jstr-1,linew)* &
573 & grid(ng)%rmask(iend+1,jstr-1)
574# endif
575# ifdef WET_DRY
576 enthalpy(iend+1,jstr-1,linew)=enthalpy(iend+1,jstr-1,linew)* &
577 & grid(ng)%rmask_wet(iend+1,jstr-1)
578# endif
579 ti(iend+1,jstr-1,linew)=enthalpy(iend+1,jstr-1,linew)/ &
580 & max(hi(iend+1,jstr-1,linew),eps)
581 IF (hi(iend+1,jstr-1,linew).LE.min_hi(ng)) THEN
582 enthalpy(iend+1,jstr-1,linew)=0.0_r8
583 ti(iend+1,jstr-1,linew)=0.0_r8
584 END IF
585 END IF
586 IF (domain(ng)%NorthWest_Corner(tile)) THEN
587 enthalpy(istr-1,jend+1,linew)=0.5_r8* &
588 & (enthalpy(istr-1,jend ,linew)+ &
589 & enthalpy(istr ,jend+1,linew))
590# ifdef MASKING
591 enthalpy(istr-1,jend+1,linew)=enthalpy(istr-1,jend+1,linew)* &
592 & grid(ng)%rmask(istr-1,jend+1)
593# endif
594# ifdef WET_DRY
595 enthalpy(istr-1,jend+1,linew)=enthalpy(istr-1,jend+1,linew)* &
596 & grid(ng)%rmask_wet(istr-1,jend+1)
597# endif
598 ti(istr-1,jend+1,linew)=enthalpy(istr-1,jend+1,linew)/ &
599 & max(hi(istr-1,jend+1,linew),eps)
600 IF (hi(istr-1,jend+1,linew).le.min_hi(ng)) THEN
601 enthalpy(istr-1,jend+1,linew)=0.0_r8
602 ti(istr-1,jend+1,linew)=0.0_r8
603 END IF
604 END IF
605 IF (domain(ng)%NorthEast_Corner(tile)) THEN
606 enthalpy(iend+1,jend+1,linew)=0.5_r8* &
607 & (enthalpy(iend+1,jend ,linew)+ &
608 & enthalpy(iend ,jend+1,linew))
609# ifdef MASKING
610 enthalpy(iend+1,jend+1,linew)=enthalpy(iend+1,jend+1,linew)* &
611 & grid(ng)%rmask(iend+1,jend+1)
612# endif
613# ifdef WET_DRY
614 enthalpy(iend+1,jend+1,linew)=enthalpy(iend+1,jend+1,linew)* &
615 & grid(ng)%rmask_wet(iend+1,jend+1)
616# endif
617 ti(iend+1,jend+1,linew)=enthalpy(iend+1,jend+1,linew)/ &
618 & max(hi(iend+1,jend+1,linew),eps)
619 IF (hi(iend+1,jend+1,linew).le.min_hi(ng)) THEN
620 enthalpy(iend+1,jend+1,linew)=0.0_r8
621 ti(iend+1,jend+1,linew)=0.0_r8
622 END IF
623 END IF
624 END IF
625!
626 RETURN
627 END SUBROUTINE ice_tibc_tile
628#endif
629 END MODULE ice_tibc_mod
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer, parameter isvice
Definition ice_mod.h:147
real(r8), dimension(:), allocatable min_hi
Definition ice_mod.h:246
integer, parameter isenth
Definition ice_mod.h:148
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, parameter istice
Definition ice_mod.h:145
integer, dimension(nices) ibice
Definition ice_mod.h:162
integer, parameter isuice
Definition ice_mod.h:146
integer, parameter ishice
Definition ice_mod.h:138
type(t_lbc), dimension(:,:,:), allocatable lbc
Definition mod_param.F:375
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
integer, parameter isouth
integer, parameter ieast
integer, parameter inorth
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