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/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/BDY – NEMO

source: branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/BDY/bdyice_lim.F90 @ 4294

Last change on this file since 4294 was 4292, checked in by cetlod, 10 years ago

dev_MERGE_2013 : 1st step of the merge, see ticket #1185

File size: 19.6 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#elif defined key_lim3
27   USE par_ice
28   USE ice             ! LIM_3 ice variables
29#endif
30   USE par_oce         ! ocean parameters
31   USE dom_oce         ! ocean space and time domain variables
32   USE dom_ice          ! sea-ice domain
33   USE sbc_oce         ! Surface boundary condition: ocean fields
34   USE bdy_oce         ! ocean open boundary conditions
35   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
36   USE in_out_manager  ! write to numout file
37   USE lib_mpp         ! distributed memory computing
38   USE lib_fortran     ! to use key_nosignedzero
39
40   IMPLICIT NONE
41   PRIVATE
42
43   PUBLIC   bdy_ice_lim    ! routine called in sbcmod
44   PUBLIC   bdy_ice_lim_dyn ! routine called in limrhg
45
46   REAL(wp) ::   epsi20 = 1.e-20_wp  ! module constants
47   !!----------------------------------------------------------------------
48   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
49   !! $Id: bdyice.F90 2715 2011-03-30 15:58:35Z rblod $
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      DO ib_bdy=1, nb_bdy
65
66         SELECT CASE( cn_ice_lim(ib_bdy) )
67         CASE('none')
68            CYCLE
69         CASE('frs')
70            CALL bdy_ice_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy )
71         CASE DEFAULT
72            CALL ctl_stop( 'bdy_ice_lim : unrecognised option for open boundaries for ice fields' )
73         END SELECT
74      ENDDO
75
76   END SUBROUTINE bdy_ice_lim
77
78   SUBROUTINE bdy_ice_frs( idx, dta, kt, ib_bdy )
79      !!------------------------------------------------------------------------------
80      !!                 ***  SUBROUTINE bdy_ice_frs  ***
81      !!                   
82      !! ** Purpose : Apply the Flow Relaxation Scheme for sea-ice fields in the case
83      !!              of unstructured open boundaries.
84      !!
85      !! Reference : Engedahl H., 1995: Use of the flow relaxation scheme in a three-
86      !!             dimensional baroclinic ocean model with realistic topography. Tellus, 365-382.
87      !!------------------------------------------------------------------------------
88      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices
89      TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data
90      INTEGER,         INTENT(in) ::   kt   ! main time-step counter
91      INTEGER,         INTENT(in) ::   ib_bdy  ! BDY set index      !!
92
93      INTEGER  ::   jb, jk, jgrd, jl   ! dummy loop indices
94      INTEGER  ::   ji, jj             ! local scalar
95      REAL(wp) ::   zwgt, zwgt1        ! local scalar
96      REAL(wp) ::   zinda, ztmelts
97      !!------------------------------------------------------------------------------
98      !
99      IF( nn_timing == 1 ) CALL timing_start('bdy_ice_frs')
100      !
101      !IF( kt /= nit000 ) THEN
102      jgrd = 1      ! Everything is at T-points here
103      !
104#if defined key_lim2
105      DO jb = 1, idx%nblen(jgrd)
106         ji    = idx%nbi(jb,jgrd)
107         jj    = idx%nbj(jb,jgrd)
108         zwgt  = idx%nbw(jb,jgrd)
109         zwgt1 = 1.e0 - idx%nbw(jb,jgrd)
110         frld (ji,jj) = ( frld (ji,jj) * zwgt1 + dta%frld (jb) * zwgt ) * tmask(ji,jj,1)     ! Leads fraction
111         hicif(ji,jj) = ( hicif(ji,jj) * zwgt1 + dta%hicif(jb) * zwgt ) * tmask(ji,jj,1)     ! Ice depth
112         hsnif(ji,jj) = ( hsnif(ji,jj) * zwgt1 + dta%hsnif(jb) * zwgt ) * tmask(ji,jj,1)     ! Snow depth
113
114!         zinda = 1.0 - MAX( 0.0_wp , SIGN ( 1.0_wp , - frld(ji,jj) ) ) ! 0 if no ice
115!         !------------------------------
116!         ! Sea ice surface temperature
117!         !------------------------------
118!         sist(ji,jj)   = zinda * 270.0 + ( 1.0 - zinda ) * tfu(ji,jj)
119!         !-----------------------------------------------
120!         ! Ice/snow temperatures and energy stored in brines
121!         !-----------------------------------------------
122!         !!! TO BE CONTIUNED (as LIM3 below) !!!
123!            zindhe = MAX( 0.e0, SIGN( 1.e0, fcor(1,jj) ) )              ! = 0 for SH, =1 for NH
124!
125!               ! Recover in situ values.
126!               zindb         = MAX( rzero, SIGN( rone, zs0a(ji,jj) - epsi06 ) )
127!               zacrith       = 1.0 - ( zindhe * acrit(1) + ( 1.0 - zindhe ) * acrit(2) )
128!               zs0a (ji,jj)  = zindb * MIN( zs0a(ji,jj), zacrith )
129!               hsnif(ji,jj)  = zindb * ( zs0sn(ji,jj) /MAX( zs0a(ji,jj), epsi16 ) )
130!               hicif(ji,jj)  = zindb * ( zs0ice(ji,jj)/MAX( zs0a(ji,jj), epsi16 ) )
131!               zindsn        = MAX( rzero, SIGN( rone, hsnif(ji,jj) - epsi06 ) )
132!               zindic        = MAX( rzero, SIGN( rone, hicif(ji,jj) - epsi03 ) )
133!               zindb         = MAX( zindsn, zindic )
134!               zs0a (ji,jj)  = zindb * zs0a(ji,jj)
135!               frld (ji,jj)  = 1.0 - zs0a(ji,jj)
136!               hsnif(ji,jj)  = zindsn * hsnif(ji,jj)
137!               hicif(ji,jj)  = zindic * hicif(ji,jj)
138!               zusvosn       = 1.0/MAX( hsnif(ji,jj) * zs0a(ji,jj), epsi16 )
139!               zusvoic       = 1.0/MAX( hicif(ji,jj) * zs0a(ji,jj), epsi16 )
140!               zignm         = MAX( rzero,  SIGN( rone, hsndif - hsnif(ji,jj) ) )
141!               zrtt          = 173.15 * rone
142!               ztsn          =          zignm   * tbif(ji,jj,1)  &
143!                              + ( 1.0 - zignm ) * MIN( MAX( zrtt, rt0_snow * zusvosn * zs0c0(ji,jj)) , tfu(ji,jj) )
144!               ztic1          = MIN( MAX( zrtt, rt0_ice * zusvoic * zs0c1(ji,jj) ) , tfu(ji,jj) )
145!               ztic2          = MIN( MAX( zrtt, rt0_ice * zusvoic * zs0c2(ji,jj) ) , tfu(ji,jj) )
146!
147!               tbif(ji,jj,1) = zindsn * ztsn  + ( 1.0 - zindsn ) * tfu(ji,jj)               
148!               tbif(ji,jj,2) = zindic * ztic1 + ( 1.0 - zindic ) * tfu(ji,jj)
149!               tbif(ji,jj,3) = zindic * ztic2 + ( 1.0 - zindic ) * tfu(ji,jj)
150!               qstoif(ji,jj) = zindb  * xlic * zs0st(ji,jj) /  MAX( zs0a(ji,jj), epsi16 )
151
152      END DO
153
154      CALL lbc_bdy_lnk( frld,  'T', 1., ib_bdy )                                         ! lateral boundary conditions
155      CALL lbc_bdy_lnk( hicif, 'T', 1., ib_bdy )
156      CALL lbc_bdy_lnk( hsnif, 'T', 1., ib_bdy )
157
158      vt_i(:,:) = hicif(:,:) * frld(:,:)
159      vt_s(:,:) = hsnif(:,:) * frld(:,:)
160      !
161#elif defined key_lim3
162
163      DO jl = 1, jpl
164         DO jb = 1, idx%nblen(jgrd)
165            ji    = idx%nbi(jb,jgrd)
166            jj    = idx%nbj(jb,jgrd)
167            zwgt  = idx%nbw(jb,jgrd)
168            zwgt1 = 1.e0 - idx%nbw(jb,jgrd)
169            a_i (ji,jj,jl) = ( a_i (ji,jj,jl) * zwgt1 + dta%a_i (jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Leads fraction
170            ht_i(ji,jj,jl) = ( ht_i(ji,jj,jl) * zwgt1 + dta%ht_i(jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Ice depth
171            ht_s(ji,jj,jl) = ( ht_s(ji,jj,jl) * zwgt1 + dta%ht_s(jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Snow depth
172         ENDDO
173         CALL lbc_bdy_lnk(  a_i(:,:,jl), 'T', 1., ib_bdy )
174         CALL lbc_bdy_lnk( ht_i(:,:,jl), 'T', 1., ib_bdy )
175         CALL lbc_bdy_lnk( ht_s(:,:,jl), 'T', 1., ib_bdy )
176
177         DO jb = 1, idx%nblen(jgrd)
178            ji    = idx%nbi(jb,jgrd)
179            jj    = idx%nbj(jb,jgrd)
180            ! clem : condition on ice thickness depends on the ice velocity
181            ! if velocity is outward (strictly), then ice thickness, volume... must be equal to adjacent values
182            IF ( u_ice(ji+1,jj  ) < 0. .AND. umask(ji-1,jj  ,1) == 0. ) THEN
183               ht_i(ji,jj,jl) = ht_i(ji+1,jj  ,jl)
184               a_i (ji,jj,jl) = a_i (ji+1,jj  ,jl)
185               ht_s(ji,jj,jl) = ht_s(ji+1,jj  ,jl)
186            ENDIF
187            IF ( u_ice(ji-1,jj  ) > 0. .AND. umask(ji+1,jj  ,1) == 0. ) THEN
188               ht_i(ji,jj,jl) = ht_i(ji-1,jj  ,jl)
189               a_i (ji,jj,jl) = a_i (ji-1,jj  ,jl)
190               ht_s(ji,jj,jl) = ht_s(ji-1,jj  ,jl)
191            ENDIF
192            IF ( v_ice(ji  ,jj+1) < 0. .AND. vmask(ji  ,jj-1,1) == 0. ) THEN
193               ht_i(ji,jj,jl) = ht_i(ji  ,jj+1,jl)
194               a_i (ji,jj,jl) = a_i (ji  ,jj+1,jl)
195               ht_s(ji,jj,jl) = ht_s(ji  ,jj+1,jl)
196            ENDIF
197            IF ( v_ice(ji  ,jj-1) > 0. .AND. vmask(ji  ,jj+1,1) == 0. ) THEN
198               ht_i(ji,jj,jl) = ht_i(ji  ,jj-1,jl)
199               a_i (ji,jj,jl) = a_i (ji  ,jj-1,jl)
200               ht_s(ji,jj,jl) = ht_s(ji  ,jj-1,jl)
201            ENDIF
202
203            zinda = 1.0 - MAX( 0.0_wp , SIGN ( 1.0_wp , - a_i(ji,jj,jl) ) ) ! 0 if no ice
204            ! --------------------
205            ! Ice and snow volumes
206            ! --------------------
207            v_i(ji,jj,jl) = ht_i(ji,jj,jl) * a_i(ji,jj,jl)
208            v_s(ji,jj,jl) = ht_s(ji,jj,jl) * a_i(ji,jj,jl)
209
210            ! -------------
211            ! Ice salinity
212            !---------------
213            sm_i(ji,jj,jl)   = zinda * bulk_sal  + ( 1.0 - zinda ) * s_i_min
214            smv_i(ji,jj,jl)  = MIN( sm_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl)
215           
216            !----------
217            ! Ice age
218            !----------
219            o_i(ji,jj,jl)    = zinda * 1.0   + ( 1.0 - zinda )
220            oa_i(ji,jj,jl)   = o_i(ji,jj,jl) * a_i(ji,jj,jl)
221           
222            !------------------------------
223            ! Sea ice surface temperature
224            !------------------------------
225            t_su(ji,jj,jl)   = zinda * 270.0 + ( 1.0 - zinda ) * t_bo(ji,jj)
226           
227            !------------------------------------
228            ! Snow temperature and heat content
229            !------------------------------------
230            DO jk = 1, nlay_s
231               t_s(ji,jj,jk,jl) = zinda * 270.00 + ( 1.0 - zinda ) * rtt
232               ! Snow energy of melting
233               e_s(ji,jj,jk,jl) = zinda * rhosn * ( cpic * ( rtt - t_s(ji,jj,jk,jl) ) + lfus )
234               ! Change dimensions
235               e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) / unit_fac
236               ! Multiply by volume, so that heat content in 10^9 Joules
237               e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * area(ji,jj) * &
238                    v_s(ji,jj,jl)  / nlay_s
239            END DO !jk
240           
241            !-----------------------------------------------
242            ! Ice salinities, temperature and heat content
243            !-----------------------------------------------
244            DO jk = 1, nlay_i
245               t_i(ji,jj,jk,jl) = zinda * 270.00 + ( 1.0 - zinda ) * rtt 
246               s_i(ji,jj,jk,jl) = zinda * bulk_sal + ( 1.0 - zinda ) * s_i_min
247               ztmelts          = - tmut * s_i(ji,jj,jk,jl) + rtt !Melting temperature in K
248               
249               ! heat content per unit volume
250               e_i(ji,jj,jk,jl) = zinda * rhoic * &
251                    (   cpic    * ( ztmelts - t_i(ji,jj,jk,jl) ) &
252                    +   lfus    * ( 1.0 - (ztmelts-rtt) / MIN((t_i(ji,jj,jk,jl)-rtt),-epsi20) ) &
253                    - rcp      * ( ztmelts - rtt ) )
254               
255               ! Correct dimensions to avoid big values
256               e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) / unit_fac 
257               
258               ! Mutliply by ice volume, and divide by number of layers to get heat content in 10^9 J
259               e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * & 
260                    area(ji,jj) * a_i(ji,jj,jl) * ht_i(ji,jj,jl) / nlay_i
261            END DO ! jk
262
263         END DO !jb
264 
265
266         CALL lbc_bdy_lnk(  a_i(:,:,jl), 'T', 1., ib_bdy )                                         ! lateral boundary conditions
267         CALL lbc_bdy_lnk( ht_i(:,:,jl), 'T', 1., ib_bdy )
268         CALL lbc_bdy_lnk( ht_s(:,:,jl), 'T', 1., ib_bdy )
269         CALL lbc_bdy_lnk(  v_i(:,:,jl), 'T', 1., ib_bdy )
270         CALL lbc_bdy_lnk(  v_s(:,:,jl), 'T', 1., ib_bdy )
271 
272         CALL lbc_bdy_lnk( smv_i(:,:,jl), 'T', 1., ib_bdy )
273         CALL lbc_bdy_lnk(  sm_i(:,:,jl), 'T', 1., ib_bdy )
274         CALL lbc_bdy_lnk(  oa_i(:,:,jl), 'T', 1., ib_bdy )
275         CALL lbc_bdy_lnk(   o_i(:,:,jl), 'T', 1., ib_bdy )
276         CALL lbc_bdy_lnk(  t_su(:,:,jl), 'T', 1., ib_bdy )
277         DO jk = 1, nlay_s
278            CALL lbc_bdy_lnk(t_s(:,:,jk,jl), 'T', 1., ib_bdy )
279            CALL lbc_bdy_lnk(e_s(:,:,jk,jl), 'T', 1., ib_bdy )
280         END DO
281         DO jk = 1, nlay_i
282            CALL lbc_bdy_lnk(t_i(:,:,jk,jl), 'T', 1., ib_bdy )
283            CALL lbc_bdy_lnk(e_i(:,:,jk,jl), 'T', 1., ib_bdy )
284         END DO
285
286      END DO !jl
287   
288#endif
289      !     
290      !ENDIF !nit000/=0
291      IF( nn_timing == 1 ) CALL timing_stop('bdy_ice_frs')
292      !
293   END SUBROUTINE bdy_ice_frs
294
295
296   SUBROUTINE bdy_ice_lim_dyn( kn, pvar1, pvar2, pvar12 )
297      !!------------------------------------------------------------------------------
298      !!                 ***  SUBROUTINE bdy_ice_lim_dyn  ***
299      !!                   
300      !! ** Purpose : Apply dynamics boundary conditions for sea-ice in the cas of unstructured open boundaries.
301      !!      kn = 1: u_ice and v_ice are equal to the value of the adjacent grid point if this latter is not ice free
302      !!              if adjacent grid point is ice free, then u_ice and v_ice are equal to ocean velocities
303      !!      kn = 2: the stress tensor is set to 0 (i.e. pvar1, pvar2, pvar12)
304      !!
305      !! 2013-06 : C. Rousset
306      !!------------------------------------------------------------------------------
307      !!
308      INTEGER,  INTENT( in    )                           ::   kn     ! set up of ice vel (kn=1) or stress tensor (kn=2)
309      REAL(wp), INTENT( inout ), DIMENSION(:,:), OPTIONAL ::   pvar1, pvar2, pvar12  ! stress tensor components (zs1,zs2,zs12)
310      INTEGER  ::   jb, jgrd   ! dummy loop indices
311      INTEGER  ::   ji, jj             ! local scalar
312      INTEGER  ::   ib_bdy ! Loop index
313      REAL(wp) ::   zmsk1, zmsk2, zflag, zinda 
314     !!------------------------------------------------------------------------------
315      !
316      IF( nn_timing == 1 ) CALL timing_start('bdy_ice_lim_dyn')
317      !
318      DO ib_bdy=1, nb_bdy
319      !
320      SELECT CASE( nn_ice_lim(ib_bdy) )
321
322      CASE(jp_none)
323         CYCLE
324
325      CASE(jp_frs)
326
327         IF ( kn == 1 ) THEN ! set up of u_ice and v_ice
328           
329            jgrd = 2      ! u velocity
330            DO jb = 1, idx_bdy(ib_bdy)%nblen(jgrd)
331               ji    = idx_bdy(ib_bdy)%nbi(jb,jgrd)
332               jj    = idx_bdy(ib_bdy)%nbj(jb,jgrd)
333               zflag = idx_bdy(ib_bdy)%flagu(jb)
334
335               IF ( ABS( zflag ) == 1. ) THEN  ! eastern and western boundaries
336                  ! one of the two zmsk is always 0 (because of zflag)
337                  zmsk1 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji+1,jj) ) ) ! 0 if no ice
338                  zmsk2 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji-1,jj) ) ) ! 0 if no ice
339                 
340                  ! u_ice = u_ice of the adjacent grid point except if this grid point is ice-free (then u_ice = u_oce)
341                  u_ice (ji,jj) = u_ice(ji+1,jj) * 0.5 * ABS( zflag + 1._wp ) * zmsk1 + &
342                     &            u_ice(ji-1,jj) * 0.5 * ABS( zflag - 1._wp ) * zmsk2 + &
343                     &            u_oce(ji  ,jj) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) )
344               ELSE                             ! everywhere else
345                  u_ice(ji,jj) = u_oce(ji,jj)
346               ENDIF
347               ! mask ice velocities
348               zinda = 1.0 - MAX( 0.0_wp, SIGN( 1.0_wp , - at_i(ji,jj) ) ) ! 0 if no ice
349               u_ice(ji,jj) = zinda * u_ice(ji,jj)
350
351            ENDDO
352           
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)
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 * ABS( zflag + 1._wp ) * zmsk1 + &
366                     &            v_ice(ji,jj-1) * 0.5 * 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               ENDIF
371               ! mask ice velocities
372               zinda = 1.0 - MAX( 0.0_wp, SIGN( 1.0_wp , - at_i(ji,jj) ) ) ! 0 if no ice
373               v_ice(ji,jj) = zinda * v_ice(ji,jj)
374
375            ENDDO
376           
377         CALL lbc_bdy_lnk( u_ice(:,:), 'U', -1., ib_bdy )
378         CALL lbc_bdy_lnk( v_ice(:,:), 'V', -1., ib_bdy )
379 
380         ELSEIF ( kn == 2 ) THEN ! set up of stress tensor (not sure it works)
381         
382            jgrd = 1      ! T grid
383            DO jb = 1, idx_bdy(ib_bdy)%nblen(jgrd)
384               ji    = idx_bdy(ib_bdy)%nbi(jb,jgrd)
385               jj    = idx_bdy(ib_bdy)%nbj(jb,jgrd)
386
387               pvar1 (ji,jj) = 0._wp
388               pvar2 (ji,jj) = 0._wp
389               pvar12(ji,jj) = 0._wp
390
391            ENDDO
392
393         ENDIF
394
395      CASE DEFAULT
396         CALL ctl_stop( 'bdy_ice_lim_dyn : unrecognised option for open boundaries for ice fields' )
397      END SELECT
398
399      ENDDO
400
401      IF( nn_timing == 1 ) CALL timing_stop('bdy_ice_lim_dyn')
402 
403    END SUBROUTINE bdy_ice_lim_dyn
404
405#else
406   !!---------------------------------------------------------------------------------
407   !!   Default option                                                    Empty module
408   !!---------------------------------------------------------------------------------
409CONTAINS
410   SUBROUTINE bdy_ice_lim( kt )      ! Empty routine
411      WRITE(*,*) 'bdy_ice_lim: You should not have seen this print! error?', kt
412   END SUBROUTINE bdy_ice_lim
413#endif
414
415   !!=================================================================================
416END MODULE bdyice_lim
Note: See TracBrowser for help on using the repository browser.