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.
limistate.F90 on Ticket #856 – Attachment – NEMO

Ticket #856: limistate.F90

File limistate.F90, 23.3 KB (added by vancop, 13 years ago)

routine

Line 
1MODULE limistate
2   !!======================================================================
3   !!                     ***  MODULE  limistate  ***
4   !!              Initialisation of diagnostics ice variables
5   !!======================================================================
6#if defined key_lim3
7   !!----------------------------------------------------------------------
8   !!   'key_lim3' :                                    LIM3 sea-ice model
9   !!----------------------------------------------------------------------
10   !!   lim_istate      :  Initialisation of diagnostics ice variables
11   !!   lim_istate_init :  initialization of ice state and namelist read
12   !!----------------------------------------------------------------------
13   !! * Modules used
14   USE phycst
15   USE oce             ! dynamics and tracers variables
16   USE dom_oce
17   USE sbc_oce         ! Surface boundary condition: ocean fields
18   USE par_ice         ! ice parameters
19   USE ice_oce         ! ice variables
20   USE eosbn2          ! equation of state
21   USE in_out_manager
22   USE dom_ice
23   USE ice
24   USE lbclnk
25
26   IMPLICIT NONE
27   PRIVATE
28
29   !! * Accessibility
30   PUBLIC lim_istate      ! routine called by lim_init.F90
31
32   !! * Module variables
33   REAL(wp) ::             & !!! ** init namelist (namiceini) **
34      ttest    = 2.0  ,    &  ! threshold water temperature for initial sea ice
35      hninn    = 0.5  ,    &  ! initial snow thickness in the north
36      hnins    = 0.1  ,    &  ! initial snow thickness in the south
37      hginn    = 2.5  ,    &  ! initial ice thickness in the north
38      hgins    = 1.0  ,    &  ! initial ice thickness in the south
39      aginn    = 0.7  ,    &  ! initial leads area in the north
40      agins    = 0.7  ,    &  ! initial leads area in the south
41      sinn     = 6.301 ,   &  ! initial salinity
42      sins     = 6.301
43
44   REAL(wp)  ::            &  ! constant values
45      zzero   = 0.0     ,  &
46      zone    = 1.0
47
48   !!----------------------------------------------------------------------
49   !!   LIM 3.0,  UCL-LOCEAN-IPSL (2008)
50   !! $Id: limistate.F90 1156 2008-06-26 16:06:45Z rblod $
51   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
52   !!----------------------------------------------------------------------
53
54CONTAINS
55
56   SUBROUTINE lim_istate
57      !!-------------------------------------------------------------------
58      !!                    ***  ROUTINE lim_istate  ***
59      !!
60      !! ** Purpose :   defined the sea-ice initial state
61      !!
62      !! ** Method  :   
63      !!                This routine will put some ice where ocean
64      !!                is at the freezing point, then fill in ice
65      !!                state variables using prescribed initial
66      !!                values in the namelist           
67      !!
68      !! ** Steps   :   
69      !!                1) Read namelist
70      !!                2) Basal temperature; ice and hemisphere masks
71      !!                3) Fill in the ice thickness distribution using gaussian
72      !!                4) Fill in space-dependent arrays for state variables
73      !!                5) Diagnostic arrays
74      !!                6) Lateral boundary conditions
75      !!
76      !! History :
77      !!   2.0  !  01-04  (C. Ethe, G. Madec)  Original code
78      !!   3.0  !  2007   (M. Vancoppenolle)   Rewrite for ice cats
79      !!   4.0  !  09-11  (M. Vancoppenolle)   Enhanced version for ice cats
80      !!--------------------------------------------------------------------
81
82      !! * Local variables
83      INTEGER  ::   ji, jj, jk, jl             ! dummy loop indices
84
85      REAL(wp) ::       &  ! temporary scalar
86         epsi06, epsi20, ztmelts
87
88      REAL(wp), DIMENSION(jpi,jpj) ::   zidto    ! ice indicator
89
90      INTEGER  :: i_hemis, i_fill, jl0 
91
92      REAL(wp) ::  &
93         ztest_1, ztest_2, ztest_3, ztest_4, ztests, zsigma, zarg, zA, zV, &
94         zA_cons, zV_cons, zconv
95
96      REAL(wp), DIMENSION(jpl,2) :: &
97         zht_i_ini, za_i_ini, zv_i_ini
98
99      REAL(wp), DIMENSION(2) :: &
100         zhm_i_ini, zat_i_ini, zvt_i_ini, &
101         zhm_s_ini, zsm_i_ini
102
103      INTEGER(wp), DIMENSION(jpi,jpj) ::   zhemis   ! hemispheric index
104      !--------------------------------------------------------------------
105
106      epsi06   = 1.0e-6
107      epsi20   = 1.0e-20
108      IF(lwp) WRITE(numout,*)
109      IF(lwp) WRITE(numout,*) 'lim_istate : Ice initialization '
110      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ '
111
112      !--------------------------------------------------------------------
113      ! 1) Read namelist
114      !--------------------------------------------------------------------
115
116      CALL lim_istate_init     !  reading the initials parameters of the ice
117
118!!gm  in lim2  the initialisation if only done if required in the namelist :
119!!gm      IF( .NOT. ln_limini ) THEN
120!!gm  this should be added in lim3 namelist...
121
122      !--------------------------------------------------------------------
123      ! 2) Basal temperature, ice mask and hemispheric index
124      !--------------------------------------------------------------------
125
126      ! Basal temperature is set to the freezing point of seawater in Celsius
127      t_bo(:,:) = tfreez( sn(:,:,1) ) * tmask(:,:,1)       ! freezing/melting point of sea water [Celcius]
128
129      DO jj = 1, jpj                                       ! ice if sst <= t-freez + ttest
130         DO ji = 1, jpi
131            IF( tn(ji,jj,1)  - t_bo(ji,jj) >= ttest ) THEN   ;   zidto(ji,jj) = 0.e0      ! no ice
132            ELSE                                             ;   zidto(ji,jj) = 1.e0      !    ice
133            ENDIF
134         END DO
135      END DO
136
137      t_bo(:,:) = t_bo(:,:) + rt0                          ! conversion to Kelvin
138
139      ! Hemispheric index
140      ! MV 2011 new initialization
141      DO jj = 1, jpj
142         DO ji = 1, jpi
143            IF( fcor(ji,jj) >= 0.e0 ) THEN   
144               zhemis(ji,jj) = 1 ! Northern hemisphere
145            ELSE
146               zhemis(ji,jj) = 2 ! Southern hemisphere
147            ENDIF
148         END DO
149      END DO
150      ! END MV 2011 new initialization
151
152      !--------------------------------------------------------------------
153      ! 3) Initialization of sea ice state variables
154      !--------------------------------------------------------------------
155
156      !-----------------------------
157      ! 3.1) Hemisphere-dependent arrays
158      !-----------------------------
159      ! assign initial thickness, concentration, snow depth and salinity to
160      ! an hemisphere-dependent array
161      zhm_i_ini(1) = hginn ; zhm_i_ini(2) = hgins  ! ice thickness
162      zat_i_ini(1) = aginn ; zat_i_ini(2) = agins  ! ice concentration
163      zvt_i_ini(:) = zhm_i_ini(:) * zat_i_ini(:)   ! ice volume
164      zhm_s_ini(1) = hninn ; zhm_s_ini(2) = hnins  ! snow depth
165      zsm_i_ini(1) = sinn  ; zsm_i_ini(2) = sins   ! bulk ice salinity
166
167      WRITE(numout,*) ' hgins, agins :' , hgins, agins
168
169      !---------------------------------------------------------------------
170      ! 3.2) Distribute ice concentration and thickness into the categories
171      !---------------------------------------------------------------------
172      ! a gaussian distribution for ice concentration is used
173      ! then we check whether the distribution fullfills
174      ! volume and area conservation, positivity and ice categories bounds
175      DO i_hemis = 1, 2
176
177      ztest_1 = 0 ; ztest_2 = 0 ; ztest_3 = 0 ; ztest_4 = 0
178
179      ! note for the great nemo engineers:
180      ! only very few of the WRITE statements are necessary for the reference version
181      ! they were one day useful, but now i personally doubt of their
182      ! potential for bringing anything useful
183
184      WRITE(numout,*)
185      WRITE(numout,*) ' ------------------------------------------ '
186      WRITE(numout,*) ' HEMISPHERE ', i_hemis
187      WRITE(numout,*) ' ------------------------------------------ '
188      WRITE(numout,*)
189      WRITE(numout,*) ' zhm_i_ini : ', zhm_i_ini(i_hemis)
190      WRITE(numout,*) ' zat_i_ini : ', zat_i_ini(i_hemis)
191      WRITE(numout,*) ' zvt_i_ini : ', zvt_i_ini(i_hemis)
192      WRITE(numout,*)
193
194      DO i_fill = jpl, 1, -1
195         IF ( ( ztest_1 + ztest_2 + ztest_3 + ztest_4 ) .NE. 4 ) THEN
196            WRITE(numout,*) ' ------------------------------------------ '
197            WRITE(numout,*) ' i_fill = ', i_fill
198            WRITE(numout,*) ' ------------------------------------------ '
199            WRITE(numout,*) ' hi_max : ', ( hi_max(jl), jl = 0, jpl )
200            !----------------------------
201            ! fill the i_fill categories
202            !----------------------------
203            ! *** 1 category to fill
204            IF ( i_fill .EQ. 1 ) THEN
205               zht_i_ini(1,i_hemis)       = zhm_i_ini(i_hemis)
206               za_i_ini(1,i_hemis)        = zat_i_ini(i_hemis)
207               zht_i_ini(2:jpl,i_hemis)   = 0.0d0
208               za_i_ini(2:jpl,i_hemis)    = 0.0d0
209               WRITE(numout,*) ' * zht_i_ini : ', ( zht_i_ini(jl,i_hemis), jl = 1, jpl )
210               WRITE(numout,*) ' * za_i_ini :  ', ( za_i_ini(jl,i_hemis), jl = 1, jpl )
211            ELSE
212
213            ! *** >1 categores to fill
214            !--- Ice thicknesses in the i_fill - 1 first categories
215               DO jl = 1, i_fill - 1
216                  zht_i_ini(jl,i_hemis)    = ( hi_max(jl) + hi_max(jl-1) ) / 2.
217               END DO
218               WRITE(numout,*) ' * zht_i_ini : ', ( zht_i_ini(jl,i_hemis), jl = 1, i_fill )
219
220            !--- jl0: most likely index where cc will be maximum
221               DO jl = 1, jpl
222                  IF ( ( zhm_i_ini(i_hemis) .GT. hi_max(jl-1) ) .AND. &
223                       ( zhm_i_ini(i_hemis) .LE. hi_max(jl)   ) ) THEN
224                     jl0 = jl
225                  ENDIF
226               END DO
227               jl0 = MIN(jl0, i_fill)
228               WRITE(numout,*) ' * jl0 = ', jl0
229
230            !--- Concentrations
231               za_i_ini(jl0,i_hemis)      = zat_i_ini(i_hemis) / SQRT(REAL(jpl))
232               DO jl = 1, i_fill - 1
233                  IF ( jl .NE. jl0 ) THEN
234                     zsigma               = zhm_i_ini(i_hemis) / 2.
235                     zarg                 = ( zht_i_ini(jl,i_hemis) - zhm_i_ini(i_hemis) ) / zsigma
236                     za_i_ini(jl,i_hemis) = za_i_ini(jl0,i_hemis) * EXP(-zarg**2)
237                  ENDIF
238               END DO
239
240               zA = 0. ! sum of the areas in the jpl categories
241               DO jl = 1, i_fill - 1
242                 zA = zA + za_i_ini(jl,i_hemis)
243               END DO
244               za_i_ini(i_fill,i_hemis)   = zat_i_ini(i_hemis) - zA ! ice conc in the last category
245               IF ( i_fill .LT. jpl ) za_i_ini(i_fill+1:jpl, i_hemis) = 0.
246         
247               WRITE(numout,*) ' * za_i_ini : ', ( za_i_ini(jl,i_hemis), jl = 1, jpl )
248
249            !--- Ice thickness in the last category
250               zV = 0. ! sum of the volumes of the N-1 categories
251               DO jl = 1, i_fill - 1
252                  zV = zV + za_i_ini(jl,i_hemis)*zht_i_ini(jl,i_hemis)
253               END DO
254               zht_i_ini(i_fill,i_hemis) = ( zvt_i_ini(i_hemis) - zV ) / za_i_ini(i_fill,i_hemis) 
255               IF ( i_fill .LT. jpl ) zht_i_ini(i_fill+1:jpl, i_hemis) = 0.
256
257               WRITE(numout,*) ' * zht_i_ini : ', ( zht_i_ini(jl,i_hemis), jl = 1, jpl )
258
259            !--- volumes
260               zv_i_ini(:,i_hemis) = za_i_ini(:,i_hemis) * zht_i_ini(:,i_hemis)
261               IF ( i_fill .LT. jpl ) zv_i_ini(i_fill+1:jpl, i_hemis) = 0.
262               WRITE(numout,*) ' * zv_i_ini : ', ( zv_i_ini(jl,i_hemis), jl = 1, jpl )
263
264            ENDIF ! i_fill
265
266            !---------------------
267            ! Compatibility tests
268            !---------------------
269            ! Test 1: area conservation
270            zA_cons = SUM(za_i_ini(:,i_hemis)) ; zconv = ABS(zat_i_ini(i_hemis) - zA_cons )
271            IF ( zconv .LT. 1.0e-6 ) THEN
272!              WRITE(numout,*) ' * TEST1 AREA CONSERVED '
273               ztest_1 = 1
274            ELSE 
275              ! this write is useful
276               WRITE(numout,*) ' * TEST1 AREA NOT CONSERVED *** zA_cons = ', zA_cons, &
277                               ' zat_i_ini = ',zat_i_ini(i_hemis) 
278               ztest_1 = 0
279            ENDIF
280
281            ! Test 2: volume conservation
282            zV_cons = SUM(zv_i_ini(:,i_hemis))
283            zconv = ABS(zvt_i_ini(i_hemis) - zV_cons)
284
285            IF ( zconv .LT. 1.0e-6 ) THEN
286!              WRITE(numout,*) ' * TEST2 VOLUME CONSERVED '
287               ztest_2 = 1
288            ELSE
289              ! this write is useful
290               WRITE(numout,*) ' * TEST2 VOLUME NOT CONSERVED *** zV_cons = ', zV_cons, &
291                            ' zvt_i_ini = ', zvt_i_ini(i_hemis)
292               ztest_2 = 0
293            ENDIF
294
295            ! Test 3: thickness of the last category is in-bounds ?
296            IF ( zht_i_ini(i_fill, i_hemis) .GT. hi_max(i_fill-1) ) THEN
297!              WRITE(numout,*) ' * TEST3 THICKNESS OF THE LAST CATEGORY OK '
298               ztest_3 = 1
299            ELSE
300               ! this write is useful
301               WRITE(numout,*) ' * TEST 3 THICKNESS OF THE LAST CATEGORY OUT OF BOUNDS *** zht_i_ini(i_fill,i_hemis) = ', &
302               zht_i_ini(i_fill,i_hemis), ' hi_max(jpl-1) = ', hi_max(i_fill-1)
303               ztest_3 = 0
304            ENDIF
305
306            ! Test 4: positivity of ice concentrations
307            ztest_4 = 1
308            DO jl = 1, jpl
309               IF ( za_i_ini(jl,i_hemis) .LT. 0.0d0 ) THEN 
310                  ! this write is useful
311                  WRITE(numout,*) ' * TEST 4 POSITIVITY NOT OK FOR CAT ', jl, ' WITH A = ', za_i_ini(jl,i_hemis)
312                  ztest_4 = 0
313               ENDIF
314            END DO
315!           IF ( ztest_4 .EQ. 1 ) WRITE(numout,*) ' * TEST 4 POSITIVITY OK '
316
317!           WRITE(numout,*), ' ztest_1 : ', ztest_1
318!           WRITE(numout,*), ' ztest_2 : ', ztest_2
319!           WRITE(numout,*), ' ztest_3 : ', ztest_3
320!           WRITE(numout,*), ' ztest_4 : ', ztest_4
321!           WRITE(numout,*), ' ztests : ', ztest_1 + ztest_2 + ztest_3 + ztest_4
322
323         ENDIF ! ztest_1 + ztest_2 + ztest_3 + ztest_4
324 
325         ztests = ztest_1 + ztest_2 + ztest_3 + ztest_4
326
327         WRITE(numout,*), ' ztests : ', ztests
328
329         IF ( ztests .NE. 4 ) THEN
330            WRITE(numout,*)
331            WRITE(numout,*), ' !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! '
332            WRITE(numout,*), ' !!!! RED ALERT                  !!! '
333            WRITE(numout,*), ' !!!! BIIIIP BIIIP BIIIIP BIIIIP !!!'
334            WRITE(numout,*), ' !!!! Something is wrong in the LIM3 initialization procedure '
335            WRITE(numout,*), ' !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! '
336            WRITE(numout,*)
337            WRITE(numout,*), ' *** ztests is not equal to 4 '
338            WRITE(numout,*), ' *** ztest_i (i=1,4) = ', ztest_1, ztest_2, ztest_3, ztest_4
339            WRITE(numout,*), ' zat_i_ini : ', zat_i_ini(i_hemis)
340            WRITE(numout,*), ' zhm_i_ini : ', zhm_i_ini(i_hemis)
341         ENDIF ! ztests .NE. 4
342
343      END DO ! i_fill
344
345      END DO ! i_hemis
346
347      !---------------------------------------------------------------------
348      ! 3.3) Space-dependent arrays for ice state variables
349      !---------------------------------------------------------------------
350
351      ! Ice concentration, thickness and volume, snow depth, ice
352      ! salinity, ice age, surface temperature
353      DO jl = 1, jpl ! loop over categories
354         DO jj = 1, jpj
355            DO ji = 1, jpi
356               a_i(ji,jj,jl)   = zidto(ji,jj) * za_i_ini (jl,zhemis(ji,jj))  ! concentration
357               ht_i(ji,jj,jl)  = zidto(ji,jj) * zht_i_ini(jl,zhemis(ji,jj))  ! ice thickness
358               ht_s(ji,jj,jl)  = zidto(ji,jj) * zhm_s_ini(zhemis(ji,jj))     ! snow depth
359               sm_i(ji,jj,jl)  = zidto(ji,jj) * zsm_i_ini(zhemis(ji,jj))     ! salinity
360               o_i(ji,jj,jl)   = zidto(ji,jj) * 1.0 + ( 1.0 - zidto(ji,jj) ) ! age
361               t_su(ji,jj,jl)  = zidto(ji,jj) * 270.0 + ( 1.0 - zidto(ji,jj) ) * t_bo(ji,jj) ! surf temp
362 
363               ! ice volume, snow volume, salt content, age content
364               v_i(ji,jj,jl)   = ht_i(ji,jj,jl) * a_i(ji,jj,jl)              ! ice volume
365               v_s(ji,jj,jl)   = ht_s(ji,jj,jl) * a_i(ji,jj,jl)              ! snow volume
366               smv_i(ji,jj,jl) = MIN( sm_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl) ! salt content
367               oa_i(ji,jj,jl)  = o_i(ji,jj,jl) * a_i(ji,jj,jl)               ! age content
368            END DO ! ji
369         END DO ! jj
370      END DO ! jl
371
372      ! Snow temperature and heat content
373      DO jk = 1, nlay_s
374         DO jl = 1, jpl ! loop over categories
375            DO jj = 1, jpj
376               DO ji = 1, jpi
377                   t_s(ji,jj,jk,jl) = zidto(ji,jj) * 270.00 + ( 1.0 - zidto(ji,jj) ) * rtt
378                   ! Snow energy of melting
379                   e_s(ji,jj,jk,jl) = zidto(ji,jj) * rhosn * ( cpic * ( rtt - t_s(ji,jj,jk,jl) ) + lfus )
380                   ! Change dimensions
381                   e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) / unit_fac
382                   ! Multiply by volume, so that heat content in 10^9 Joules
383                   e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * area(ji,jj) * v_s(ji,jj,jl) / nlay_s
384               END DO ! ji
385            END DO ! jj
386         END DO ! jl
387      END DO ! jk
388
389      ! Ice salinity, temperature and heat content
390      DO jk = 1, nlay_i
391         DO jl = 1, jpl ! loop over categories
392            DO jj = 1, jpj
393               DO ji = 1, jpi
394                   t_i(ji,jj,jk,jl) = zidto(ji,jj)*270.00 + ( 1.0 - zidto(ji,jj) ) * rtt 
395                   s_i(ji,jj,jk,jl) = zidto(ji,jj) * zsm_i_ini(zhemis(ji,jj)) + ( 1.0 - zidto(ji,jj) ) * 0.1
396                   ztmelts          = - tmut * s_i(ji,jj,jk,jl) + rtt !Melting temperature in K
397
398                   ! heat content per unit volume
399                   e_i(ji,jj,jk,jl) = zidto(ji,jj) * rhoic * &
400                      (   cpic    * ( ztmelts - t_i(ji,jj,jk,jl) ) &
401                      +   lfus    * ( 1.0 - (ztmelts-rtt) / MIN((t_i(ji,jj,jk,jl)-rtt),-epsi20) ) &
402                      - rcp      * ( ztmelts - rtt ) &
403                      )
404
405                   ! Correct dimensions to avoid big values
406                   e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) / unit_fac 
407
408                   ! Mutliply by ice volume, and divide by number of layers
409                   ! to get heat content in 10^9 J
410                   e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * area(ji,jj) * v_i(ji,jj,jl) / nlay_i
411               END DO ! ji
412            END DO ! jj
413         END DO ! jl
414      END DO ! jk
415
416      !--------------------------------------------------------------------
417      ! 4) Global ice variables for output diagnostics                    |
418      !--------------------------------------------------------------------
419      fsbbq (:,:)     = 0.e0
420      u_ice (:,:)     = 0.e0
421      v_ice (:,:)     = 0.e0
422      stress1_i(:,:)  = 0.0
423      stress2_i(:,:)  = 0.0
424      stress12_i(:,:) = 0.0
425
426# if defined key_coupled
427      albege(:,:)   = 0.8 * tms(:,:)
428# endif
429
430      !--------------------------------------------------------------------
431      ! 5) Moments for advection
432      !--------------------------------------------------------------------
433
434      sxice (:,:,:)  = 0.e0   ;   sxsn (:,:,:)  = 0.e0   ;   sxa  (:,:,:)  = 0.e0
435      syice (:,:,:)  = 0.e0   ;   sysn (:,:,:)  = 0.e0   ;   sya  (:,:,:)  = 0.e0
436      sxxice(:,:,:)  = 0.e0   ;   sxxsn(:,:,:)  = 0.e0   ;   sxxa (:,:,:)  = 0.e0
437      syyice(:,:,:)  = 0.e0   ;   syysn(:,:,:)  = 0.e0   ;   syya (:,:,:)  = 0.e0
438      sxyice(:,:,:)  = 0.e0   ;   sxysn(:,:,:)  = 0.e0   ;   sxya (:,:,:)  = 0.e0
439
440      sxc0  (:,:,:)  = 0.e0   ;   sxe  (:,:,:,:)= 0.e0   
441      syc0  (:,:,:)  = 0.e0   ;   sye  (:,:,:,:)= 0.e0   
442      sxxc0 (:,:,:)  = 0.e0   ;   sxxe (:,:,:,:)= 0.e0   
443      syyc0 (:,:,:)  = 0.e0   ;   syye (:,:,:,:)= 0.e0   
444      sxyc0 (:,:,:)  = 0.e0   ;   sxye (:,:,:,:)= 0.e0   
445
446      sxsal  (:,:,:)  = 0.e0
447      sysal  (:,:,:)  = 0.e0
448      sxxsal (:,:,:)  = 0.e0
449      syysal (:,:,:)  = 0.e0
450      sxysal (:,:,:)  = 0.e0
451
452      !--------------------------------------------------------------------
453      ! 6) Lateral boundary conditions                                    |
454      !--------------------------------------------------------------------
455
456      DO jl = 1, jpl
457
458         CALL lbc_lnk( a_i(:,:,jl)  , 'T', 1. )
459         CALL lbc_lnk( v_i(:,:,jl)  , 'T', 1. )
460         CALL lbc_lnk( v_s(:,:,jl)  , 'T', 1. )
461         CALL lbc_lnk( smv_i(:,:,jl), 'T', 1. )
462         CALL lbc_lnk( oa_i(:,:,jl) , 'T', 1. )
463
464         CALL lbc_lnk( ht_i(:,:,jl) , 'T', 1. )
465         CALL lbc_lnk( ht_s(:,:,jl) , 'T', 1. )
466         CALL lbc_lnk( sm_i(:,:,jl) , 'T', 1. )
467         CALL lbc_lnk( o_i(:,:,jl)  , 'T', 1. )
468         CALL lbc_lnk( t_su(:,:,jl) , 'T', 1. )
469         DO jk = 1, nlay_s
470            CALL lbc_lnk(t_s(:,:,jk,jl), 'T', 1. )
471            CALL lbc_lnk(e_s(:,:,jk,jl), 'T', 1. )
472         END DO
473         DO jk = 1, nlay_i
474            CALL lbc_lnk(t_i(:,:,jk,jl), 'T', 1. )
475            CALL lbc_lnk(e_i(:,:,jk,jl), 'T', 1. )
476         END DO
477
478         a_i (:,:,jl) = tms(:,:) * a_i(:,:,jl)
479
480      END DO
481
482      CALL lbc_lnk( at_i , 'T', 1. )
483      at_i(:,:) = tms(:,:) * at_i(:,:)                       ! put 0 over land
484
485      CALL lbc_lnk( fsbbq  , 'T', 1. )
486
487   END SUBROUTINE lim_istate
488
489   SUBROUTINE lim_istate_init
490      !!-------------------------------------------------------------------
491      !!                   ***  ROUTINE lim_istate_init  ***
492      !!       
493      !! ** Purpose : Definition of initial state of the ice
494      !!
495      !! ** Method : Read the namiceini namelist and check the parameter
496      !!       values called at the first timestep (nit000)
497      !!
498      !! ** input :
499      !!        Namelist namiceini
500      !!
501      !! history :
502      !!  8.5  ! 03-08 (C. Ethe) original code
503      !!  8.5  ! 07-11 (M. Vancoppenolle) rewritten initialization
504      !!-----------------------------------------------------------------------------
505      NAMELIST/namiceini/ ttest, hninn, hnins, hginn, hgins, aginn, agins, sinn, sins
506      !!-----------------------------------------------------------------------------
507
508      ! Define the initial parameters
509      ! -------------------------
510
511      ! Read Namelist namiceini
512      REWIND ( numnam_ice )
513      READ   ( numnam_ice , namiceini )
514      IF(lwp) THEN
515         WRITE(numout,*)
516         WRITE(numout,*) 'lim_istate_init : ice parameters inititialisation '
517         WRITE(numout,*) '~~~~~~~~~~~~~~~'
518         WRITE(numout,*) '   threshold water temp. for initial sea-ice    ttest      = ', ttest
519         WRITE(numout,*) '   initial snow thickness in the north          hninn      = ', hninn
520         WRITE(numout,*) '   initial snow thickness in the south          hnins      = ', hnins 
521         WRITE(numout,*) '   initial ice thickness  in the north          hginn      = ', hginn
522         WRITE(numout,*) '   initial ice thickness  in the south          hgins      = ', hgins
523         WRITE(numout,*) '   initial ice concentr.  in the north          aginn      = ', aginn
524         WRITE(numout,*) '   initial ice concentr.  in the north          agins      = ', agins
525         WRITE(numout,*) '   initial  ice salinity  in the north          sinn       = ', sinn
526         WRITE(numout,*) '   initial  ice salinity  in the south          sins       = ', sins
527      ENDIF
528
529   END SUBROUTINE lim_istate_init
530
531#else
532   !!----------------------------------------------------------------------
533   !!   Default option :         Empty module          NO LIM sea-ice model
534   !!----------------------------------------------------------------------
535CONTAINS
536   SUBROUTINE lim_istate          ! Empty routine
537   END SUBROUTINE lim_istate
538#endif
539
540   !!======================================================================
541END MODULE limistate