ROMS
Loading...
Searching...
No Matches
ice_vbc.F
Go to the documentation of this file.
1#include "cppdefs.h"
3
4#ifdef ICE_MODEL
5!
6!git $Id$
7!================================================== Hernan G. Arango ===
8! Copyright (c) 2002-2025 The ROMS Group Paul Budgell !
9! Licensed under a MIT/X style license Katherine Hedstrom !
10! See License_ROMS.md Scott M. Durski !
11!=======================================================================
12! !
13! This routine sets the ice-water and ice-air stresses for the ice !
14! momentum equations. !
15! !
16!=======================================================================
17!
18 USE mod_param
19 USE mod_coupling
20 USE mod_forces
21 USE mod_grid
22 USE mod_ice
23# ifdef LMD_SKPP
24 USE mod_mixing
25# endif
26 USE mod_ocean
27 USE mod_scalars
28!
29 USE bc_2d_mod
30# ifdef DISTRIBUTE
32# endif
33!
34 implicit none
35!
36 PRIVATE
37 PUBLIC :: ice_vbc
38!
39 CONTAINS
40!
41!***********************************************************************
42 SUBROUTINE ice_vbc (ng, tile, model)
43!***********************************************************************
44!
45 USE mod_stepping
46!
47! Imported variable declarations.
48!
49 integer, intent(in) :: ng, tile, model
50!
51! Local variable declarations.
52!
53 character (len=*), parameter :: MyFile = &
54 & __FILE__
55!
56# include "tile.h"
57!
58# ifdef PROFILE
59 CALL wclock_on (ng, model, 42, __line__, myfile)
60# endif
61 CALL ice_vbc_tile (ng, tile, model, &
62 & lbi, ubi, lbj, ubj, &
63 & imins, imaxs, jmins, jmaxs, &
64 & liold(ng), liuol(ng), &
65# ifdef ICE_SHOREFAST
66 & grid(ng) % h, &
67# endif
68# ifdef WET_DRY
69 & grid(ng) % umask_wet, &
70 & grid(ng) % vmask_wet, &
71# endif
72 & grid(ng) % z_r, &
73 & grid(ng) % z_w, &
74# ifdef ICESHELF
75 & grid(ng) % zice, &
76# endif
77 & coupling(ng) % Zt_avg1, &
78 & ocean(ng) % rho, &
79 & ocean(ng) % u, &
80 & ocean(ng) % v, &
81# ifdef ICE_BULK_FLUXES
82 & forces(ng) % sustr_ai, &
83 & forces(ng) % svstr_ai, &
84 & forces(ng) % sustr_ao, &
85 & forces(ng) % svstr_ao, &
86# endif
87 & forces(ng) % sustr, &
88 & forces(ng) % svstr, &
89 & ice(ng) % Fi, &
90 & ice(ng) % Si)
91# ifdef PROFILE
92 CALL wclock_off (ng, model, 42, __line__, myfile)
93# endif
94!
95 RETURN
96 END SUBROUTINE ice_vbc
97!
98!***********************************************************************
99 SUBROUTINE ice_vbc_tile (ng, tile, model, &
100 & LBi, UBi, LBj, UBj, &
101 & IminS, ImaxS, JminS, JmaxS, &
102 & liold, liuol, &
103# ifdef ICE_SHOREFAST
104 & h, &
105# endif
106# ifdef WET_DRY
107 & umask_wet, vmask_wet, &
108# endif
109 & z_r, z_w, &
110# ifdef ICESHELF
111 & zice, &
112# endif
113 & Zt_avg1, &
114 & rho, u, v, &
115# ifdef ICE_BULK_FLUXES
116 & sustr_ai, svstr_ai, &
117 & sustr_ao, svstr_ao, &
118# endif
119 & sustr, svstr, &
120 & Fi, Si)
121!***********************************************************************
122!
123! Imported variable declarations.
124!
125 integer, intent(in) :: ng, tile, model
126 integer, intent(in) :: LBi, UBi, LBj, UBj
127 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
128 integer, intent(in) :: liold, liuol
129!
130# ifdef ASSUMED_SHAPE
131# ifdef ICE_SHOREFAST
132 real(r8), intent(in) :: h(LBi:,LBj:)
133# endif
134# ifdef WET_DRY
135 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
136 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
137# endif
138 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
139 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
140# ifdef ICESHELF
141 real(r8), intent(in) :: zice(LBi:,LBj:)
142# endif
143 real(r8), intent(in) :: Zt_avg1(LBi:,LBj:)
144 real(r8), intent(in) :: rho(LBi:,LBj:,:)
145 real(r8), intent(in) :: u(LBi:,LBj:,:,:)
146 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
147# ifdef ICE_BULK_FLUXES
148 real(r8), intent(in) :: sustr_ao(LBi:,LBj:)
149 real(r8), intent(in) :: svstr_ao(LBi:,LBj:)
150 real(r8), intent(in) :: sustr_ai(LBi:,LBj:)
151 real(r8), intent(in) :: svstr_ai(LBi:,LBj:)
152# endif
153 real(r8), intent(in) :: Si(LBi:,LBj:,:,:)
154 real(r8), intent(inout) :: sustr(LBi:,LBj:)
155 real(r8), intent(inout) :: svstr(LBi:,LBj:)
156 real(r8), intent(inout) :: Fi(LBi:,LBj:,:)
157# else
158# ifdef ICE_SHOREFAST
159 real(r8), intent(in) :: h(LBi:UBi,LBj:UBj)
160# endif
161# ifdef WET_DRY
162 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
163 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
164# endif
165 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
166 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
167# ifdef ICESHELF
168 real(r8), intent(in) :: zice(LBi:UBi,LBj:UBj)
169# endif
170 real(r8), intent(in) :: Zt_avg1(LBi:UBi,LBj:UBj)
171 real(r8), intent(in) :: rho(LBi:UBi,LBj:UBj,N(ng))
172 real(r8), intent(in) :: u(LBi:UBi,LBj:UBj,N(ng),2)
173 real(r8), intent(in) :: v(LBi:UBi,LBj:UBj,N(ng),2)
174# ifdef ICE_BULK_FLUXES
175 real(r8), intent(in) :: sustr_ao(LBi:UBi,LBj:UBj)
176 real(r8), intent(in) :: svstr_ao(LBi:UBi,LBj:UBj)
177 real(r8), intent(in) :: sustr_ai(LBi:UBi,LBj:UBj)
178 real(r8), intent(in) :: svstr_ai(LBi:UBi,LBj:UBj)
179# endif
180 real(r8), intent(in) :: Si(LBi:UBi,LBj:UBj,2,nIceS)
181 real(r8), intent(inout) :: sustr(LBi:UBi,LBj:UBj)
182 real(r8), intent(inout) :: svstr(LBi:UBi,LBj:UBj)
183 real(r8), intent(inout) :: Fi(LBi:UBi,LBj:UBj,nIceF)
184# endif
185!
186! Local variable declarations.
187!
188 logical :: IceCavity
189!
190 integer :: i, j
191!
192 real(r8) :: aix, aiy, cff, chux, chuy, chuax, chuay, dztop
193 real(r8) :: hix, hiy, rhoO, spd, thic, zdz0, z0
194 real(r8) :: tauiwu, tauiwv
195# ifdef ICE_SHOREFAST
196 real(r8) :: clear, hh
197# endif
198!
199 real(r8), parameter :: kappa = 0.4_r8
200 real(r8), parameter :: z0ii = 0.01_r8
201 real(r8), parameter :: eps = 1.0e-20
202!
203 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: spdiw
204 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: chuiw
205 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: chuai
206 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: utauiw
207# ifndef ICE_MK
208 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: uwind
209 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: vwind
210 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: wind_speed
211 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: windu
212 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: windv
213# endif
214
215# include "set_bounds.h"
216
217# ifdef ICE_MK
218!
219!-----------------------------------------------------------------------
220! Compute ice-water and ice-air stresses according with the Mellor-
221! Kantha ice thermodynamics formulation.
222!-----------------------------------------------------------------------
223!
224 DO j=jstr-1,jend
225 DO i=istr-1,iend
226 rhoo=rho0+rho(i,j,n(ng))
227 spd=fi(i,j,iciovs)
228# ifdef ICE_SHOREFAST
229 spd=max(spd, 0.10_r8)
230 clear=h(i,j)+ &
231 & min(zt_avg1(i,j), 0.0_r8)-0.9*si(i,j,liold,ishice)
232 clear=max(clear, 0.001_r8)
233 IF (clear.lt.5.0_r8) THEN
234 spd=0.2_r8*clear*spd ! 0.2 = 1/5
235 END IF
236# else
237 spd=max(spd, 0.15_r8)
238# endif
239# ifdef ICE_MOM_BULK
240 utauiw(i,j)=spd
241 chuiw(i,j)=cd_io(ng)*spd
242# else
243 thic=si(i,j,liold,ishice)/ &
244 & max(si(i,j,liold,isaice), max(min_ai(ng), eps))
245!
246! This parameterization for z0 seems to be designed to vary the
247! ice-ocean drag for ice between half a meter and 5 m thick (SMD).
248!
249 z0=max(z0ii*thic, 0.01_r8)
250 z0=min(z0, 0.1_r8)
251 dztop=z_w(i,j,n(ng))-z_r(i,j,n(ng))
252 zdz0=dztop/z0
253!
254! This formulation for the ice-ocean momentum transfer coeficient,
255! chuiw, leads to extremely high values along the Alaskan coast where
256! the ice is particularly thick and the water is particularly shallow.
257! Therfore, the 'upper limit' is changed from 3 to 6 to cap the ice
258! water drag at a more 'reasonable' value in shallow water (SMD).
259!
260 IF (zdz0.lt.6.0_r8) zdz0=6.0_r8
261 utauiw(i,j)=sqrt(fi(i,j,iciomt)*spd)
262 utauiw(i,j)=max(utauiw(i,j), 1.0e-04_r8)
263 chuiw(i,j)=kappa*utauiw(i,j)/log(zdz0)
264!
265! Add an ice thickness dependent factor for the wind stress similar
266! to what is below for the othe formulation option (SMD).
267!
268 chuai(i,j)=1.0_r8-cos(1.0_r8*pi*min((thic+0.05_r8), 1.0_r8))
269# ifdef ICE_SHOREFAST
270 hh=h(i,j)+min(zt_avg1(i,j),0.0_r8)
271 clear=hh-0.9_r8*si(i,j,liold,ishice)
272 clear=max(clear, 0.0_r8)
273 IF (clear.lt.5.0_r8) &
274 & chuiw(i,j)=(max(clear-1.0_r8, 0.0_r8)*0.25_r8)*chuiw(i,j)
275# endif
276# endif
277 END DO
278 END DO
279!
280! Compute surface air-sea U-stress.
281!
282 DO j=jstr,jend
283 DO i=istrp,iend
284 rhoo=1000.0_r8+0.5_r8*(rho(i,j,n(ng))+rho(i-1,j,n(ng)))
285 aix =0.5_r8*(si(i,j,liold,isaice)+si(i-1,j,liold,isaice))
286 chux=0.5_r8*(chuiw(i,j)+chuiw(i-1,j))
287!
288! Add thickness dependent wind ice drag (SMD).
289!
290# ifdef ICE_MOM_BULK
291 fi(i,j,icaius)=0.5_r8*aix* &
292 & (sustr_ai(i,j)+sustr_ai(i-1,j))/icerho(ng)
293# else
294 chuax=0.5_r8*(chuai(i,j)+chuai(i-1,j))
295 fi(i,j,icaius)=0.5_r8*aix*chuax* &
296 & (sustr_ai(i,j)+sustr_ai(i-1,j))/icerho(ng)
297# endif
298# ifdef WET_DRY
299 fi(i,j,icaius)=fi(i,j,icaius)*umask_wet(i,j)
300# endif
301# ifdef ICE_BULK_FLUXES
302!
303# ifdef ICESHELF
304 icecavity=zice(i,j).ne.0.0_r8
305# else
306 icecavity=.false.
307# endif
308 IF (.not.icecavity) THEN
309 sustr(i,j)=aix*chux*(si(i,j,liuol,isuice)-fi(i,j,icuavg))+ &
310 & (1.0_r8-aix)*sustr_ao(i,j)
311# ifdef WET_DRY
312 sustr(i,j)=sustr(i,j)*umask_wet(i,j)
313# endif
314 END IF
315# endif
316 END DO
317 END DO
318!
319! Compute surface air-ice interface V-stress.
320!
321 DO j=jstrp,jend
322 DO i=istr,iend
323 rhoo=1000.0_r8+0.5_r8*(rho(i,j,n(ng))+rho(i,j-1,n(ng)))
324 aiy =0.5_r8*(si(i,j,liold,isaice)+si(i,j-1,liold,isaice))
325 chuy=0.5_r8*(chuiw(i,j)+chuiw(i,j-1))
326!
327! Add thickness dependent wind ice drag (SMD).
328!
329# ifdef ICE_MOM_BULK
330 fi(i,j,icaivs)=0.5_r8*aiy* &
331 & (svstr_ai(i,j)+svstr_ai(i,j-1))/icerho(ng)
332
333# else
334 chuay=0.5_r8*(chuai(i,j)+chuai(i,j-1))
335 fi(i,j,icaivs)=0.5_r8*aiy*chuay* &
336 & (svstr_ai(i,j)+svstr_ai(i,j-1))/icerho(ng)
337# endif
338# ifdef WET_DRY
339 fi(i,j,icaivs)=fi(i,j,icaivs)*vmask_wet(i,j)
340# endif
341
342# ifdef ICE_BULK_FLUXES
343!
344# ifdef ICESHELF
345 icecavity=zice(i,j).ne.0.0_r8
346# else
347 icecavity=.false.
348# endif
349 IF (.not.icecavity) THEN
350 svstr(i,j)=aiy*chuy*(si(i,j,liuol,isvice)-fi(i,j,icvavg))+ &
351 & (1.0_r8-aiy)*svstr_ao(i,j)
352# ifdef WET_DRY
353 svstr(i,j)=svstr(i,j)*vmask_wet(i,j)
354# endif
355 END IF
356# endif
357 END DO
358 END DO
359# else
360!
361!-----------------------------------------------------------------------
362! Alternate formulation to the Mellor-Kantha ice thermodynamics.
363!-----------------------------------------------------------------------
364!
365 DO j=jstr,jend
366 DO i=istrp,iend
367 windu(i,j)=0.5_r8*(uwind(i-1,j)+uwind(i,j))
368 END DO
369 END DO
370 DO j=jstrp,jend
371 DO i=istr,iend
372 windv(i,j)=0.5_r8*(vwind(i,j-1)+vwind(i,j))
373 END DO
374 END DO
375!
376 DO j=jstr,jend
377 DO i=istr,iend
378 cff=(si(i,j,liold,isaice)+0.02)*rho0*cd_io(ng)/icerho(ng)
379 spdiw(i,j)=cff*spd_iw(i,j)
380 spdiw(i,j)=max(spdiw(i,j), cff*0.1_r8)
381 wind_speed(i,j)=0.5*sqrt((windu(i,j)+windu(i+1,j))**2 + &
382 & (windv(i,j)+windv(i,j+1))**2)
383# ifdef ICE_SHOREFAST
384 clear=h(i,j)+min(zt_avg1(i,j), 0.0_r8)-0.9*hi(i,j,liold)
385 clear=max(clear, 0.001_r8)
386 IF (clear.lt.5.0_r8) THEN
387 spd=0.2_r8*clear*spd
388 spdiw(i,j)=0.2_r8*clear*spdiw(i,j) ! 0.2 = 1/5
389 END IF
390# endif
391 END DO
392 END DO
393!
394! Compute surface air-ice interface U-stress.
395!
396 DO j=jstr,jend
397 DO i=istrp,iend
398 rhoo=1000.0_r8+0.5_r8*(rho(i,j,n(ng))+rho(i-1,j,n(ng)))
399 aix=0.5_r8*(si(i,j,liold,isaice)+si(i-1,j,liold,isaice))
400 hix=0.5_r8*(si(i,j,liold,ishice)+si(i-1,j,liold,ishice))
401 spd=0.5_r8*(wind_speed(i,j)+wind_speed(i-1,j))
402 fi(i,j,icaius)=aix*airrho(ng)* &
403 & (1.0_r8*cd_ai(ng)* &
404 & (1.0_r8-cos(1.0_r8*pi* &
405 & min((hix/(aix+0.02_r8)+0.01_r8), &
406 & 1.0_r8))))* &
407 & spd*windu(i,j)/icerho(ng)
408
409# ifdef ICE_BULK_FLUXES
410!
411# ifdef ICESHELF
412 icecavity=zice(i,j).ne.0.0_r8
413# else
414 icecavity=.false.
415# endif
416 IF (.not.icecavity) THEN
417 tauiwu=0.5_r8*(spdiw(i,j)+spdiw(i-1,j))* &
418 & (si(i,j,liuol,isuice)-fi(i,j,icuavg))*icerho(ng)/rhoo
419 sustr(i,j)=tauiwu + &
420 & (1.0_r8-aix)* &
421 & 0.5_r8*(sustr_ao(i,j)+sustr_ao(i-1,j))
422# ifdef WET_DRY
423 sustr(i,j)=sustr(i,j)*umask_wet(i,j)
424# endif
425 END IF
426# endif
427 END DO
428 END DO
429!
430! Compute surface air-ice interface V-stress.
431!
432 DO j=jstrp,jend
433 DO i=istr,iend
434 rhoo=1000.0_r8+0.5_r8*(rho(i,j,n(ng))+rho(i,j-1,n(ng)))
435 aiy=0.5_r8*(si(i,j,liold,isaice)+si(i,j-1,liold,isaice))
436 hiy=0.5_r8*(si(i,j,liold,ishice)+si(i,j-1,liold,ishice))
437 spd=0.5_r8*(wind_speed(i,j)+wind_speed(i,j-1))
438 fi(i,j,icaivs)=aiy*airrho(ng)* &
439 & (0.5_r8*cd_ai(ng)* &
440 & (1.0_r8-cos(2.0_r8*pi* &
441 & min((hiy/(aiy+0.02_r8)+0.1_r8), &
442 & 0.5_r8))))* &
443 & spd*windv(i,j)/icerho(ng)
444
445# ifdef ICE_BULK_FLUXES
446!
447# ifdef ICESHELF
448 icecavity=zice(i,j).ne.0.0_r8
449# else
450 icecavity=.false.
451# endif
452 IF (.not.icecavity) THEN
453 tauiwv=0.5_r8*(spdiw(i,j)+spdiw(i,j-1))* &
454 & (si(i,j,liuol,isvice)-fi(i,j,icvavg))*icerho(ng)/rho0
455 svstr(i,j)=tauiwv+ &
456 & (1.0_r8-aiy)* &
457 & 0.5_r8*(svstr_ao(i,j)+svstr_ao(i,j-1))
458# ifdef WET_DRY
459 svstr(i,j)=svstr(i,j)*vmask_wet(i,j)
460# endif
461 END IF
462# endif
463 END DO
464 END DO
465# endif
466!
467# ifdef ICE_MK
468! Load sea surface elevation (isHsse), ice-ocean momentum transfer
469! coefficient (icIOmt), ice-ocean friction velocity (icIOfv).
470# else
471! Load sea surface elevation (isHsse).
472# endif
473!
474 DO j=jstr,jend
475 DO i=istr,iend
476 fi(i,j,ichsse)=zt_avg1(i,j)
477# ifdef ICE_MK
478 fi(i,j,iciomt)=chuiw(i,j)
479 fi(i,j,iciofv)=utauiw(i,j)
480# endif
481 END DO
482 END DO
483!
484! Apply boundary conditions.
485!
486 CALL bc_r2d_tile (ng, tile, &
487 & lbi, ubi, lbj, ubj, &
488 & fi(:,:,ichsse))
489
490 CALL bc_u2d_tile (ng, tile, &
491 & lbi, ubi, lbj, ubj, &
492 & fi(:,:,icaius))
493
494 CALL bc_v2d_tile (ng, tile, &
495 & lbi, ubi, lbj, ubj, &
496 & fi(:,:,icaivs))
497
498# ifdef ICE_MK
499 CALL bc_r2d_tile (ng, tile, &
500 & lbi, ubi, lbj, ubj, &
501 & fi(:,:,iciomt))
502
503 CALL bc_r2d_tile (ng, tile, &
504 & lbi, ubi, lbj, ubj, &
505 & fi(:,:,iciofv))
506# endif
507
508# ifdef ICE_BULK_FLUXES
509 CALL bc_u2d_tile (ng, tile, &
510 & lbi, ubi, lbj, ubj, &
511 & sustr)
512
513 CALL bc_v2d_tile (ng, tile, &
514 & lbi, ubi, lbj, ubj, &
515 & svstr)
516# endif
517
518# ifdef DISTRIBUTE
519!
520 CALL mp_exchange2d (ng, tile, inlm, 3, &
521 & lbi, ubi, lbj, ubj, &
522 & nghostpoints, ewperiodic(ng), nsperiodic(ng), &
523 & fi(:,:,ichsse), &
524 & fi(:,:,icaius), &
525 & fi(:,:,icaivs))
526
527# ifdef ICE_MK
528 CALL mp_exchange2d (ng, tile, inlm, 2, &
529 & lbi, ubi, lbj, ubj, &
530 & nghostpoints, ewperiodic(ng), nsperiodic(ng), &
531 & fi(:,:,iciomt), &
532 & fi(:,:,iciofv))
533# endif
534
535# ifdef ICE_BULK_FLUXES
536 CALL mp_exchange2d (ng, tile, inlm, 2, &
537 & lbi, ubi, lbj, ubj, &
538 & nghostpoints, ewperiodic(ng), nsperiodic(ng), &
539 & sustr, svstr)
540# endif
541# endif
542!
543 RETURN
544 END SUBROUTINE ice_vbc_tile
545!
546#endif
547 END MODULE ice_vbc_mod
subroutine bc_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition bc_2d.F:44
subroutine bc_v2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition bc_2d.F:345
subroutine bc_u2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition bc_2d.F:167
type(t_coupling), dimension(:), allocatable coupling
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer, parameter isvice
Definition ice_mod.h:147
integer, parameter ichsse
Definition ice_mod.h:174
integer, parameter icuavg
Definition ice_mod.h:187
real(r8), dimension(:), allocatable airrho
Definition ice_mod.h:222
real(r8), dimension(:), allocatable icerho
Definition ice_mod.h:223
type(t_ice), dimension(:), allocatable ice
Definition ice_mod.h:283
integer, parameter icvavg
Definition ice_mod.h:188
integer, parameter iciomt
Definition ice_mod.h:177
real(r8), dimension(:), allocatable min_ai
Definition ice_mod.h:241
integer, parameter iciovs
Definition ice_mod.h:178
integer, parameter iciofv
Definition ice_mod.h:175
integer, parameter icaius
Definition ice_mod.h:171
integer, parameter isaice
Definition ice_mod.h:137
integer, parameter isuice
Definition ice_mod.h:146
real(r8), dimension(:), allocatable cd_io
Definition ice_mod.h:230
integer, parameter ishice
Definition ice_mod.h:138
integer, parameter icaivs
Definition ice_mod.h:172
real(r8), dimension(:), allocatable cd_ai
Definition ice_mod.h:229
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer, parameter inlm
Definition mod_param.F:662
integer nghostpoints
Definition mod_param.F:710
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
real(dp) rho0
real(dp), parameter pi
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
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