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 in branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/LIM_SRC_3 – NEMO

source: branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/LIM_SRC_3/limistate.F90 @ 4317

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