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 trunk/NEMOGCM/NEMO/LIM_SRC_3 – NEMO

source: trunk/NEMOGCM/NEMO/LIM_SRC_3/limistate.F90 @ 3319

Last change on this file since 3319 was 3294, checked in by rblod, 12 years ago

Merge of 3.4beta into the trunk

  • Property svn:keywords set to Id
File size: 25.2 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   !!----------------------------------------------------------------------
9#if defined key_lim3
10   !!----------------------------------------------------------------------
11   !!   'key_lim3' :                                    LIM3 sea-ice model
12   !!----------------------------------------------------------------------
13   !!   lim_istate      :  Initialisation of diagnostics ice variables
14   !!   lim_istate_init :  initialization of ice state and namelist read
15   !!----------------------------------------------------------------------
16   USE phycst           ! physical constant
17   USE oce              ! dynamics and tracers variables
18   USE dom_oce          ! ocean domain
19   USE sbc_oce          ! Surface boundary condition: ocean fields
20   USE eosbn2           ! equation of state
21   USE ice              ! sea-ice variables
22   USE par_ice          ! ice parameters
23   USE dom_ice          ! sea-ice domain
24   USE in_out_manager   ! I/O manager
25   USE lbclnk           ! lateral boundary condition - MPP exchanges
26   USE lib_mpp          ! MPP library
27   USE wrk_nemo         ! work arrays
28
29   IMPLICIT NONE
30   PRIVATE
31
32   PUBLIC   lim_istate      ! routine called by lim_init.F90
33
34   !                                  !!** init namelist (namiceini) **
35   REAL(wp) ::   ttest    = 2.0_wp     ! threshold water temperature for initial sea ice
36   REAL(wp) ::   hninn    = 0.5_wp     ! initial snow thickness in the north
37   REAL(wp) ::   hginn_u  = 2.5_wp     ! initial ice thickness in the north
38   REAL(wp) ::   aginn_u  = 0.7_wp     ! initial leads area in the north
39   REAL(wp) ::   hginn_d  = 5.0_wp     ! initial ice thickness in the north
40   REAL(wp) ::   aginn_d  = 0.25_wp    ! initial leads area in the north
41   REAL(wp) ::   hnins    = 0.1_wp     ! initial snow thickness in the south
42   REAL(wp) ::   hgins_u  = 1.0_wp     ! initial ice thickness in the south
43   REAL(wp) ::   agins_u  = 0.7_wp     ! initial leads area in the south
44   REAL(wp) ::   hgins_d  = 2.0_wp     ! initial ice thickness in the south
45   REAL(wp) ::   agins_d  = 0.2_wp     ! initial leads area in the south
46   REAL(wp) ::   sinn     = 6.301_wp   ! initial salinity
47   REAL(wp) ::   sins     = 6.301_wp   !
48
49   !!----------------------------------------------------------------------
50   !! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2011)
51   !! $Id$
52   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
53   !!----------------------------------------------------------------------
54CONTAINS
55
56   SUBROUTINE lim_istate
57      !!-------------------------------------------------------------------
58      !!                    ***  ROUTINE lim_istate  ***
59      !!
60      !! ** Purpose :   defined the sea-ice initial state
61      !!
62      !! ** Method  :   restart from a state defined in a binary file
63      !!                or from arbitrary sea-ice conditions
64      !!-------------------------------------------------------------------
65      INTEGER  ::   ji, jj, jk, jl             ! dummy loop indices
66      REAL(wp) ::   zeps6, zeps, ztmelts, epsi06   ! local scalars
67      REAL(wp) ::   zvol, zare, zh, zh1, zh2, zh3, zan, zbn, zas, zbs 
68      REAL(wp), POINTER, DIMENSION(:)   ::   zgfactorn, zhin 
69      REAL(wp), POINTER, DIMENSION(:)   ::   zgfactors, zhis
70      REAL(wp), POINTER, DIMENSION(:,:) ::   zidto      ! ice indicator
71      !--------------------------------------------------------------------
72
73      CALL wrk_alloc( jpm, zgfactorn, zgfactors, zhin, zhis )
74      CALL wrk_alloc( jpi, jpj, zidto )
75
76      !--------------------------------------------------------------------
77      ! 1) Preliminary things
78      !--------------------------------------------------------------------
79      epsi06 = 1.e-6_wp
80
81      CALL lim_istate_init     !  reading the initials parameters of the ice
82
83!!gm  in lim2  the initialisation if only done if required in the namelist :
84!!gm      IF( .NOT. ln_limini ) THEN
85!!gm  this should be added in lim3 namelist...
86
87      !--------------------------------------------------------------------
88      ! 2) Ice initialization (hi,hs,frld,t_su,sm_i,t_i,t_s)              |
89      !--------------------------------------------------------------------
90
91      IF(lwp) WRITE(numout,*)
92      IF(lwp) WRITE(numout,*) 'lim_istate : Ice initialization '
93      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ '
94
95      t_bo(:,:) = tfreez( tsn(:,:,1,jp_sal) ) * tmask(:,:,1)       ! freezing/melting point of sea water [Celcius]
96
97      DO jj = 1, jpj                                       ! ice if sst <= t-freez + ttest
98         DO ji = 1, jpi
99            IF( tsn(ji,jj,1,jp_tem)  - t_bo(ji,jj) >= ttest ) THEN   ;   zidto(ji,jj) = 0.e0      ! no ice
100            ELSE                                                     ;   zidto(ji,jj) = 1.e0      !    ice
101            ENDIF
102         END DO
103      END DO
104
105      t_bo(:,:) = t_bo(:,:) + rt0                          ! t_bo converted from Celsius to Kelvin (rt0 over land)
106
107      ! constants for heat contents
108      zeps   = 1.e-20_wp
109      zeps6  = 1.e-06_wp
110
111      ! zgfactor for initial ice distribution
112      zgfactorn(:) = 0._wp
113      zgfactors(:) = 0._wp
114
115      ! first ice type
116      DO jl = ice_cat_bounds(1,1), ice_cat_bounds(1,2)
117         zhin (1)     = ( hi_max(jl-1) + hi_max(jl) ) * 0.5_wp
118         zgfactorn(1) = zgfactorn(1) + exp(-(zhin(1)-hginn_u)*(zhin(1)-hginn_u) * 0.5_wp )
119         zhis (1)     = ( hi_max(jl-1) + hi_max(jl) ) * 0.5_wp
120         zgfactors(1) = zgfactors(1) + exp(-(zhis(1)-hgins_u)*(zhis(1)-hgins_u) * 0.5_wp )
121      END DO ! jl
122      zgfactorn(1) = aginn_u / zgfactorn(1)
123      zgfactors(1) = agins_u / zgfactors(1)
124
125      ! -------------
126      ! new distribution, polynom of second order, conserving area and volume
127      zh1 = 0._wp
128      zh2 = 0._wp
129      zh3 = 0._wp
130      DO jl = 1, jpl
131         zh = ( hi_max(jl-1) + hi_max(jl) ) * 0.5_wp
132         zh1 = zh1 + zh
133         zh2 = zh2 + zh * zh
134         zh3 = zh3 + zh * zh * zh
135      END DO
136      IF(lwp) WRITE(numout,*) ' zh1 : ', zh1
137      IF(lwp) WRITE(numout,*) ' zh2 : ', zh2
138      IF(lwp) WRITE(numout,*) ' zh3 : ', zh3
139
140      zvol = aginn_u * hginn_u
141      zare = aginn_u
142      IF( jpl >= 2 ) THEN
143         zbn = ( zvol*zh2 - zare*zh3 ) / ( zh2*zh2 - zh1*zh3)
144         zan = ( zare - zbn*zh1 ) / zh2
145      ENDIF
146
147      IF(lwp) WRITE(numout,*) ' zvol: ', zvol
148      IF(lwp) WRITE(numout,*) ' zare: ', zare
149      IF(lwp) WRITE(numout,*) ' zbn : ', zbn 
150      IF(lwp) WRITE(numout,*) ' zan : ', zan 
151
152      zvol = agins_u * hgins_u
153      zare = agins_u
154      IF( jpl >= 2 ) THEN
155         zbs = ( zvol*zh2 - zare*zh3 ) / ( zh2*zh2 - zh1*zh3)
156         zas = ( zare - zbs*zh1 ) / zh2
157      ENDIF
158
159      IF(lwp) WRITE(numout,*) ' zvol: ', zvol
160      IF(lwp) WRITE(numout,*) ' zare: ', zare
161      IF(lwp) WRITE(numout,*) ' zbn : ', zbn 
162      IF(lwp) WRITE(numout,*) ' zan : ', zan 
163
164      !end of new lines
165      ! -------------
166!!!
167      ! retour a LIMA_MEC
168      !     ! second ice type
169      !     zdummy  = hi_max(ice_cat_bounds(2,1)-1)
170      !     hi_max(ice_cat_bounds(2,1)-1) = 0.0
171
172      !     ! here to change !!!!
173      !     jm = 2
174      !     DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2)
175      !        zhin (2)     = ( hi_max(jl-1) + hi_max(jl) ) / 2.0
176      !        zhin (2)     = ( hi_max_typ(jl-ice_cat_bounds(2,1),jm    ) + &
177      !                         hi_max_typ(jl-ice_cat_bounds(2,1) + 1,jm)   ) / 2.0
178      !        zgfactorn(2) = zgfactorn(2) + exp(-(zhin(2)-hginn_d)*(zhin(2)-hginn_d)/2.0)
179      !        zhis (2)     = ( hi_max(jl-1) + hi_max(jl) ) / 2.0
180      !        zhis (2)     = ( hi_max_typ(jl-ice_cat_bounds(2,1),jm    ) + &
181      !                         hi_max_typ(jl-ice_cat_bounds(2,1) + 1,jm)   ) / 2.0
182      !        zgfactors(2) = zgfactors(2) + exp(-(zhis(2)-hgins_d)*(zhis(2)-hgins_d)/2.0)
183      !     END DO ! jl
184      !     zgfactorn(2) = aginn_d / zgfactorn(2)
185      !     zgfactors(2) = agins_d / zgfactors(2)
186      !     hi_max(ice_cat_bounds(2,1)-1) = zdummy
187      ! END retour a LIMA_MEC
188!!!
189
190!!gm  optimisation :  loop over the ice categories inside the ji, jj loop !!!
191
192      DO jj = 1, jpj
193         DO ji = 1, jpi
194
195            !--- Northern hemisphere
196            !----------------------------------------------------------------
197            IF( fcor(ji,jj) >= 0._wp ) THEN   
198
199               !-----------------------
200               ! Ice area / thickness
201               !-----------------------
202
203               IF ( jpl .EQ. 1) THEN ! one category
204
205                  DO jl = ice_cat_bounds(1,1), ice_cat_bounds(1,2) ! loop over ice thickness categories
206                     a_i(ji,jj,jl)    = zidto(ji,jj) * aginn_u
207                     ht_i(ji,jj,jl)   = zidto(ji,jj) * hginn_u
208                     v_i(ji,jj,jl)    = ht_i(ji,jj,jl)*a_i(ji,jj,jl)
209                  END DO
210
211               ELSE ! several categories
212
213                  DO jl = ice_cat_bounds(1,1), ice_cat_bounds(1,2) ! loop over ice thickness categories
214                     zhin(1)          = ( hi_max(jl-1) + hi_max(jl) ) / 2.0
215                     a_i(ji,jj,jl)    = zidto(ji,jj) * MAX( zgfactorn(1) * exp(-(zhin(1)-hginn_u)* & 
216                        (zhin(1)-hginn_u)/2.0) , epsi06)
217                     ! new line
218                     a_i(ji,jj,jl)    = zidto(ji,jj) * ( zan * zhin(1) * zhin(1) + zbn * zhin(1) )
219                     ht_i(ji,jj,jl)   = zidto(ji,jj) * zhin(1) 
220                     v_i(ji,jj,jl)    = ht_i(ji,jj,jl)*a_i(ji,jj,jl)
221                  END DO
222
223               ENDIF
224
225
226!!!
227               ! retour a LIMA_MEC
228               !              !ridged ice
229               !              zdummy  = hi_max(ice_cat_bounds(2,1)-1)
230               !              hi_max(ice_cat_bounds(2,1)-1) = 0.0
231               !              DO jl = ice_cat_bounds(2,1), ice_cat_bounds(2,2) ! loop over ice thickness categories
232               !                 zhin(2)          = ( hi_max(jl-1) + hi_max(jl) ) / 2.0
233               !                 a_i(ji,jj,jl)    = zidto(ji,jj) * MAX( zgfactorn(2) * exp(-(zhin(2)-hginn_d)* &
234               !                                    (zhin(2)-hginn_d)/2.0) , epsi06)
235               !                 ht_i(ji,jj,jl)   = zidto(ji,jj) * zhin(2)
236               !                 v_i(ji,jj,jl)    = ht_i(ji,jj,jl)*a_i(ji,jj,jl)
237               !              END DO
238               !              hi_max(ice_cat_bounds(2,1)-1) = zdummy
239
240               !              !rafted ice
241               !              jl = 6
242               !              a_i(ji,jj,jl)       = 0.0
243               !              ht_i(ji,jj,jl)      = 0.0
244               !              v_i(ji,jj,jl)       = 0.0
245               ! END retour a LIMA_MEC
246!!!
247
248               DO jl = 1, jpl
249
250                  !-------------
251                  ! Snow depth
252                  !-------------
253                  ht_s(ji,jj,jl)   = zidto(ji,jj) * hninn
254                  v_s(ji,jj,jl)    = ht_s(ji,jj,jl)*a_i(ji,jj,jl)
255
256                  !---------------
257                  ! Ice salinity
258                  !---------------
259                  sm_i(ji,jj,jl)   = zidto(ji,jj) * sinn  + ( 1.0 - zidto(ji,jj) ) * 0.1
260                  smv_i(ji,jj,jl)  = MIN( sm_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl)
261
262                  !----------
263                  ! Ice age
264                  !----------
265                  o_i(ji,jj,jl)    = zidto(ji,jj) * 1.0   + ( 1.0 - zidto(ji,jj) )
266                  oa_i(ji,jj,jl)   = o_i(ji,jj,jl) * a_i(ji,jj,jl)
267
268                  !------------------------------
269                  ! Sea ice surface temperature
270                  !------------------------------
271
272                  t_su(ji,jj,jl)   = zidto(ji,jj) * 270.0 + ( 1.0 - zidto(ji,jj) ) * t_bo(ji,jj)
273
274                  !------------------------------------
275                  ! Snow temperature and heat content
276                  !------------------------------------
277
278                  DO jk = 1, nlay_s
279                     t_s(ji,jj,jk,jl) = zidto(ji,jj) * 270.00 + ( 1.0 - zidto(ji,jj) ) * rtt
280                     ! Snow energy of melting
281                     e_s(ji,jj,jk,jl) = zidto(ji,jj) * rhosn * ( cpic * ( rtt - t_s(ji,jj,jk,jl) ) + lfus )
282                     ! Change dimensions
283                     e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) / unit_fac
284                     ! Multiply by volume, so that heat content in 10^9 Joules
285                     e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * area(ji,jj) * &
286                        v_s(ji,jj,jl)  / nlay_s
287                  END DO !jk
288
289                  !-----------------------------------------------
290                  ! Ice salinities, temperature and heat content
291                  !-----------------------------------------------
292
293                  DO jk = 1, nlay_i
294                     t_i(ji,jj,jk,jl) = zidto(ji,jj)*270.00 + ( 1.0 - zidto(ji,jj) ) * rtt 
295                     s_i(ji,jj,jk,jl) = zidto(ji,jj) * sinn + ( 1.0 - zidto(ji,jj) ) * 0.1
296                     ztmelts          = - tmut * s_i(ji,jj,jk,jl) + rtt !Melting temperature in K
297
298                     ! heat content per unit volume
299                     e_i(ji,jj,jk,jl) = zidto(ji,jj) * rhoic * &
300                        (   cpic    * ( ztmelts - t_i(ji,jj,jk,jl) ) &
301                        +   lfus    * ( 1.0 - (ztmelts-rtt) / MIN((t_i(ji,jj,jk,jl)-rtt),-zeps) ) &
302                        - rcp      * ( ztmelts - rtt ) &
303                        )
304
305                     ! Correct dimensions to avoid big values
306                     e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) / unit_fac 
307
308                     ! Mutliply by ice volume, and divide by number of layers to get heat content in 10^9 J
309                     e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * & 
310                        area(ji,jj) * a_i(ji,jj,jl) * ht_i(ji,jj,jl) / &
311                        nlay_i
312                  END DO ! jk
313
314               END DO ! jl
315
316            ELSE ! on fcor
317
318               !--- Southern hemisphere
319               !----------------------------------------------------------------
320
321               !-----------------------
322               ! Ice area / thickness
323               !-----------------------
324
325               IF ( jpl .EQ. 1) THEN ! one category
326
327                  DO jl = ice_cat_bounds(1,1), ice_cat_bounds(1,2) ! loop over ice thickness categories
328                     a_i(ji,jj,jl)    = zidto(ji,jj) * agins_u
329                     ht_i(ji,jj,jl)   = zidto(ji,jj) * hgins_u
330                     v_i(ji,jj,jl)    = ht_i(ji,jj,jl)*a_i(ji,jj,jl)
331                  END DO
332
333               ELSE ! several categories
334
335                  !level ice
336                  DO jl = ice_cat_bounds(1,1), ice_cat_bounds(1,2) !over thickness categories
337
338                     zhis(1)       = ( hi_max(jl-1) + hi_max(jl) ) / 2.0
339                     a_i(ji,jj,jl) = zidto(ji,jj) * MAX( zgfactors(1) * exp(-(zhis(1)-hgins_u) * & 
340                        (zhis(1)-hgins_u)/2.0) , epsi06 )
341                     ! new line square distribution volume conserving
342                     a_i(ji,jj,jl)    = zidto(ji,jj) * ( zas * zhis(1) * zhis(1) + zbs * zhis(1) )
343                     ht_i(ji,jj,jl)   = zidto(ji,jj) * zhis(1) 
344                     v_i(ji,jj,jl)    = ht_i(ji,jj,jl)*a_i(ji,jj,jl)
345
346                  END DO ! jl
347
348               ENDIF
349
350!!!
351               ! retour a LIMA_MEC
352               !              !ridged ice
353               !              zdummy  = hi_max(ice_cat_bounds(2,1)-1)
354               !              hi_max(ice_cat_bounds(2,1)-1) = 0.0
355               !              DO jl = ice_cat_bounds(2,1), ice_cat_bounds(2,2) !over thickness categories
356               !                 zhis(2)       = ( hi_max(jl-1) + hi_max(jl) ) / 2.0
357               !                 a_i(ji,jj,jl) = zidto(ji,jj)*MAX( zgfactors(2)   &
358               !                    &          * exp(-(zhis(2)-hgins_d)*(zhis(2)-hgins_d)/2.0), epsi06 )
359               !                 ht_i(ji,jj,jl)   = zidto(ji,jj) * zhis(2)
360               !                 v_i(ji,jj,jl)    = ht_i(ji,jj,jl)*a_i(ji,jj,jl)
361               !              END DO
362               !              hi_max(ice_cat_bounds(2,1)-1) = zdummy
363
364               !              !rafted ice
365               !              jl = 6
366               !              a_i(ji,jj,jl)       = 0.0
367               !              ht_i(ji,jj,jl)      = 0.0
368               !              v_i(ji,jj,jl)       = 0.0
369               ! END retour a LIMA_MEC
370!!!
371
372               DO jl = 1, jpl !over thickness categories
373
374                  !---------------
375                  ! Snow depth
376                  !---------------
377
378                  ht_s(ji,jj,jl)   = zidto(ji,jj) * hnins
379                  v_s(ji,jj,jl)    = ht_s(ji,jj,jl)*a_i(ji,jj,jl)
380
381                  !---------------
382                  ! Ice salinity
383                  !---------------
384
385                  sm_i(ji,jj,jl)   = zidto(ji,jj) * sins  + ( 1.0 - zidto(ji,jj) ) * 0.1
386                  smv_i(ji,jj,jl)  = MIN( sm_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl)
387
388                  !----------
389                  ! Ice age
390                  !----------
391
392                  o_i(ji,jj,jl)    = zidto(ji,jj) * 1.0   + ( 1.0 - zidto(ji,jj) )
393                  oa_i(ji,jj,jl)   = o_i(ji,jj,jl) * a_i(ji,jj,jl)
394
395                  !------------------------------
396                  ! Sea ice surface temperature
397                  !------------------------------
398
399                  t_su(ji,jj,jl)   = zidto(ji,jj) * 270.0 + ( 1.0 - zidto(ji,jj) ) * t_bo(ji,jj)
400
401                  !----------------------------------
402                  ! Snow temperature / heat content
403                  !----------------------------------
404
405                  DO jk = 1, nlay_s
406                     t_s(ji,jj,jk,jl) = zidto(ji,jj) * 270.00 + ( 1.0 - zidto(ji,jj) ) * rtt
407                     ! Snow energy of melting
408                     e_s(ji,jj,jk,jl) = zidto(ji,jj) * rhosn * ( cpic * ( rtt - t_s(ji,jj,jk,jl) ) + lfus )
409                     ! Change dimensions
410                     e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) / unit_fac
411                     ! Multiply by volume, so that heat content in 10^9 Joules
412                     e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * area(ji,jj) * &
413                        v_s(ji,jj,jl)  / nlay_s
414                  END DO
415
416                  !---------------------------------------------
417                  ! Ice temperature, salinity and heat content
418                  !---------------------------------------------
419
420                  DO jk = 1, nlay_i
421                     t_i(ji,jj,jk,jl) = zidto(ji,jj)*270.00 + ( 1.0 - zidto(ji,jj) ) * rtt 
422                     s_i(ji,jj,jk,jl) = zidto(ji,jj) * sins + ( 1.0 - zidto(ji,jj) ) * 0.1
423                     ztmelts          = - tmut * s_i(ji,jj,jk,jl) + rtt !Melting temperature in K
424
425                     ! heat content per unit volume
426                     e_i(ji,jj,jk,jl) = zidto(ji,jj) * rhoic * &
427                        (   cpic    * ( ztmelts - t_i(ji,jj,jk,jl) ) &
428                        +   lfus  * ( 1.0 - (ztmelts-rtt) / MIN((t_i(ji,jj,jk,jl)-rtt),-zeps) ) &
429                        - rcp      * ( ztmelts - rtt ) &
430                        )
431
432                     ! Correct dimensions to avoid big values
433                     e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) / unit_fac 
434
435                     ! Mutliply by ice volume, and divide by number of layers to get heat content in 10^9 J
436                     e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * & 
437                        area(ji,jj) * a_i(ji,jj,jl) * ht_i(ji,jj,jl) / &
438                        nlay_i
439                  END DO !jk
440
441               END DO ! jl
442
443            ENDIF ! on fcor
444
445         END DO
446      END DO
447
448      !--------------------------------------------------------------------
449      ! 3) Global ice variables for output diagnostics                    |
450      !--------------------------------------------------------------------
451
452      fsbbq (:,:)     = 0.e0
453      u_ice (:,:)     = 0.e0
454      v_ice (:,:)     = 0.e0
455      stress1_i(:,:)  = 0.0
456      stress2_i(:,:)  = 0.0
457      stress12_i(:,:) = 0.0
458
459      !--------------------------------------------------------------------
460      ! 4) Moments for advection
461      !--------------------------------------------------------------------
462
463      sxice (:,:,:)  = 0.e0   ;   sxsn (:,:,:)  = 0.e0   ;   sxa  (:,:,:)  = 0.e0
464      syice (:,:,:)  = 0.e0   ;   sysn (:,:,:)  = 0.e0   ;   sya  (:,:,:)  = 0.e0
465      sxxice(:,:,:)  = 0.e0   ;   sxxsn(:,:,:)  = 0.e0   ;   sxxa (:,:,:)  = 0.e0
466      syyice(:,:,:)  = 0.e0   ;   syysn(:,:,:)  = 0.e0   ;   syya (:,:,:)  = 0.e0
467      sxyice(:,:,:)  = 0.e0   ;   sxysn(:,:,:)  = 0.e0   ;   sxya (:,:,:)  = 0.e0
468
469      sxc0  (:,:,:)  = 0.e0   ;   sxe  (:,:,:,:)= 0.e0   
470      syc0  (:,:,:)  = 0.e0   ;   sye  (:,:,:,:)= 0.e0   
471      sxxc0 (:,:,:)  = 0.e0   ;   sxxe (:,:,:,:)= 0.e0   
472      syyc0 (:,:,:)  = 0.e0   ;   syye (:,:,:,:)= 0.e0   
473      sxyc0 (:,:,:)  = 0.e0   ;   sxye (:,:,:,:)= 0.e0   
474
475      sxsal  (:,:,:)  = 0.e0
476      sysal  (:,:,:)  = 0.e0
477      sxxsal (:,:,:)  = 0.e0
478      syysal (:,:,:)  = 0.e0
479      sxysal (:,:,:)  = 0.e0
480
481      !--------------------------------------------------------------------
482      ! 5) Lateral boundary conditions                                    |
483      !--------------------------------------------------------------------
484
485      DO jl = 1, jpl
486         CALL lbc_lnk( a_i(:,:,jl)  , 'T', 1. )
487         CALL lbc_lnk( v_i(:,:,jl)  , 'T', 1. )
488         CALL lbc_lnk( v_s(:,:,jl)  , 'T', 1. )
489         CALL lbc_lnk( smv_i(:,:,jl), 'T', 1. )
490         CALL lbc_lnk( oa_i(:,:,jl) , 'T', 1. )
491         !
492         CALL lbc_lnk( ht_i(:,:,jl) , 'T', 1. )
493         CALL lbc_lnk( ht_s(:,:,jl) , 'T', 1. )
494         CALL lbc_lnk( sm_i(:,:,jl) , 'T', 1. )
495         CALL lbc_lnk( o_i(:,:,jl)  , 'T', 1. )
496         CALL lbc_lnk( t_su(:,:,jl) , 'T', 1. )
497         DO jk = 1, nlay_s
498            CALL lbc_lnk(t_s(:,:,jk,jl), 'T', 1. )
499            CALL lbc_lnk(e_s(:,:,jk,jl), 'T', 1. )
500         END DO
501         DO jk = 1, nlay_i
502            CALL lbc_lnk(t_i(:,:,jk,jl), 'T', 1. )
503            CALL lbc_lnk(e_i(:,:,jk,jl), 'T', 1. )
504         END DO
505         !
506         a_i(:,:,jl) = tms(:,:) * a_i(:,:,jl)
507      END DO
508
509      CALL lbc_lnk( at_i , 'T', 1. )
510      at_i(:,:) = tms(:,:) * at_i(:,:)                       ! put 0 over land
511      !
512      CALL lbc_lnk( fsbbq  , 'T', 1. )
513      !
514      CALL wrk_dealloc( jpm, zgfactorn, zgfactors, zhin, zhis )
515      CALL wrk_dealloc( jpi, jpj, zidto )
516      !
517   END SUBROUTINE lim_istate
518
519
520   SUBROUTINE lim_istate_init
521      !!-------------------------------------------------------------------
522      !!                   ***  ROUTINE lim_istate_init  ***
523      !!       
524      !! ** Purpose : Definition of initial state of the ice
525      !!
526      !! ** Method :   Read the namiceini namelist and check the parameter
527      !!             values called at the first timestep (nit000)
528      !!
529      !! ** input  :   namelist namiceini
530      !!-----------------------------------------------------------------------------
531      NAMELIST/namiceini/ ttest, hninn, hginn_u, aginn_u, hginn_d, aginn_d, hnins,   &
532         &                hgins_u, agins_u, hgins_d, agins_d, sinn, sins
533      !!-----------------------------------------------------------------------------
534      !
535      REWIND ( numnam_ice )               ! Read Namelist namiceini
536      READ   ( numnam_ice , namiceini )
537      !
538      IF(lwp) THEN                        ! control print
539         WRITE(numout,*)
540         WRITE(numout,*) 'lim_istate_init : ice parameters inititialisation '
541         WRITE(numout,*) '~~~~~~~~~~~~~~~'
542         WRITE(numout,*) '   threshold water temp. for initial sea-ice    ttest      = ', ttest
543         WRITE(numout,*) '   initial snow thickness in the north          hninn      = ', hninn
544         WRITE(numout,*) '   initial undef ice thickness in the north     hginn_u    = ', hginn_u
545         WRITE(numout,*) '   initial undef ice concentr. in the north     aginn_u    = ', aginn_u         
546         WRITE(numout,*) '   initial  def  ice thickness in the north     hginn_d    = ', hginn_d
547         WRITE(numout,*) '   initial  def  ice concentr. in the north     aginn_d    = ', aginn_d         
548         WRITE(numout,*) '   initial snow thickness in the south          hnins      = ', hnins 
549         WRITE(numout,*) '   initial undef ice thickness in the north     hgins_u    = ', hgins_u
550         WRITE(numout,*) '   initial undef ice concentr. in the north     agins_u    = ', agins_u         
551         WRITE(numout,*) '   initial  def  ice thickness in the north     hgins_d    = ', hgins_d
552         WRITE(numout,*) '   initial  def  ice concentr. in the north     agins_d    = ', agins_d         
553         WRITE(numout,*) '   initial  ice salinity       in the north     sinn       = ', sinn
554         WRITE(numout,*) '   initial  ice salinity       in the south     sins       = ', sins
555      ENDIF
556      !
557   END SUBROUTINE lim_istate_init
558
559#else
560   !!----------------------------------------------------------------------
561   !!   Default option :         Empty module          NO LIM sea-ice model
562   !!----------------------------------------------------------------------
563CONTAINS
564   SUBROUTINE lim_istate          ! Empty routine
565   END SUBROUTINE lim_istate
566#endif
567
568   !!======================================================================
569END MODULE limistate
Note: See TracBrowser for help on using the repository browser.