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

Context Navigation

source: branches/2016/dev_CNRS_2016/NEMOGCM/NEMO/LIM_SRC_3/limistate.F90 @ 7280

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

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

Download in other formats: