New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
bdyice_lim.F90 in branches/2015/dev_r5803_NOC_WAD/NEMOGCM/NEMO/OPA_SRC/BDY – NEMO

source: branches/2015/dev_r5803_NOC_WAD/NEMOGCM/NEMO/OPA_SRC/BDY/bdyice_lim.F90 @ 5870

Last change on this file since 5870 was 5870, checked in by acc, 9 years ago

Branch 2015/dev_r5803_NOC_WAD. Merge in trunk changes from 5803 to 5869 in preparation for merge. Also tidied and reorganised some wetting and drying code. Renamed wadlmt.F90 to wetdry.F90. Wetting drying code changes restricted to domzgr.F90, domvvl.F90 nemogcm.F90 sshwzv.F90, dynspg_ts.F90, wetdry.F90 and dynhpg.F90. Code passes full SETTE tests with ln_wd=.false.. Still awaiting test case for checking with ln_wd=.false.

  • Property svn:keywords set to Id
File size: 18.8 KB
Line 
1MODULE bdyice_lim
2   !!======================================================================
3   !!                       ***  MODULE  bdyice_lim  ***
4   !! Unstructured Open Boundary Cond. :  Open boundary conditions for sea-ice (LIM2 and LIM3)
5   !!======================================================================
6   !!  History :  3.3  !  2010-09 (D. Storkey)  Original code
7   !!             3.4  !  2011    (D. Storkey)  rewrite in preparation for OBC-BDY merge
8   !!              -   !  2012-01 (C. Rousset)  add lim3 and remove useless jk loop
9   !!----------------------------------------------------------------------
10#if defined   key_bdy   &&  ( defined key_lim2 || defined key_lim3 )
11   !!----------------------------------------------------------------------
12   !!   'key_bdy'            and                 Unstructured Open Boundary Conditions
13   !!   'key_lim2'                                                 LIM-2 sea ice model
14   !!   'key_lim3'                                                 LIM-3 sea ice model
15   !!----------------------------------------------------------------------
16   !!   bdy_ice_lim        : Application of open boundaries to ice
17   !!   bdy_ice_frs        : Application of Flow Relaxation Scheme
18   !!----------------------------------------------------------------------
19   USE timing          ! Timing
20   USE phycst          ! physical constant
21   USE eosbn2          ! equation of state
22   USE oce             ! ocean dynamics and tracers variables
23#if defined key_lim2
24   USE par_ice_2
25   USE ice_2           ! LIM_2 ice variables
26   USE dom_ice_2       ! sea-ice domain
27#elif defined key_lim3
28   USE ice             ! LIM_3 ice variables
29   USE dom_ice         ! sea-ice domain
30   USE limvar
31#endif
32   USE par_oce         ! ocean parameters
33   USE dom_oce         ! ocean space and time domain variables
34   USE sbc_oce         ! Surface boundary condition: ocean fields
35   USE bdy_oce         ! ocean open boundary conditions
36   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
37   USE in_out_manager  ! write to numout file
38   USE lib_mpp         ! distributed memory computing
39   USE lib_fortran     ! to use key_nosignedzero
40
41   IMPLICIT NONE
42   PRIVATE
43
44   PUBLIC   bdy_ice_lim     ! routine called in sbcmod
45   PUBLIC   bdy_ice_lim_dyn ! routine called in limrhg
46
47   !!----------------------------------------------------------------------
48   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
49   !! $Id$
50   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
51   !!----------------------------------------------------------------------
52CONTAINS
53
54   SUBROUTINE bdy_ice_lim( kt )
55      !!----------------------------------------------------------------------
56      !!                  ***  SUBROUTINE bdy_ice_lim  ***
57      !!
58      !! ** Purpose : - Apply open boundary conditions for ice (LIM2 and LIM3)
59      !!
60      !!----------------------------------------------------------------------
61      INTEGER, INTENT( in ) ::   kt   ! Main time step counter
62      !
63      INTEGER ::   ib_bdy   ! Loop index
64      !!----------------------------------------------------------------------
65      !
66#if defined key_lim3
67      CALL lim_var_glo2eqv
68#endif
69      !
70      DO ib_bdy=1, nb_bdy
71         !
72         SELECT CASE( cn_ice_lim(ib_bdy) )
73         CASE('none')
74            CYCLE
75         CASE('frs')
76            CALL bdy_ice_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy )
77         CASE DEFAULT
78            CALL ctl_stop( 'bdy_ice_lim : unrecognised option for open boundaries for ice fields' )
79         END SELECT
80         !
81      END DO
82      !
83#if defined key_lim3
84      CALL lim_var_zapsmall
85      CALL lim_var_agg(1)
86#endif
87      !
88   END SUBROUTINE bdy_ice_lim
89
90
91   SUBROUTINE bdy_ice_frs( idx, dta, kt, ib_bdy )
92      !!------------------------------------------------------------------------------
93      !!                 ***  SUBROUTINE bdy_ice_frs  ***
94      !!                   
95      !! ** Purpose : Apply the Flow Relaxation Scheme for sea-ice fields in the case
96      !!              of unstructured open boundaries.
97      !!
98      !! Reference : Engedahl H., 1995: Use of the flow relaxation scheme in a three-
99      !!             dimensional baroclinic ocean model with realistic topography. Tellus, 365-382.
100      !!------------------------------------------------------------------------------
101      TYPE(OBC_INDEX), INTENT(in) ::   idx     ! OBC indices
102      TYPE(OBC_DATA),  INTENT(in) ::   dta     ! OBC external data
103      INTEGER,         INTENT(in) ::   kt      ! main time-step counter
104      INTEGER,         INTENT(in) ::   ib_bdy  ! BDY set index
105      !
106      INTEGER  ::   jpbound            ! 0 = incoming ice
107      !                                ! 1 = outgoing ice
108      INTEGER  ::   jb, jk, jgrd, jl   ! dummy loop indices
109      INTEGER  ::   ji, jj, ii, ij     ! local scalar
110      REAL(wp) ::   zwgt, zwgt1        ! local scalar
111      REAL(wp) ::   ztmelts, zdh
112#if  defined key_lim2 && ! defined key_lim2_vp && defined key_agrif
113     USE ice_2, vt_s => hsnm
114     USE ice_2, vt_i => hicm
115#endif
116      !!------------------------------------------------------------------------------
117      !
118      IF( nn_timing == 1 )   CALL timing_start('bdy_ice_frs')
119      !
120      jgrd = 1      ! Everything is at T-points here
121      !
122#if defined key_lim2
123      DO jb = 1, idx%nblen(jgrd)
124         ji    = idx%nbi(jb,jgrd)
125         jj    = idx%nbj(jb,jgrd)
126         zwgt  = idx%nbw(jb,jgrd)
127         zwgt1 = 1.e0 - idx%nbw(jb,jgrd)
128         frld (ji,jj) = ( frld (ji,jj) * zwgt1 + dta%frld (jb) * zwgt ) * tmask(ji,jj,1)     ! Leads fraction
129         hicif(ji,jj) = ( hicif(ji,jj) * zwgt1 + dta%hicif(jb) * zwgt ) * tmask(ji,jj,1)     ! Ice depth
130         hsnif(ji,jj) = ( hsnif(ji,jj) * zwgt1 + dta%hsnif(jb) * zwgt ) * tmask(ji,jj,1)     ! Snow depth
131      END DO
132
133      CALL lbc_bdy_lnk( frld,  'T', 1., ib_bdy )                                         ! lateral boundary conditions
134      CALL lbc_bdy_lnk( hicif, 'T', 1., ib_bdy )
135      CALL lbc_bdy_lnk( hsnif, 'T', 1., ib_bdy )
136
137      vt_i(:,:) = hicif(:,:) * frld(:,:)
138      vt_s(:,:) = hsnif(:,:) * frld(:,:)
139      !
140#elif defined key_lim3
141
142      DO jl = 1, jpl
143         DO jb = 1, idx%nblen(jgrd)
144            ji    = idx%nbi(jb,jgrd)
145            jj    = idx%nbj(jb,jgrd)
146            zwgt  = idx%nbw(jb,jgrd)
147            zwgt1 = 1.e0 - idx%nbw(jb,jgrd)
148            a_i (ji,jj,jl) = ( a_i (ji,jj,jl) * zwgt1 + dta%a_i (jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Leads fraction
149            ht_i(ji,jj,jl) = ( ht_i(ji,jj,jl) * zwgt1 + dta%ht_i(jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Ice depth
150            ht_s(ji,jj,jl) = ( ht_s(ji,jj,jl) * zwgt1 + dta%ht_s(jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Snow depth
151
152            ! -----------------
153            ! Pathological case
154            ! -----------------
155            ! In case a) snow load would be in excess or b) ice is coming into a warmer environment that would lead to
156            ! very large transformation from snow to ice (see limthd_dh.F90)
157
158            ! Then, a) transfer the snow excess into the ice (different from limthd_dh)
159            zdh = MAX( 0._wp, ( rhosn * ht_s(ji,jj,jl) + ( rhoic - rau0 ) * ht_i(ji,jj,jl) ) * r1_rau0 )
160            ! Or, b) transfer all the snow into ice (if incoming ice is likely to melt as it comes into a warmer environment)
161            !zdh = MAX( 0._wp, ht_s(ji,jj,jl) * rhosn / rhoic )
162
163            ! recompute ht_i, ht_s
164            ht_i(ji,jj,jl) = MIN( hi_max(jl), ht_i(ji,jj,jl) + zdh )
165            ht_s(ji,jj,jl) = MAX( 0._wp, ht_s(ji,jj,jl) - zdh * rhoic / rhosn ) 
166
167         ENDDO
168         CALL lbc_bdy_lnk(  a_i(:,:,jl), 'T', 1., ib_bdy )
169         CALL lbc_bdy_lnk( ht_i(:,:,jl), 'T', 1., ib_bdy )
170         CALL lbc_bdy_lnk( ht_s(:,:,jl), 'T', 1., ib_bdy )
171      ENDDO
172      ! retrieve at_i
173      at_i(:,:) = 0._wp
174      DO jl = 1, jpl
175         at_i(:,:) = a_i(:,:,jl) + at_i(:,:)
176      END DO
177
178      DO jl = 1, jpl
179         DO jb = 1, idx%nblen(jgrd)
180            ji    = idx%nbi(jb,jgrd)
181            jj    = idx%nbj(jb,jgrd)
182
183            ! condition on ice thickness depends on the ice velocity
184            ! if velocity is outward (strictly), then ice thickness, volume... must be equal to adjacent values
185            jpbound = 0   ;   ii = ji   ;   ij = jj
186            !
187            IF( u_ice(ji+1,jj  ) < 0. .AND. umask(ji-1,jj  ,1) == 0. ) jpbound = 1; ii = ji+1; ij = jj
188            IF( u_ice(ji-1,jj  ) > 0. .AND. umask(ji+1,jj  ,1) == 0. ) jpbound = 1; ii = ji-1; ij = jj
189            IF( v_ice(ji  ,jj+1) < 0. .AND. vmask(ji  ,jj-1,1) == 0. ) jpbound = 1; ii = ji  ; ij = jj+1
190            IF( v_ice(ji  ,jj-1) > 0. .AND. vmask(ji  ,jj+1,1) == 0. ) jpbound = 1; ii = ji  ; ij = jj-1
191            !
192            IF( nn_ice_lim_dta(ib_bdy) == 0 ) jpbound = 0; ii = ji; ij = jj   ! case ice boundaries = initial conditions
193            !                                                                 !      do not make state variables dependent on velocity
194            !
195            rswitch = MAX( 0.0_wp , SIGN ( 1.0_wp , at_i(ii,ij) - 0.01 ) ) ! 0 if no ice
196            !
197            ! concentration and thickness
198            a_i (ji,jj,jl) = a_i (ii,ij,jl) * rswitch
199            ht_i(ji,jj,jl) = ht_i(ii,ij,jl) * rswitch
200            ht_s(ji,jj,jl) = ht_s(ii,ij,jl) * rswitch
201            !
202            ! Ice and snow volumes
203            v_i(ji,jj,jl) = ht_i(ji,jj,jl) * a_i(ji,jj,jl)
204            v_s(ji,jj,jl) = ht_s(ji,jj,jl) * a_i(ji,jj,jl)
205            !
206            SELECT CASE( jpbound )
207            !
208            CASE( 0 )   ! velocity is inward
209               !
210               ! Ice salinity, age, temperature
211               sm_i(ji,jj,jl)   = rswitch * rn_ice_sal(ib_bdy)  + ( 1.0 - rswitch ) * rn_simin
212               oa_i(ji,jj,jl)   = rswitch * rn_ice_age(ib_bdy) * a_i(ji,jj,jl)
213               t_su(ji,jj,jl)   = rswitch * rn_ice_tem(ib_bdy)  + ( 1.0 - rswitch ) * rn_ice_tem(ib_bdy)
214               DO jk = 1, nlay_s
215                  t_s(ji,jj,jk,jl) = rswitch * rn_ice_tem(ib_bdy) + ( 1.0 - rswitch ) * rt0
216               END DO
217               DO jk = 1, nlay_i
218                  t_i(ji,jj,jk,jl) = rswitch * rn_ice_tem(ib_bdy) + ( 1.0 - rswitch ) * rt0 
219                  s_i(ji,jj,jk,jl) = rswitch * rn_ice_sal(ib_bdy) + ( 1.0 - rswitch ) * rn_simin
220               END DO
221               !
222            CASE( 1 )   ! velocity is outward
223               !
224               ! Ice salinity, age, temperature
225               sm_i(ji,jj,jl)   = rswitch * sm_i(ii,ij,jl)  + ( 1.0 - rswitch ) * rn_simin
226               oa_i(ji,jj,jl)   = rswitch * oa_i(ii,ij,jl)
227               t_su(ji,jj,jl)   = rswitch * t_su(ii,ij,jl)  + ( 1.0 - rswitch ) * rt0
228               DO jk = 1, nlay_s
229                  t_s(ji,jj,jk,jl) = rswitch * t_s(ii,ij,jk,jl) + ( 1.0 - rswitch ) * rt0
230               END DO
231               DO jk = 1, nlay_i
232                  t_i(ji,jj,jk,jl) = rswitch * t_i(ii,ij,jk,jl) + ( 1.0 - rswitch ) * rt0
233                  s_i(ji,jj,jk,jl) = rswitch * s_i(ii,ij,jk,jl) + ( 1.0 - rswitch ) * rn_simin
234               END DO
235               !
236            END SELECT
237            !
238            IF( nn_icesal == 1 ) THEN     ! constant salinity : overwrite rn_ice_sal
239               sm_i(ji,jj  ,jl) = rn_icesal
240               s_i (ji,jj,:,jl) = rn_icesal
241            ENDIF
242            !
243            ! contents
244            smv_i(ji,jj,jl)  = MIN( sm_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl)
245            DO jk = 1, nlay_s
246               ! Snow energy of melting
247               e_s(ji,jj,jk,jl) = rswitch * rhosn * ( cpic * ( rt0 - t_s(ji,jj,jk,jl) ) + lfus )
248               ! Multiply by volume, so that heat content in J/m2
249               e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * v_s(ji,jj,jl) * r1_nlay_s
250            END DO
251            DO jk = 1, nlay_i
252               ztmelts          = - tmut * s_i(ji,jj,jk,jl) + rt0 !Melting temperature in K                 
253               ! heat content per unit volume
254               e_i(ji,jj,jk,jl) = rswitch * rhoic * &
255                  (   cpic    * ( ztmelts - t_i(ji,jj,jk,jl) ) &
256                  +   lfus    * ( 1.0 - (ztmelts-rt0) / MIN((t_i(ji,jj,jk,jl)-rt0),-epsi20) ) &
257                  - rcp      * ( ztmelts - rt0 ) )
258               ! Mutliply by ice volume, and divide by number of layers to get heat content in J/m2
259               e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * a_i(ji,jj,jl) * ht_i(ji,jj,jl) * r1_nlay_i
260            END DO
261            !
262         END DO
263         !
264         CALL lbc_bdy_lnk(  a_i(:,:,jl), 'T', 1., ib_bdy )
265         CALL lbc_bdy_lnk( ht_i(:,:,jl), 'T', 1., ib_bdy )
266         CALL lbc_bdy_lnk( ht_s(:,:,jl), 'T', 1., ib_bdy )
267         CALL lbc_bdy_lnk(  v_i(:,:,jl), 'T', 1., ib_bdy )
268         CALL lbc_bdy_lnk(  v_s(:,:,jl), 'T', 1., ib_bdy )
269         !
270         CALL lbc_bdy_lnk( smv_i(:,:,jl), 'T', 1., ib_bdy )
271         CALL lbc_bdy_lnk(  sm_i(:,:,jl), 'T', 1., ib_bdy )
272         CALL lbc_bdy_lnk(  oa_i(:,:,jl), 'T', 1., ib_bdy )
273         CALL lbc_bdy_lnk(  t_su(:,:,jl), 'T', 1., ib_bdy )
274         DO jk = 1, nlay_s
275            CALL lbc_bdy_lnk(t_s(:,:,jk,jl), 'T', 1., ib_bdy )
276            CALL lbc_bdy_lnk(e_s(:,:,jk,jl), 'T', 1., ib_bdy )
277         END DO
278         DO jk = 1, nlay_i
279            CALL lbc_bdy_lnk(t_i(:,:,jk,jl), 'T', 1., ib_bdy )
280            CALL lbc_bdy_lnk(e_i(:,:,jk,jl), 'T', 1., ib_bdy )
281         END DO
282         !
283      END DO !jl
284      !
285#endif
286      !     
287      IF( nn_timing == 1 )   CALL timing_stop('bdy_ice_frs')
288      !
289   END SUBROUTINE bdy_ice_frs
290
291
292   SUBROUTINE bdy_ice_lim_dyn( cd_type )
293      !!------------------------------------------------------------------------------
294      !!                 ***  SUBROUTINE bdy_ice_lim_dyn  ***
295      !!                   
296      !! ** Purpose : Apply dynamics boundary conditions for sea-ice in the cas of unstructured open boundaries.
297      !!              u_ice and v_ice are equal to the value of the adjacent grid point if this latter is not ice free
298      !!              if adjacent grid point is ice free, then u_ice and v_ice are equal to ocean velocities
299      !!
300      !! 2013-06 : C. Rousset
301      !!------------------------------------------------------------------------------
302      CHARACTER(len=1), INTENT(in)  ::   cd_type   ! nature of velocity grid-points
303      !
304      INTEGER  ::   jb, jgrd           ! dummy loop indices
305      INTEGER  ::   ji, jj             ! local scalar
306      INTEGER  ::   ib_bdy             ! Loop index
307      REAL(wp) ::   zmsk1, zmsk2, zflag
308      !!------------------------------------------------------------------------------
309      !
310      IF( nn_timing == 1 ) CALL timing_start('bdy_ice_lim_dyn')
311      !
312      DO ib_bdy=1, nb_bdy
313         !
314         SELECT CASE( cn_ice_lim(ib_bdy) )
315         !
316         CASE('none')
317            CYCLE
318            !
319         CASE('frs')
320            !
321            IF( nn_ice_lim_dta(ib_bdy) == 0 ) CYCLE            ! case ice boundaries = initial conditions
322            !                                                  !      do not change ice velocity (it is only computed by rheology)
323            SELECT CASE ( cd_type )
324            !     
325            CASE ( 'U' ) 
326               jgrd = 2      ! u velocity
327               DO jb = 1, idx_bdy(ib_bdy)%nblen(jgrd)
328                  ji    = idx_bdy(ib_bdy)%nbi(jb,jgrd)
329                  jj    = idx_bdy(ib_bdy)%nbj(jb,jgrd)
330                  zflag = idx_bdy(ib_bdy)%flagu(jb,jgrd)
331                  !
332                  IF ( ABS( zflag ) == 1. ) THEN  ! eastern and western boundaries
333                     ! one of the two zmsk is always 0 (because of zflag)
334                     zmsk1 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji+1,jj) ) ) ! 0 if no ice
335                     zmsk2 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji-1,jj) ) ) ! 0 if no ice
336                     
337                     ! u_ice = u_ice of the adjacent grid point except if this grid point is ice-free (then u_ice = u_oce)
338                     u_ice (ji,jj) = u_ice(ji+1,jj) * 0.5_wp * ABS( zflag + 1._wp ) * zmsk1 + &
339                        &            u_ice(ji-1,jj) * 0.5_wp * ABS( zflag - 1._wp ) * zmsk2 + &
340                        &            u_oce(ji  ,jj) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) )
341                  ELSE                             ! everywhere else
342                     !u_ice(ji,jj) = u_oce(ji,jj)
343                     u_ice(ji,jj) = 0._wp
344                  ENDIF
345                  ! mask ice velocities
346                  rswitch = MAX( 0.0_wp , SIGN ( 1.0_wp , at_i(ji,jj) - 0.01_wp ) ) ! 0 if no ice
347                  u_ice(ji,jj) = rswitch * u_ice(ji,jj)
348                  !
349               END DO
350               CALL lbc_bdy_lnk( u_ice(:,:), 'U', -1., ib_bdy )
351               !
352            CASE ( 'V' )
353               jgrd = 3      ! v velocity
354               DO jb = 1, idx_bdy(ib_bdy)%nblen(jgrd)
355                  ji    = idx_bdy(ib_bdy)%nbi(jb,jgrd)
356                  jj    = idx_bdy(ib_bdy)%nbj(jb,jgrd)
357                  zflag = idx_bdy(ib_bdy)%flagv(jb,jgrd)
358                  !
359                  IF ( ABS( zflag ) == 1. ) THEN  ! northern and southern boundaries
360                     ! one of the two zmsk is always 0 (because of zflag)
361                     zmsk1 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji,jj+1) ) ) ! 0 if no ice
362                     zmsk2 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji,jj-1) ) ) ! 0 if no ice
363                     
364                     ! u_ice = u_ice of the adjacent grid point except if this grid point is ice-free (then u_ice = u_oce)
365                     v_ice (ji,jj) = v_ice(ji,jj+1) * 0.5_wp * ABS( zflag + 1._wp ) * zmsk1 + &
366                        &            v_ice(ji,jj-1) * 0.5_wp * ABS( zflag - 1._wp ) * zmsk2 + &
367                        &            v_oce(ji,jj  ) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) )
368                  ELSE                             ! everywhere else
369                     !v_ice(ji,jj) = v_oce(ji,jj)
370                     v_ice(ji,jj) = 0._wp
371                  ENDIF
372                  ! mask ice velocities
373                  rswitch = MAX( 0.0_wp , SIGN ( 1.0_wp , at_i(ji,jj) - 0.01 ) ) ! 0 if no ice
374                  v_ice(ji,jj) = rswitch * v_ice(ji,jj)
375                  !
376               END DO
377               CALL lbc_bdy_lnk( v_ice(:,:), 'V', -1., ib_bdy )
378               !
379            END SELECT
380            !
381         CASE DEFAULT
382            CALL ctl_stop( 'bdy_ice_lim_dyn : unrecognised option for open boundaries for ice fields' )
383         END SELECT
384         !
385      END DO
386      !
387      IF( nn_timing == 1 ) CALL timing_stop('bdy_ice_lim_dyn')
388      !
389    END SUBROUTINE bdy_ice_lim_dyn
390
391#else
392   !!---------------------------------------------------------------------------------
393   !!   Default option                                                    Empty module
394   !!---------------------------------------------------------------------------------
395CONTAINS
396   SUBROUTINE bdy_ice_lim( kt )      ! Empty routine
397      WRITE(*,*) 'bdy_ice_lim: You should not have seen this print! error?', kt
398   END SUBROUTINE bdy_ice_lim
399#endif
400
401   !!=================================================================================
402END MODULE bdyice_lim
Note: See TracBrowser for help on using the repository browser.