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

Context Navigation

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

Last change on this file since 4737 was 4688, checked in by clem, 10 years ago
new version of LIM3 with perfect conservation of heat, see ticket #1352
Property svn:keywords set to `Id`
File size: 24.2 KB

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

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

Download in other formats: