ROMS
Loading...
Searching...
No Matches
ice_smolar.h
Go to the documentation of this file.
1 MODULE ice_advect_mod
2!
3!git $Id
4!=======================================================================
5! Copyright (c) 2002-2025 The ROMS Group !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md W. Paul Budgell !
8!================================================== Hernan G. Arango ===
9! !
10! This subroutine performs advection of ice scalars using the !
11! Smolarkiewicz second-order upwind scheme. !
12! !
13! Reference: !
14! !
15! Smolarkiewicz, P.K. and W.W. Grabowski, 1990: The multidimensional !
16! Positive definite advection transport algorithm: Nonoscillatory !
17! option, J. Comp. Phys., 86, 355-375. !
18! !
19!=======================================================================
20!
21 USE mod_param
22 USE mod_boundary
23 USE mod_forces
24 USE mod_grid
25 USE mod_ice
26 USE mod_ocean
27 USE mod_scalars
28!
29 USE exchange_2d_mod, ONLY : exchange_r2d_tile
30#ifdef DISTRIBUTE
31 USE mp_exchange_mod, ONLY : mp_exchange2d
32#endif
33 USE ice_bc2d_mod, ONLY : ice_bc2d_tile
34 USE ice_tibc_mod, ONLY : ice_tibc_tile
35!
36 implicit none
37!
38 PUBLIC :: ice_advect
39 PRIVATE :: ice_advect_tile
40 PRIVATE :: ice_mpdata_tile
41!
42 CONTAINS
43!
44!***********************************************************************
45 SUBROUTINE ice_advect (ng, tile, model)
46!***********************************************************************
47!
48 USE mod_stepping
49!
50! Imported variable declarations.
51!
52 integer, intent(in) :: ng, tile, model
53!
54! Local variable declarations.
55!
56 character (len=*), parameter :: myfile = &
57 & __FILE__
58!
59#include "tile.h"
60!
61#ifdef PROFILE
62 CALL wclock_on (ng, model, 42, __line__, myfile)
63#endif
64 CALL ice_advect_tile (ng, tile, model, &
65 & lbi, ubi, lbj, ubj, &
66 & imins, imaxs, jmins, jmaxs, &
67 & linew(ng), liold(ng), liunw(ng))
68#ifdef PROFILE
69 CALL wclock_off (ng, model, 42, __line__, myfile)
70#endif
71!
72 RETURN
73 END SUBROUTINE ice_advect
74!
75!***********************************************************************
76 SUBROUTINE ice_advect_tile (ng, tile, model, &
77 & LBi, UBi, LBj, UBj, &
78 & IminS, ImaxS, JminS, JmaxS, &
79 & linew, liold, liunw)
80!***********************************************************************
81!
82! Imported variable declarations.
83!
84 integer, intent(in) :: ng, tile, model
85 integer, intent(in) :: lbi, ubi, lbj, ubj
86 integer, intent(in) :: imins, imaxs, jmins, jmaxs
87 integer, intent(in) :: linew, liold, liunw
88
89!
90! Local variable declarations.
91!
92 integer :: i, j
93
94#include "set_bounds.h"
95!
96!-----------------------------------------------------------------------
97! Advect the ice concentation, isAice.
98!-----------------------------------------------------------------------
99!
100 CALL ice_mpdata_tile (ng, tile, model, &
101 & lbi, ubi, lbj, ubj, &
102 & imins, imaxs, jmins, jmaxs, &
103 & linew, liold, liunw, &
104#ifdef MASKING
105 & grid(ng) % rmask, &
106 & grid(ng) % umask, &
107 & grid(ng) % vmask, &
108#endif
109#ifdef WET_DRY
110 & grid(ng) % rmask_wet, &
111 & grid(ng) % umask_wet, &
112 & grid(ng) % vmask_wet, &
113#endif
114#ifdef ICESHELF
115 & grid(ng) % zice, &
116#endif
117#ifndef ICE_UPWIND
118 & grid(ng) % pm, &
119 & grid(ng) % pn, &
120#endif
121 & grid(ng) % on_u, &
122 & grid(ng) % om_v, &
123 & grid(ng) % omn, &
124 & ice(ng) % Si(:,:,:,isuice), &
125 & ice(ng) % Si(:,:,:,isvice), &
126 & ice(ng) % Si(:,:,:,isaice))
127!
128 CALL ice_bc2d_tile (ng, tile, model, isaice, &
129 & lbi, ubi, lbj, ubj, &
130 & imins, imaxs, jmins, jmaxs, &
131 & liold, linew, &
132 & ice(ng) % Si(:,:,:,isuice), &
133 & ice(ng) % Si(:,:,:,isvice), &
134 & ice(ng) % Si(:,:,:,isaice), &
135 & lbc(:,ibice(isaice),ng))
136!
137!-----------------------------------------------------------------------
138! Advect the ice thickness, isHice.
139!-----------------------------------------------------------------------
140!
141 CALL ice_mpdata_tile (ng, tile, model, &
142 & lbi, ubi, lbj, ubj, &
143 & imins, imaxs, jmins, jmaxs, &
144 & linew, liold, liunw, &
145#ifdef MASKING
146 & grid(ng) % rmask, &
147 & grid(ng) % umask, &
148 & grid(ng) % vmask, &
149#endif
150#ifdef WET_DRY
151 & grid(ng) % rmask_wet, &
152 & grid(ng) % umask_wet, &
153 & grid(ng) % vmask_wet, &
154#endif
155#ifdef ICESHELF
156 & grid(ng) % zice, &
157#endif
158#ifndef ICE_UPWIND
159 & grid(ng) % pm, &
160 & grid(ng) % pn, &
161#endif
162 & grid(ng) % on_u, &
163 & grid(ng) % om_v, &
164 & grid(ng) % omn, &
165 & ice(ng) % Si(:,:,:,isuice), &
166 & ice(ng) % Si(:,:,:,isvice), &
167 & ice(ng) % Si(:,:,:,ishice))
168!
169 CALL ice_bc2d_tile (ng, tile, model, ishice, &
170 & lbi, ubi, lbj, ubj, &
171 & imins, imaxs, jmins, jmaxs, &
172 & liold, linew, &
173 & ice(ng) % Si(:,:,:,isuice), &
174 & ice(ng) % Si(:,:,:,isvice), &
175 & ice(ng) % Si(:,:,:,ishice), &
176 & lbc(:,ibice(ishice),ng))
177!
178! Compute rate of ice divergence (m3/s).
179!
180 DO j=jstrt,jendt
181 DO i=istrt,iendt
182 ice(ng)%Fi(i,j,icwdiv)=(ice(ng)%Si(i,j,linew,ishice)- &
183 & ice(ng)%Si(i,j,liold,ishice))/dt(ng)
184 END DO
185 END DO
186
187#ifdef ICE_THERMO
188!
189!-----------------------------------------------------------------------
190! Advect the snow thickness, isHsno.
191!-----------------------------------------------------------------------
192!
193 CALL ice_mpdata_tile (ng, tile, model, &
194 & lbi, ubi, lbj, ubj, &
195 & imins, imaxs, jmins, jmaxs, &
196 & linew, liold, liunw, &
197# ifdef MASKING
198 & grid(ng) % rmask, &
199 & grid(ng) % umask, &
200 & grid(ng) % vmask, &
201# endif
202# ifdef WET_DRY
203 & grid(ng) % rmask_wet, &
204 & grid(ng) % umask_wet, &
205 & grid(ng) % vmask_wet, &
206# endif
207# ifdef ICESHELF
208 & grid(ng) % zice, &
209# endif
210# ifndef ICE_UPWIND
211 & grid(ng) % pm, &
212 & grid(ng) % pn, &
213# endif
214 & grid(ng) % on_u, &
215 & grid(ng) % om_v, &
216 & grid(ng) % omn, &
217 & ice(ng) % Si(:,:,:,isuice), &
218 & ice(ng) % Si(:,:,:,isvice), &
219 & ice(ng) % Si(:,:,:,ishsno))
220!
221 CALL ice_bc2d_tile (ng, tile, model, ishsno, &
222 & lbi, ubi, lbj, ubj, &
223 & imins, imaxs, jmins, jmaxs, &
224 & liold, linew, &
225 & ice(ng) % Si(:,:,:,isuice), &
226 & ice(ng) % Si(:,:,:,isvice), &
227 & ice(ng) % Si(:,:,:,ishsno), &
228 & lbc(:,ibice(ishsno),ng))
229!
230!-----------------------------------------------------------------------
231! Advect the surface melt water thickness, isHmel.
232!-----------------------------------------------------------------------
233!
234 CALL ice_mpdata_tile (ng, tile, model, &
235 & lbi, ubi, lbj, ubj, &
236 & imins, imaxs, jmins, jmaxs, &
237 & linew, liold, liunw, &
238# ifdef MASKING
239 & grid(ng) % rmask, &
240 & grid(ng) % umask, &
241 & grid(ng) % vmask, &
242# endif
243# ifdef WET_DRY
244 & grid(ng) % rmask_wet, &
245 & grid(ng) % umask_wet, &
246 & grid(ng) % vmask_wet, &
247# endif
248# ifdef ICESHELF
249 & grid(ng) % zice, &
250# endif
251# ifndef ICE_UPWIND
252 & grid(ng) % pm, &
253 & grid(ng) % pn, &
254# endif
255 & grid(ng) % on_u, &
256 & grid(ng) % om_v, &
257 & grid(ng) % omn, &
258 & ice(ng) % Si(:,:,:,isuice), &
259 & ice(ng) % Si(:,:,:,isvice), &
260 & ice(ng) % Si(:,:,:,ishmel))
261!
262 CALL ice_bc2d_tile (ng, tile, model, ishmel, &
263 & lbi, ubi, lbj, ubj, &
264 & imins, imaxs, jmins, jmaxs, &
265 & liold, linew, &
266 & ice(ng) % Si(:,:,:,isuice), &
267 & ice(ng) % Si(:,:,:,isvice), &
268 & ice(ng) % Si(:,:,:,ishmel), &
269 & lbc(:,ibice(ishmel),ng))
270!
271!-----------------------------------------------------------------------
272! Advect the interior ice temperature, isTice.
273!-----------------------------------------------------------------------
274!
275 CALL ice_mpdata_tile (ng, tile, model, &
276 & lbi, ubi, lbj, ubj, &
277 & imins, imaxs, jmins, jmaxs, &
278 & linew, liold, liunw, &
279# ifdef MASKING
280 & grid(ng) % rmask, &
281 & grid(ng) % umask, &
282 & grid(ng) % vmask, &
283# endif
284# ifdef WET_DRY
285 & grid(ng) % rmask_wet, &
286 & grid(ng) % umask_wet, &
287 & grid(ng) % vmask_wet, &
288# endif
289# ifdef ICESHELF
290 & grid(ng) % zice, &
291# endif
292# ifndef ICE_UPWIND
293 & grid(ng) % pm, &
294 & grid(ng) % pn, &
295# endif
296 & grid(ng) % on_u, &
297 & grid(ng) % om_v, &
298 & grid(ng) % omn, &
299 & ice(ng) % Si(:,:,:,isuice), &
300 & ice(ng) % Si(:,:,:,isvice), &
301 & ice(ng) % Si(:,:,:,isenth))
302!
303 DO j=jstrt,jendt
304 DO i=istrt,iendt
305 ice(ng)%Si(i,j,linew,istice)=ice(ng)%Si(i,j,linew,isenth)/ &
306 & max(ice(ng)%Si(i,j,linew,ishice),&
307 & 1.0e-6_r8)
308 IF (ice(ng)%Si(i,j,linew,ishice).le.min_hi(ng)) THEN
309 ice(ng)%Si(i,j,linew,isenth)=0.0_r8
310 ice(ng)%Si(i,j,linew,istice)=0.0_r8
311 END IF
312 END DO
313 END DO
314!
315 CALL ice_tibc_tile (ng, tile, model, &
316 & lbi, ubi, lbj, ubj, &
317 & liold, linew, &
318 & ice(ng) % Si(:,:,:,isuice), &
319 & ice(ng) % Si(:,:,:,isvice), &
320 & ice(ng) % Si(:,:,:,ishice), &
321 & ice(ng) % Si(:,:,:,istice), &
322 & ice(ng) % Si(:,:,:,isenth))
323!
324!-----------------------------------------------------------------------
325! Advect thickness associated with age of ice, isHage. Then, compute
326! ice age (isIage).
327!-----------------------------------------------------------------------
328!
329 DO j=jstrt,jendt
330 DO i=istrt,iendt
331 ice(ng)%Si(i,j,liold,ishage)=ice(ng)%Si(i,j,liold,ishice)* &
332 & ice(ng)%Si(i,j,liold,isiage)
333 ice(ng)%Si(i,j,linew,ishage)=ice(ng)%Si(i,j,linew,ishice)* &
334 & ice(ng)%Si(i,j,linew,isiage)
335 IF (ice(ng)%Si(i,j,liold,ishice).le.min_hi(ng)) THEN
336 ice(ng)%Si(i,j,liold,ishage)=0.0_r8
337 ice(ng)%Si(i,j,liold,isiage)=0.0_r8
338 END IF
339 END DO
340 END DO
341!
342 CALL ice_mpdata_tile (ng, tile, model, &
343 & lbi, ubi, lbj, ubj, &
344 & imins, imaxs, jmins, jmaxs, &
345 & linew, liold, liunw, &
346# ifdef MASKING
347 & grid(ng) % rmask, &
348 & grid(ng) % umask, &
349 & grid(ng) % vmask, &
350# endif
351# ifdef WET_DRY
352 & grid(ng) % rmask_wet, &
353 & grid(ng) % umask_wet, &
354 & grid(ng) % vmask_wet, &
355# endif
356# ifdef ICESHELF
357 & grid(ng) % zice, &
358# endif
359# ifndef ICE_UPWIND
360 & grid(ng) % pm, &
361 & grid(ng) % pn, &
362# endif
363 & grid(ng) % on_u, &
364 & grid(ng) % om_v, &
365 & grid(ng) % omn, &
366 & ice(ng) % Si(:,:,:,isuice), &
367 & ice(ng) % Si(:,:,:,isvice), &
368 & ice(ng) % Si(:,:,:,ishage))
369!
370 DO j=jstrt,jendt
371 DO i=istrt,iendt
372 ice(ng)%Si(i,j,linew,isiage)=ice(ng)%Si(i,j,linew,ishage)/ &
373 & max(ice(ng)%Si(i,j,linew,ishice),&
374 & 1.0e-6_r8)
375 IF (ice(ng)%Si(i,j,linew,ishice).le.min_hi(ng)) THEN
376 ice(ng)%Si(i,j,linew,ishage)=0.0_r8
377 ice(ng)%Si(i,j,linew,isiage)=0.0_r8
378 END IF
379 END DO
380 END DO
381!
382 CALL ice_bc2d_tile (ng, tile, model, isiage, &
383 & lbi, ubi, lbj, ubj, &
384 & imins, imaxs, jmins, jmaxs, &
385 & liold, linew, &
386 & ice(ng) % Si(:,:,:,isuice), &
387 & ice(ng) % Si(:,:,:,isvice), &
388 & ice(ng) % Si(:,:,:,isiage), &
389 & lbc(:,ibice(isiage),ng))
390#endif
391
392#if defined ICE_THERMO && defined ICE_BIO
393!
394!-----------------------------------------------------------------------
395! Advect the ice algae concentration, IcePhL.
396!-----------------------------------------------------------------------
397!
398 CALL ice_mpdata_tile (ng, tile, model, &
399 & lbi, ubi, lbj, ubj, &
400 & imins, imaxs, jmins, jmaxs, &
401 & linew, liold, liunw, &
402# ifdef MASKING
403 & grid(ng) % rmask, &
404 & grid(ng) % rmask, &
405 & grid(ng) % vmask, &
406# endif
407# ifdef WET_DRY
408 & grid(ng) % rmask_wet, &
409 & grid(ng) % umask_wet, &
410 & grid(ng) % vmask_wet, &
411# endif
412# ifdef ICESHELF
413 & grid(ng) % zice, &
414# endif
415# ifndef ICE_UPWIND
416 & grid(ng) % pm, &
417 & grid(ng) % pn, &
418# endif
419 & grid(ng) % on_u, &
420 & grid(ng) % om_v, &
421 & grid(ng) % omn, &
422 & ice(ng) % Si(:,:,:,isuice), &
423 & ice(ng) % Si(:,:,:,isvice), &
424 & ice(ng) % Si(:,:,:,isiphy))
425!
426! Need to change this to ice_bc2d_tile calls.
427!
428 CALL ice_bc2d_tile (ng, tile, model, isiphy, &
429 & lbi, ubi, lbj, ubj, &
430 & liold, linew, &
431 & ice(ng) % Si(:,:,:,isuice), &
432 & ice(ng) % Si(:,:,:,isvice), &
433 & ice(ng) % Si(:,:,:,isiphy), &
434 & lbc(:,ibice(isiphy),ng))
435!
436!-----------------------------------------------------------------------
437! Advect the ice nitrate concentration, isINO3.
438!-----------------------------------------------------------------------
439!
440 CALL ice_mpdata_tile (ng, tile, model, &
441 & lbi, ubi, lbj, ubj, &
442 & imins, imaxs, jmins, jmaxs, &
443 & linew, liold, liunw, &
444# ifdef MASKING
445 & grid(ng) % rmask, &
446 & grid(ng) % umask, &
447 & grid(ng) % vmask, &
448# endif
449# ifdef WET_DRY
450 & grid(ng) % rmask_wet, &
451 & grid(ng) % umask_wet, &
452 & grid(ng) % vmask_wet, &
453# endif
454# ifdef ICESHELF
455 & grid(ng) % zice, &
456# endif
457# ifndef ICE_UPWIND
458 & grid(ng) % pm, &
459 & grid(ng) % pn, &
460# endif
461 & grid(ng) % on_u, &
462 & grid(ng) % om_v, &
463 & grid(ng) % omn, &
464 & ice(ng) % Si(:,:,:,isuice), &
465 & ice(ng) % Si(:,:,:,isvice), &
466 & ice(ng) % Si(:,:,:,isino3))
467!
468 CALL ice_bc2d_tile (ng, tile, model, isino3, &
469 & lbi, ubi, lbj, ubj, &
470 & liold, linew, &
471 & ice(ng) % Si(:,:,:,isuice), &
472 & ice(ng) % Si(:,:,:,isvice), &
473 & ice(ng) % Si(:,:,:,isino3), &
474 & lbc(:,ibice(isino3),ng))
475!
476!-----------------------------------------------------------------------
477! Advect the ice ammonia concentration, isINH4.
478!-----------------------------------------------------------------------
479!
480 CALL ice_mpdata_tile (ng, tile, model, &
481 & lbi, ubi, lbj, ubj, &
482 & imins, imaxs, jmins, jmaxs, &
483 & linew, liold, liunw, &
484# ifdef MASKING
485 & grid(ng) % rmask, &
486 & grid(ng) % umask, &
487 & grid(ng) % vmask, &
488# endif
489# ifdef WET_DRY
490 & grid(ng) % rmask_wet, &
491 & grid(ng) % umask_wet, &
492 & grid(ng) % vmask_wet, &
493# endif
494# ifdef ICESHELF
495 & grid(ng) % zice, &
496# endif
497# ifndef ICE_UPWIND
498 & grid(ng) % pm, &
499 & grid(ng) % pn, &
500# endif
501 & grid(ng) % on_u, &
502 & grid(ng) % om_v, &
503 & grid(ng) % omn, &
504 & ice(ng) % Si(:,:,:,isuice), &
505 & ice(ng) % Si(:,:,:,isvice), &
506 & ice(ng) % Si(:,:,:,isinh4))
507!
508 CALL ice_bc2d_tile (ng, tile, model, isinh4, &
509 & lbi, ubi, lbj, ubj, &
510 & liold, linew, &
511 & ice(ng) % Si(:,:,:,isuice), &
512 & ice(ng) % Si(:,:,:,isvice), &
513 & ice(ng) % Si(:,:,:,isinh4), &
514 & lbc(:,ibice(isinh4),ng))
515#endif
516!
517!-----------------------------------------------------------------------
518! Exchange boundary information.
519!-----------------------------------------------------------------------
520!
521 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
522 CALL exchange_r2d_tile (ng, tile, &
523 & lbi, ubi, lbj, ubj, &
524 & ice(ng)%Si(:,:,linew,isaice))
525
526 CALL exchange_r2d_tile (ng, tile, &
527 & lbi, ubi, lbj, ubj, &
528 & ice(ng)%Si(:,:,linew,ishice))
529
530#ifdef ICE_THERMO
531 CALL exchange_r2d_tile (ng, tile, &
532 & lbi, ubi, lbj, ubj, &
533 & ice(ng)%Si(:,:,linew,ishsno))
534
535 CALL exchange_r2d_tile (ng, tile, &
536 & lbi, ubi, lbj, ubj, &
537 & ice(ng)%Si(:,:,linew,ishmel))
538
539 CALL exchange_r2d_tile (ng, tile, &
540 & lbi, ubi, lbj, ubj, &
541 & ice(ng)%Si(:,:,linew,istice))
542
543 CALL exchange_r2d_tile (ng, tile, &
544 & lbi, ubi, lbj, ubj, &
545 & ice(ng)%Si(:,:,linew,isenth))
546
547 CALL exchange_r2d_tile (ng, tile, &
548 & lbi, ubi, lbj, ubj, &
549 & ice(ng)%Si(:,:,linew,isiage))
550
551 CALL exchange_r2d_tile (ng, tile, &
552 & lbi, ubi, lbj, ubj, &
553 & ice(ng)%Si(:,:,linew,ishage))
554
555# ifdef ICE_BIO
556 CALL exchange_r2d_tile (ng, tile, &
557 & lbi, ubi, lbj, ubj, &
558 & ice(ng)%Si(:,:,linew,isiphy))
559
560 CALL exchange_r2d_tile (ng, tile, &
561 & lbi, ubi, lbj, ubj, &
562 & ice(ng)%Si(:,:,linew,isino3))
563
564 CALL exchange_r2d_tile (ng, tile, &
565 & lbi, ubi, lbj, ubj, &
566 & ice(ng)%Si(:,:,linew,isinh4))
567# endif
568#endif
569 END IF
570
571#ifdef DISTRIBUTE
572!
573 CALL mp_exchange2d (ng, tile, model, 2, &
574 & lbi, ubi, lbj, ubj, &
575 & nghostpoints, ewperiodic(ng), nsperiodic(ng), &
576 & ice(ng)%Si(:,:,linew,isaice), &
577 & ice(ng)%Si(:,:,linew,ishice))
578
579# ifdef ICE_THERMO
580 CALL mp_exchange2d (ng, tile, model, 2, &
581 & lbi, ubi, lbj, ubj, &
582 & nghostpoints, ewperiodic(ng), nsperiodic(ng), &
583 & ice(ng)%Si(:,:,linew,ishsno), &
584 & ice(ng)%Si(:,:,linew,ishmel))
585
586 CALL mp_exchange2d (ng, tile, model, 4, &
587 & lbi, ubi, lbj, ubj, &
588 & nghostpoints, ewperiodic(ng), nsperiodic(ng), &
589 & ice(ng)%Si(:,:,linew,istice), &
590 & ice(ng)%Si(:,:,linew,isenth), &
591 & ice(ng)%Si(:,:,linew,isiage), &
592 & ice(ng)%Si(:,:,linew,ishage))
593
594# if defined ICE_BIO
595 CALL mp_exchange2d (ng, tile, model, 3, &
596 & lbi, ubi, lbj, ubj, &
597 & nghostpoints, ewperiodic(ng), nsperiodic(ng), &
598 & ice(ng)%Si(:,:,linew,isiphy), &
599 & ice(ng)%Si(:,:,linew,isino3), &
600 & ice(ng)%Si(:,:,linew,isinh4))
601# endif
602# endif
603#endif
604!
605 RETURN
606 END SUBROUTINE ice_advect_tile
607!
608!***********************************************************************
609 SUBROUTINE ice_mpdata_tile (ng, tile, model, &
610 & LBi, UBi, LBj, UBj, &
611 & IminS, ImaxS, JminS, JmaxS, &
612 & linew, liold, liunw, &
613#ifdef MASKING
614 & rmask, umask, vmask, &
615#endif
616#ifdef WET_DRY
617 & rmask_wet, umask_wet, vmask_wet, &
618#endif
619#ifdef ICESHELF
620 & zice, &
621#endif
622#ifndef ICE_UPWIND
623 & pm, pn, &
624#endif
625 & on_u, om_v, omn, &
626 & ui, vi, field)
627!***********************************************************************
628!
629! Imported variable declarations.
630!
631 integer, intent(in) :: ng, tile, model
632 integer, intent(in) :: lbi, ubi, lbj, ubj
633 integer, intent(in) :: imins, imaxs, jmins, jmaxs
634 integer, intent(in) :: linew, liold, liunw
635!
636#ifdef ASSUMED_SHAPE
637# ifdef MASKING
638 real(r8), intent(in) :: rmask(lbi:,lbj:)
639 real(r8), intent(in) :: umask(lbi:,lbj:)
640 real(r8), intent(in) :: vmask(lbi:,lbj:)
641# endif
642# ifdef WET_DRY
643 real(r8), intent(in) :: rmask_wet(lbi:,lbj:)
644 real(r8), intent(in) :: umask_wet(lbi:,lbj:)
645 real(r8), intent(in) :: vmask_wet(lbi:,lbj:)
646# endif
647# ifdef ICESHELF
648 real(r8), intent(in) :: zice(lbi:,lbj:)
649# endif
650# ifndef ICE_UPWIND
651 real(r8), intent(in) :: pm(lbi:,lbj:)
652 real(r8), intent(in) :: pn(lbi:,lbj:)
653# endif
654 real(r8), intent(in) :: on_u(lbi:,lbj:)
655 real(r8), intent(in) :: om_v(lbi:,lbj:)
656 real(r8), intent(in) :: omn(lbi:,lbj:)
657 real(r8), intent(in) :: ui(lbi:,lbj:,:)
658 real(r8), intent(in) :: vi(lbi:,lbj:,:)
659 real(r8), intent(inout) :: field(lbi:,lbj:,:)
660#else
661# ifdef MASKING
662 real(r8), intent(in) :: rmask(lbi:ubi,lbj:ubj)
663 real(r8), intent(in) :: umask(lbi:ubi,lbj:ubj)
664 real(r8), intent(in) :: vmask(lbi:ubi,lbj:ubj)
665# endif
666# ifdef WET_DRY
667 real(r8), intent(in) :: rmask_wet(lbi:ubi,lbj:ubj)
668 real(r8), intent(in) :: umask_wet(lbi:ubi,lbj:ubj)
669 real(r8), intent(in) :: vmask_wet(lbi:ubi,lbj:ubj)
670# endif
671# ifdef ICESHELF
672 real(r8), intent(in) :: zice(lbi:ubi,lbj:ubj)
673# endif
674# ifndef ICE_UPWIND
675 real(r8), intent(in) :: pm(lbi:ubi,lbj:ubj)
676 real(r8), intent(in) :: pn(lbi:ubi,lbj:ubj)
677# endif
678 real(r8), intent(in) :: on_u(lbi:ubi,lbj:ubj)
679 real(r8), intent(in) :: om_v(lbi:ubi,lbj:ubj)
680 real(r8), intent(in) :: omn(lbi:ubi,lbj:ubj)
681 real(r8), intent(in) :: ui(lbi:ubi,lbj:ubj,2)
682 real(r8), intent(in) :: vi(lbi:ubi,lbj:ubj,2)
683 real(r8), intent(inout) :: field(lbi:ubi,lbj:ubj,2)
684#endif
685!
686! Local variable definitions
687!
688 integer :: imin, imax, jmin, jmax
689 integer :: i, j
690!
691 real(r8) :: cu_crss, cu
692 real(r8) :: cff1, cff2, rateu, ratev, rateyiu, ratexiv
693 real(r8) :: uspeed, vspeed
694!
695 real(r8), parameter :: epsil = 1.0e-15_r8
696 real(r8), parameter :: add = 3.0e+3_r8
697!
698 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ar
699 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: aflxu
700 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: aflxv
701 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: aif
702 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: fx
703 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: fe
704
705#include "set_bounds.h"
706!
707!-----------------------------------------------------------------------
708! Horizontal advection of generic ice field.
709!-----------------------------------------------------------------------
710!
711#ifndef ICE_UPWIND
712 IF (ewperiodic(ng)) THEN
713 imin=istr-1
714 imax=iend+1
715 ELSE
716 imin=max(istr-1,1)
717 imax=min(iend+1,lm(ng))
718 END IF
719 IF (nsperiodic(ng)) THEN
720 jmin=jstr-1
721 jmax=jend+1
722 ELSE
723 jmin=max(jstr-1,1)
724 jmax=min(jend+1,mm(ng))
725 END IF
726#else
727 imin=istr
728 imax=iend
729 jmin=jstr
730 jmax=jend
731#endif
732!
733! Upstream fluxes. (HGA: the advection is not in flux form?)
734!
735 DO j=jmin,jmax
736 DO i=imin,imax+1
737 cff1=max(0.0_r8,ui(i,j,liunw))
738 cff2=min(0.0_r8,ui(i,j,liunw))
739 aflxu(i,j)=on_u(i,j)* &
740 & (cff1*field(i-1,j,liold)+ &
741 & cff2*field(i ,j,liold))
742 END DO
743 END DO
744 DO j=jmin,jmax+1
745 DO i=imin,imax
746 cff1=max(0.0_r8,vi(i,j,liunw))
747 cff2=min(0.0_r8,vi(i,j,liunw))
748 aflxv(i,j)=om_v(i,j)* &
749 & (cff1*field(i,j-1,liold)+ &
750 & cff2*field(i,j ,liold))
751 END DO
752 END DO
753!
754! Step 1: first guess advection operator.
755!
756 DO j=jmin,jmax
757 DO i=imin,imax
758 ar(i,j)=1.0_r8/omn(i,j)
759 aif(i,j)=(field(i,j,liold)- &
760 & dtice(ng)*(aflxu(i+1,j)-aflxu(i,j)+ &
761 & aflxv(i,j+1)-aflxv(i,j))*ar(i,j))
762#ifdef MASKING
763 aif(i,j)=aif(i,j)*rmask(i,j)
764#endif
765#ifdef WET_DRY
766 aif(i,j)=aif(i,j)*rmask_wet(i,j)
767#endif
768#ifdef ICESHELF
769 IF (zice(i,j).ne.0.0_r8) aif(i,j)=0.0_r8
770#endif
771 END DO
772 END DO
773!
774! Set values at the open boundaries.
775!
776 IF (.not.ewperiodic(ng)) THEN
777 IF (domain(ng)%Western_Edge(tile)) THEN
778 DO j=jmin,jmax
779 aif(istr-1,j)=aif(istr,j)
780 END DO
781 END IF
782 IF (domain(ng)%Eastern_Edge(tile)) THEN
783 DO j=jmin,jmax
784 aif(iend+1,j)=aif(iend,j)
785 END DO
786 END IF
787 END IF
788!
789 IF (.not.nsperiodic(ng)) THEN
790 IF (domain(ng)%Southern_Edge(tile)) THEN
791 DO i=imin,imax
792 aif(i,jstr-1)=aif(i,jstr)
793 END DO
794 END IF
795 IF (domain(ng)%Northern_Edge(tile)) THEN
796 DO i=imin,imax
797 aif(i,jend+1)=aif(i,jend)
798 END DO
799 END IF
800 IF (.not.ewperiodic(ng)) THEN
801 IF (domain(ng)%SouthWest_Corner(tile)) THEN
802 aif(istr-1,jstr-1)=aif(istr,jstr)
803 END IF
804 IF (domain(ng)%NorthWest_Corner(tile)) THEN
805 aif(istr-1,jend+1)=aif(istr,jend)
806 END IF
807 IF (domain(ng)%SouthEast_Corner(tile)) THEN
808 aif(iend+1,jstr-1)=aif(iend,jstr)
809 END IF
810 IF (domain(ng)%NorthEast_Corner(tile)) THEN
811 aif(iend+1,jend+1)=aif(iend,jend)
812 END IF
813 END IF
814 END IF
815
816#ifdef MASKING
817!
818! Apply land/sea mask.
819!
820 DO j=jmin,jmax
821 DO i=imin,imax
822 aif(i,j)=aif(i,j)*rmask(i,j)
823# ifdef WET_DRY
824 aif(i,j)=aif(i,j)*rmask_wet(i,j)
825# endif
826 END DO
827 END DO
828#endif
829#ifdef ICESHELF
830 DO j=jmin,jmax
831 DO i=imin,jmax
832 IF (zice(i,j).ne.0.0_r8) THEN
833 aif(i,j)=0.0_r8
834 END IF
835 END DO
836 END DO
837#endif
838
839#ifndef ICE_UPWIND
840!
841! Step 2: Compute antidiffusive velocities .
842!
843! This is needed to avoid touching "aif" under land mask.
844! Note that only aif(i,j) and aif(i-1,j) are allowed to appear
845! explicitly in the code segment below. This is okay because if
846! either of them is masked, then "ui" is zero at that point,
847! and therefore no additional masking is required.
848!
849 DO j=jstr,jend+1
850 DO i=istr,iend+1
851 fe(i,j)=0.5* &
852# ifdef MASKING
853 & vmask(i,j)* &
854# endif
855# ifdef WET_DRY
856 & vmask_wet(i,j)* &
857# endif
858 & (aif(i,j)-aif(i,j-1))
859!
860 fx(i,j)=0.5* &
861# ifdef MASKING
862 & umask(i,j)* &
863# endif
864# ifdef WET_DRY
865 & umask_wet(i,j)* &
866# endif
867 & (aif(i,j)-aif(i-1,j))
868 END DO
869 END DO
870!
871 DO j=jstr,jend
872 DO i=istr,iend+1
873 rateu=(aif(i,j)-aif(i-1,j))/ &
874 & max(epsil, aif(i,j)+aif(i-1,j))
875
876 rateyiu=(fe(i,j+1)+fe(i,j)+fe(i-1,j+1)+fe(i-1,j))/ &
877 & (max(epsil, aif(i ,j)+fe(i ,j+1)-fe(i ,j)+ &
878 & aif(i-1,j)+fe(i-1,j+1)-fe(i-1,j)))
879
880 cu=0.5*dtice(ng)*(pm(i,j)+pm(i-1,j))*ui(i,j,liunw)
881
882 cu_crss=0.5_r8*dtice(ng)* &
883 & 0.0625_r8*(pn(i-1,j+1)+pn(i,j+1)+ &
884 & pn(i-1,j-1)+pn(i,j-1))* &
885 & (vi(i-1,j+1,liunw)+vi(i,j+1,liunw)+ &
886 & vi(i-1,j ,liunw)+vi(i,j ,liunw))
887
888 uspeed=rateu*(abs(ui(i,j,liunw))-cu*ui(i,j,liunw))- &
889 & rateyiu*cu_crss*ui(i,j,liunw)
890
891 cff1=max(0.0_r8,uspeed)
892 cff2=min(0.0_r8,uspeed)
893 aflxu(i,j)=on_u(i,j)*(cff1*aif(i-1,j)+ &
894 & cff2*aif(i ,j))
895 END DO
896 END DO
897!
898 DO j=jstr,jend+1
899 DO i=istr,iend
900 ratev=(aif(i,j)-aif(i,j-1))/ &
901 & max(epsil, aif(i,j)+aif(i,j-1))
902
903 ratexiv=(fx(i+1,j)+fx(i,j) +fx(i+1,j-1)+fx(i,j-1))/ &
904 & (max(epsil, aif(i,j )+fx(i+1,j )-fx(i,j )+ &
905 & aif(i,j-1)+fx(i+1,j-1)-fx(i,j-1)))
906
907 cu=0.5*dtice(ng)*(pn(i,j)+pn(i,j-1))*vi(i,j,liunw)
908
909 cu_crss=0.5_r8*dtice(ng)* &
910 & 0.0625_r8*(pm(i+1,j)+pm(i+1,j-1)+ &
911 & pm(i-1,j)+pm(i-1,j-1))* &
912 & (ui(i,j ,liunw)+ui(i+1,j ,liunw)+ &
913 & ui(i,j-1,liunw)+ui(i+1,j-1,liunw))
914
915 vspeed=ratev*(abs(vi(i,j,liunw))-cu*vi(i,j,liunw))- &
916 & ratexiv*cu_crss*vi(i,j,liunw)
917
918 cff1=max(0.0_r8,vspeed)
919 cff2=min(0.0_r8,vspeed)
920 aflxv(i,j)=om_v(i,j)*(cff1*aif(i,j-1)+ &
921 & cff2*aif(i,j ))
922 END DO
923 END DO
924!
925! Correct advection operator by substraction the antidifusive fluxes.
926!
927 DO j=jstr,jend
928 DO i=istr,iend
929 aif(i,j)=aif(i,j)- &
930 & dtice(ng)*pm(i,j)*pn(i,j)* &
931 & (aflxu(i+1,j)-aflxu(i,j)+ &
932 & aflxv(i,j+1)-aflxv(i,j))
933# ifdef MASKING
934 aif(i,j)=aif(i,j)*rmask(i,j)
935# endif
936# ifdef WET_DRY
937 aif(i,j)=aif(i,j)*rmask_wet(i,j)
938# endif
939# ifdef ICESHELF
940 IF (zice(i,j).ne.0.0_r8) aif(i,j)=0.
941# endif
942 END DO
943 END DO
944#endif
945!
946! Load advected solution.
947!
948 DO j=jstr,jend
949 DO i=istr,iend
950 field(i,j,linew)=aif(i,j)
951 END DO
952 END DO
953!
954 RETURN
955 END SUBROUTINE ice_mpdata_tile
956
957 END MODULE ice_advect_mod
exchange_2d_mod
Definition exchange_2d.F:2
exchange_2d_mod::exchange_r2d_tile
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition exchange_2d.F:253
ice_advect_mod
Definition ice_smolar.h:1
ice_advect_mod::ice_advect_tile
subroutine, private ice_advect_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, linew, liold, liunw)
Definition ice_smolar.h:80
ice_advect_mod::ice_mpdata_tile
subroutine, private ice_mpdata_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, linew, liold, liunw, rmask, umask, vmask, rmask_wet, umask_wet, vmask_wet, zice, pm, pn, on_u, om_v, omn, ui, vi, field)
Definition ice_smolar.h:627
ice_advect_mod::ice_advect
subroutine, public ice_advect(ng, tile, model)
Definition ice_smolar.h:46
ice_bc2d_mod
Definition ice_bc2d.F:2
ice_tibc_mod
Definition ice_tibc.F:2
mod_boundary
Definition mod_boundary.F:2
mod_forces
Definition mod_forces.F:2
mod_grid
Definition mod_grid.F:2
mod_grid::grid
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
mod_ice
Definition ice_mod.h:1
mod_ice::isvice
integer, parameter isvice
Definition ice_mod.h:147
mod_ice::dtice
integer, dimension(:), allocatable dtice
Definition ice_mod.h:217
mod_ice::min_hi
real(r8), dimension(:), allocatable min_hi
Definition ice_mod.h:246
mod_ice::ishage
integer, parameter ishage
Definition ice_mod.h:149
mod_ice::isino3
integer, parameter isino3
Definition ice_mod.h:154
mod_ice::icwdiv
integer, parameter icwdiv
Definition ice_mod.h:189
mod_ice::isenth
integer, parameter isenth
Definition ice_mod.h:148
mod_ice::ishsno
integer, parameter ishsno
Definition ice_mod.h:140
mod_ice::ice
type(t_ice), dimension(:), allocatable ice
Definition ice_mod.h:283
mod_ice::istice
integer, parameter istice
Definition ice_mod.h:145
mod_ice::isiphy
integer, parameter isiphy
Definition ice_mod.h:153
mod_ice::ibice
integer, dimension(nices) ibice
Definition ice_mod.h:162
mod_ice::ishmel
integer, parameter ishmel
Definition ice_mod.h:139
mod_ice::isiage
integer, parameter isiage
Definition ice_mod.h:141
mod_ice::isaice
integer, parameter isaice
Definition ice_mod.h:137
mod_ice::isinh4
integer, parameter isinh4
Definition ice_mod.h:155
mod_ice::isuice
integer, parameter isuice
Definition ice_mod.h:146
mod_ice::ishice
integer, parameter ishice
Definition ice_mod.h:138
mod_ocean
Definition mod_ocean.F:2
mod_param
Definition mod_param.F:2
mod_param::nghostpoints
integer nghostpoints
Definition mod_param.F:710
mod_param::lbc
type(t_lbc), dimension(:,:,:), allocatable lbc
Definition mod_param.F:375
mod_param::domain
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
mod_param::lm
integer, dimension(:), allocatable lm
Definition mod_param.F:455
mod_param::mm
integer, dimension(:), allocatable mm
Definition mod_param.F:456
mod_scalars
Definition mod_scalars.F:2
mod_scalars::dt
real(dp), dimension(:), allocatable dt
Definition mod_scalars.F:269
mod_scalars::ewperiodic
logical, dimension(:), allocatable ewperiodic
Definition mod_scalars.F:510
mod_scalars::nsperiodic
logical, dimension(:), allocatable nsperiodic
Definition mod_scalars.F:511
mod_stepping
Definition mod_stepping.F:2
mp_exchange_mod
Definition mp_exchange.F:2
mp_exchange_mod::mp_exchange2d
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
Definition mp_exchange.F:294
wclock_off
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
wclock_on
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3