Context Navigation

limistate.F90 @ 6140

Last change on this file since 6140 was 6140, checked in by timgraham, 8 years ago

Merge of branches/2015/dev_merge_2015 back into trunk. Merge excludes NEMOGCM/TOOLS/OBSTOOLS/ for now due to issues with the change of file type. Will sort these manually with further commits.

Branch merged as follows:
In the working copy of branch ran:
svn merge svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/trunk@HEAD
Small conflicts due to bug fixes applied to trunk since the dev_merge_2015 was copied. Bug fixes were applied to the branch as well so these were easy to resolve.
Branch committed at this stage

In working copy run:
svn switch svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/trunk
to switch working copy

Run:
svn merge --reintegrate svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/branches/2015/dev_merge_2015
to merge the branch into the trunk and then commit - no conflicts at this stage.

Property svn:keywords set to Id

File size: 28.6 KB

Line
1	MODULE 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	!! - ! 2014 (C. Rousset) add N/S initializations
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_oce ! ocean parameters
25	USE dom_ice ! sea-ice domain
26	USE in_out_manager ! I/O manager
27	USE lib_mpp ! MPP library
28	USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)
29	USE wrk_nemo ! work arrays
30	USE fldread ! read input fields
31	USE iom
32
33	IMPLICIT NONE
34	PRIVATE
35
36	PUBLIC lim_istate ! routine called by lim_init.F90
37
38	! !! init namelist (namiceini)
39	REAL(wp) :: rn_thres_sst ! threshold water temperature for initial sea ice
40	REAL(wp) :: rn_hts_ini_n ! initial snow thickness in the north
41	REAL(wp) :: rn_hts_ini_s ! initial snow thickness in the south
42	REAL(wp) :: rn_hti_ini_n ! initial ice thickness in the north
43	REAL(wp) :: rn_hti_ini_s ! initial ice thickness in the south
44	REAL(wp) :: rn_ati_ini_n ! initial leads area in the north
45	REAL(wp) :: rn_ati_ini_s ! initial leads area in the south
46	REAL(wp) :: rn_smi_ini_n ! initial salinity
47	REAL(wp) :: rn_smi_ini_s ! initial salinity
48	REAL(wp) :: rn_tmi_ini_n ! initial temperature
49	REAL(wp) :: rn_tmi_ini_s ! initial temperature
50
51	INTEGER , PARAMETER :: jpfldi = 6 ! maximum number of files to read
52	INTEGER , PARAMETER :: jp_hti = 1 ! index of ice thickness (m) at T-point
53	INTEGER , PARAMETER :: jp_hts = 2 ! index of snow thicknes (m) at T-point
54	INTEGER , PARAMETER :: jp_ati = 3 ! index of ice fraction (%) at T-point
55	INTEGER , PARAMETER :: jp_tsu = 4 ! index of ice surface temp (K) at T-point
56	INTEGER , PARAMETER :: jp_tmi = 5 ! index of ice temp at T-point
57	INTEGER , PARAMETER :: jp_smi = 6 ! index of ice sali at T-point
58	TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: si ! structure of input fields (file informations, fields read)
59
60	LOGICAL :: ln_iceini ! initialization or not
61	LOGICAL :: ln_iceini_file ! Ice initialization state from 2D netcdf file
62	!!----------------------------------------------------------------------
63	!! LIM 3.0, UCL-LOCEAN-IPSL (2008)
64	!! $Id$
65	!! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
66	!!----------------------------------------------------------------------
67	CONTAINS
68
69	SUBROUTINE lim_istate
70	!!-------------------------------------------------------------------
71	!! * ROUTINE lim_istate *
72	!!
73	!! ** Purpose : defined the sea-ice initial state
74	!!
75	!! ** Method :
76	!! This routine will put some ice where ocean
77	!! is at the freezing point, then fill in ice
78	!! state variables using prescribed initial
79	!! values in the namelist
80	!!
81	!! ** Steps :
82	!! 1) Read namelist
83	!! 2) Basal temperature; ice and hemisphere masks
84	!! 3) Fill in the ice thickness distribution using gaussian
85	!! 4) Fill in space-dependent arrays for state variables
86	!! 5) Diagnostic arrays
87	!! 6) Lateral boundary conditions
88	!!
89	!! ** Notes : o_i, t_su, t_s, t_i, s_i must be filled everywhere, even
90	!! where there is no ice (clem: I do not know why, is it mandatory?)
91	!!
92	!! History :
93	!! 2.0 ! 01-04 (C. Ethe, G. Madec) Original code
94	!! 3.0 ! 2007 (M. Vancoppenolle) Rewrite for ice cats
95	!! 4.0 ! 09-11 (M. Vancoppenolle) Enhanced version for ice cats
96	!!--------------------------------------------------------------------
97
98	!! * Local variables
99	INTEGER :: ji, jj, jk, jl ! dummy loop indices
100	REAL(wp) :: ztmelts, zdh
101	INTEGER :: i_hemis, i_fill, jl0
102	REAL(wp) :: ztest_1, ztest_2, ztest_3, ztest_4, ztests, zsigma, zarg, zA, zV, zA_cons, zV_cons, zconv
103	REAL(wp), POINTER, DIMENSION(:,:) :: zswitch ! ice indicator
104	REAL(wp), POINTER, DIMENSION(:,:) :: zht_i_ini, zat_i_ini, zvt_i_ini !data from namelist or nc file
105	REAL(wp), POINTER, DIMENSION(:,:) :: zts_u_ini, zht_s_ini, zsm_i_ini, ztm_i_ini !data from namelist or nc file
106	REAL(wp), POINTER, DIMENSION(:,:,:) :: zh_i_ini, za_i_ini, zv_i_ini !data by cattegories to fill
107	!--------------------------------------------------------------------
108
109	CALL wrk_alloc( jpi, jpj, jpl, zh_i_ini, za_i_ini, zv_i_ini )
110	CALL wrk_alloc( jpi, jpj, zht_i_ini, zat_i_ini, zvt_i_ini, zts_u_ini, zht_s_ini, zsm_i_ini, ztm_i_ini )
111	CALL wrk_alloc( jpi, jpj, zswitch )
112
113	IF(lwp) WRITE(numout,*)
114	IF(lwp) WRITE(numout,*) 'lim_istate : Ice initialization '
115	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ '
116
117	!--------------------------------------------------------------------
118	! 1) Read namelist
119	!--------------------------------------------------------------------
120
121	CALL lim_istate_init ! reading the initials parameters of the ice
122
123	! surface temperature
124	DO jl = 1, jpl ! loop over categories
125	t_su (:,:,jl) = rt0 * tmask(:,:,1)
126	tn_ice(:,:,jl) = rt0 * tmask(:,:,1)
127	END DO
128
129	! basal temperature (considered at freezing point)
130	CALL eos_fzp( sss_m(:,:), t_bo(:,:) )
131	t_bo(:,:) = ( t_bo(:,:) + rt0 ) * tmask(:,:,1)
132
133
134	IF( ln_iceini ) THEN
135
136	!--------------------------------------------------------------------
137	! 2) Basal temperature, ice mask and hemispheric index
138	!--------------------------------------------------------------------
139
140	DO jj = 1, jpj ! ice if sst <= t-freez + ttest
141	DO ji = 1, jpi
142	IF( ( sst_m(ji,jj) - ( t_bo(ji,jj) - rt0 ) ) * tmask(ji,jj,1) >= rn_thres_sst ) THEN
143	zswitch(ji,jj) = 0._wp * tmask(ji,jj,1) ! no ice
144	ELSE
145	zswitch(ji,jj) = 1._wp * tmask(ji,jj,1) ! ice
146	ENDIF
147	END DO
148	END DO
149
150	!--------------------------------------------------------------------
151	! 3) Initialization of sea ice state variables
152	!--------------------------------------------------------------------
153	IF( ln_iceini_file )THEN
154
155	zht_i_ini(:,:) = si(jp_hti)%fnow(:,:,1)
156	zht_s_ini(:,:) = si(jp_hts)%fnow(:,:,1)
157	zat_i_ini(:,:) = si(jp_ati)%fnow(:,:,1)
158	zts_u_ini(:,:) = si(jp_tsu)%fnow(:,:,1)
159	ztm_i_ini(:,:) = si(jp_tmi)%fnow(:,:,1)
160	zsm_i_ini(:,:) = si(jp_smi)%fnow(:,:,1)
161
162	ELSE ! ln_iceini_file = F
163
164	!-----------------------------
165	! 3.1) Hemisphere-dependent arrays
166	!-----------------------------
167	! assign initial thickness, concentration, snow depth and salinity to an hemisphere-dependent array
168	DO jj = 1, jpj
169	DO ji = 1, jpi
170	IF( fcor(ji,jj) >= 0._wp ) THEN
171	zht_i_ini(ji,jj) = rn_hti_ini_n
172	zht_s_ini(ji,jj) = rn_hts_ini_n
173	zat_i_ini(ji,jj) = rn_ati_ini_n
174	zts_u_ini(ji,jj) = rn_tmi_ini_n
175	zsm_i_ini(ji,jj) = rn_smi_ini_n
176	ztm_i_ini(ji,jj) = rn_tmi_ini_n
177	ELSE
178	zht_i_ini(ji,jj) = rn_hti_ini_s
179	zht_s_ini(ji,jj) = rn_hts_ini_s
180	zat_i_ini(ji,jj) = rn_ati_ini_s
181	zts_u_ini(ji,jj) = rn_tmi_ini_s
182	zsm_i_ini(ji,jj) = rn_smi_ini_s
183	ztm_i_ini(ji,jj) = rn_tmi_ini_s
184	ENDIF
185	END DO
186	END DO
187
188	ENDIF ! ln_iceini_file
189
190	zvt_i_ini(:,:) = zht_i_ini(:,:) * zat_i_ini(:,:) ! ice volume
191
192	!---------------------------------------------------------------------
193	! 3.2) Distribute ice concentration and thickness into the categories
194	!---------------------------------------------------------------------
195	! a gaussian distribution for ice concentration is used
196	! then we check whether the distribution fullfills
197	! volume and area conservation, positivity and ice categories bounds
198	zh_i_ini(:,:,:) = 0._wp
199	za_i_ini(:,:,:) = 0._wp
200	zv_i_ini(:,:,:) = 0._wp
201
202	DO jj = 1, jpj
203	DO ji = 1, jpi
204
205	IF( zat_i_ini(ji,jj) > 0._wp .AND. zht_i_ini(ji,jj) > 0._wp )THEN
206
207	ztest_1 = 0 ; ztest_2 = 0 ; ztest_3 = 0 ; ztest_4 = 0
208	! ztests = 0
209
210	DO i_fill = jpl, 1, -1
211
212	! IF( ztests .NE. 4 ) THEN
213	IF ( ( ztest_1 + ztest_2 + ztest_3 + ztest_4 ) .NE. 4 ) THEN
214	!----------------------------
215	! fill the i_fill categories
216	!----------------------------
217	! *** 1 category to fill
218	IF ( i_fill .EQ. 1 ) THEN
219	zh_i_ini(ji,jj, 1) = zht_i_ini(ji,jj)
220	za_i_ini(ji,jj, 1) = zat_i_ini(ji,jj)
221	zh_i_ini(ji,jj,2:jpl) = 0._wp
222	za_i_ini(ji,jj,2:jpl) = 0._wp
223	ELSE
224
225	! *** >1 categores to fill
226	!--- Ice thicknesses in the i_fill - 1 first categories
227	DO jl = 1, i_fill - 1
228	zh_i_ini(ji,jj,jl) = hi_mean(jl)
229	END DO
230
231	!--- jl0: most likely index where cc will be maximum
232	DO jl = 1, jpl
233	IF ( ( zht_i_ini(ji,jj) > hi_max(jl-1) ) .AND. &
234	& ( zht_i_ini(ji,jj) <= hi_max(jl) ) ) THEN
235	jl0 = jl
236	ENDIF
237	END DO
238	jl0 = MIN(jl0, i_fill)
239
240	!--- Concentrations
241	za_i_ini(ji,jj,jl0) = zat_i_ini(ji,jj) / SQRT(REAL(jpl))
242	DO jl = 1, i_fill - 1
243	IF( jl .NE. jl0 )THEN
244	zsigma = 0.5 * zht_i_ini(ji,jj)
245	zarg = ( zh_i_ini(ji,jj,jl) - zht_i_ini(ji,jj) ) / zsigma
246	za_i_ini(ji,jj,jl) = za_i_ini(ji,jj,jl0) * EXP(-zarg**2)
247	ENDIF
248	END DO
249
250	zA = 0. ! sum of the areas in the jpl categories
251	DO jl = 1, i_fill - 1
252	zA = zA + za_i_ini(ji,jj,jl)
253	END DO
254	za_i_ini(ji,jj,i_fill) = zat_i_ini(ji,jj) - zA ! ice conc in the last category
255	IF ( i_fill .LT. jpl ) za_i_ini(ji,jj,i_fill+1:jpl) = 0._wp
256
257	!--- Ice thickness in the last category
258	zV = 0. ! sum of the volumes of the N-1 categories
259	DO jl = 1, i_fill - 1
260	zV = zV + za_i_ini(ji,jj,jl)*zh_i_ini(ji,jj,jl)
261	END DO
262	zh_i_ini(ji,jj,i_fill) = ( zvt_i_ini(ji,jj) - zV ) / za_i_ini(ji,jj,i_fill)
263	IF ( i_fill .LT. jpl ) zh_i_ini(ji,jj,i_fill+1:jpl) = 0._wp
264
265	!--- volumes
266	zv_i_ini(ji,jj,:) = za_i_ini(ji,jj,:) * zh_i_ini(ji,jj,:)
267	IF ( i_fill .LT. jpl ) zv_i_ini(ji,jj,i_fill+1:jpl) = 0._wp
268
269	ENDIF ! i_fill
270
271	!---------------------
272	! Compatibility tests
273	!---------------------
274	! Test 1: area conservation
275	zA_cons = SUM(za_i_ini(ji,jj,:)) ; zconv = ABS(zat_i_ini(ji,jj) - zA_cons )
276	IF ( zconv .LT. 1.0e-6 ) THEN
277	ztest_1 = 1
278	ELSE
279	!this write is useful
280	IF(lwp) WRITE(numout,) ' TEST1 AREA NOT CONSERVED *** zA_cons = ', zA_cons,' zat_i_ini = ',zat_i_ini(ji,jj)
281	ztest_1 = 0
282	ENDIF
283
284	! Test 2: volume conservation
285	zV_cons = SUM(zv_i_ini(ji,jj,:))
286	zconv = ABS(zvt_i_ini(ji,jj) - zV_cons)
287
288	IF( zconv .LT. 1.0e-6 ) THEN
289	ztest_2 = 1
290	ELSE
291	!this write is useful
292	IF(lwp) WRITE(numout,) ' TEST2 VOLUME NOT CONSERVED *** zV_cons = ', zV_cons, &
293	' zvt_i_ini = ', zvt_i_ini(ji,jj)
294	ztest_2 = 0
295	ENDIF
296
297	! Test 3: thickness of the last category is in-bounds ?
298	IF ( zh_i_ini(ji,jj,i_fill) > hi_max(i_fill-1) ) THEN
299	ztest_3 = 1
300	ELSE
301	! this write is useful
302	IF(lwp) WRITE(numout,) ' TEST 3 THICKNESS OF THE LAST CATEGORY OUT OF BOUNDS *** zh_i_ini(ji,jj,i_fill) = ', &
303	zh_i_ini(ji,jj,i_fill), ' hi_max(jpl-1) = ', hi_max(i_fill-1)
304	IF(lwp) WRITE(numout,*) ' ji,jj,i_fill ',ji,jj,i_fill
305	IF(lwp) WRITE(numout,*) 'zht_i_ini ',zht_i_ini(ji,jj)
306	ztest_3 = 0
307	ENDIF
308
309	! Test 4: positivity of ice concentrations
310	ztest_4 = 1
311	DO jl = 1, jpl
312	IF ( za_i_ini(ji,jj,jl) .LT. 0._wp ) THEN
313	! this write is useful
314	IF(lwp) WRITE(numout,) ' TEST 4 POSITIVITY NOT OK FOR CAT ', jl, ' WITH A = ', za_i_ini(ji,jj,jl)
315	ztest_4 = 0
316	ENDIF
317	END DO
318
319	ENDIF ! ztest_1 + ztest_2 + ztest_3 + ztest_4
320
321	ztests = ztest_1 + ztest_2 + ztest_3 + ztest_4
322
323	END DO ! i_fill
324
325	IF(lwp) THEN
326	WRITE(numout,*) ' ztests : ', ztests
327	IF( ztests .NE. 4 )THEN
328	WRITE(numout,*)
329	WRITE(numout,*) ' !!!! ALERT !!! '
330	WRITE(numout,*) ' !!!! Something is wrong in the LIM3 initialization procedure '
331	WRITE(numout,*)
332	WRITE(numout,) ' ** ztests is not equal to 4 '
333	WRITE(numout,) ' ** ztest_i (i=1,4) = ', ztest_1, ztest_2, ztest_3, ztest_4
334	WRITE(numout,*) ' zat_i_ini : ', zat_i_ini(ji,jj)
335	WRITE(numout,*) ' zht_i_ini : ', zht_i_ini(ji,jj)
336	ENDIF ! ztests .NE. 4
337	ENDIF
338
339	ENDIF ! zat_i_ini(ji,jj) > 0._wp .AND. zhm_i_ini(ji,jj) > 0._wp
340
341	ENDDO
342	ENDDO
343
344	!---------------------------------------------------------------------
345	! 3.3) Space-dependent arrays for ice state variables
346	!---------------------------------------------------------------------
347
348	! Ice concentration, thickness and volume, ice salinity, ice age, surface temperature
349	DO jl = 1, jpl ! loop over categories
350	DO jj = 1, jpj
351	DO ji = 1, jpi
352	a_i(ji,jj,jl) = zswitch(ji,jj) * za_i_ini(ji,jj,jl) ! concentration
353	ht_i(ji,jj,jl) = zswitch(ji,jj) * zh_i_ini(ji,jj,jl) ! ice thickness
354	sm_i(ji,jj,jl) = zswitch(ji,jj) * zsm_i_ini(ji,jj) ! salinity
355	o_i(ji,jj,jl) = zswitch(ji,jj) * 1._wp ! age (1 day)
356	t_su(ji,jj,jl) = zswitch(ji,jj) * zts_u_ini(ji,jj) + ( 1._wp - zswitch(ji,jj) ) * rt0 ! surf temp
357
358	IF( zht_i_ini(ji,jj) > 0._wp )THEN
359	ht_s(ji,jj,jl)= ht_i(ji,jj,jl) * ( zht_s_ini(ji,jj) / zht_i_ini(ji,jj) ) ! snow depth
360	ELSE
361	ht_s(ji,jj,jl)= 0._wp
362	ENDIF
363
364	! This case below should not be used if (ht_s/ht_i) is ok in namelist
365	! In case snow load is in excess that would lead to transformation from snow to ice
366	! Then, transfer the snow excess into the ice (different from limthd_dh)
367	zdh = MAX( 0._wp, ( rhosn * ht_s(ji,jj,jl) + ( rhoic - rau0 ) * ht_i(ji,jj,jl) ) * r1_rau0 )
368	! recompute ht_i, ht_s avoiding out of bounds values
369	ht_i(ji,jj,jl) = MIN( hi_max(jl), ht_i(ji,jj,jl) + zdh )
370	ht_s(ji,jj,jl) = MAX( 0._wp, ht_s(ji,jj,jl) - zdh * rhoic * r1_rhosn )
371
372	! ice volume, salt content, age content
373	v_i(ji,jj,jl) = ht_i(ji,jj,jl) * a_i(ji,jj,jl) ! ice volume
374	v_s(ji,jj,jl) = ht_s(ji,jj,jl) * a_i(ji,jj,jl) ! snow volume
375	smv_i(ji,jj,jl) = MIN( sm_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl) ! salt content
376	oa_i(ji,jj,jl) = o_i(ji,jj,jl) * a_i(ji,jj,jl) ! age content
377	END DO
378	END DO
379	END DO
380
381	! Snow temperature and heat content
382	DO jk = 1, nlay_s
383	DO jl = 1, jpl ! loop over categories
384	DO jj = 1, jpj
385	DO ji = 1, jpi
386	t_s(ji,jj,jk,jl) = zswitch(ji,jj) * ztm_i_ini(ji,jj) + ( 1._wp - zswitch(ji,jj) ) * rt0
387	! Snow energy of melting
388	e_s(ji,jj,jk,jl) = zswitch(ji,jj) * rhosn * ( cpic * ( rt0 - t_s(ji,jj,jk,jl) ) + lfus )
389
390	! Mutliply by volume, and divide by number of layers to get heat content in J/m2
391	e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * v_s(ji,jj,jl) * r1_nlay_s
392	END DO
393	END DO
394	END DO
395	END DO
396
397	! Ice salinity, temperature and heat content
398	DO jk = 1, nlay_i
399	DO jl = 1, jpl ! loop over categories
400	DO jj = 1, jpj
401	DO ji = 1, jpi
402	t_i(ji,jj,jk,jl) = zswitch(ji,jj) * ztm_i_ini(ji,jj) + ( 1._wp - zswitch(ji,jj) ) * rt0
403	s_i(ji,jj,jk,jl) = zswitch(ji,jj) * zsm_i_ini(ji,jj) + ( 1._wp - zswitch(ji,jj) ) * rn_simin
404	ztmelts = - tmut * s_i(ji,jj,jk,jl) + rt0 !Melting temperature in K
405
406	! heat content per unit volume
407	e_i(ji,jj,jk,jl) = zswitch(ji,jj) * rhoic * ( cpic * ( ztmelts - t_i(ji,jj,jk,jl) ) &
408	+ lfus * ( 1._wp - (ztmelts-rt0) / MIN((t_i(ji,jj,jk,jl)-rt0),-epsi20) ) &
409	- rcp * ( ztmelts - rt0 ) )
410
411	! Mutliply by ice volume, and divide by number of layers to get heat content in J/m2
412	e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * v_i(ji,jj,jl) * r1_nlay_i
413	END DO
414	END DO
415	END DO
416	END DO
417
418	tn_ice (:,:,:) = t_su (:,:,:)
419
420	ELSE ! if ln_iceini=false
421	a_i (:,:,:) = 0._wp
422	v_i (:,:,:) = 0._wp
423	v_s (:,:,:) = 0._wp
424	smv_i(:,:,:) = 0._wp
425	oa_i (:,:,:) = 0._wp
426	ht_i (:,:,:) = 0._wp
427	ht_s (:,:,:) = 0._wp
428	sm_i (:,:,:) = 0._wp
429	o_i (:,:,:) = 0._wp
430
431	e_i(:,:,:,:) = 0._wp
432	e_s(:,:,:,:) = 0._wp
433
434	DO jl = 1, jpl
435	DO jk = 1, nlay_i
436	t_i(:,:,jk,jl) = rt0 * tmask(:,:,1)
437	END DO
438	DO jk = 1, nlay_s
439	t_s(:,:,jk,jl) = rt0 * tmask(:,:,1)
440	END DO
441	END DO
442
443	ENDIF ! ln_iceini
444
445	at_i (:,:) = 0.0_wp
446	DO jl = 1, jpl
447	at_i (:,:) = at_i (:,:) + a_i (:,:,jl)
448	END DO
449	!
450	!--------------------------------------------------------------------
451	! 4) Global ice variables for output diagnostics \|
452	!--------------------------------------------------------------------
453	u_ice (:,:) = 0._wp
454	v_ice (:,:) = 0._wp
455	stress1_i(:,:) = 0._wp
456	stress2_i(:,:) = 0._wp
457	stress12_i(:,:) = 0._wp
458
459	!--------------------------------------------------------------------
460	! 5) Moments for advection
461	!--------------------------------------------------------------------
462
463	sxopw (:,:) = 0._wp
464	syopw (:,:) = 0._wp
465	sxxopw(:,:) = 0._wp
466	syyopw(:,:) = 0._wp
467	sxyopw(:,:) = 0._wp
468
469	sxice (:,:,:) = 0._wp ; sxsn (:,:,:) = 0._wp ; sxa (:,:,:) = 0._wp
470	syice (:,:,:) = 0._wp ; sysn (:,:,:) = 0._wp ; sya (:,:,:) = 0._wp
471	sxxice(:,:,:) = 0._wp ; sxxsn(:,:,:) = 0._wp ; sxxa (:,:,:) = 0._wp
472	syyice(:,:,:) = 0._wp ; syysn(:,:,:) = 0._wp ; syya (:,:,:) = 0._wp
473	sxyice(:,:,:) = 0._wp ; sxysn(:,:,:) = 0._wp ; sxya (:,:,:) = 0._wp
474
475	sxc0 (:,:,:) = 0._wp ; sxe (:,:,:,:)= 0._wp
476	syc0 (:,:,:) = 0._wp ; sye (:,:,:,:)= 0._wp
477	sxxc0 (:,:,:) = 0._wp ; sxxe (:,:,:,:)= 0._wp
478	syyc0 (:,:,:) = 0._wp ; syye (:,:,:,:)= 0._wp
479	sxyc0 (:,:,:) = 0._wp ; sxye (:,:,:,:)= 0._wp
480
481	sxsal (:,:,:) = 0._wp
482	sysal (:,:,:) = 0._wp
483	sxxsal (:,:,:) = 0._wp
484	syysal (:,:,:) = 0._wp
485	sxysal (:,:,:) = 0._wp
486
487	sxage (:,:,:) = 0._wp
488	syage (:,:,:) = 0._wp
489	sxxage (:,:,:) = 0._wp
490	syyage (:,:,:) = 0._wp
491	sxyage (:,:,:) = 0._wp
492
493
494	CALL wrk_dealloc( jpi, jpj, jpl, zh_i_ini, za_i_ini, zv_i_ini )
495	CALL wrk_dealloc( jpi, jpj, zht_i_ini, zat_i_ini, zvt_i_ini, zts_u_ini, zht_s_ini, zsm_i_ini, ztm_i_ini )
496	CALL wrk_dealloc( jpi, jpj, zswitch )
497
498	END SUBROUTINE lim_istate
499
500	SUBROUTINE lim_istate_init
501	!!-------------------------------------------------------------------
502	!! * ROUTINE lim_istate_init *
503	!!
504	!! ** Purpose : Definition of initial state of the ice
505	!!
506	!! ** Method : Read the namiceini namelist and check the parameter
507	!! values called at the first timestep (nit000)
508	!!
509	!! ** input :
510	!! Namelist namiceini
511	!!
512	!! history :
513	!! 8.5 ! 03-08 (C. Ethe) original code
514	!! 8.5 ! 07-11 (M. Vancoppenolle) rewritten initialization
515	!!-----------------------------------------------------------------------------
516	!
517	INTEGER :: ios,ierr,inum_ice ! Local integer output status for namelist read
518	INTEGER :: ji,jj
519	INTEGER :: ifpr, ierror
520	!
521	CHARACTER(len=100) :: cn_dir ! Root directory for location of ice files
522	TYPE(FLD_N) :: sn_hti, sn_hts, sn_ati, sn_tsu, sn_tmi, sn_smi
523	TYPE(FLD_N), DIMENSION(jpfldi) :: slf_i ! array of namelist informations on the fields to read
524	!
525	NAMELIST/namiceini/ ln_iceini, ln_iceini_file, rn_thres_sst, rn_hts_ini_n, rn_hts_ini_s, &
526	& rn_hti_ini_n, rn_hti_ini_s, rn_ati_ini_n, rn_ati_ini_s, rn_smi_ini_n, &
527	& rn_smi_ini_s, rn_tmi_ini_n, rn_tmi_ini_s, &
528	& sn_hti, sn_hts, sn_ati, sn_tsu, sn_tmi, sn_smi, cn_dir
529	!!-----------------------------------------------------------------------------
530	!
531	REWIND( numnam_ice_ref ) ! Namelist namiceini in reference namelist : Ice initial state
532	READ ( numnam_ice_ref, namiceini, IOSTAT = ios, ERR = 901)
533	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namiceini in reference namelist', lwp )
534
535	REWIND( numnam_ice_cfg ) ! Namelist namiceini in configuration namelist : Ice initial state
536	READ ( numnam_ice_cfg, namiceini, IOSTAT = ios, ERR = 902 )
537	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namiceini in configuration namelist', lwp )
538	IF(lwm) WRITE ( numoni, namiceini )
539
540	slf_i(jp_hti) = sn_hti ; slf_i(jp_hts) = sn_hts
541	slf_i(jp_ati) = sn_ati ; slf_i(jp_tsu) = sn_tsu
542	slf_i(jp_tmi) = sn_tmi ; slf_i(jp_smi) = sn_smi
543
544	! Define the initial parameters
545	! -------------------------
546
547	IF(lwp) THEN
548	WRITE(numout,*)
549	WRITE(numout,*) 'lim_istate_init : ice parameters inititialisation '
550	WRITE(numout,*) '~~~~~~~~~~~~~~~'
551	WRITE(numout,*) ' initialization with ice (T) or not (F) ln_iceini = ', ln_iceini
552	WRITE(numout,*) ' ice initialization from a netcdf file ln_iceini_file = ', ln_iceini_file
553	WRITE(numout,*) ' threshold water temp. for initial sea-ice rn_thres_sst = ', rn_thres_sst
554	WRITE(numout,*) ' initial snow thickness in the north rn_hts_ini_n = ', rn_hts_ini_n
555	WRITE(numout,*) ' initial snow thickness in the south rn_hts_ini_s = ', rn_hts_ini_s
556	WRITE(numout,*) ' initial ice thickness in the north rn_hti_ini_n = ', rn_hti_ini_n
557	WRITE(numout,*) ' initial ice thickness in the south rn_hti_ini_s = ', rn_hti_ini_s
558	WRITE(numout,*) ' initial ice concentr. in the north rn_ati_ini_n = ', rn_ati_ini_n
559	WRITE(numout,*) ' initial ice concentr. in the north rn_ati_ini_s = ', rn_ati_ini_s
560	WRITE(numout,*) ' initial ice salinity in the north rn_smi_ini_n = ', rn_smi_ini_n
561	WRITE(numout,*) ' initial ice salinity in the south rn_smi_ini_s = ', rn_smi_ini_s
562	WRITE(numout,*) ' initial ice/snw temp in the north rn_tmi_ini_n = ', rn_tmi_ini_n
563	WRITE(numout,*) ' initial ice/snw temp in the south rn_tmi_ini_s = ', rn_tmi_ini_s
564	ENDIF
565
566	IF( ln_iceini_file ) THEN ! Ice initialization using input file
567	!
568	! set si structure
569	ALLOCATE( si(jpfldi), STAT=ierror )
570	IF( ierror > 0 ) THEN
571	CALL ctl_stop( 'Ice_ini in limistate: unable to allocate si structure' ) ; RETURN
572	ENDIF
573
574	DO ifpr = 1, jpfldi
575	ALLOCATE( si(ifpr)%fnow(jpi,jpj,1) )
576	ALLOCATE( si(ifpr)%fdta(jpi,jpj,1,2) )
577	END DO
578
579	! fill si with slf_i and control print
580	CALL fld_fill( si, slf_i, cn_dir, 'lim_istate', 'lim istate ini', 'numnam_ice' )
581
582	CALL fld_read( nit000, 1, si ) ! input fields provided at the current time-step
583
584	ENDIF
585
586	END SUBROUTINE lim_istate_init
587
588	#else
589	!!----------------------------------------------------------------------
590	!! Default option : Empty module NO LIM sea-ice model
591	!!----------------------------------------------------------------------
592	CONTAINS
593	SUBROUTINE lim_istate ! Empty routine
594	END SUBROUTINE lim_istate
595	#endif
596
597	!!======================================================================
598	END MODULE limistate

Note: See TracBrowser for help on using the repository browser.

New URL for NEMO forge! http://forge.nemo-ocean.eu

Context Navigation

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

Download in other formats: