Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Normal
Revision Log

sisvat_gsn.f90 @ 5428

Last change on this file since 5428 was 5246, checked in by abarral, 2 months ago
Convert fixed-form to free-form sources .F -> .{f,F}90 (WIP: some .F remain, will be handled in subsequent commits)
File size: 42.0 KB

Rev	Line
[3792]	1
[5246]	2	subroutine SISVAT_GSn
[3792]	3
[5246]	4	! +------------------------------------------------------------------------+
	5	! \| MAR SISVAT_GSn 20-09-2003 MAR \|
	6	! \| SubRoutine SISVAT_GSn simulates SNOW Metamorphism \|
	7	! +------------------------------------------------------------------------+
	8	! \| \|
	9	! \| PARAMETERS: knonv: Total Number of columns = \|
	10	! \| ^^^^^^^^^^ = Total Number of continental grid boxes \|
	11	! \| X Number of Mosaic Cell per grid box \|
	12	! \| \|
	13	! \| INPUT / isnoSV = total Nb of Ice/Snow Layers \|
	14	! \| OUTPUT: iiceSV = total Nb of Ice Layers \|
	15	! \| ^^^^^^ istoSV = 0,...,5 : Snow History (see istdSV data) \|
	16	! \| \|
	17	! \| INPUT: TsisSV : Soil/Ice Temperatures (layers -nsol,-nsol+1,..,0)\|
	18	! \| ^^^^^ & Snow Temperatures (layers 1,2,...,nsno) [K] \|
	19	! \| ro__SV : Soil/Snow Volumic Mass [kg/m3] \|
	20	! \| eta_SV : Soil/Snow Water Content [m3/m3] \|
	21	! \| slopSV : Surface Slope [-] \|
	22	! \| dzsnSV : Snow Layer Thickness [m] \|
	23	! \| dt__SV2 : Time Step [s] \|
	24	! \| \|
	25	! \| INPUT / G1snSV : Dendricity (<0) or Sphericity (>0) of Snow Layer \|
	26	! \| OUTPUT: G2snSV : Sphericity (>0) or Size of Snow Layer \|
	27	! \| ^^^^^^ \|
	28	! \| \|
	29	! \| Formalisme adopte pour la Representation des Grains: \|
	30	! \| Formalism for the Representation of Grains: \|
	31	! \| ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ \|
	32	! \| \|
	33	! \| 1 - -1 Neige Fraiche \|
	34	! \| / \ \| ------------- \|
	35	! \| / \ \| Dendricite decrite par Dendricite \|
	36	! \| / \ \| Dendricity et Sphericite \|
	37	! \| / \ \| \|
	38	! \| 2---------3 - 0 described by Dendricity \|
	39	! \| and Sphericity \|
	40	! \| \|---------\| \|
	41	! \| 0 1 \|
	42	! \| Sphericite \|
	43	! \| Sphericity \|
	44	! \| \|
	45	! \| 4---------5 - \|
	46	! \| \| \| \| \|
	47	! \| \| \| \| Diametre (1/10eme de mm) (ou Taille) \|
	48	! \| \| \| \| Diameter (1/10th of mm) (or Size ) \|
	49	! \| \| \| \| \|
	50	! \| \| \| \| Neige non dendritique \|
	51	! \| 6---------7 - --------------------- \|
	52	! \| decrite par Sphericite \|
	53	! \| et Taille \|
	54	! \| described by Sphericity \|
	55	! \| and Size \|
	56	! \| \|
	57	! \| Les Variables du Modele: \|
	58	! \| Model Variables: \|
	59	! \| ^^^^^^^^^^^^^^^^^^^^^^^^ \|
	60	! \| Cas Dendritique Cas non Dendritique \|
	61	! \| \|
	62	! \| G1snSV : Dendricite G1snSV : Sphericite \|
	63	! \| G2snSV : Sphericite G2snSV : Taille (1/10e mm) \|
	64	! \| Size \|
	65	! \| \|
	66	! \| Cas Dendritique/ Dendritic Case \|
	67	! \| ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ \|
	68	! \| Dendricite(Dendricity) G1snSV \|
	69	! \| varie de -G1_dSV (-99 par defaut / etoile) a 0 \|
	70	! \| division par -G1_dSV pour obtenir des valeurs entre 1 et 0 \|
	71	! \| varies from -G1_dSV (default -99 / fresh snow) to 0 \|
	72	! \| division by -G1_dSV to obtain values between 1 and 0 \|
	73	! \| \|
	74	! \| Sphericite(Sphericity) G2snSV \|
	75	! \| varie de 0 (cas completement anguleux) \|
	76	! \| a G1_dSV (99 par defaut, cas spherique) \|
	77	! \| division par G1_dSV pour obtenir des valeurs entre 0 et 1 \|
	78	! \| varies from 0 (full faceted) to G1_dSV \|
	79	! \| \|
	80	! \| Cas non Dendritique / non Dendritic Case \|
	81	! \| ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ \|
	82	! \| Sphericite(Sphericity) G1snSV \|
	83	! \| varie de 0 (cas completement anguleux) \|
	84	! \| a G1_dSV (99 par defaut, cas spherique) \|
	85	! \| division par G1_dSV pour obtenir des valeurs entre 0 et 1 \|
	86	! \| varies from 0 (full faceted) to G1_dSV \|
	87	! \| \|
	88	! \| Taille (Size) G2snSV \|
	89	! \| superieure a ADSdSV (.4 mm) et ne fait que croitre \|
	90	! \| greater than ADSdSV (.4 mm) always increases \|
	91	! \| \|
	92	! \| Exemples: Points caracteristiques des Figures ci-dessus \|
	93	! \| ^^^^^^^^^ \|
	94	! \| \|
	95	! \| G1snSV G2snSV dendricite sphericite taille \|
	96	! \| dendricity sphericity size \|
	97	! \| ------------------------------------------------------------------ \|
	98	! \| [1/10 mm] \|
	99	! \| 1 -G1_dSV sph3SN 1 0.5 \|
	100	! \| 2 0 0 0 0 \|
	101	! \| 3 0 G1_dSV 0 1 \|
	102	! \| 4 0 ADSdSV 0 4. \|
	103	! \| 5 G1_dSV ADSdSV-vsphe1 1 3. \|
	104	! \| 6 0 -- 0 -- \|
	105	! \| 7 G1_dSV -- 1 -- \|
	106	! \| \|
	107	! \| par defaut: G1_dSV=99. \|
	108	! \| sph3SN=50. \|
	109	! \| ADSdSV= 4. \|
	110	! \| vsphe1=1. \|
	111	! \| \|
	112	! \| Methode: \|
	113	! \| ^^^^^^^^ \|
	114	! \| 1. Evolution Types de Grains selon Lois de Brun et al. (1992): \|
	115	! \| Grain metamorphism according to Brun et al. (1992): \|
	116	! \| Plusieurs Cas sont a distiguer / the different Cases are: \|
	117	! \| 1.1 Metamorphose Neige humide / wet Snow \|
	118	! \| 1.2 Metamorphose Neige seche / dry Snow \|
	119	! \| 1.2.1 Gradient faible / low Temperature Gradient \|
	120	! \| 1.2.2 Gradient moyen / moderate Temperature Gradient \|
	121	! \| 1.2.3 Gradient fort / high Temperature Gradient \|
	122	! \| Dans chaque Cas on separe Neige Dendritique et non Dendritique \|
	123	! \| le Passage Dendritique -> non Dendritique \|
	124	! \| se fait lorsque G1snSV devient > 0 \|
	125	! \| the Case of Dentritic or non Dendritic Snow is treated separately \|
	126	! \| the Limit Dentritic -> non Dendritic is reached when G1snSV > 0 \|
	127	! \| \|
	128	! \| 2. Tassement: Loi de Viscosite adaptee selon le Type de Grains \|
	129	! \| Snow Settling: Viscosity depends on the Grain Type \|
	130	! \| \|
	131	! \| 3. Update Variables historiques (cas non dendritique seulement) \|
	132	! \| nhSNow defaut \|
	133	! \| 0 Cas normal \|
	134	! \| istdSV(1) 1 Grains anguleux / faceted cristal \|
	135	! \| istdSV(2) 2 Grains ayant ete en presence d eau liquide \|
	136	! \| mais n'ayant pas eu de caractere anguleux / \|
	137	! \| liquid water and no faceted cristals before \|
	138	! \| istdSV(3) 3 Grains ayant ete en presence d eau liquide \|
	139	! \| ayant eu auparavant un caractere anguleux / \|
	140	! \| liquid water and faceted cristals before \|
	141	! \| \|
	142	! \| REFER. : Brun et al. 1989, J. Glaciol 35 pp. 333--342 \|
	143	! \| ^^^^^^^^ Brun et al. 1992, J. Glaciol 38 pp. 13-- 22 \|
	144	! \| (CROCUS Model, adapted to MAR at CEN by H.Gallee) \|
	145	! \| \|
	146	! \| REFER. : Marbouty, D. 1980, J. Glaciol 26 pp. xxx--xxx \|
	147	! \| ^^^^^^^^ (CROCUS Model, adapted to MAR at CEN by H.Gallee) \|
	148	! \| (for angular shapes) \|
	149	! \| \|
	150	! \| Preprocessing Option: SISVAT IO (not always a standard preprocess.) \|
	151	! \| ^^^^^^^^^^^^^^^^^^^^^ ^^^^^^^^^ \|
	152	! \| FILE \| CONTENT \|
	153	! \| ~~~~~~~~~~~~~~~~~~~~~+~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ \|
	154	! \| # SISVAT_GSn.vp \| #vp: OUTPUT/Verification: Snow Properties \|
	155	! \| \| unit 47, SubRoutines SISVAT_zSn, _GSn \|
	156	! \| # stdout \| #wp: OUTPUT/Verification: Snow Properties \|
	157	! \| \| unit 6, SubRoutine SISVAT_GSn \|
	158	! \| \|
	159	! +------------------------------------------------------------------------+
[3792]	160
[5246]	161
	162
	163
	164	! +--Global Variables
	165	! + ================
	166
	167	use VARphy
	168	use VAR_SV
	169	use VARdSV
	170	use VAR0SV
	171	use VARxSV
	172	use VARtSV
	173
	174
	175	IMPLICIT NONE
	176
	177
	178
	179	! +--INPUT/OUTPUT
	180	! + ------------
	181
	182
	183	! +--OUTPUT
	184	! + ------
	185
	186	integer :: dt__SV2
	187
	188
	189	! +--Local Variables
	190	! + ================
	191
	192	logical :: vector !
	193	integer :: ikl !
	194	integer :: isn ,isnp !
	195	integer :: istoOK !
	196	real :: G1_bak,G2_bak ! Old Values of G1, G2
	197	real :: ro_dry(knonv, nsno) ! Dry Density [g/cm3]
	198	real :: etaSno(knonv, nsno) ! Liquid Water Content [g/cm2]
	199	real :: SnMass(knonv) ! Snow Mass [kg/m2]
	200	real :: dTsndz ! Temperature Gradient
	201	real :: sWater ! Water Content [%]
	202	real :: exp1Wa !
	203	real :: dDENDR ! Dendricity Increment
	204	real :: DENDRn ! Normalized Dendricity
	205	real :: SPHERn ! Normalized Sphericity
	206	real :: Wet_OK ! Wet Metamorphism Switch
	207	real :: OK__DE !
	208	real :: OK__wd ! New G*, from wet Dendritic
	209	real :: G1__wd ! New G1, from wet Dendritic
	210	real :: G2__wd ! New G2, from wet Dendritic
	211	real :: OKlowT !
	212	real :: facVap !
	213	real :: OK_ldd !
	214	real :: G1_ldd !
	215	real :: G2_ldd !
	216	real :: DiamGx !
	217	real :: DiamOK !
	218	real :: No_Big !
	219	real :: dSPHER !
	220	real :: SPHER0 !
	221	real :: SPHbig !
	222	real :: G1_lds !
	223	real :: OK_mdT !
	224	real :: OKmidT !
	225	real :: OKhigT !
	226	real :: OK_mdd !
	227	real :: G1_mdd !
	228	real :: G2_mdd !
	229	real :: G1_mds !
	230	real :: OK_hdd !
	231	real :: G1_hdd !
	232	real :: G2_hdd !
	233	real :: OK_hds !
	234	real :: G1_hds !
	235	real :: T1__OK,T2__OK !
	236	real :: T3_xOK,T3__OK,T3_nOK !
	237	real :: ro1_OK,ro2_OK !
	238	real :: dT1_OK,dT2_OK,dT3xOK,dT3_OK !
	239	real :: dT4xOK,dT4_OK,dT4nOK,AngSno !
	240	real :: G2_hds,SphrOK,HISupd !
	241	real :: H1a_OK,H1b_OK,H1__OK !
	242	real :: H23aOK,H23bOK,H23_OK !
	243	real :: H2__OK,H3__OK !
	244	real :: H45_OK,H4__OK,H5__OK !
	245	real :: ViscSn,OK_Liq,OK_Ang,OKxLiq !
	246	real :: dSnMas,dzsnew,rosnew,rosmax,smb_old,smb_new
	247	real :: zn_old,zn_new
	248
	249	real :: epsi5 ! Alpha ev67 single precision
	250	real :: vdiam1 ! Small Grains Min.Diam.[.0001m]
	251	real :: vdiam2 ! Spher.Variat.Max Diam. [mm]
	252	real :: vdiam3 ! Min.Diam.\|Limit Spher. [mm]
	253	real :: vdiam4 ! Min.Diam.\|Viscosity Change
	254	real :: vsphe1 ! Max Sphericity
	255	real :: vsphe2 ! Low T Metamorphism Coeff.
	256	real :: vsphe3 ! Max.Sphericity (history=1)
	257	real :: vsphe4 ! Min.Sphericity=>history=1
	258	real :: vtang1,vtang2,vtang3,vtang4 ! Temperature Contribution
	259	real :: vtang5,vtang6,vtang7,vtang8 !
	260	real :: vtang9,vtanga,vtangb,vtangc !
	261	real :: vrang1,vrang2 ! Density Contribution
	262	real :: vgang1,vgang2,vgang3,vgang4 ! Grad(T) Contribution
	263	real :: vgang5,vgang6,vgang7,vgang8 !
	264	real :: vgang9,vganga,vgangb,vgangc !
	265	real :: vgran6 ! Max.Sphericity for Settling
	266	real :: vtelv1 ! Threshold \| history = 2, 3
	267	real :: vvap1 ! Vapor Pressure Coefficient
	268	real :: vvap2 ! Vapor Pressure Exponent
	269	real :: vgrat1 ! Boundary weak/mid grad(T)
	270	real :: vgrat2 ! Boundary mid/strong grad(T)
	271	real :: vfi ! PHI, strong grad(T)
	272	real :: vvisc1,vvisc2,vvisc3,vvisc4 ! Viscosity Coefficients
	273	real :: vvisc5,vvisc6,vvisc7 ! id., wet Snow
	274	real :: rovisc ! Wet Snow Density Influence
	275	real :: vdz3 ! Maximum Layer Densification
	276	real :: OK__ws ! New G2
	277	real :: G1__ws ! New G1, from wet Spheric
	278	real :: G2__ws ! New G2, from wet Spheric
	279	real :: husi_0,husi_1,husi_2,husi_3 ! Constants for New G2
	280	real :: vtail1,vtail2 ! Constants for New G2
	281	real :: frac_j ! Time Step [Day]
	282
	283	real :: vdent1 ! Wet Snow Metamorphism
	284	integer :: nvdent1 ! (Coefficients for
	285	integer :: nvdent2 ! Dendricity)
	286
	287	! +--Snow Properties: IO
	288	! + ~~~~~~~~~~~~~~~~~~~
	289	! #vp real G_curr(18),Gcases(18)
	290	! #vp common /GSnLOC/ Gcases
	291	! #wp real D__MAX
	292	! #wp common /GSnMAX/ D__MAX
	293
	294
	295	! +--DATA
	296	! + ====
	297
	298	data vector/.true./ ! Vectorization Switch
	299	data vdent1/ 0.5e8/ ! Wet Snow Metamorphism
	300	!XF tuned for Greenland (2.e8=old value)
	301	data nvdent1/ 3 / ! (Coefficients for
	302	data nvdent2/16 / ! Dendricity)
	303
	304	data husi_0 /20. / ! 10 * 2
	305	data husi_1 / 0.23873 / ! (3/4) /pi
	306	data husi_2 / 4.18880 / ! (4/3) *pi
	307	data husi_3 / 0.33333 / ! 1/3
	308	data vtail1 / 1.28e-08/ ! Wet Metamorphism
	309	data vtail2 / 4.22e-10/ ! (NON Dendritic / Spheric)
	310
	311	data epsi5 / 1.0e-5 / !
	312
	313	data vdiam1 / 4.0 / ! Small Grains Min.Diameter
	314
	315	data vdiam2 / 0.5 / ! Spher.Variat.Max Diam.[mm]
	316	data vdiam3 / 3.0 / ! Min.Diam.\|Limit Spher.[mm]
	317	data vdiam4 / 2.0 / ! Min.Diam.\|Viscosity Change
	318
	319	data vsphe1 / 1.0 / ! Max Sphericity
	320	data vsphe2 / 1.0e9 / ! Low T Metamorphism Coeff.
	321	data vsphe3 / 0.5 / ! Max.Sphericity (history=1)
	322	data vsphe4 / 0.1 / ! Min.Sphericity=>history=1
	323
	324	data vgran6 / 51. / ! Max.Sphericity for Settling
	325	data vtelv1 / 5.e-1 / ! Threshold \| history = 2, 3
	326
	327	data vvap1 /-6.e3 / ! Vapor Pressure Coefficient
	328	data vvap2 / 0.4 / ! Vapor Pressure Exponent
	329
	330	data vgrat1 /0.05 / ! Boundary weak/mid grad(T)
	331	data vgrat2 /0.15 / ! Boundary mid/strong grad(T)
	332	data vfi /0.09 / ! PHI, strong grad(T)
	333
	334	data vvisc1 / 0.70 / ! Viscosity Coefficients
	335	data vvisc2 / 1.11e5 / !
	336	data vvisc3 /23.00 / !
	337	data vvisc4 / 0.10 / !
	338	data vvisc5 / 1.00 / ! id., wet Snow
	339	data vvisc6 / 2.00 / !
	340	data vvisc7 /10.00 / !
	341	data rovisc / 0.25 / ! Wet Snow Density Influence
	342	data vdz3 / 0.30 / ! Maximum Layer Densification
	343
	344
	345	! +--DATA (Coefficient Fonction fort Gradient Marbouty)
	346	! + --------------------------------------------------
	347
	348	data vtang1 /40.0/ ! Temperature Contribution
	349	data vtang2 / 6.0/ !
	350	data vtang3 /22.0/ !
	351	data vtang4 / 0.7/ !
	352	data vtang5 / 0.3/ !
	353	data vtang6 / 6.0/ !
	354	data vtang7 / 1.0/ !
	355	data vtang8 / 0.8/ !
	356	data vtang9 /16.0/ !
	357	data vtanga / 0.2/ !
	358	data vtangb / 0.2/ !
	359	data vtangc /18.0/ !
	360
	361	data vrang1 / 0.40/ ! Density Contribution
	362	data vrang2 / 0.15/ !
	363
	364	data vgang1 / 0.70/ ! Grad(T) Contribution
	365	data vgang2 / 0.25/ !
	366	data vgang3 / 0.40/ !
	367	data vgang4 / 0.50/ !
	368	data vgang5 / 0.10/ !
	369	data vgang6 / 0.15/ !
	370	data vgang7 / 0.10/ !
	371	data vgang8 / 0.55/ !
	372	data vgang9 / 0.65/ !
	373	data vganga / 0.20/ !
	374	data vgangb / 0.85/ !
	375	data vgangc / 0.15/ !
	376
	377	! #wp data D__MAX / 4.00/ !
	378
	379
	380	! +-- 1. Metamorphoses dans les Strates
	381	! + Metamorphism
	382	! + ==============================
	383
	384	dt__SV2= dt__SV
	385	frac_j = dt__SV2 / 86400. ! Time Step [Day]
	386
	387	zn4_SV = 0
	388
	389
	390	! +-- 1.1 Initialisation: teneur en eau liquide et gradient de temperature
	391	! + ------------------ liquid water content and temperature gradient
	392
	393	DO ikl=1,knonv
	394	DO isn=1,isnoSV(ikl)
	395
	396	ro_dry(ikl,isn) = 1.e-3 *ro__SV(ikl,isn) & ! Dry Density
	397	*(1. -eta_SV(ikl,isn)) ! [g/cm3]
	398	etaSno(ikl,isn) = 1.e-1 *dzsnSV(ikl,isn) & ! Liquid Water
	399	* ro__SV(ikl,isn) & ! Content [g/cm2]
	400	* eta_SV(ikl,isn) !
	401	END DO
	402	END DO
	403
	404	!!$OMP PARALLEL DO default(firstprivate)
	405	!!$OMP.shared (/xSISVAT_I/,/xSISVAT_R/,/SoR0SV/,/SoI0SV/,/Sn_dSV/)
	406	DO ikl=1,knonv
	407	DO isn=1,isnoSV(ikl)
	408	isnp = min(isn+1,isnoSV(ikl))
	409
	410	dTsndz = abs( (TsisSV(ikl,isnp)-TsisSV(ikl,isn-1)) *2.e-2 &
	411	/max(((dzsnSV(ikl,isnp)+dzsnSV(ikl,isn) ) &
	412	*( isnp - isn) &
	413	+(dzsnSV(ikl,isn )+dzsnSV(ikl,isn-1))),epsi))
	414	! +... Factor 1.d-2 for Conversion K/m --> K/cm
	415
	416
	417	! +-- 1.2 Metamorphose humide
	418	! + Wet Snow Metamorphism
	419	! + ---------------------
	420
	421	Wet_OK = max(zero,sign(unun,eta_SV(ikl,isn)-epsi))
	422
	423
	424	! +-- Vitesse de diminution de la dendricite
	425	! + Rate of the dendricity decrease
	426	! + ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	427	sWater=1.d-1ro__SV(ikl,isn)eta_SV(ikl,isn) &
	428	/max(epsi,ro_dry(ikl,isn))
	429	! +... sWater:Water Content [%]
	430	! + 1.d-1= 1.d2(1->%) * 1.d-3(ro__SV*eta_SV:kg/m3->g/cm3)
	431
	432	exp1Wa= sWater**nvdent1
	433	dDENDR=max(exp1Wa/nvdent2,vdent1*exp(vvap1/TfSnow))
	434
	435	! +-- 1.2.1 Cas dendritique/dendritic Case
	436	! + ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	437	OK__wd=max(zero, & !
	438	sign(unun,-G1snSV(ikl,isn) & !
	439	-epsi )) !
	440
	441	DENDRn=-G1snSV(ikl,isn)/G1_dSV ! Normalized Dendricity (+)
	442	SPHERn= G2snSV(ikl,isn)/G1_dSV ! Normalized Sphericity
	443	DENDRn= DENDRn -dDENDR *frac_j ! New Dendricity (+)
	444	SPHERn= SPHERn +dDENDR *frac_j ! New Sphericity
	445
	446	OK__DE=max(zero, & ! IF 1.,
	447	sign(unun, DENDRn & ! NO change
	448	-epsi )) ! Dendr. -> Spheric
	449
	450	G1__wd=OK__DE * ( -DENDRn*G1_dSV) & ! Dendritic
	451	+(1.-OK__DE)* min(G1_dSV,SPHERn*G1_dSV) ! Dendr. -> Spheric
	452	G2__wd=OK__DE * min(G1_dSV,SPHERn*G1_dSV) & ! Spheric
	453	+(1.-OK__DE)*(ADSdSV-min(SPHERn,vsphe1)) ! Spher. -> Size
	454
	455	! +-- 1.2.2 Cas non dendritique non completement spherique
	456	! + Evolution de la Sphericite seulement.
	457	! + Non dendritic and not completely spheric Case
	458	! + Evolution of Sphericity only (not size)
	459	! + ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	460	OK__ws=max(zero, & !
	461	sign(unun, G1_dSV & !
	462	-epsi5 & !
	463	-G1snSV(ikl,isn))) !
	464
	465	SPHERn= G1snSV(ikl,isn)/G1_dSV
	466	SPHERn= SPHERn +dDENDR *frac_j
	467	G1__ws= min(G1_dSV,SPHERn*G1_dSV)
	468
	469	! +-- 1.2.3 Cas non dendritique et spherique / non dendritic and spheric
	470	! + Evolution de la Taille seulement / Evolution of Size only
	471	! + ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	472	G2__ws = husi_0 &
	473	*( husi_1 &
	474	(husi_2 ( G2snSV(ikl,isn)/husi_0)**3 &
	475	+(vtail1 +vtail2 exp1Wa )dt__SV2)) &
	476	** husi_3
	477
	478
	479	! +-- 1.3 Metamorposes seches / Dry Metamorphism
	480	! + --------------------------------------
	481
	482
	483	! +-- 1.3.1 Calcul Metamorphose faible/low Gradient (0.00-0.05 deg/cm)
	484	! + ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	485	OKlowT=max(zero, & !
	486	sign(unun, vgrat1 & !
	487	-dTsndz )) !
	488
	489	facVap=exp(vvap1/TsisSV(ikl,isn))
	490
	491	! +-- 1.3.1.1 Cas dendritique / dendritic Case
	492
	493	OK_ldd=max(zero, & !
	494	sign(unun,-G1snSV(ikl,isn) & !
	495	-epsi )) !
	496
	497	DENDRn=-G1snSV(ikl,isn) /G1_dSV
	498	SPHERn= G2snSV(ikl,isn) /G1_dSV
	499	DENDRn= DENDRn-vdent1facVapfrac_j
	500	SPHERn= SPHERn+vsphe2facVapfrac_j
	501
	502	OK__DE=max(zero, & ! IF 1.,
	503	sign(unun, DENDRn & ! NO change
	504	-epsi )) ! Dendr. -> Spheric
	505
	506	G1_ldd= OK__DE * ( -DENDRn*G1_dSV) & ! Dendritic
	507	+(1.-OK__DE)* min(G1_dSV,SPHERn*G1_dSV) ! Dendr. -> Spheric
	508	G2_ldd= OK__DE * min(G1_dSV,SPHERn*G1_dSV) & ! Spheric
	509	+(1.-OK__DE)*(ADSdSV-min(SPHERn,vsphe1)) ! Spher. -> Size
	510
	511	! +-- 1.3.1.2 Cas non dendritique / non dendritic Case
	512
	513	SPHERn=G1snSV(ikl,isn)/G1_dSV
	514	DiamGx=G2snSV(ikl,isn)*0.1
	515
	516	istoOK=min( abs(istoSV(ikl,isn)- &
	517	istdSV(1) ),1) ! zero if istoSV = 1
	518	DiamOK=max(zero, sign(unun,vdiam2-DiamGx))
	519	No_Big= istoOK+DiamOK
	520	No_Big=min(No_Big,unun)
	521
	522	dSPHER= vsphe2facVapfrac_j !
	523	SPHER0= SPHERn+dSPHER ! small grains
	524	SPHbig= SPHERn+dSPHER & ! big grains
	525	*exp(min(zero,vdiam3-G2snSV(ikl,isn))) ! (history = 2 or 3)
	526	SPHbig= min(vsphe3,SPHbig) ! limited sphericity
	527	SPHERn= No_Big * SPHER0 &
	528	+ (1.-No_Big)* SPHbig
	529
	530	G1_lds= min(G1_dSV,SPHERn*G1_dSV)
	531
	532	! +-- 1.3.2 Calcul Metamorphose Gradient Moyen/Moderate (0.05-0.15)
	533	! + ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	534	OK_mdT=max(zero, & !
	535	sign(unun, vgrat2 & !
	536	-dTsndz)) !
	537	OKmidT= OK_mdT *(1.-OKlowT) !
	538	OKhigT= (1. -OK_mdT) *(1.-OKlowT) !
	539
	540	facVap=vdent1*exp(vvap1/TsisSV(ikl,isn)) &
	541	* (1.e2 dTsndz)*vvap2
	542
	543	! +-- 1.3.2.1 cas dendritique / dendritic case.
	544
	545	OK_mdd=max(zero, & !
	546	sign(unun,-G1snSV(ikl,isn) & !
	547	-epsi )) !
	548
	549	DENDRn=-G1snSV(ikl,isn)/G1_dSV
	550	SPHERn= G2snSV(ikl,isn)/G1_dSV
	551	DENDRn= DENDRn - facVap*frac_j
	552	SPHERn= SPHERn - facVap*frac_j
	553
	554	OK__DE=max(zero, & ! IF 1.,
	555	sign(unun, DENDRn & ! NO change
	556	-epsi )) ! Dendr. -> Spheric
	557
	558	G1_mdd= OK__DE * ( -DENDRn*G1_dSV) & ! Dendritic
	559	+(1.-OK__DE)* max(zero ,SPHERn*G1_dSV) ! Dendr. -> Spheric
	560	G2_mdd= OK__DE * max(zero ,SPHERn*G1_dSV) & ! Spheric
	561	+(1.-OK__DE)*(ADSdSV-max(SPHERn,zero )) ! Spher. -> Size
	562
	563	! +-- 1.3.2.2 Cas non dendritique / non dendritic Case
	564
	565	SPHERn=G1snSV(ikl,isn)/G1_dSV
	566	SPHERn= SPHERn-facVap*frac_j
	567	G1_mds=max(zero,SPHERn*G1_dSV)
	568
	569	! +-- 1.3.3 Calcul Metamorphose fort / high Gradient
	570	! + ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	571	facVap=vdent1*exp(vvap1/TsisSV(ikl,isn)) &
	572	* (1.e2 dTsndz)*vvap2
	573
	574	! +-- 1.3.3.1 Cas dendritique / dendritic Case
	575
	576	OK_hdd=max(zero, & !
	577	sign(unun,-G1snSV(ikl,isn) & !
	578	-epsi )) !
	579
	580	DENDRn=-G1snSV(ikl,isn)/G1_dSV !
	581	SPHERn= G2snSV(ikl,isn)/G1_dSV !
	582	DENDRn= DENDRn - facVap*frac_j !
	583	SPHERn= SPHERn - facVap*frac_j ! Non dendritic
	584	! + ! and angular
	585	OK__DE=max(zero, & ! IF 1.,
	586	sign(unun, DENDRn & ! NO change
	587	-epsi )) ! Dendr. -> Spheric
	588
	589	G1_hdd= OK__DE * ( -DENDRn*G1_dSV) & ! Dendritic
	590	+(1.-OK__DE)* max(zero ,SPHERn*G1_dSV) ! Dendr. -> Spheric
	591	G2_hdd= OK__DE * max(zero ,SPHERn*G1_dSV) & ! Spheric
	592	+(1.-OK__DE)*(ADSdSV-max(SPHERn,zero )) ! Spher. -> Size
	593
	594	! +-- 1.3.3.2 Cas non dendritique non completement anguleux.
	595	! + non dendritic and spericity gt. 0
	596
	597	OK_hds=max(zero, & !
	598	sign(unun, G1snSV(ikl,isn) & !
	599	-epsi )) !
	600
	601	SPHERn= G1snSV(ikl,isn)/G1_dSV
	602	SPHERn= SPHERn - facVap*frac_j
	603	G1_hds= max(zero,SPHERn*G1_dSV)
	604
	605	! +-- 1.3.3.3 Cas non dendritique et anguleux
	606	! + dendritic and spericity = 0.
	607
	608	T1__OK = max(zero,sign(unun,TsisSV(ikl,isn)-TfSnow+vtang1))
	609	T2__OK = max(zero,sign(unun,TsisSV(ikl,isn)-TfSnow+vtang2))
	610	T3_xOK = max(zero,sign(unun,TsisSV(ikl,isn)-TfSnow+vtang3))
	611	T3__OK = T3_xOK * (1. - T2__OK)
	612	T3_nOK = (1. - T3_xOK) * (1. - T2__OK)
	613	ro1_OK = max(zero,sign(unun,vrang1-ro_dry(ikl,isn)))
	614	ro2_OK = max(zero,sign(unun,ro_dry(ikl,isn)-vrang2))
	615	dT1_OK = max(zero,sign(unun,vgang1-dTsndz ))
	616	dT2_OK = max(zero,sign(unun,vgang2-dTsndz ))
	617	dT3xOK = max(zero,sign(unun,vgang3-dTsndz ))
	618	dT3_OK = dT3xOK * (1. - dT2_OK)
	619	dT4xOK = max(zero,sign(unun,vgang4-dTsndz ))
	620	dT4_OK = dT4xOK * (1. - dT3_OK) &
	621	* (1. - dT2_OK)
	622	dT4nOK = (1. - dT4xOK) * (1. - dT3_OK) &
	623	* (1. - dT2_OK)
	624
	625	! +-- Influence de la Temperature /Temperature Influence
	626	! + ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	627	AngSno = &
	628	T1__OK & ! 11
	629	(T2__OK(vtang4+vtang5*(TfSnow -TsisSV(ikl,isn)) & ! 12
	630	/vtang6) & !
	631	+T3__OK(vtang7-vtang8(TfSnow-vtang2-TsisSV(ikl,isn)) & ! 13
	632	/vtang9) & !
	633	+T3_nOK(vtanga-vtangb(TfSnow-vtang3-TsisSV(ikl,isn)) & ! 14
	634	/vtangc)) & !
	635
	636	! +-- Influence de la Masse Volumique /Density Influence
	637	! + ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	638	* ro1_OK &
	639	( ro2_OK(1. - (ro_dry(ikl,isn)-vrang2) & !
	640	/(vrang1-vrang2)) & !
	641	+1.-ro2_OK ) & !
	642
	643	! +-- Influence du Gradient de Temperature /Temperature Gradient Influence
	644	! + ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	645	( dT1_OK(dT2_OKvgang5(dTsndz-vgang6) & ! 15
	646	/(vgang2-vgang6) & !
	647	+dT3_OK*vgang7 & ! 16
	648	+dT4_OK*vgang9 & ! 17
	649	+dT4nOK*vgangb ) & ! 18
	650	+1.-dT1_OK ) & !
	651	+ ro1_OK &
	652	* dT1_OK(dT3_OKvgang8*(dTsndz-vgang2) &
	653	/(vgang3-vgang2) &
	654	+dT4_OKvganga(dTsndz-vgang3) &
	655	/(vgang4-vgang3) &
	656	+dT4nOKvgangc(dTsndz-vgang4) &
	657	/(vgang1-vgang4))
	658
	659	G2_hds = G2snSV(ikl,isn) + 1.d2 AngSnovfi *frac_j
	660
	661
	662	! +--New Properties
	663	! + --------------
	664
	665	G1_bak = G1snSV(ikl,isn)
	666	G2_bak = G2snSV(ikl,isn)
	667
	668	G1snSV(ikl,isn) = Wet_OK * ( OK__wd *G1__wd & ! 1
	669	+(1.-OK__wd)* OK__ws *G1__ws & ! 2
	670	+(1.-OK__wd)(1.-OK__ws)G1_bak) & ! 3
	671	+(1. - Wet_OK) & !
	672	( OKlowT ( OK_ldd *G1_ldd & ! 4
	673	+(1.-OK_ldd) *G1_lds) & ! 5
	674	+ OKmidT ( OK_mdd G1_mdd & ! 6
	675	+(1.-OK_mdd) *G1_mds) & ! 7
	676	+ OKhigT ( OK_hdd G1_hdd & ! 8
	677	+(1.-OK_hdd)* OK_hds *G1_hds & ! 9
	678	+(1.-OK_hdd)(1.-OK_hds)G1_bak)) ! 10
	679
	680	!XF
	681	if(G1snSV(ikl,isn)<0.1) &
	682	G2_hds = G2snSV(ikl,isn) + 1.d1 AngSnovfi *frac_j
	683	!XF
	684
	685
	686	G2snSV(ikl,isn) = Wet_OK * ( OK__wd *G2__wd & ! 1
	687	+(1.-OK__wd)* OK__ws *G2_bak & ! 2
	688	+(1.-OK__wd)(1.-OK__ws)G2__ws) & ! 3
	689	+(1. - Wet_OK) & !
	690	( OKlowT ( OK_ldd *G2_ldd & ! 4
	691	+(1.-OK_ldd) *G2_bak) & ! 5
	692	+ OKmidT ( OK_mdd G2_mdd & ! 6
	693	+(1.-OK_mdd) *G2_bak) & ! 7
	694	+ OKhigT ( OK_hdd G2_hdd & ! 8
	695	+(1.-OK_hdd)* OK_hds *G2_bak & ! 9
	696	+(1.-OK_hdd)(1.-OK_hds)G2_hds)) ! 10
	697
	698	! +--Snow Properties: IO Set Up
	699	! + ~~~~~~~~~~~~~~~~~~~~~~~~~~
	700	! #vp G_curr( 1) = Wet_OK * OK__wd
	701	! #vp G_curr( 2) = Wet_OK (1.-OK__wd) OK__ws
	702	! #vp G_curr( 3) = Wet_OK (1.-OK__wd)(1.-OK__ws)
	703	! #vp G_curr( 4) = (1.-Wet_OK)* OKlowT * OK_ldd
	704	! #vp G_curr( 5) = (1.-Wet_OK)* OKlowT *(1.-OK_ldd)
	705	! #vp G_curr( 6) = (1.-Wet_OK)* OKmidT * OK_mdd
	706	! #vp G_curr( 7) = (1.-Wet_OK)* OKmidT *(1.-OK_mdd)
	707	! #vp G_curr( 8) = (1.-Wet_OK)* OKhigT * OK_hdd
	708	! #vp G_curr( 9) = (1.-Wet_OK)* OKhigT (1.-OK_hdd) OK_hds
	709	! #vp G_curr(10) = (1.-Wet_OK)* OKhigT (1.-OK_hdd)(1.-OK_hds)
	710	! #vp G_curr(11) = T1__OK * G_curr(10)
	711	! #vp G_curr(12) = T2__OK * G_curr(10)
	712	! #vp G_curr(13) = T3__OK * G_curr(10)
	713	! #vp G_curr(14) = T3_nOK * G_curr(10)
	714	! #vp G_curr(15) = ro1_OK* dT1_OK * dT2_OK * G_curr(10)
	715	! #vp G_curr(16) = ro1_OK* dT1_OK * dT3_OK * G_curr(10)
	716	! #vp G_curr(17) = ro1_OK* dT1_OK * dT4_OK * G_curr(10)
	717	! #vp G_curr(18) = ro1_OK* dT1_OK * dT4nOK * G_curr(10)
	718
	719	! #vp Gcases( 1) = max(Gcases( 1),G_curr( 1))
	720	! #vp Gcases( 2) = max(Gcases( 2),G_curr( 2))
	721	! #vp Gcases( 3) = max(Gcases( 3),G_curr( 3))
	722	! #vp Gcases( 4) = max(Gcases( 4),G_curr( 4))
	723	! #vp Gcases( 5) = max(Gcases( 5),G_curr( 5))
	724	! #vp Gcases( 6) = max(Gcases( 6),G_curr( 6))
	725	! #vp Gcases( 7) = max(Gcases( 7),G_curr( 7))
	726	! #vp Gcases( 8) = max(Gcases( 8),G_curr( 8))
	727	! #vp Gcases( 9) = max(Gcases( 9),G_curr( 9))
	728	! #vp Gcases(10) = max(Gcases(10),G_curr(10))
	729	! #vp Gcases(11) = max(Gcases(11),G_curr(11))
	730	! #vp Gcases(12) = max(Gcases(12),G_curr(12))
	731	! #vp Gcases(13) = max(Gcases(13),G_curr(13))
	732	! #vp Gcases(14) = max(Gcases(14),G_curr(14))
	733	! #vp Gcases(15) = max(Gcases(15),G_curr(15))
	734	! #vp Gcases(16) = max(Gcases(16),G_curr(16))
	735	! #vp Gcases(17) = max(Gcases(17),G_curr(17))
	736	! #vp Gcases(18) = max(Gcases(18),G_curr(18))
	737
	738	! +--Snow Properties: IO
	739	! + ~~~~~~~~~~~~~~~~~~~
	740	! #vp IF (isn .le. isnoSV(ikl))
	741	! #vp. write(47,471)isn ,isnoSV(ikl) ,
	742	! #vp. TsisSV(ikl,isn),ro__SV(ikl,isn),eta_SV(ikl,isn),
	743	! #vp. G1_bak ,G2_bak ,istoSV(ikl,isn),
	744	! #vp. dTsndz,
	745	! #vp. ( k ,k=1,18),
	746	! #vp. (G_curr(k),k=1,18),
	747	! #vp. (Gcases(k),k=1,18),
	748	! #vp. Wet_OK,OK__wd,G1__wd,G2__wd,
	749	! #vp. 1.-OK__wd,OK__ws,G1__ws,1.-OK__ws,G2__ws,
	750	! #vp. 1.-Wet_OK,OKlowT,OK_ldd,G1_ldd, G2_ldd,
	751	! #vp. 1.-OK_ldd,G1_lds,
	752	! #vp. OKmidT,OK_mdd,G1_mdd, G1_mdd,
	753	! #vp. 1.-OK_mdd,G1_mds,
	754	! #vp. OKhigT,OK_hdd,G1_hdd, G2_hdd,
	755	! #vp. 1.-OK_hdd,OK_hds, G1_hds,
	756	! #vp. 1.-OK_hds,G2_hds,
	757	! #vp. G1snSV(ikl,isn),
	758	! #vp. G2snSV(ikl,isn)
	759
	760	END DO
	761	END DO
	762	!!$OMP END PARALLEL DO
	763
	764	! +-- 2. Mise a Jour Variables Historiques (Cas non dendritique)
	765	! + Update of the historical Variables
	766	! + =======================================================
	767
	768	IF (vector) THEN
	769	!XF
	770	DO ikl=1,knonv
	771	DO isn=1,isnoSV(ikl)
	772	SphrOK = max(zero,sign(unun, G1snSV(ikl,isn)))
	773	H1a_OK = max(zero,sign(unun,vsphe4-G1snSV(ikl,isn)))
	774	H1b_OK = 1 - min(1 , istoSV(ikl,isn))
	775	H1__OK = H1a_OK*H1b_OK
	776	H23aOK = max(zero,sign(unun,vsphe4-G1_dSV &
	777	+G1snSV(ikl,isn)))
	778	H23bOK = max(zero,sign(unun,etaSno(ikl,isn) &
	779	/max(epsi,dzsnSV(ikl,isn)) &
	780	-vtelv1 ))
	781	H23_OK = H23aOK*H23bOK
	782	H2__OK = 1 - min(1 , istoSV(ikl,isn))
	783	H3__OK = 1 - min(1 , abs(istoSV(ikl,isn)-istdSV(1)))
	784	H45_OK = max(zero,sign(unun,TfSnow-TsisSV(ikl,isn)+epsi))
	785	H4__OK = 1 - min(1 , abs(istoSV(ikl,isn)-istdSV(2)))
	786	H5__OK = 1 - min(1 , abs(istoSV(ikl,isn)-istdSV(3)))
	787
	788	HISupd = &
	789	SphrOK(H1__OK istdSV(1) &
	790	+(1.-H1__OK)* H23_OK (H2__OKistdSV(2) &
	791	+H3__OK*istdSV(3)) &
	792	+(1.-H1__OK)(1.-H23_OK) H45_OK(H4__OKistdSV(4) &
	793	+H5__OK*istdSV(5)))
	794	istoSV(ikl,isn) = HISupd + &
	795	(1.-min(unun,HISupd)) *istoSV(ikl,isn)
	796	END DO
	797	END DO
	798	ELSE
	799
	800
	801	! +-- 2. Mise a Jour Variables Historiques (Cas non dendritique)
	802	! + Update of the historical Variables
	803	! + =======================================================
	804
	805	DO ikl=1,knonv
	806	DO isn=iiceSV(ikl),isnoSV(ikl)
	807	IF (G1snSV(ikl,isn).ge.0.) THEN
	808	IF(G1snSV(ikl,isn).lt.vsphe4.and.istoSV(ikl,isn).eq.0) THEN
	809	istoSV(ikl,isn)=istdSV(1)
	810	ELSEIF(G1_dSV-G1snSV(ikl,isn) .lt.vsphe4.and. &
	811	etaSno(ikl,isn)/dzsnSV(ikl,isn).gt.vtelv1) THEN
	812	IF (istoSV(ikl,isn).eq.0) &
	813	istoSV(ikl,isn)= istdSV(2)
	814	IF (istoSV(ikl,isn).eq.istdSV(1)) &
	815	istoSV(ikl,isn)= istdSV(3)
	816	ELSEIF(TsisSV(ikl,isn).lt.TfSnow) THEN
	817	IF (istoSV(ikl,isn).eq.istdSV(2)) &
	818	istoSV(ikl,isn)= istdSV(4)
	819	IF (istoSV(ikl,isn).eq.istdSV(3)) &
	820	istoSV(ikl,isn)= istdSV(5)
	821	END IF
[3792]	822	END IF
[5246]	823	END DO
	824	END DO
	825	END IF
[3792]	826
[5246]	827
	828	! +-- 3. Tassement mecanique /mechanical Settlement
	829	! + ==========================================
	830
	831	DO ikl=1,knonv
	832	SnMass(ikl) = 0.
	833	END DO
	834	!XF
	835	DO ikl=1,knonv
	836
	837	smb_old = 0.
	838	zn_old = 0
	839	DO isn = 1, isnoSV(ikl)
	840	smb_old = smb_old + dzsnSV(ikl,isn) *ro__SV(ikl,isn)
	841	zn_old = zn_old + dzsnSV(ikl,isn)
	842	ENDDO
	843
	844	DO isn=isnoSV(ikl),1,-1
	845	dSnMas = 100.dzsnSV(ikl,isn)ro_dry(ikl,isn)
	846	SnMass(ikl)= SnMass(ikl)+0.5*dSnMas
	847	ViscSn = vvisc1 *vvisc2 &
	848	exp(vvisc3 ro_dry(ikl,isn) &
	849	+vvisc4*abs(TfSnow-TsisSV(ikl,isn))) &
	850	*ro_dry(ikl,isn)/rovisc
	851
	852	! +-- Changement de Viscosite si Teneur en Eau liquide
	853	! + Change of the Viscosity if liquid Water Content
	854	! + ------------------------------------------------
	855
	856	OK_Liq = max(zero,sign(unun,etaSno(ikl,isn)-epsi))
	857	OK_Ang = max(zero,sign(unun,vgran6-G1snSV(ikl,isn))) &
	858	*(1-min(1 , abs(istoSV(ikl,isn)-istdSV(1))))
	859	! #wp IF (G1snSV(ikl,isn).gt.0..AND.G1snSV(ikl,isn).lt.vsphe4
	860	! #wp. .AND.istoSV(ikl,isn).eq. 0)
	861	! #wp. THEN
	862	! #wp write(6,*) ikl,isn,' G1,G2,hist,OK_Ang ',
	863	! #wp. G1snSV(ikl,isn), G2snSV(ikl,isn),istoSV(ikl,isn),OK_Ang
	864	! #wp stop "Grains anguleux mal d?finis"
	865	! #wp END IF
	866	OKxLiq = max(zero,sign(unun,vtelv1-etaSno(ikl,isn) &
	867	/max(epsi,dzsnSV(ikl,isn)))) &
	868	* max(0 ,sign(1 ,istoSV(ikl,isn) &
	869	-istdSV(1) ))
	870	ViscSn = &
	871	ViscSn( OK_Liq/(vvisc5+vvisc6etaSno(ikl,isn) &
	872	/max(epsi,dzsnSV(ikl,isn))) &
	873	+(1.-OK_Liq) ) &
	874	( OK_Angexp(min(ADSdSV,G2snSV(ikl,isn)-vdiam4)) &
	875	+(1.-OK_Ang) ) &
	876	( OKxLiq vvisc7 &
	877	+(1.-OKxLiq) )
	878
	879
	880	! +-- Calcul nouvelle Epaisseur / new Thickness
	881	! + -----------------------------------------
	882
	883	dzsnew = &
	884	dzsnSV(ikl,isn) &
	885	*max(vdz3, &
	886	(unun-dt__SV2max(SnMass(ikl)cos(slopSV(ikl)),unun) &
	887	/max(ViscSn ,epsi)))
	888	rosnew = ro__SV(ikl,isn) *dzsnSV(ikl,isn) &
	889	/max(1e-10,dzsnew)
	890	rosmax = 1. /( (1. -eta_SV(ikl,isn)) /ro_Ice &
	891	+ eta_SV(ikl,isn) /ro_Wat)
	892	rosnew = min(rosnew ,rosmax)
	893	dzsnew = dzsnSV(ikl,isn) *ro__SV(ikl,isn) &
	894	/max(1e-10,rosnew)
	895	ro__SV(ikl,isn)= rosnew
	896	dzsnSV(ikl,isn)= dzsnew
	897	ro_dry(ikl,isn)= ro__SV(ikl,isn)(1.-eta_SV(ikl,isn))1.e-3
	898	! +... ro_dry: Dry Density (g/cm3)
	899	! +
	900	SnMass(ikl) = SnMass(ikl)+dSnMas*0.5
	901	END DO
	902
	903	smb_new = 0.
	904	DO isn = 1, isnoSV(ikl)
	905	smb_new = smb_new + dzsnSV(ikl,isn) *ro__SV(ikl,isn)
	906	ENDDO
	907
	908	isn=1
	909	if (dzsnSV(ikl,isn)>0.and.ro__SV(ikl,isn)>0) then
	910	dzsnSV(ikl,isn) = dzsnSV(ikl,isn) +0.9999*(smb_old-smb_new) &
	911	/ ro__SV(ikl,isn)
	912	endif
	913
	914	zn_new = 0
	915	DO isn = 1, isnoSV(ikl)
	916	zn_new = zn_new + dzsnSV(ikl,isn)
	917	ENDDO
	918	zn4_SV(ikl) = zn4_SV(ikl) + (zn_new - zn_old)
	919
	920	END DO
	921
	922
	923
	924	return
	925	end subroutine sisvat_gsn

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

Download in other formats:

Original Format