Context Navigation

sisvat_qsn.F @ 4743

Last change on this file since 4743 was 3900, checked in by evignon, 3 years ago
Commission de la nouvelle interface LMDZ-SISVAT la prochaine commission consistera a supprimer l'ancien repertoire sisvat et a faire un peu de nettoyage.
File size: 27.8 KB

Rev	Line
[3792]	1
	2
	3	subroutine SISVAT_qSn
	4	. (
	5	! #e1. EqSn_0,EqSn_1,EqSn_d
	6	! #m1. ,SIsubl,SImelt,SIrnof
	7	. )
	8
	9	C +------------------------------------------------------------------------+
	10	C \| MAR SISVAT_qSn Fri 29-Jul-2011 MAR \|
	11	C \| SubRoutine SISVAT_qSn updates the Snow Water Content \|
	12	C +------------------------------------------------------------------------+
	13	C \| \|
	14	C \| PARAMETERS: knonv: Total Number of columns = \|
	15	C \| ^^^^^^^^^^ = Total Number of continental grid boxes \|
	16	C \| X Number of Mosaic Cell per grid box \|
	17	C \| \|
	18	C \| INPUT: isnoSV = total Nb of Ice/Snow Layers \|
	19	C \| ^^^^^ \|
	20	C \| \|
	21	C \| INPUT: TaT_SV : SBL Top Temperature [K] \|
	22	C \| ^^^^^ dt__SV : Time Step [s] \|
	23	C \| \|
	24	C \| INPUT / drr_SV : Rain Intensity [kg/m2/s] \|
	25	C \| OUTPUT: dzsnSV : Snow Layer Thickness [m] \|
	26	C \| ^^^^^^ eta_SV : Snow Water Content [m3/m3] \|
	27	C \| ro__SV : Snow/Soil Volumic Mass [kg/m3] \|
	28	C \| TsisSV : Soil/Ice Temperatures (layers -nsol,-nsol+1,..,0)\|
	29	C \| & Snow Temperatures (layers 1,2,...,nsno) [K] \|
	30	C \| \|
	31	C \| OUTPUT: SWS_SV : Surficial Water Status \|
	32	C \| ^^^^^^ \|
	33	C \| EExcsv : Snow Energy in Excess, initial Forcing [J/m2] \|
	34	C \| EqSn_d : Snow Energy in Excess, remaining [J/m2] \|
	35	C \| EqSn_0 : Snow Energy, before Phase Change [J/m2] \|
	36	C \| EqSn_1 : Snow Energy, after Phase Change [J/m2] \|
	37	C \| SIsubl : Snow sublimed/deposed Mass [mm w.e.] \|
	38	C \| SImelt : Snow Melted Mass [mm w.e.] \|
	39	C \| SIrnof : Surficial Water + Run OFF Change [mm w.e.] \|
	40	C \| \|
	41	C \| Internal Variables: \|
	42	C \| ^^^^^^^^^^^^^^^^^^ \|
	43	C \| \|
	44	C \| # OPTIONS: #E0: IO for Verification: Energy Budget \|
	45	C \| # ^^^^^^^ \|
	46	C \| # #su: IO for Verification: Slush Diagnostic \|
	47	C \| \|
	48	C \| \|
	49	C +------------------------------------------------------------------------+
	50
	51
	52
	53
	54	C +--Global Variables
	55	C + ================
	56
	57	use VARphy
	58	use VAR_SV
	59	use VARdSV
	60	use VAR0SV
	61	use VARxSV
	62	use VARySV
[3900]	63	use surface_data, only: is_ok_slush,opt_runoff_ac
[3792]	64
[3900]	65
[3792]	66	IMPLICIT NONE
	67
	68
	69	! Energy Budget
	70	! ~~~~~~~~~~~~~~~~~~~~~~
	71	! #e1 real EqSn_d(knonv) ! Energy in Excess, initial
	72	! #e1 real EqSn_0(knonv) ! Snow Energy, befor Phase Change
	73	! #vm real EqSn01(knonv) ! Snow Energy, after Phase Change
	74	! #vm real EqSn02(knonv) ! Snow Energy, after Phase Change
	75	! .AND. Last Melting
	76	! #e1 real EqSn_1(knonv) ! Snow Energy, after Phase Change
	77	! .AND. Mass Redistr.
	78	! Snow/Ice (Mass) Budget
	79	! ~~~~~~~~~~~~~~~~~~~~~~
	80	! #m1 real SIsubl(knonv) ! Snow Deposed Mass
	81	! #m1 real SImelt(knonv) ! Snow Melted Mass
	82	! #m1 real SIrnof(knonv) ! Local Surficial Water + Run OFF
	83
	84
	85	C +--Internal Variables
	86	C + ==================
	87
	88	integer ikl ,isn !
	89	integer nh ! Non erodible Snow: up.lay.Index
	90	integer LayrOK ! 1 (0) if In(Above) Snow Pack
	91	integer k_face ! 1 (0) if Crystal(no) faceted
	92	integer LastOK ! 1 ==> 1! Snow Layer
	93	integer NOLayr ! 1 Layer Update
	94	integer noSnow(knonv) ! Nb of Layers Updater
	95	integer kSlush ! Slush Switch
	96	real dTSnow ! Temperature [C]
	97	real EExdum(knonv) ! Energy in Excess when no Snow
	98	real OKmelt ! 1 (0) if (no) Melting
	99	real EnMelt ! Energy in excess, for Melting
	100	real SnHLat ! Energy consumed in Melting
	101	real AdEnrg,B_Enrg ! Additional Energy from Vapor
	102	real dzVap0,dzVap1 ! Vaporized Thickness [m]
	103	real dzMelt(knonv) ! Melted Thickness [m]
	104	real rosDry ! Snow volumic Mass if no Water in
	105	real PorVol ! Pore volume
	106	real PClose ! Pore Hole Close OFF Switch
	107	real SGDiam ! Snow Grain Diameter
	108	real SGDmax ! Max. Snow Grain Diameter
	109	real rWater ! Retained Water [kg/m2]
	110	real drrNEW ! New available Water [kg/m2]
	111	real rdzNEW ! Snow Mass [kg/m2]
	112	real rdzsno ! Snow Mass [kg/m2]
	113	real EnFrez ! Energy Release in Freezing
	114	real WaFrez ! Water consumed in Melting
	115	real RapdOK ! 1. ==> Snow melts rapidly
	116	real ThinOK ! 1. ==> Snow Layer is thin
	117	real dzepsi ! Minim. Snow Layer Thickness (!)
	118	real dz_Min ! Minim. Snow Layer Thickness
	119	real z_Melt ! Last (thin) Layer Melting
	120	real rusnew ! Surficial Water Thickness [mm]
	121	real zWater ! Max Slush Water Thickness [mm]
	122	real zSlush ! Slush Water Thickness [mm]
	123	real ro_new ! New Snow/ice Density [kg/m3]
	124	real zc,zt ! Non erod.Snow Thickness[mm w.e.]
	125	real rusnSV0(knonv)
	126	real Tsave
	127
	128	C +--OUTPUT of SISVAT Trace Statistics (see assignation in PHY_SISVAT)
	129	C + ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	130	integer isnnew,isinew,isnUpD,isnitr
	131
	132	! OUTPUT in SISVAT at specified i,j,k,n (see assignation in PHY_SISVAT)
	133	! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	134	! #wx integer iSV_v1,jSV_v1,nSV_v1,kSV_v1,lSV_v1
	135	! #wx common/SISVAT_EV/ iSV_v1,jSV_v1,nSV_v1,kSV_v1,lSV_v1
	136
	137	C +--Energy and Mass Budget
	138	C + ~~~~~~~~~~~~~~~~~~~~~~
	139	! #vm real WqSn_0(knonv) ! Snow Water+Forcing Initial
	140	! #vm real WqSn_1(knonv) ! Snow Water+Forcing, Final
	141	! #vm logical emopen ! IO Switch
	142	! #vm common/Se_qSn_L/emopen !
	143	! #vm integer no_err !
	144	! #vm common/Se_qSn_I/no_err !
	145	! #vm real hourer,timeer !
	146	! #vm common/Se_qSn_R/timeer !
	147
	148	C +--Slush Diagnostic: IO
	149	C + ~~~~~~~~~~~~~~~~~~~~
	150	! #vu logical su_opn ! IO Switch
	151	! #vu common/SI_qSn_L/su_opn !
	152
	153
	154	C +--DATA
	155	C + ====
	156
	157	data dzepsi/0.0001/ ! Minim. Snow Layer Thickness (!)
	158	c #?? data dz_Min/0.005/ ! Minim. Snow Layer Thickness
	159	c ... Warning: Too high for Col de Porte: precludes 1st snow (layer) apparition
	160	data dz_Min/2.5e-3/ ! Minim. Snow Layer Thickness
	161	data SGDmax/0.003/ ! Maxim. Snow Grain Diameter [m]
	162	! (Rowe et al. 1995, JGR p.16268)
	163
	164	C +--Energy Budget (IN)
	165	C + ==================
	166
	167	! #e1 DO ikl=1,knonv
	168	! #e1 EqSn_0(ikl) = 0.
	169	! #e1 END DO
	170	! #e1 DO isn=nsno,1,-1
	171	! #e1 DO ikl=1,knonv
	172	! #e1 EqSn_0(ikl) = EqSn_0(ikl) + ro__SV(ikl,isn) *dzsnSV(ikl,isn)
	173	! #e1. (Cn_dSV (TsisSV(ikl,isn) -TfSnow )
	174	! #e1. -Lf_H2O *(1. -eta_SV(ikl,isn)))
	175	! #e1 END DO
	176	! #e1 END DO
	177
	178
	179	C +--Water Budget (IN)
	180	C + ==================
	181
	182	! #vm DO ikl=1,knonv
	183	! #vm WqSn_0(ikl) = drr_SV(ikl) * dt__SV
	184	! #vm. +rusnSV(ikl)
	185	! #vm END DO
	186	! #vm DO isn=nsno,1,-1
	187	! #vm DO ikl=1,knonv
	188	! #vm WqSn_0(ikl) = WqSn_0(ikl) + ro__SV(ikl,isn) *dzsnSV(ikl,isn)
	189	! #vm END DO
	190	! #vm END DO
	191
	192
	193	C +--Snow Melt Budget
	194	C + ================
	195
	196	! #m1 DO ikl=1,knonv
	197	! #m1 SImelt(ikl) = 0.
	198	! #m1 SIrnof(ikl) = rusnSV(ikl) + RnofSV(ikl) * dt__SV
	199	! #m1 END DO
	200
	201
	202	C +--Initialization
	203	C + ==============
	204
	205	DO ikl=1,knonv
	206	noSnow(ikl) = 0 ! Nb of Layers Updater
	207	ispiSV(ikl) = 0 ! Pore Hole Close OFF Index
	208	! (assumed to be the Top of
	209	! the surimposed Ice Layer)
	210	zn5_SV(ikl) = 0.
	211	rusnSV0(ikl) = 0.
	212
	213	END DO
	214
	215
	216	C +--Melting/Freezing Energy
	217	C + =======================
	218
	219	C +...REMARK: Snow liquid Water Temperature assumed = TfSnow
	220	C + ^^^^^^
	221	DO ikl=1,knonv
	222	EExdum(ikl) = drr_SV(ikl) * C__Wat *(TaT_SV(ikl)-TfSnow)
	223	. * dt__SV
	224	EExcsv(ikl) = EExdum(ikl) * min(1,isnoSV(ikl)) ! Snow exists
	225	EExdum(ikl) = EExdum(ikl) - EExcsv(ikl) !
	226	! #e1 EqSn_d(ikl) = EExcsv(ikl) !
	227	END DO
	228
	229
	230	C +--Surficial Water Status
	231	C + ----------------------
	232
	233	DO ikl=1,knonv
	234	SWS_SV(ikl) = max(zero,sign(unun,TfSnow
	235	. -TsisSV(ikl,isnoSV(ikl))))
	236	END DO
	237
	238	DO ikl=1,knonv
[3900]	239
	240	DO isn=min(nsno,isnoSV(ikl)+1),1,-1
[3792]	241	! EV DO isn=nsno,1,-1
	242	C +--Energy, store Previous Content
	243	C + ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
	244	dTSnow = TsisSV(ikl,isn) - TfSnow
	245	EExcsv(ikl) = EExcsv(ikl)
	246	. + ro__SV(ikl,isn) * Cn_dSV * dTSnow
	247	. * dzsnSV(ikl,isn)
	248	TsisSV(ikl,isn) = TfSnow
	249
	250	C +--Water, store Previous Content
	251	C + ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
	252	drr_SV(ikl) = drr_SV(ikl)
	253	. + ro__SV(ikl,isn) * eta_SV(ikl,isn)
	254	. * dzsnSV(ikl,isn)
	255	. / dt__SV
	256	ro__SV(ikl,isn) =
	257	. ro__SV(ikl,isn) *(1. - eta_SV(ikl,isn))
	258	eta_SV(ikl,isn) = 0.
	259
	260
	261	C +--Melting if EExcsv > 0
	262	C + ======================
	263
	264	EnMelt = max(zero, EExcsv(ikl) )
	265
	266	C +--Energy Consumption
	267	C + ^^^^^^^^^^^^^^^^^^
	268	SnHLat = ro__SV(ikl,isn) * Lf_H2O
	269	dzMelt(ikl) = EnMelt / max(SnHLat, epsi )
	270	noSnow(ikl) = noSnow(ikl)
	271	. + max(zero ,sign(unun,dzMelt(ikl) !
	272	. -dzsnSV(ikl ,isn))) ! 1 if full Melt
	273	. *min(1 , max(0 ,1+isnoSV(ikl)-isn)) ! 1 in the Pack
	274	dzMelt(ikl) =
	275	. min(dzsnSV(ikl, isn),dzMelt(ikl))
	276	dzsnSV(ikl,isn) =
	277	. dzsnSV(ikl,isn) -dzMelt(ikl)
	278	zn5_SV(ikl) = zn5_SV(ikl) +dzMelt(ikl)
	279	EExcsv(ikl) = EExcsv(ikl) -dzMelt(ikl)*SnHLat
	280	wem_SV(ikl) = wem_SV(ikl) -dzMelt(ikl)*ro__SV(ikl,isn)
	281
	282	C +--Water Production
	283	C + ^^^^^^^^^^^^^^^^^
	284	drr_SV(ikl) = drr_SV(ikl)
	285	. + ro__SV(ikl,isn) * dzMelt(ikl)/dt__SV
	286	! #m1 SImelt(ikl) = SImelt(ikl)
	287	! #m1. + ro__SV(ikl,isn) * dzMelt(ikl)
	288	OKmelt =max(zero,sign(unun,drr_SV(ikl)-epsi))
	289
	290	C +--Snow History
	291	C + ^^^^^^^^^^^^
	292	k_face = min( istoSV(ikl,isn),istdSV(1)) ! = 1 if
	293	. *max(0,2-istoSV(ikl,isn) ) ! faceted
	294	istoSV(ikl,isn) = !
	295	. (1.-OKmelt) * istoSV(ikl,isn) !
	296	. + OKmelt ((1-k_face) istdSV(2) !
	297	. + k_face * istdSV(3) ) !
	298
	299
	300	C +--Freezing if EExcsv < 0
	301	C + ======================
	302
	303	rdzsno = ro__SV(ikl,isn) * dzsnSV(ikl ,isn)
	304	LayrOK = min( 1, max(0 , isnoSV(ikl)-isn+1))
	305	EnFrez = min(zero, EExcsv(ikl))
	306	WaFrez = -( EnFrez * LayrOK / Lf_H2O)
	307	drrNEW = max(zero,drr_SV(ikl) - WaFrez / dt__SV)
	308	WaFrez = ( drr_SV(ikl) - drrNEW)* dt__SV
	309	drr_SV(ikl) = drrNEW
	310	EExcsv(ikl) = EExcsv(ikl) + WaFrez * Lf_H2O
	311	EnFrez = min(zero,EExcsv(ikl)) * LayrOK
	312	rdzNEW = WaFrez + rdzsno
	313	ro__SV(ikl,isn) = rdzNEW /max(epsi, dzsnSV(ikl,isn))
	314	TsisSV(ikl,isn) = TfSnow
	315	. + EnFrez /(Cn_dSV *max(epsi, rdzNEW) )
	316	EExcsv(ikl) = EExcsv(ikl) - EnFrez
	317	wer_SV(ikl) = WaFrez
	318	. + wer_SV(ikl)
	319
	320
	321
	322	C +--Snow Water Content
	323	C + ==================
	324
	325	C +--Percolation Velocity
	326	C + ^^^^^^^^^^^^^^^^^^^^
	327	c #PW SGDiam = 1.6d-4
	328	c #PW. + 1.1d-13 (ro__SV(ikl,isn)ro__SV(ikl,isn)
	329	c #PW. ro__SV(ikl,isn)ro__SV(ikl,isn))
	330
	331	C +--Pore Volume [-]
	332	C + ^^^^^^^^^^^^^^^^^
	333	rosDry =(1. - eta_SV(ikl,isn))* ro__SV(ikl,isn) !
	334	PorVol = 1. - rosDry / ro_Ice !
	335	PorVol = max(PorVol , zero ) !
	336
	337	C +--Water Retention
	338	C + ^^^^^^^^^^^^^^^^
	339	rWater = ws0dSV * PorVol * ro_Wat * dzsnSV(ikl,isn)
	340	drrNEW = max(zero,drr_SV(ikl) - rWater /dt__SV)
	341	rWater = ( drr_SV(ikl) - drrNEW)*dt__SV
	342	drr_SV(ikl) = drrNEW
	343	rdzNEW = rWater
	344	. + rosDry * dzsnSV(ikl,isn)
	345	eta_SV(ikl,isn) = rWater / max(epsi,rdzNEW)
	346	ro__SV(ikl,isn) = rdzNEW / max(epsi,dzsnSV(ikl,isn))
	347
	348	C +--Pore Hole Close OFF
	349	C + ^^^^^^^^^^^^^^^^^^^
	350	PClose = max(zero,
	351	. sign(unun,ro__SV(ikl,isn)
	352	. -roCdSV ))
	353	ispiSV(ikl) = ispiSV(ikl) *(1.-PClose)
	354	. + max(ispiSV(ikl),isn) * Pclose
	355	PClose = max(0 , ! Water under SuPer.Ice
	356	. min (1 ,ispiSV(ikl) ! contributes to
	357	. -isn )) ! Surficial Water
	358
	359	cXF
	360	if(ro__SV(ikl,isn) >= roCdSV.and.ro__SV(ikl,1)<900)
	361	. PClose = min(0.50,PClose *
	362	. (1.-(ro_ice-ro__SV(ikl,isn))/(ro_ice-roCdSV)))
	363
	364	PClose = max(0.,min(1.,PClose))
	365
	366	if(isn==1) then
	367	PClose = 1
	368	ispiSV(ikl)= max(ispiSV(ikl),1)
	369	endif
	370
	371	if(drr_SV(ikl) >0 .and.TsisSV(ikl,isn)>273.14) then
	372	if((ro__SV(ikl,isn)>900.and.ro__SV(ikl,isn)<920).or.
	373	. ro__SV(ikl,isn)>950) then
	374	dzsnSV(ikl,isn) = dzsnSV(ikl,isn)*ro__SV(ikl,isn)/ro_ice
	375	ro__SV(ikl,isn) = ro_ice
	376	PClose = 1
	377	endif
	378	endif
	379
	380	c if (isn>1.and.isn<nsno .and.
	381	c . ro__SV(ikl,isn-1)>900 .and.
	382	c . ro__SV(ikl,isn) >roCdSV .and.
	383	c . ro__SV(ikl,isn) <900 .and.
	384	c . TsisSV(ikl,isn) >273.14 .and.
	385	c . TsisSV(ikl,isn+1)<273.15 .and.
	386	c . drr_SV(ikl) >0) then
	387	c TsisSV(ikl,isn)=273.14
	388	c PClose = 1
	389	c endif
	390
	391	cXF
	392	rusnSV(ikl) = rusnSV(ikl)
	393	. + drr_SV(ikl) dt__SV PClose
	394	rusnSV0(ikl)= rusnSV0(ikl)
	395	. + drr_SV(ikl) dt__SV PClose
	396	drr_SV(ikl) = drr_SV(ikl) *(1.-PClose)
	397
	398	END DO
	399
	400	END DO
	401
	402
	403	C +--Remove Zero-Thickness Layers
	404	C + ============================
	405
	406	1000 CONTINUE
	407	isnitr = 0
	408	DO ikl=1,knonv
	409	isnUpD = 0
	410	isinew = 0
	411	cXF
	412
	413
	414	DO isn=1,min(nsno-1,isnoSV(ikl))
	415	isnnew =(unun-max(zero ,sign(unun,dzsnSV(ikl,isn)-dzepsi)))
	416	. * max(0 , min(1 ,isnoSV(ikl) +1 -isn ))
	417	isnUpD = max(isnUpD, isnnew)
	418	isnitr = max(isnitr, isnnew)
	419	isinew = isnisnUpD max(0, 1-isinew) ! LowerMost 0-Layer
	420	. +isinew ! Index
	421	dzsnSV(ikl,isn) = dzsnSV(ikl,isn+isnnew)
	422	ro__SV(ikl,isn) = ro__SV(ikl,isn+isnnew)
	423	TsisSV(ikl,isn) = TsisSV(ikl,isn+isnnew)
	424	eta_SV(ikl,isn) = eta_SV(ikl,isn+isnnew)
	425	G1snSV(ikl,isn) = G1snSV(ikl,isn+isnnew)
	426	G2snSV(ikl,isn) = G2snSV(ikl,isn+isnnew)
	427	dzsnSV(ikl,isn+isnnew) =(1-isnnew)*dzsnSV(ikl,isn+isnnew)
	428	ro__SV(ikl,isn+isnnew) =(1-isnnew)*ro__SV(ikl,isn+isnnew)
	429	eta_SV(ikl,isn+isnnew) =(1-isnnew)*eta_SV(ikl,isn+isnnew)
	430	G1snSV(ikl,isn+isnnew) =(1-isnnew)*G1snSV(ikl,isn+isnnew)
	431	G2snSV(ikl,isn+isnnew) =(1-isnnew)*G2snSV(ikl,isn+isnnew)
	432
	433	END DO
	434	isnoSV(ikl) = isnoSV(ikl)-isnUpD ! Nb of Snow Layer
	435	ispiSV(ikl) = ispiSV(ikl) ! Nb of SuperI Layer
	436	. -isnUpD *max(0,min(ispiSV(ikl)-isinew,1)) ! Update if I=0
	437
	438	END DO
	439
	440	IF (isnitr.GT.0) GO TO 1000
	441
	442
	443	C +--New upper Limit of the non erodible Snow (istoSV .GT. 1)
	444	C + ========================================
	445
	446	DO ikl=1,knonv
	447	nh = 0
	448	cXF
	449	DO isn= isnoSV(ikl),1,-1
	450	nh = nh + isn* min(istoSV(ikl,isn)-1,1)*max(0,1-nh)
	451	ENDDO
	452	zc = 0.
	453	zt = 0.
	454	cXF
	455	DO isn=1,isnoSV(ikl)
	456	zc = zc + dzsnSV(ikl,isn) *ro__SV(ikl,isn)
	457	. * max(0,min(1,nh+1-isn))
	458	zt = zt + dzsnSV(ikl,isn) *ro__SV(ikl,isn)
	459	END DO
	460	zWE_SV(ikl) = zt
	461	zWEcSV(ikl) = min(zWEcSV(ikl),zt)
	462	zWEcSV(ikl) = max(zWEcSV(ikl),zc)
	463	END DO
	464
	465
	466	C +--Energy Budget (OUT)
	467	C + ===================
	468
	469	! #vm DO ikl=1,knonv
	470	! #vm EqSn01(ikl) =-EqSn_0(ikl)
	471	! #vm. -EExcsv(ikl)
	472	! #vm END DO
	473	! #vm DO isn=nsno,1,-1
	474	! #vm DO ikl=1,knonv
	475	! #vm EqSn01(ikl) = EqSn01(ikl) + ro__SV(ikl,isn) *dzsnSV(ikl,isn)
	476	! #vm. (Cn_dSV (TsisSV(ikl,isn) -TfSnow )
	477	! #vm. -Lf_H2O *(1. -eta_SV(ikl,isn)))
	478	! #vm END DO
	479	! #vm END DO
	480
	481
	482	C +--"Negative Heat" from supercooled rain
	483	C + ------------------------------------
	484
	485	DO ikl=1,knonv
	486	EExcsv(ikl) = EExcsv(ikl) + EExdum(ikl)
	487
	488
	489	C +--Surficial Water Run OFF
	490	C + -----------------------
	491
	492	rusnew = rusnSV(ikl) * SWf_SV(ikl)
	493
[3900]	494	if(isnoSV(ikl)<=1 .OR. opt_runoff_ac) rusnew = 0.
[3792]	495	!if(ivgtSV(ikl)>=1) rusnew = 0.
	496
	497	c #EU rusnew = 0.
[3900]	498	c #AC rusnew = 0.
	499
[3792]	500	RnofSV(ikl) = RnofSV(ikl)
	501	. +(rusnSV(ikl) - rusnew ) / dt__SV
	502	RuofSV(ikl,1) = RuofSV(ikl,1)
	503	. +(rusnSV(ikl) - rusnew ) / dt__SV
	504	RuofSV(ikl,4) = RuofSV(ikl,4)
	505	. +(rusnSV0(ikl) ) / dt__SV
	506	rusnSV(ikl) = rusnew
	507	END DO
	508
	509
	510	C +--Percolation down the Continental Ice Pack
	511	C + -----------------------------------------
	512
	513	DO ikl=1,knonv
	514	drr_SV(ikl) = drr_SV(ikl) + rusnSV(ikl)
	515	. * (1-min(1,ispiSV(ikl)))/ dt__SV
	516	rusnSV(ikl) = rusnSV(ikl)
	517	. * min(1,ispiSV(ikl))
	518	END DO
	519
	520	cXF removal of too thin snowlayers if TT> 275.15 + bug if TT>> 273.15
	521	DO ikl=1,knonv
	522	zt=0.
	523	DO isn=1,isnoSV(ikl)
	524	zt=zt+dzsnSV(ikl,isn)
	525	ENDDO
	526
	527	if(zt<0.005+(TaT_SV(ikl)-TfSnow)/1000..and.
	528	. isnoSV(ikl) >0 .and.
	529	. TaT_SV(ikl) >=TfSnow .and.
	530	. istoSV(ikl,isnoSV(ikl)) >1 ) then
	531	DO isn=1,isnoSV(ikl)
	532	drr_SV(ikl) = drr_SV(ikl)
	533	. + dzsnSV(ikl,isn)*ro__SV(ikl,isn) /dt__SV
	534	dzsnSV(ikl,isn)= 0.
	535
	536	ENDDO
	537	isnoSV(ikl) = 0
	538	endif
	539	ENDDO
	540
[3900]	541	C +--Slush Formation (Activated. CAUTION: ADD RunOff Possibility before Activation)
[3792]	542	C + --------------- ^^^^^^^ ^^^
	543
[3900]	544	IF (is_ok_slush) THEN
	545
	546	DO ikl=1,knonv
	547	DO isn=1,isnoSV(ikl)
	548	kSlush = min(1,max(0,isn+1-ispiSV(ikl))) ! Slush Switch
[3792]	549
	550	C +--Available Additional Pore Volume [-]
	551	C + ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
[3900]	552	PorVol = 1. - ro__SV(ikl,isn) ! [--]
	553	. *(1. - eta_SV(ikl,isn))/ ro_Ice !
	554	. - eta_SV(ikl,isn) !
	555	. *ro__SV(ikl,isn) / ro_Wat !
	556	PorVol = max(PorVol , zero ) !
	557	zWater = dzsnSV(ikl,isn) * PorVol * 1000. ! [mm] OR [kg/m2]
	558	. * (1. -SWS_SV(ikl) ! 0 <=> freezing
	559	. *(1 -min(1,iabs(isn-isnoSV(ikl))))) ! 1 <=> isn=isnoSV
	560	zSlush = min(rusnSV(ikl) , zWater) ! [mm] OR [kg/m2]
	561	ro_new =(dzsnSV(ikl,isn) * ro__SV(ikl,isn) !
	562	. +zSlush ) !
	563	. / max(dzsnSV(ikl,isn) , epsi ) !
	564	if(ro_new<ro_Ice+20) then ! MAX 940kg/m3 !
	565	rusnSV(ikl) = rusnSV(ikl) - zSlush ! [mm] OR [kg/m2]
	566	RuofSV(ikl,4)= max(0.,RuofSV(ikl,4) - zSlush/dt__SV)
	567	eta_SV(ikl,isn) =(ro_new - ro__SV(ikl,isn) !
	568	. *(1. - eta_SV(ikl,isn))) !
	569	. / max (ro_new , epsi ) !
	570	ro__SV(ikl,isn) = ro_new !
	571	endif
	572	END DO
	573	END DO
	574	END IF
[3792]	575
	576	C +--Impact of the Sublimation/Deposition on the Surface Mass Balance
	577	C + ================================================================
	578
	579	DO ikl=1,knonv
	580	isn = isnoSV(ikl)
	581	dzVap0 = dt__SV
	582	. * HLs_sv(ikl) * min(isn , 1 )
	583	. /(Lx_H2O(ikl) * max(ro__SV(ikl,isn) , epsi))
	584	NOLayr=min(zero,sign(unun,dzsnSV(ikl,isn) + dzVap0))
	585	dzVap1=min(zero, dzsnSV(ikl,isn) + dzVap0)
	586
	587
	588	C +--Additional Energy
	589	C + -----------------
	590
	591	c #VH AdEnrg = dzVap0 * ro__SV(ikl,isnoSV(ikl)) ! Water Vapor
	592	c #VH. C__Wat (TsisSV(ikl,isnoSV(ikl)) -TfSnow) ! Sensible Heat
	593
	594	c #aH B_Enrg =(Cn_dSV *(TsisSV(ikl,isn) -TfSnow )
	595	c #aH. -Lf_H2O *(1. -eta_SV(ikl,isn)))
	596	c #aH. /(1. + dzVap0 /max(epsi,dzsnSV(ikl,isn)))
	597	c #aH eta_SV(ikl,isn) =
	598	c #aH. max(zero,unun +(B_Enrg
	599	c #aH. -(TsisSV(ikl,isn) -TfSnow)*Cn_dSV)
	600	c #aH. /Lf_H2O )
	601	c #aH TsisSV(ikl,isn) = ( B_Enrg
	602	c #aH. +(1. -eta_SV(ikl,isn))
	603	c #aH. *Lf_H2O )
	604	c #aH. / Cn_dSV
	605	c #aH. + TfSnow
	606
	607	! #e1 STOP "PLEASE add Energy (#aH) from deposition/sublimation"
	608
	609
	610	C +--Update of the upper Snow layer Thickness
	611	C + ----------------------------------------
	612
	613	dzsnSV(ikl,isn) =
	614	. max(zero, dzsnSV(ikl,isnoSV(ikl)) + dzVap0)
	615	isnoSV(ikl) = isnoSV(ikl) + NOLayr
	616	isn = isnoSV(ikl)
	617	dzsnSV(ikl,isn) = dzsnSV(ikl,isn) + dzVap1
	618	wes_SV(ikl) = ro__SV(ikl,isn) * dzVap0
	619
	620	END DO
	621
	622
	623	C +--Energy Budget (OUT)
	624	C + ===================
	625
	626	! #vm DO ikl=1,knonv
	627	! #vm EqSn02(ikl) =-EqSn_0(ikl)
	628	! #vm. -EExcsv(ikl)
	629	! #vm END DO
	630	! #vm DO isn=nsno,1,-1
	631	! #vm DO ikl=1,knonv
	632	! #vm EqSn02(ikl) = EqSn01(ikl) + ro__SV(ikl,isn) *dzsnSV(ikl,isn)
	633	! #vm. (Cn_dSV (TsisSV(ikl,isn) -TfSnow )
	634	! #vm. -Lf_H2O *(1. -eta_SV(ikl,isn)))
	635	! #vm END DO
	636	! #vm END DO
	637
	638
	639	C +--Snow/I Budget
	640	C + -------------
	641
	642	! #m1 DO ikl=1,knonv
	643	! #m1 SIsubl(ikl) = dt__SVHLs_sv(ikl)min(isnoSV(ikl),1)
	644	! #m1. /Lx_H2O(ikl)
	645	! #m1 SIrnof(ikl) = rusnSV(ikl) + RnofSV(ikl) * dt__SV
	646	! #m1. - SIrnof(ikl)
	647	! #m1 END DO
	648
	649
	650	C +--Anticipated Disappearance of a rapidly Melting too thin Last Snow Layer
	651	C + =======================================================================
	652
	653	DO ikl=1,knonv
	654	LastOK = min(1 , max(0 ,iiceSV(ikl)-isnoSV(ikl)+2)
	655	. *min(1 ,isnoSV(ikl)-iiceSV(ikl))
	656	. +min(1 ,isnoSV(ikl)) )
	657	RapdOK = max(zero,sign(unun,dzMelt(ikl)-epsi ))
	658	ThinOK = max(zero,sign(unun,dz_Min -dzsnSV(ikl,1)))
	659	z_Melt = LastOK RapdOKThinOK
	660	noSnow(ikl) = noSnow(ikl) + z_Melt
	661	z_Melt = z_Melt *dzsnSV(ikl,1)
	662	dzsnSV(ikl,1) = dzsnSV(ikl,1) - z_Melt
	663	EExcsv(ikl) = EExcsv(ikl) - z_Melt *ro__SV(ikl,1)
	664	. (1. -eta_SV(ikl,1))Lf_H2O
	665
	666	C +--Water Production
	667	C + ^^^^^^^^^^^^^^^^^
	668	drr_SV(ikl) = drr_SV(ikl)
	669	. + ro__SV(ikl,1) * z_Melt /dt__SV
	670	END DO
	671
	672
	673	C +--Update Nb of Layers
	674	C + ===================
	675
	676	DO ikl=1,knonv
	677	isnoSV(ikl) = isnoSV(ikl)
	678	. * min(1,iabs(isnoSV(ikl)-noSnow(ikl)))
	679	END DO
	680
	681
	682	! Energy Budget (OUT)
	683	! ===================
	684
	685	! #e1 DO ikl=1,knonv
	686	! #e1 EqSn_1(ikl) = 0.
	687	! #e1 END DO
	688	! #e1 DO isn=nsno,1,-1
	689	! #e1 DO ikl=1,knonv
	690	! #e1 EqSn_1(ikl) = EqSn_1(ikl) + ro__SV(ikl,isn) *dzsnSV(ikl,isn)
	691	! #e1. (Cn_dSV (TsisSV(ikl,isn) -TfSnow )
	692	! #e1. -Lf_H2O *(1. -eta_SV(ikl,isn)))
	693	! #e1 END DO
	694	! #e1 END DO
	695
	696
	697	C +--Water Budget (OUT)
	698	C + ===================
	699
	700	! #vm DO ikl=1,knonv
	701	! #vm WqSn_0(ikl) = WqSn_0(ikl)
	702	! #vm. + HLs_sv(ikl) * dt__SV
	703	! #vm. *min(isnoSV(ikl),1) / Lx_H2O(ikl)
	704	! #vm WqSn_1(ikl) = drr_SV(ikl) * dt__SV
	705	! #vm. + rusnSV(ikl)
	706	! #vm. + RnofSV(ikl) * dt__SV
	707	! #vm END DO
	708	! #vm DO isn=nsno,1,-1
	709	! #vm DO ikl=1,knonv
	710	! #vm WqSn_1(ikl) = WqSn_1(ikl)
	711	! #vm. + ro__SV(ikl,isn)* dzsnSV(ikl,isn)
	712	! #vm END DO
	713	! #vm END DO
	714
	715
	716	return
	717	end

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

Download in other formats:

Original Format