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

source: branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limistate.F90 @ 3625

Last change on this file since 3625 was 3625, checked in by acc, 11 years ago

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