Context Navigation

limthd_dh.F90 @ 4990

Last change on this file since 4990 was 4990, checked in by timgraham, 9 years ago

Merged branches/2014/dev_MERGE_2014 back onto the trunk as follows:

In the working copy of branch ran:
svn merge svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/trunk@HEAD
1 conflict in LIM_SRC_3/limdiahsb.F90
Resolved by keeping the version from dev_MERGE_2014 branch
and commited at r4989

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/2014/dev_MERGE_2014
to merge the branch into the trunk - no conflicts at this stage.

Property svn:keywords set to Id

File size: 36.2 KB

Rev	Line
[825]	1	MODULE limthd_dh
[1572]	2	!!======================================================================
	3	!! * MODULE limthd_dh *
	4	!! LIM-3 : thermodynamic growth and decay of the ice
	5	!!======================================================================
	6	!! History : LIM ! 2003-05 (M. Vancoppenolle) Original code in 1D
	7	!! ! 2005-06 (M. Vancoppenolle) 3D version
[4688]	8	!! 3.2 ! 2009-07 (M. Vancoppenolle, Y. Aksenov, G. Madec) bug correction in wfx_snw & wfx_ice
[3625]	9	!! 3.4 ! 2011-02 (G. Madec) dynamical allocation
	10	!! 3.5 ! 2012-10 (G. Madec & co) salt flux + bug fixes
[1572]	11	!!----------------------------------------------------------------------
[825]	12	#if defined key_lim3
[834]	13	!!----------------------------------------------------------------------
	14	!! 'key_lim3' LIM3 sea-ice model
	15	!!----------------------------------------------------------------------
[3625]	16	!! lim_thd_dh : vertical accr./abl. and lateral ablation of sea ice
[825]	17	!!----------------------------------------------------------------------
[3625]	18	USE par_oce ! ocean parameters
	19	USE phycst ! physical constants (OCE directory)
	20	USE sbc_oce ! Surface boundary condition: ocean fields
	21	USE ice ! LIM variables
	22	USE par_ice ! LIM parameters
	23	USE thd_ice ! LIM thermodynamics
	24	USE in_out_manager ! I/O manager
	25	USE lib_mpp ! MPP library
	26	USE wrk_nemo ! work arrays
	27	USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)
[4688]	28
[825]	29	IMPLICIT NONE
	30	PRIVATE
	31
[1572]	32	PUBLIC lim_thd_dh ! called by lim_thd
[825]	33
	34	!!----------------------------------------------------------------------
[4161]	35	!! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2010)
[1156]	36	!! $Id$
[2715]	37	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[825]	38	!!----------------------------------------------------------------------
	39	CONTAINS
	40
[4688]	41	SUBROUTINE lim_thd_dh( kideb, kiut )
[921]	42	!!------------------------------------------------------------------
	43	!! * ROUTINE lim_thd_dh *
	44	!!
[1572]	45	!! ** Purpose : determines variations of ice and snow thicknesses.
[921]	46	!!
[1572]	47	!! ** Method : Ice/Snow surface melting arises from imbalance in surface fluxes
	48	!! Bottom accretion/ablation arises from flux budget
	49	!! Snow thickness can increase by precipitation and decrease by sublimation
	50	!! If snow load excesses Archmiede limit, snow-ice is formed by
	51	!! the flooding of sea-water in the snow
[921]	52	!!
[1572]	53	!! 1) Compute available flux of heat for surface ablation
	54	!! 2) Compute snow and sea ice enthalpies
	55	!! 3) Surface ablation and sublimation
	56	!! 4) Bottom accretion/ablation
	57	!! 5) Case of Total ablation
	58	!! 6) Snow ice formation
[921]	59	!!
[1572]	60	!! References : Bitz and Lipscomb, 1999, J. Geophys. Res.
	61	!! Fichefet T. and M. Maqueda 1997, J. Geophys. Res., 102(C6), 12609-12646
	62	!! Vancoppenolle, Fichefet and Bitz, 2005, Geophys. Res. Let.
	63	!! Vancoppenolle et al.,2009, Ocean Modelling
[921]	64	!!------------------------------------------------------------------
[1572]	65	INTEGER , INTENT(in) :: kideb, kiut ! Start/End point on which the the computation is applied
	66	!!
	67	INTEGER :: ji , jk ! dummy loop indices
[4161]	68	INTEGER :: ii, ij ! 2D corresponding indices to ji
[1572]	69	INTEGER :: iter
[825]	70
[4688]	71	REAL(wp) :: ztmelts ! local scalar
	72	REAL(wp) :: zdh, zfdum !
[1572]	73	REAL(wp) :: zfracs ! fractionation coefficient for bottom salt entrapment
	74	REAL(wp) :: zcoeff ! dummy argument for snowfall partitioning over ice and leads
[4688]	75	REAL(wp) :: zs_snic ! snow-ice salinity
[1572]	76	REAL(wp) :: zswi1 ! switch for computation of bottom salinity
	77	REAL(wp) :: zswi12 ! switch for computation of bottom salinity
	78	REAL(wp) :: zswi2 ! switch for computation of bottom salinity
	79	REAL(wp) :: zgrr ! bottom growth rate
[4688]	80	REAL(wp) :: zt_i_new ! bottom formation temperature
	81
	82	REAL(wp) :: zQm ! enthalpy exchanged with the ocean (J/m2), >0 towards the ocean
	83	REAL(wp) :: zEi ! specific enthalpy of sea ice (J/kg)
	84	REAL(wp) :: zEw ! specific enthalpy of exchanged water (J/kg)
	85	REAL(wp) :: zdE ! specific enthalpy difference (J/kg)
	86	REAL(wp) :: zfmdt ! exchange mass flux x time step (J/m2), >0 towards the ocean
	87	REAL(wp) :: zsstK ! SST in Kelvin
	88
[3294]	89	REAL(wp), POINTER, DIMENSION(:) :: zh_s ! snow layer thickness
[4688]	90	REAL(wp), POINTER, DIMENSION(:) :: zqprec ! energy of fallen snow (J.m-3)
	91	REAL(wp), POINTER, DIMENSION(:) :: zq_su ! heat for surface ablation (J.m-2)
	92	REAL(wp), POINTER, DIMENSION(:) :: zq_bo ! heat for bottom ablation (J.m-2)
	93	REAL(wp), POINTER, DIMENSION(:) :: zq_1cat ! corrected heat in case 1-cat and hmelt>15cm (J.m-2)
	94	REAL(wp), POINTER, DIMENSION(:) :: zq_rema ! remaining heat at the end of the routine (J.m-2)
	95	REAL(wp), POINTER, DIMENSION(:) :: zf_tt ! Heat budget to determine melting or freezing(W.m-2)
	96	INTEGER , POINTER, DIMENSION(:) :: icount ! number of layers vanished by melting
[3294]	97
[3625]	98	REAL(wp), POINTER, DIMENSION(:) :: zdh_s_mel ! snow melt
	99	REAL(wp), POINTER, DIMENSION(:) :: zdh_s_pre ! snow precipitation
	100	REAL(wp), POINTER, DIMENSION(:) :: zdh_s_sub ! snow sublimation
[3294]	101
	102	REAL(wp), POINTER, DIMENSION(:,:) :: zdeltah
[4688]	103	REAL(wp), POINTER, DIMENSION(:,:) :: zh_i ! ice layer thickness
[3294]	104
[4688]	105	REAL(wp), POINTER, DIMENSION(:) :: zqh_i ! total ice heat content (J.m-2)
	106	REAL(wp), POINTER, DIMENSION(:) :: zqh_s ! total snow heat content (J.m-2)
	107	REAL(wp), POINTER, DIMENSION(:) :: zq_s ! total snow enthalpy (J.m-3)
[3294]	108
[4161]	109	! mass and salt flux (clem)
[4990]	110	REAL(wp) :: zdvres, zswitch_sal
[4161]	111
[3294]	112	! Heat conservation
[4688]	113	INTEGER :: num_iter_max
	114
[1572]	115	!!------------------------------------------------------------------
[825]	116
[4688]	117	! Discriminate between varying salinity (num_sal=2) and prescribed cases (other values)
	118	SELECT CASE( num_sal ) ! varying salinity or not
	119	CASE( 1, 3, 4 ) ; zswitch_sal = 0 ! prescribed salinity profile
	120	CASE( 2 ) ; zswitch_sal = 1 ! varying salinity profile
	121	END SELECT
[825]	122
[4688]	123	CALL wrk_alloc( jpij, zh_s, zqprec, zq_su, zq_bo, zf_tt, zq_1cat, zq_rema )
	124	CALL wrk_alloc( jpij, zdh_s_mel, zdh_s_pre, zdh_s_sub, zqh_i, zqh_s, zq_s )
[4873]	125	CALL wrk_alloc( jpij, nlay_i+1, zdeltah, zh_i )
[4688]	126	CALL wrk_alloc( jpij, icount )
[4161]	127
[4688]	128	dh_i_surf (:) = 0._wp ; dh_i_bott (:) = 0._wp ; dh_snowice(:) = 0._wp
	129	dsm_i_se_1d(:) = 0._wp ; dsm_i_si_1d(:) = 0._wp
	130
	131	zqprec (:) = 0._wp ; zq_su (:) = 0._wp ; zq_bo (:) = 0._wp ; zf_tt (:) = 0._wp
	132	zq_1cat(:) = 0._wp ; zq_rema(:) = 0._wp
[2715]	133
[4688]	134	zh_s (:) = 0._wp
	135	zdh_s_pre(:) = 0._wp
	136	zdh_s_mel(:) = 0._wp
	137	zdh_s_sub(:) = 0._wp
	138	zqh_s (:) = 0._wp
	139	zqh_i (:) = 0._wp
[4161]	140
[4688]	141	zh_i (:,:) = 0._wp
	142	zdeltah (:,:) = 0._wp
	143	icount (:) = 0
	144
	145	! initialize layer thicknesses and enthalpies
	146	h_i_old (:,0:nlay_i+1) = 0._wp
	147	qh_i_old(:,0:nlay_i+1) = 0._wp
	148	DO jk = 1, nlay_i
	149	DO ji = kideb, kiut
[4872]	150	h_i_old (ji,jk) = ht_i_1d(ji) / REAL( nlay_i )
	151	qh_i_old(ji,jk) = q_i_1d(ji,jk) * h_i_old(ji,jk)
[4688]	152	ENDDO
	153	ENDDO
[921]	154	!
	155	!------------------------------------------------------------------------------!
[4688]	156	! 1) Calculate available heat for surface and bottom ablation !
[921]	157	!------------------------------------------------------------------------------!
	158	!
[2715]	159	DO ji = kideb, kiut
[4990]	160	rswitch = 1._wp - MAX( 0._wp , SIGN( 1._wp , - ht_s_1d(ji) ) )
	161	ztmelts = rswitch * rtt + ( 1._wp - rswitch ) * rtt
[825]	162
[4990]	163	zfdum = qns_ice_1d(ji) + ( 1._wp - i0(ji) ) * qsr_ice_1d(ji) - fc_su(ji)
	164	zf_tt(ji) = fc_bo_i(ji) + fhtur_1d(ji) + fhld_1d(ji)
[4688]	165
[4872]	166	zq_su (ji) = MAX( 0._wp, zfdum * rdt_ice ) * MAX( 0._wp , SIGN( 1._wp, t_su_1d(ji) - ztmelts ) )
[4688]	167	zq_bo (ji) = MAX( 0._wp, zf_tt(ji) * rdt_ice )
	168	END DO
	169
[921]	170	!
	171	!------------------------------------------------------------------------------!
[4688]	172	! If snow temperature is above freezing point, then snow melts
	173	! (should not happen but sometimes it does)
[921]	174	!------------------------------------------------------------------------------!
[4688]	175	DO ji = kideb, kiut
[4872]	176	IF( t_s_1d(ji,1) > rtt ) THEN !!! Internal melting
[4688]	177	! Contribution to heat flux to the ocean [W.m-2], < 0
[4872]	178	hfx_res_1d(ji) = hfx_res_1d(ji) + q_s_1d(ji,1) * ht_s_1d(ji) * a_i_1d(ji) * r1_rdtice
[4688]	179	! Contribution to mass flux
[4872]	180	wfx_snw_1d(ji) = wfx_snw_1d(ji) + rhosn * ht_s_1d(ji) * a_i_1d(ji) * r1_rdtice
[4688]	181	! updates
[4872]	182	ht_s_1d(ji) = 0._wp
	183	q_s_1d (ji,1) = 0._wp
	184	t_s_1d (ji,1) = rtt
[4688]	185	END IF
	186	END DO
	187
	188	!------------------------------------------------------------!
	189	! 2) Computing layer thicknesses and enthalpies. !
	190	!------------------------------------------------------------!
[921]	191	!
[4688]	192	DO ji = kideb, kiut
[4872]	193	zh_s(ji) = ht_s_1d(ji) / REAL( nlay_s )
[825]	194	END DO
[2715]	195	!
[825]	196	DO jk = 1, nlay_s
[2715]	197	DO ji = kideb, kiut
[4872]	198	zqh_s(ji) = zqh_s(ji) + q_s_1d(ji,jk) * zh_s(ji)
[825]	199	END DO
	200	END DO
[2715]	201	!
[825]	202	DO jk = 1, nlay_i
[2715]	203	DO ji = kideb, kiut
[4872]	204	zh_i(ji,jk) = ht_i_1d(ji) / REAL( nlay_i )
	205	zqh_i(ji) = zqh_i(ji) + q_i_1d(ji,jk) * zh_i(ji,jk)
[825]	206	END DO
	207	END DO
[921]	208	!
	209	!------------------------------------------------------------------------------\|
	210	! 3) Surface ablation and sublimation \|
	211	!------------------------------------------------------------------------------\|
	212	!
[834]	213	!-------------------------
	214	! 3.1 Snow precips / melt
	215	!-------------------------
[825]	216	! Snow accumulation in one thermodynamic time step
	217	! snowfall is partitionned between leads and ice
	218	! if snow fall was uniform, a fraction (1-at_i) would fall into leads
	219	! but because of the winds, more snow falls on leads than on sea ice
	220	! and a greater fraction (1-at_i)^beta of the total mass of snow
[834]	221	! (beta < 1) falls in leads.
[825]	222	! In reality, beta depends on wind speed,
	223	! and should decrease with increasing wind speed but here, it is
[834]	224	! considered as a constant. an average value is 0.66
[825]	225	! Martin Vancoppenolle, December 2006
	226
	227	DO ji = kideb, kiut
[4688]	228	!-----------
	229	! Snow fall
	230	!-----------
	231	! thickness change
[4872]	232	zcoeff = ( 1._wp - ( 1._wp - at_i_1d(ji) )**betas ) / at_i_1d(ji)
[825]	233	zdh_s_pre(ji) = zcoeff * sprecip_1d(ji) * rdt_ice / rhosn
[4688]	234	! enthalpy of the precip (>0, J.m-3) (tatm_ice is now in K)
	235	zqprec (ji) = rhosn * ( cpic * ( rtt - MIN( tatm_ice_1d(ji), rt0_snow) ) + lfus )
	236	IF( sprecip_1d(ji) == 0._wp ) zqprec(ji) = 0._wp
	237	! heat flux from snow precip (>0, W.m-2)
[4872]	238	hfx_spr_1d(ji) = hfx_spr_1d(ji) + zdh_s_pre(ji) * a_i_1d(ji) * zqprec(ji) * r1_rdtice
[4688]	239	! mass flux, <0
[4872]	240	wfx_spr_1d(ji) = wfx_spr_1d(ji) - rhosn * a_i_1d(ji) * zdh_s_pre(ji) * r1_rdtice
[4688]	241	! update thickness
[4872]	242	ht_s_1d (ji) = MAX( 0._wp , ht_s_1d(ji) + zdh_s_pre(ji) )
[825]	243
[4688]	244	!---------------------
	245	! Melt of falling snow
	246	!---------------------
	247	! thickness change
	248	IF( zdh_s_pre(ji) > 0._wp ) THEN
[4990]	249	rswitch = 1._wp - MAX( 0._wp , SIGN( 1._wp , - zqprec(ji) + epsi20 ) )
	250	zdh_s_mel (ji) = - rswitch * zq_su(ji) / MAX( zqprec(ji) , epsi20 )
[4688]	251	zdh_s_mel (ji) = MAX( - zdh_s_pre(ji), zdh_s_mel(ji) ) ! bound melting
	252	! heat used to melt snow (W.m-2, >0)
[4872]	253	hfx_snw_1d(ji) = hfx_snw_1d(ji) - zdh_s_mel(ji) * a_i_1d(ji) * zqprec(ji) * r1_rdtice
[4688]	254	! snow melting only = water into the ocean (then without snow precip), >0
[4872]	255	wfx_snw_1d(ji) = wfx_snw_1d(ji) - rhosn * a_i_1d(ji) * zdh_s_mel(ji) * r1_rdtice
[4688]	256
	257	! updates available heat + thickness
	258	zq_su (ji) = MAX( 0._wp , zq_su (ji) + zdh_s_mel(ji) * zqprec(ji) )
[4872]	259	ht_s_1d(ji) = MAX( 0._wp , ht_s_1d(ji) + zdh_s_mel(ji) )
	260	zh_s (ji) = ht_s_1d(ji) / REAL( nlay_s )
[4688]	261
	262	ENDIF
[825]	263	END DO
	264
[4688]	265	! If heat still available, then melt more snow
	266	zdeltah(:,:) = 0._wp ! important
[825]	267	DO jk = 1, nlay_s
	268	DO ji = kideb, kiut
[4688]	269	! thickness change
[4990]	270	rswitch = 1._wp - MAX( 0._wp, SIGN( 1._wp, - ht_s_1d(ji) ) )
	271	rswitch = rswitch * ( 1._wp - MAX( 0._wp, SIGN( 1._wp, - q_s_1d(ji,jk) + epsi20 ) ) )
	272	zdeltah (ji,jk) = - rswitch * zq_su(ji) / MAX( q_s_1d(ji,jk), epsi20 )
[4688]	273	zdeltah (ji,jk) = MAX( zdeltah(ji,jk) , - zh_s(ji) ) ! bound melting
	274	zdh_s_mel(ji) = zdh_s_mel(ji) + zdeltah(ji,jk)
	275	! heat used to melt snow(W.m-2, >0)
[4872]	276	hfx_snw_1d(ji) = hfx_snw_1d(ji) - zdeltah(ji,jk) * a_i_1d(ji) * q_s_1d(ji,jk) * r1_rdtice
[4688]	277	! snow melting only = water into the ocean (then without snow precip)
[4872]	278	wfx_snw_1d(ji) = wfx_snw_1d(ji) - rhosn * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice
[825]	279
[4688]	280	! updates available heat + thickness
[4872]	281	zq_su (ji) = MAX( 0._wp , zq_su (ji) + zdeltah(ji,jk) * q_s_1d(ji,jk) )
	282	ht_s_1d(ji) = MAX( 0._wp , ht_s_1d(ji) + zdeltah(ji,jk) )
[825]	283
[1572]	284	END DO
	285	END DO
[825]	286
[4688]	287	!----------------------
	288	! 3.2 Snow sublimation
	289	!----------------------
	290	! qla_ice is always >=0 (upwards), heat goes to the atmosphere, therefore snow sublimates
	291	! clem comment: not counted in mass exchange in limsbc since this is an exchange with atm. (not ocean)
	292	! clem comment: ice should also sublimate
	293	IF( lk_cpl ) THEN
	294	! coupled mode: sublimation already included in emp_ice (to do in limsbc_ice)
	295	zdh_s_sub(:) = 0._wp
	296	ELSE
	297	! forced mode: snow thickness change due to sublimation
[921]	298	DO ji = kideb, kiut
[4872]	299	zdh_s_sub(ji) = MAX( - ht_s_1d(ji) , - parsub * qla_ice_1d(ji) / ( rhosn * lsub ) * rdt_ice )
[4688]	300	! Heat flux by sublimation [W.m-2], < 0
	301	! sublimate first snow that had fallen, then pre-existing snow
	302	zcoeff = ( MAX( zdh_s_sub(ji), - MAX( 0._wp, zdh_s_pre(ji) + zdh_s_mel(ji) ) ) * zqprec(ji) + &
[4872]	303	& ( zdh_s_sub(ji) - MAX( zdh_s_sub(ji), - MAX( 0._wp, zdh_s_pre(ji) + zdh_s_mel(ji) ) ) ) * q_s_1d(ji,1) ) &
	304	& * a_i_1d(ji) * r1_rdtice
[4688]	305	hfx_sub_1d(ji) = hfx_sub_1d(ji) + zcoeff
	306	! Mass flux by sublimation
[4872]	307	wfx_sub_1d(ji) = wfx_sub_1d(ji) - rhosn * a_i_1d(ji) * zdh_s_sub(ji) * r1_rdtice
[4688]	308	! new snow thickness
[4872]	309	ht_s_1d(ji) = MAX( 0._wp , ht_s_1d(ji) + zdh_s_sub(ji) )
[1572]	310	END DO
	311	ENDIF
[825]	312
[4688]	313	! --- Update snow diags --- !
[825]	314	DO ji = kideb, kiut
[4688]	315	dh_s_tot(ji) = zdh_s_mel(ji) + zdh_s_pre(ji) + zdh_s_sub(ji)
[4872]	316	zh_s(ji) = ht_s_1d(ji) / REAL( nlay_s )
[4688]	317	END DO ! ji
[825]	318
[4688]	319	!-------------------------------------------
	320	! 3.3 Update temperature, energy
	321	!-------------------------------------------
	322	! new temp and enthalpy of the snow (remaining snow precip + remaining pre-existing snow)
	323	zq_s(:) = 0._wp
[825]	324	DO jk = 1, nlay_s
	325	DO ji = kideb,kiut
[4990]	326	rswitch = MAX( 0._wp , SIGN( 1._wp, - ht_s_1d(ji) + epsi20 ) )
	327	q_s_1d(ji,jk) = ( 1._wp - rswitch ) / MAX( ht_s_1d(ji), epsi20 ) * &
[4688]	328	& ( ( MAX( 0._wp, dh_s_tot(ji) ) ) * zqprec(ji) + &
[4872]	329	& ( - MAX( 0._wp, dh_s_tot(ji) ) + ht_s_1d(ji) ) * rhosn * ( cpic * ( rtt - t_s_1d(ji,jk) ) + lfus ) )
	330	zq_s(ji) = zq_s(ji) + q_s_1d(ji,jk)
[825]	331	END DO
	332	END DO
	333
[4688]	334	!--------------------------
	335	! 3.4 Surface ice ablation
	336	!--------------------------
	337	zdeltah(:,:) = 0._wp ! important
	338	DO jk = 1, nlay_i
	339	DO ji = kideb, kiut
[4872]	340	zEi = - q_i_1d(ji,jk) / rhoic ! Specific enthalpy of layer k [J/kg, <0]
[4688]	341
[4872]	342	ztmelts = - tmut * s_i_1d(ji,jk) + rtt ! Melting point of layer k [K]
[4688]	343
	344	zEw = rcp * ( ztmelts - rt0 ) ! Specific enthalpy of resulting meltwater [J/kg, <0]
	345
	346	zdE = zEi - zEw ! Specific enthalpy difference < 0
	347
	348	zfmdt = - zq_su(ji) / zdE ! Mass flux to the ocean [kg/m2, >0]
	349
	350	zdeltah(ji,jk) = - zfmdt / rhoic ! Melt of layer jk [m, <0]
	351
	352	zdeltah(ji,jk) = MIN( 0._wp , MAX( zdeltah(ji,jk) , - zh_i(ji,jk) ) ) ! Melt of layer jk cannot exceed the layer thickness [m, <0]
	353
	354	zq_su(ji) = MAX( 0._wp , zq_su(ji) - zdeltah(ji,jk) * rhoic * zdE ) ! update available heat
	355
	356	dh_i_surf(ji) = dh_i_surf(ji) + zdeltah(ji,jk) ! Cumulate surface melt
	357
	358	zfmdt = - rhoic * zdeltah(ji,jk) ! Recompute mass flux [kg/m2, >0]
	359
	360	zQm = zfmdt * zEw ! Energy of the melt water sent to the ocean [J/m2, <0]
	361
[4872]	362	! Contribution to salt flux (clem: using sm_i_1d and not s_i_1d(jk) is ok)
	363	sfx_sum_1d(ji) = sfx_sum_1d(ji) - sm_i_1d(ji) * a_i_1d(ji) * zdeltah(ji,jk) * rhoic * r1_rdtice
[4688]	364
	365	! Contribution to heat flux [W.m-2], < 0
[4872]	366	hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_1d(ji) * zEw * r1_rdtice
[4688]	367
	368	! Total heat flux used in this process [W.m-2], > 0
[4872]	369	hfx_sum_1d(ji) = hfx_sum_1d(ji) - zfmdt * a_i_1d(ji) * zdE * r1_rdtice
[4688]	370
	371	! Contribution to mass flux
[4872]	372	wfx_sum_1d(ji) = wfx_sum_1d(ji) - rhoic * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice
[4688]	373
	374	! record which layers have disappeared (for bottom melting)
	375	! => icount=0 : no layer has vanished
	376	! => icount=5 : 5 layers have vanished
[4990]	377	rswitch = MAX( 0._wp , SIGN( 1._wp , - ( zh_i(ji,jk) + zdeltah(ji,jk) ) ) )
	378	icount(ji) = icount(ji) + NINT( rswitch )
[4688]	379	zh_i(ji,jk) = MAX( 0._wp , zh_i(ji,jk) + zdeltah(ji,jk) )
	380
	381	! update heat content (J.m-2) and layer thickness
[4872]	382	qh_i_old(ji,jk) = qh_i_old(ji,jk) + zdeltah(ji,jk) * q_i_1d(ji,jk)
[4688]	383	h_i_old (ji,jk) = h_i_old (ji,jk) + zdeltah(ji,jk)
[825]	384	END DO
[921]	385	END DO
[4688]	386	! update ice thickness
	387	DO ji = kideb, kiut
[4872]	388	ht_i_1d(ji) = MAX( 0._wp , ht_i_1d(ji) + dh_i_surf(ji) )
[4688]	389	END DO
[825]	390
[921]	391	!
	392	!------------------------------------------------------------------------------!
	393	! 4) Basal growth / melt !
	394	!------------------------------------------------------------------------------!
	395	!
[4688]	396	!------------------
	397	! 4.1 Basal growth
	398	!------------------
	399	! Basal growth is driven by heat imbalance at the ice-ocean interface,
	400	! between the inner conductive flux (fc_bo_i), from the open water heat flux
	401	! (fhld) and the turbulent ocean flux (fhtur).
	402	! fc_bo_i is positive downwards. fhtur and fhld are positive to the ice
[825]	403
[4688]	404	! If salinity varies in time, an iterative procedure is required, because
	405	! the involved quantities are inter-dependent.
	406	! Basal growth (dh_i_bott) depends upon new ice specific enthalpy (zEi),
	407	! which depends on forming ice salinity (s_i_new), which depends on dh/dt (dh_i_bott)
	408	! -> need for an iterative procedure, which converges quickly
	409
	410	IF ( num_sal == 2 ) THEN
	411	num_iter_max = 5
	412	ELSE
	413	num_iter_max = 1
[1572]	414	ENDIF
[825]	415
[4688]	416	!clem debug. Just to be sure that enthalpy at nlay_i+1 is null
	417	DO ji = kideb, kiut
[4872]	418	q_i_1d(ji,nlay_i+1) = 0._wp
[4688]	419	END DO
[825]	420
[4688]	421	! Iterative procedure
	422	DO iter = 1, num_iter_max
	423	DO ji = kideb, kiut
	424	IF( zf_tt(ji) < 0._wp ) THEN
[825]	425
[4688]	426	! New bottom ice salinity (Cox & Weeks, JGR88 )
	427	!--- zswi1 if dh/dt < 2.0e-8
	428	!--- zswi12 if 2.0e-8 < dh/dt < 3.6e-7
	429	!--- zswi2 if dh/dt > 3.6e-7
	430	zgrr = MIN( 1.0e-3, MAX ( dh_i_bott(ji) * r1_rdtice , epsi10 ) )
	431	zswi2 = MAX( 0._wp , SIGN( 1._wp , zgrr - 3.6e-7 ) )
	432	zswi12 = MAX( 0._wp , SIGN( 1._wp , zgrr - 2.0e-8 ) ) * ( 1.0 - zswi2 )
	433	zswi1 = 1. - zswi2 * zswi12
	434	zfracs = MIN ( zswi1 * 0.12 + zswi12 * ( 0.8925 + 0.0568 * LOG( 100.0 * zgrr ) ) &
	435	& + zswi2 * 0.26 / ( 0.26 + 0.74 * EXP ( - 724300.0 * zgrr ) ) , 0.5 )
[825]	436
[4688]	437	ii = MOD( npb(ji) - 1, jpi ) + 1 ; ij = ( npb(ji) - 1 ) / jpi + 1
[825]	438
[4688]	439	s_i_new(ji) = zswitch_sal * zfracs * sss_m(ii,ij) & ! New ice salinity
[4872]	440	+ ( 1. - zswitch_sal ) * sm_i_1d(ji)
[4688]	441	! New ice growth
	442	ztmelts = - tmut * s_i_new(ji) + rtt ! New ice melting point (K)
[825]	443
[4872]	444	zt_i_new = zswitch_sal * t_bo_1d(ji) + ( 1. - zswitch_sal) * t_i_1d(ji, nlay_i)
[4688]	445
	446	zEi = cpic * ( zt_i_new - ztmelts ) & ! Specific enthalpy of forming ice (J/kg, <0)
	447	& - lfus * ( 1.0 - ( ztmelts - rtt ) / ( zt_i_new - rtt ) ) &
	448	& + rcp * ( ztmelts-rtt )
	449
[4872]	450	zEw = rcp * ( t_bo_1d(ji) - rt0 ) ! Specific enthalpy of seawater (J/kg, < 0)
[4688]	451
	452	zdE = zEi - zEw ! Specific enthalpy difference (J/kg, <0)
	453
	454	dh_i_bott(ji) = rdt_ice * MAX( 0._wp , zf_tt(ji) / ( zdE * rhoic ) )
	455
[4872]	456	q_i_1d(ji,nlay_i+1) = -zEi * rhoic ! New ice energy of melting (J/m3, >0)
[4688]	457
	458	ENDIF ! fc_bo_i
	459	END DO ! ji
	460	END DO ! iter
	461
	462	! Contribution to Energy and Salt Fluxes
[825]	463	DO ji = kideb, kiut
[4688]	464	IF( zf_tt(ji) < 0._wp ) THEN
	465	! New ice growth
	466
	467	zfmdt = - rhoic * dh_i_bott(ji) ! Mass flux x time step (kg/m2, < 0)
	468
	469	ztmelts = - tmut * s_i_new(ji) + rtt ! New ice melting point (K)
	470
[4872]	471	zt_i_new = zswitch_sal * t_bo_1d(ji) + ( 1. - zswitch_sal) * t_i_1d(ji, nlay_i)
[4688]	472
	473	zEi = cpic * ( zt_i_new - ztmelts ) & ! Specific enthalpy of forming ice (J/kg, <0)
	474	& - lfus * ( 1.0 - ( ztmelts - rtt ) / ( zt_i_new - rtt ) ) &
	475	& + rcp * ( ztmelts-rtt )
	476
[4872]	477	zEw = rcp * ( t_bo_1d(ji) - rt0 ) ! Specific enthalpy of seawater (J/kg, < 0)
[4688]	478
	479	zdE = zEi - zEw ! Specific enthalpy difference (J/kg, <0)
	480
	481	! Contribution to heat flux to the ocean [W.m-2], >0
[4872]	482	hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_1d(ji) * zEw * r1_rdtice
[4688]	483
	484	! Total heat flux used in this process [W.m-2], <0
[4872]	485	hfx_bog_1d(ji) = hfx_bog_1d(ji) - zfmdt * a_i_1d(ji) * zdE * r1_rdtice
[4688]	486
	487	! Contribution to salt flux, <0
[4872]	488	sfx_bog_1d(ji) = sfx_bog_1d(ji) + s_i_new(ji) * a_i_1d(ji) * zfmdt * r1_rdtice
[4688]	489
	490	! Contribution to mass flux, <0
[4872]	491	wfx_bog_1d(ji) = wfx_bog_1d(ji) - rhoic * a_i_1d(ji) * dh_i_bott(ji) * r1_rdtice
[4688]	492
	493	! update heat content (J.m-2) and layer thickness
[4872]	494	qh_i_old(ji,nlay_i+1) = qh_i_old(ji,nlay_i+1) + dh_i_bott(ji) * q_i_1d(ji,nlay_i+1)
[4688]	495	h_i_old (ji,nlay_i+1) = h_i_old (ji,nlay_i+1) + dh_i_bott(ji)
[825]	496	ENDIF
	497	END DO
	498
[4688]	499	!----------------
	500	! 4.2 Basal melt
	501	!----------------
	502	zdeltah(:,:) = 0._wp ! important
[825]	503	DO jk = nlay_i, 1, -1
	504	DO ji = kideb, kiut
[4688]	505	IF( zf_tt(ji) >= 0._wp .AND. jk > icount(ji) ) THEN ! do not calculate where layer has already disappeared from surface melting
[825]	506
[4872]	507	ztmelts = - tmut * s_i_1d(ji,jk) + rtt ! Melting point of layer jk (K)
[825]	508
[4872]	509	IF( t_i_1d(ji,jk) >= ztmelts ) THEN !!! Internal melting
[825]	510
[4872]	511	zEi = - q_i_1d(ji,jk) / rhoic ! Specific enthalpy of melting ice (J/kg, <0)
[825]	512
[4872]	513	!!zEw = rcp * ( t_i_1d(ji,jk) - rtt ) ! Specific enthalpy of meltwater at T = t_i_1d (J/kg, <0)
[4161]	514
[4688]	515	zdE = 0._wp ! Specific enthalpy difference (J/kg, <0)
	516	! set up at 0 since no energy is needed to melt water...(it is already melted)
[4161]	517
[4688]	518	zdeltah (ji,jk) = MIN( 0._wp , - zh_i(ji,jk) ) ! internal melting occurs when the internal temperature is above freezing
	519	! this should normally not happen, but sometimes, heat diffusion leads to this
[825]	520
[4688]	521	dh_i_bott (ji) = dh_i_bott(ji) + zdeltah(ji,jk)
[825]	522
[4688]	523	zfmdt = - zdeltah(ji,jk) * rhoic ! Mass flux x time step > 0
[825]	524
[4688]	525	! Contribution to heat flux to the ocean [W.m-2], <0 (ice enthalpy zEi is "sent" to the ocean)
[4872]	526	hfx_res_1d(ji) = hfx_res_1d(ji) + zfmdt * a_i_1d(ji) * zEi * r1_rdtice
[825]	527
[4872]	528	! Contribution to salt flux (clem: using sm_i_1d and not s_i_1d(jk) is ok)
	529	sfx_res_1d(ji) = sfx_res_1d(ji) - sm_i_1d(ji) * a_i_1d(ji) * zdeltah(ji,jk) * rhoic * r1_rdtice
[4688]	530
	531	! Contribution to mass flux
[4872]	532	wfx_res_1d(ji) = wfx_res_1d(ji) - rhoic * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice
[825]	533
[4688]	534	! update heat content (J.m-2) and layer thickness
[4872]	535	qh_i_old(ji,jk) = qh_i_old(ji,jk) + zdeltah(ji,jk) * q_i_1d(ji,jk)
[4688]	536	h_i_old (ji,jk) = h_i_old (ji,jk) + zdeltah(ji,jk)
[825]	537
[4688]	538	ELSE !!! Basal melting
[825]	539
[4872]	540	zEi = - q_i_1d(ji,jk) / rhoic ! Specific enthalpy of melting ice (J/kg, <0)
[825]	541
[4688]	542	zEw = rcp * ( ztmelts - rtt )! Specific enthalpy of meltwater (J/kg, <0)
[825]	543
[4688]	544	zdE = zEi - zEw ! Specific enthalpy difference (J/kg, <0)
	545
	546	zfmdt = - zq_bo(ji) / zdE ! Mass flux x time step (kg/m2, >0)
	547
	548	zdeltah(ji,jk) = - zfmdt / rhoic ! Gross thickness change
	549
	550	zdeltah(ji,jk) = MIN( 0._wp , MAX( zdeltah(ji,jk), - zh_i(ji,jk) ) ) ! bound thickness change
	551
	552	zq_bo(ji) = MAX( 0._wp , zq_bo(ji) - zdeltah(ji,jk) * rhoic * zdE ) ! update available heat. MAX is necessary for roundup errors
	553
	554	dh_i_bott(ji) = dh_i_bott(ji) + zdeltah(ji,jk) ! Update basal melt
	555
	556	zfmdt = - zdeltah(ji,jk) * rhoic ! Mass flux x time step > 0
	557
	558	zQm = zfmdt * zEw ! Heat exchanged with ocean
	559
	560	! Contribution to heat flux to the ocean [W.m-2], <0
[4872]	561	hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_1d(ji) * zEw * r1_rdtice
[4688]	562
[4872]	563	! Contribution to salt flux (clem: using sm_i_1d and not s_i_1d(jk) is ok)
	564	sfx_bom_1d(ji) = sfx_bom_1d(ji) - sm_i_1d(ji) * a_i_1d(ji) * zdeltah(ji,jk) * rhoic * r1_rdtice
[4688]	565
	566	! Total heat flux used in this process [W.m-2], >0
[4872]	567	hfx_bom_1d(ji) = hfx_bom_1d(ji) - zfmdt * a_i_1d(ji) * zdE * r1_rdtice
[4688]	568
	569	! Contribution to mass flux
[4872]	570	wfx_bom_1d(ji) = wfx_bom_1d(ji) - rhoic * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice
[4688]	571
	572	! update heat content (J.m-2) and layer thickness
[4872]	573	qh_i_old(ji,jk) = qh_i_old(ji,jk) + zdeltah(ji,jk) * q_i_1d(ji,jk)
[4688]	574	h_i_old (ji,jk) = h_i_old (ji,jk) + zdeltah(ji,jk)
	575	ENDIF
	576
	577	ENDIF
	578	END DO ! ji
	579	END DO ! jk
	580
	581	!------------------------------------------------------------------------------!
	582	! Excessive ablation in a 1-category model
	583	! in a 1-category sea ice model, bottom ablation must not exceed hmelt (-0.15)
	584	!------------------------------------------------------------------------------!
	585	! ??? keep ???
	586	! clem bug: I think this should be included above, so we would not have to
	587	! track heat/salt/mass fluxes backwards
	588	! IF( jpl == 1 ) THEN
	589	! DO ji = kideb, kiut
	590	! IF( zf_tt(ji) >= 0._wp ) THEN
	591	! zdh = MAX( hmelt , dh_i_bott(ji) )
	592	! zdvres = zdh - dh_i_bott(ji) ! >=0
	593	! dh_i_bott(ji) = zdh
	594	!
	595	! ! excessive energy is sent to lateral ablation
[4990]	596	! rswitch = MAX( 0._wp, SIGN( 1._wp , 1._wp - at_i_1d(ji) - epsi20 ) )
	597	! zq_1cat(ji) = rswitch * rhoic * lfus * at_i_1d(ji) / MAX( 1._wp - at_i_1d(ji) , epsi20 ) * zdvres ! J.m-2 >=0
[4688]	598	!
	599	! ! correct salt and mass fluxes
[4872]	600	! sfx_bom_1d(ji) = sfx_bom_1d(ji) - sm_i_1d(ji) * a_i_1d(ji) * zdvres * rhoic * r1_rdtice ! this is only a raw approximation
	601	! wfx_bom_1d(ji) = wfx_bom_1d(ji) - rhoic * a_i_1d(ji) * zdvres * r1_rdtice
[4688]	602	! ENDIF
	603	! END DO
	604	! ENDIF
	605
[825]	606	!-------------------------------------------
[4688]	607	! Update temperature, energy
[825]	608	!-------------------------------------------
	609	DO ji = kideb, kiut
[4872]	610	ht_i_1d(ji) = MAX( 0._wp , ht_i_1d(ji) + dh_i_bott(ji) )
[4688]	611	END DO
[825]	612
[4688]	613	!-------------------------------------------
	614	! 5. What to do with remaining energy
	615	!-------------------------------------------
	616	! If heat still available for melting and snow remains, then melt more snow
	617	!-------------------------------------------
	618	zdeltah(:,:) = 0._wp ! important
	619	DO ji = kideb, kiut
	620	zq_rema(ji) = zq_su(ji) + zq_bo(ji)
[4872]	621	! zindh = 1._wp - MAX( 0._wp, SIGN( 1._wp, - ht_s_1d(ji) ) ) ! =1 if snow
[4688]	622	! zindq = 1._wp - MAX( 0._wp, SIGN( 1._wp, - zq_s(ji) + epsi20 ) )
	623	! zdeltah (ji,1) = - zindh * zindq * zq_rema(ji) / MAX( zq_s(ji), epsi20 )
[4872]	624	! zdeltah (ji,1) = MIN( 0._wp , MAX( zdeltah(ji,1) , - ht_s_1d(ji) ) ) ! bound melting
[4688]	625	! zdh_s_mel(ji) = zdh_s_mel(ji) + zdeltah(ji,1)
	626	! dh_s_tot (ji) = dh_s_tot(ji) + zdeltah(ji,1)
[4872]	627	! ht_s_1d (ji) = ht_s_1d(ji) + zdeltah(ji,1)
[4688]	628	!
	629	! zq_rema(ji) = zq_rema(ji) + zdeltah(ji,1) * zq_s(ji) ! update available heat (J.m-2)
	630	! ! heat used to melt snow
[4872]	631	! hfx_snw_1d(ji) = hfx_snw_1d(ji) - zdeltah(ji,1) * a_i_1d(ji) * zq_s(ji) * r1_rdtice ! W.m-2 (>0)
[4688]	632	! ! Contribution to mass flux
[4872]	633	! wfx_snw_1d(ji) = wfx_snw_1d(ji) - rhosn * a_i_1d(ji) * zdeltah(ji,1) * r1_rdtice
[4688]	634	!
	635	ii = MOD( npb(ji) - 1, jpi ) + 1 ; ij = ( npb(ji) - 1 ) / jpi + 1
	636	! Remaining heat flux (W.m-2) is sent to the ocean heat budget
[4872]	637	hfx_out(ii,ij) = hfx_out(ii,ij) + ( zq_1cat(ji) + zq_rema(ji) * a_i_1d(ji) ) * r1_rdtice
[825]	638
[4688]	639	IF( ln_nicep .AND. zq_rema(ji) < 0. .AND. lwp ) WRITE(numout,*) 'ALERTE zq_rema <0 = ', zq_rema(ji)
	640	END DO
	641
[921]	642	!
	643	!------------------------------------------------------------------------------\|
	644	! 6) Snow-Ice formation \|
	645	!------------------------------------------------------------------------------\|
[1572]	646	! When snow load excesses Archimede's limit, snow-ice interface goes down under sea-level,
	647	! flooding of seawater transforms snow into ice dh_snowice is positive for the ice
[825]	648	DO ji = kideb, kiut
[1572]	649	!
[4872]	650	dh_snowice(ji) = MAX( 0._wp , ( rhosn * ht_s_1d(ji) + (rhoic-rau0) * ht_i_1d(ji) ) / ( rhosn+rau0-rhoic ) )
[825]	651
[4872]	652	ht_i_1d(ji) = ht_i_1d(ji) + dh_snowice(ji)
	653	ht_s_1d(ji) = ht_s_1d(ji) - dh_snowice(ji)
[825]	654
[4688]	655	! Salinity of snow ice
	656	ii = MOD( npb(ji) - 1, jpi ) + 1 ; ij = ( npb(ji) - 1 ) / jpi + 1
[4872]	657	zs_snic = zswitch_sal * sss_m(ii,ij) * ( rhoic - rhosn ) / rhoic + ( 1. - zswitch_sal ) * sm_i_1d(ji)
[825]	658
	659	! entrapment during snow ice formation
[4688]	660	! new salinity difference stored (to be used in limthd_ent.F90)
[4161]	661	IF ( num_sal == 2 ) THEN
[4990]	662	rswitch = MAX( 0._wp , SIGN( 1._wp , ht_i_1d(ji) - epsi10 ) )
[4161]	663	! salinity dif due to snow-ice formation
[4990]	664	dsm_i_si_1d(ji) = ( zs_snic - sm_i_1d(ji) ) * dh_snowice(ji) / MAX( ht_i_1d(ji), epsi10 ) * rswitch
[4161]	665	! salinity dif due to bottom growth
[4688]	666	IF ( zf_tt(ji) < 0._wp ) THEN
[4990]	667	dsm_i_se_1d(ji) = ( s_i_new(ji) - sm_i_1d(ji) ) * dh_i_bott(ji) / MAX( ht_i_1d(ji), epsi10 ) * rswitch
[4161]	668	ENDIF
	669	ENDIF
[825]	670
[4688]	671	! Contribution to energy flux to the ocean [J/m2], >0 (if sst<0)
	672	ii = MOD( npb(ji) - 1, jpi ) + 1 ; ij = ( npb(ji) - 1 ) / jpi + 1
	673	zfmdt = ( rhosn - rhoic ) * MAX( dh_snowice(ji), 0._wp ) ! <0
	674	zsstK = sst_m(ii,ij) + rt0
	675	zEw = rcp * ( zsstK - rt0 )
	676	zQm = zfmdt * zEw
	677
	678	! Contribution to heat flux
[4872]	679	hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_1d(ji) * zEw * r1_rdtice
[825]	680
[4688]	681	! Contribution to salt flux
[4872]	682	sfx_sni_1d(ji) = sfx_sni_1d(ji) + sss_m(ii,ij) * a_i_1d(ji) * zfmdt * r1_rdtice
[4688]	683
	684	! Contribution to mass flux
	685	! All snow is thrown in the ocean, and seawater is taken to replace the volume
[4872]	686	wfx_sni_1d(ji) = wfx_sni_1d(ji) - a_i_1d(ji) * dh_snowice(ji) * rhoic * r1_rdtice
	687	wfx_snw_1d(ji) = wfx_snw_1d(ji) + a_i_1d(ji) * dh_snowice(ji) * rhosn * r1_rdtice
[4688]	688
	689	! update heat content (J.m-2) and layer thickness
[4872]	690	qh_i_old(ji,0) = qh_i_old(ji,0) + dh_snowice(ji) * q_s_1d(ji,1) + zfmdt * zEw
[4688]	691	h_i_old (ji,0) = h_i_old (ji,0) + dh_snowice(ji)
	692
	693	! Total ablation (to debug)
[4872]	694	IF( ht_i_1d(ji) <= 0._wp ) a_i_1d(ji) = 0._wp
[825]	695
[4688]	696	END DO !ji
[4161]	697
[2715]	698	!
[4688]	699	!-------------------------------------------
	700	! Update temperature, energy
	701	!-------------------------------------------
	702	!clem bug: we should take snow into account here
	703	DO ji = kideb, kiut
[4990]	704	rswitch = 1.0 - MAX( 0._wp , SIGN( 1._wp , - ht_i_1d(ji) ) )
	705	t_su_1d(ji) = rswitch * t_su_1d(ji) + ( 1.0 - rswitch ) * rtt
[4688]	706	END DO ! ji
	707
	708	DO jk = 1, nlay_s
	709	DO ji = kideb,kiut
	710	! mask enthalpy
[4990]	711	rswitch = MAX( 0._wp , SIGN( 1._wp, - ht_s_1d(ji) ) )
	712	q_s_1d(ji,jk) = ( 1.0 - rswitch ) * q_s_1d(ji,jk)
[4872]	713	! recalculate t_s_1d from q_s_1d
[4990]	714	t_s_1d(ji,jk) = rtt + ( 1._wp - rswitch ) * ( - q_s_1d(ji,jk) / ( rhosn * cpic ) + lfus / cpic )
[4688]	715	END DO
	716	END DO
	717
	718	CALL wrk_dealloc( jpij, zh_s, zqprec, zq_su, zq_bo, zf_tt, zq_1cat, zq_rema )
	719	CALL wrk_dealloc( jpij, zdh_s_mel, zdh_s_pre, zdh_s_sub, zqh_i, zqh_s, zq_s )
[4873]	720	CALL wrk_dealloc( jpij, nlay_i+1, zdeltah, zh_i )
[4688]	721	CALL wrk_dealloc( jpij, icount )
[2715]	722	!
[4161]	723	!
[921]	724	END SUBROUTINE lim_thd_dh
[1572]	725
[825]	726	#else
[1572]	727	!!----------------------------------------------------------------------
	728	!! Default option NO LIM3 sea-ice model
	729	!!----------------------------------------------------------------------
[825]	730	CONTAINS
	731	SUBROUTINE lim_thd_dh ! Empty routine
	732	END SUBROUTINE lim_thd_dh
	733	#endif
[1572]	734
	735	!!======================================================================
[921]	736	END MODULE limthd_dh

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/limthd_dh.F90 @ 4990

Download in other formats: