subroutine SISVAT_zSn C +------------------------------------------------------------------------+ C | MAR SISVAT_zSn 12-07-2019 MAR | C | SubRoutine SISVAT_zSn manages the Snow Pack vertical Discretization | C | | C +------------------------------------------------------------------------+ C | | C | PARAMETERS: knonv: Total Number of columns = | C | ^^^^^^^^^^ = Total Number of continental grid boxes | C | X Number of Mosaic Cell per grid box | C | | C | INPUT / NLaysv = New Snow Layer Switch | C | OUTPUT: isnoSV = total Nb of Ice/Snow Layers | C | ^^^^^^ ispiSV = 0,...,nsno: Uppermost Superimposed Ice Layer | C | iiceSV = total Nb of Ice Layers | C | istoSV = 0,...,5 : Snow History (see istdSV data) | C | | C | INPUT / TsisSV : Soil/Ice Temperatures (layers -nsol,-nsol+1,..,0)| C | OUTPUT: & Snow Temperatures (layers 1,2,...,nsno) [K] | C | ^^^^^^ ro__SV : Soil/Snow Volumic Mass [kg/m3] | C | eta_SV : Soil/Snow Water Content [m3/m3] | C | dzsnSV : Snow Layer Thickness [m] | C | G1snSV : Dendricity (<0) or Sphericity (>0) of Snow Layer | C | G2snSV : Sphericity (>0) or Size of Snow Layer | C | agsnSV : Snow Age [day] | C | | C | METHOD: 1) Agregate the thinest Snow Layer | C | ^^^^^^ if a new Snow Layer has been precipitated (NLaysv = 1) | C | 2) Divide a too thick Snow Layer except | C | if the maximum Number of Layer is reached | C | in this case forces NLay_s = 1 | C | 3) Agregate the thinest Snow Layer | C | in order to divide a too thick Snow Layer | C | at next Time Step when NLay_s = 1 | C | | C +------------------------------------------------------------------------+ C +--Global Variables C + ================ use VARphy use VAR_SV use VARdSV use VAR0SV use VARxSV use VARySV IMPLICIT NONE C +--Internal Variables C + ================== integer ikl ,isn ,i ! integer NLay_s(knonv) ! Split Snow Layer Switch integer isagr1(knonv) ! 1st Layer History integer isagr2(knonv) ! 2nd Layer History integer LstLay ! 0 ====> isnoSV = 1 integer isno_n ! Snow Normal.Profile integer iice_n ! Ice Normal.Profile integer iiceOK ! Ice Switch integer icemix ! 0 ====> Agregated Snow+Ice=Snow C + ! 1 Ice integer isn1 (knonv) ! 1st layer to stagger real staggr ! stagger Switch real WEagre(knonv) ! Snow Water Equivalent Thickness real dzthin(knonv) ! Thickness of the thinest layer real OKthin ! Swich ON a new thinest layer real dz_dif ! difference from ideal discret. real thickL ! Thick Layer Indicator real OK_ICE ! Swich ON uppermost Ice Layer real Agrege(knonv) ! 1. when Agregation constrained real dzepsi ! Min Single Snw Layer Thickness real dzxmin ! Min Acceptable Layer Thickness real dz_min ! Min Layer Thickness real dz_max ! Max Layer Thickness real dzagr1(knonv) ! 1st Layer Thickness real dzagr2(knonv) ! 2nd Layer Thickness real T_agr1(knonv) ! 1st Layer Temperature real T_agr2(knonv) ! 2nd Layer Temperature real roagr1(knonv) ! 1st Layer Density real roagr2(knonv) ! 2nd Layer Density real etagr1(knonv) ! 1st Layer Water Content real etagr2(knonv) ! 2nd Layer Water Content real G1agr1(knonv) ! 1st Layer Dendricity/Spher. real G1agr2(knonv) ! 2nd Layer Dendricity/Spher. real G2agr1(knonv) ! 1st Layer Sphericity/Size real G2agr2(knonv) ! 2nd Layer Sphericity/Size real agagr1(knonv) ! 1st Layer Age real agagr2(knonv) ! 2nd Layer Age C +--DATA C + ==== data icemix / 0 / ! 0 ====> Agregated Snow+Ice=Snow data dzepsi / 0.0020/ ! Min single Layer Thickness data dzxmin / 0.0025/ ! Min accept.Layer Thickness c #EU data dz_min / 0.0050/ ! Min Local Layer Thickness < SMn data dz_min / 0.0040/ ! Min Local Layer Thickness < SMn data dz_max / 0.0300/ ! Min Gener. Layer Thickness C + CAUTION: dz_max > dz_min*2 is required ! Otherwise re-agregation is C + ! activated after splitting C +--Constrains Agregation of too thin Layers C + ================================================= C +--Search the thinest non-zero Layer C + ---------------------------------- DO ikl=1,knonv if(isnoSV(ikl)<=2) dz_min=max(0.0050,dz_min) dzepsi=0.0015 if(ro__SV(ikl,isnoSV(ikl))>920) dzepsi=0.0020 dzthin(ikl) = 0. ! Arbitrary unrealistic END DO ! Layer Thickness cXF DO ikl=1,knonv DO isn=1,isnoSV(ikl)-3 ! no agregation of 3 first snowlayers ! XF 04/07/2019 isno_n = isnoSV(ikl)-isn+1 ! Snow Normal.Profile iice_n = iiceSV(ikl)-isn ! Ice Normal.Profile iiceOK = min(1,max(0,iice_n +1)) ! Ice Switch ! #vz dz_ref(isn) = ! ! #vz. dz_min *((1-iiceOK)*isno_n*isno_n ! Theoretical Profile ! #vz. + iiceOK * 2**iice_n) ! ! #vz. /max(1,isnoSV(ikl)) ! dz_dif = max(zero, ! Actual Profile . dz_min ! . *((1-iiceOK)*isno_n*isno_n ! Theoretical Profile . + iiceOK *2. **iice_n) ! . - dzsnSV(ikl, isn) ) ! Actual Profile ! #vz dzwdif(isn) = dz_dif ! OKthin = max(zero, ! . sign(unun, ! . dz_dif-dzthin(ikl))) ! 1.=> New thinest Lay. . * max(0, ! 1 => .le. isnoSV . min(1, ! 1 => isn is in the . isnoSV(ikl)-isn +1 )) ! Snow Pack . * min(unun, ! ! ! 1st additional Condition to accept OKthin . max(zero, ! combination . sign(unun,G1snSV(ikl, isn ) ! G1 with same . *G1snSV(ikl,max(1,isn-1)))) ! sign => OK ! ! 2nd additional Condition to accept OKthin . + max(zero, ! G1>0 . sign(unun,G1snSV(ikl, isn ))) ! =>OK ! ! 3rd additional Condition to accept OKthin . + max(zero, ! dz too small . sign(unun,dzxmin ! =>OK . -dzsnSV(ikl, isn ))))! i_thin(ikl) = (1. - OKthin) * i_thin(ikl) ! Update thinest Lay. . + OKthin * isn ! Index dzthin(ikl) = (1. - OKthin) * dzthin(ikl) ! . + OKthin * dz_dif ! END DO END DO DO ikl=1,knonv DO isn=1,isnoSV(ikl) OKthin = max(zero, ! . sign(unun, ! . dz_min ! . -dzsnSV(ikl,isn))) ! . * max(zero, ! ON if dz > 0 . sign(unun, ! . dzsnSV(ikl,isn)-epsi)) ! . *min(1,max(0, ! Multiple Snow Lay. . min (1, ! Switch = 1 . isnoSV(ikl) ! if isno > iice + 1 . -iiceSV(ikl)-1)) ! C + ! . +int(max(zero, ! . sign(unun, ! . dzepsi ! Minimum accepted for . -dzsnSV(ikl,isn)))) ! 1 Snow Layer over Ice . *int(max(zero, ! ON if dz > 0 . sign(unun, ! . dzsnSV(ikl,isn)-epsi)))! . *(1 -min (abs(isnoSV(ikl) ! Switch = 1 . -iiceSV(ikl)-1),1)) ! if isno = iice + 1 C + ! . +max(0, ! Ice . min (1, ! Switch . iiceSV(ikl)+1-isn))) ! . *min(unun, ! . max(zero, ! combination . sign(unun,G1snSV(ikl, isn ) ! G1>0 + G1<0 . *G1snSV(ikl,max(1,isn-1)))) ! NO . + max(zero, ! . sign(unun,G1snSV(ikl, isn ))) ! . + max(zero, ! . sign(unun,dzxmin ! . -dzsnSV(ikl, isn ))))! i_thin(ikl) = (1. - OKthin) * i_thin(ikl) ! Update thinest Lay. . + OKthin * isn ! Index END DO END DO C + *************** call SISVAT_zCr C + *************** C +--Assign the 2 Layers to agregate C + ------------------------------- DO ikl=1,knonv isn = i_thin(ikl) if(LIndsv(ikl)>0) isn=min(nsno-1,isn) ! cXF isagr1(ikl) = istoSV(ikl,isn) isagr2(ikl) = istoSV(ikl,isn+LIndsv(ikl)) dzagr1(ikl) = dzsnSV(ikl,isn) dzagr2(ikl) = dzsnSV(ikl,isn+LIndsv(ikl)) T_agr1(ikl) = TsisSV(ikl,isn) T_agr2(ikl) = TsisSV(ikl,isn+LIndsv(ikl)) roagr1(ikl) = ro__SV(ikl,isn) roagr2(ikl) = ro__SV(ikl,isn+LIndsv(ikl)) etagr1(ikl) = eta_SV(ikl,isn) etagr2(ikl) = eta_SV(ikl,isn+LIndsv(ikl)) G1agr1(ikl) = G1snSV(ikl,isn) G1agr2(ikl) = G1snSV(ikl,isn+LIndsv(ikl)) G2agr1(ikl) = G2snSV(ikl,isn) G2agr2(ikl) = G2snSV(ikl,isn+LIndsv(ikl)) agagr1(ikl) = agsnSV(ikl,isn) agagr2(ikl) = agsnSV(ikl,isn+LIndsv(ikl)) LstLay = min(1,max( 0,isnoSV(ikl) -1)) ! 0 if single Layer isnoSV(ikl) = isnoSV(ikl) ! decrement isnoSV . -(1-LstLay)* max(zero, ! if downmost Layer . sign(unun,eps_21 ! < 1.e-21 m . -dzsnSV(ikl,1))) ! isnoSV(ikl) = max( 0, isnoSV(ikl) ) ! Agrege(ikl) = max(zero, ! . sign(unun,dz_min ! No Agregation . -dzagr1(ikl) )) ! if too thick Layer . *LstLay ! if a single Layer . * min( max(0 ,isnoSV(ikl)+1 ! if Agregation . -i_thin(ikl) ! with a Layer . -LIndsv(ikl) ),1) ! above the Pack WEagre(ikl) = 0. END DO DO ikl=1,knonv DO isn=1,isnoSV(ikl) WEagre(ikl) = WEagre(ikl) + ro__SV(ikl,isn)*dzsnSV(ikl,isn) . *min(1,max(0,i_thin(ikl)+1-isn)) ENDDO ENDDO C +--Agregates C + --------- C + *************** call SISVAT_zAg . (isagr1,isagr2,WEagre . ,dzagr1,dzagr2,T_agr1,T_agr2 . ,roagr1,roagr2,etagr1,etagr2 . ,G1agr1,G1agr2,G2agr1,G2agr2 . ,agagr1,agagr2,Agrege . ) C + *************** C +--Rearranges the Layers C + --------------------- C +--New (agregated) Snow layer C + ^^^^^^^^^^^^^^^^^^^^^^^^^^ DO ikl=1,knonv isn = i_thin(ikl) isn = min(isn,isn+LIndsv(ikl)) isnoSV(ikl) = max(0.,isnoSV(ikl) -Agrege(ikl)) iiceSV(ikl) = iiceSV(ikl) . -max(0,sign(1,iiceSV(ikl) -isn +icemix)) . *Agrege(ikl) . *max(0,sign(1,iiceSV(ikl) -1 )) istoSV(ikl,isn) = (1.-Agrege(ikl))*istoSV(ikl,isn) . + Agrege(ikl) *isagr1(ikl) dzsnSV(ikl,isn) = (1.-Agrege(ikl))*dzsnSV(ikl,isn) . + Agrege(ikl) *dzagr1(ikl) TsisSV(ikl,isn) = (1.-Agrege(ikl))*TsisSV(ikl,isn) . + Agrege(ikl) *T_agr1(ikl) ro__SV(ikl,isn) = (1.-Agrege(ikl))*ro__SV(ikl,isn) . + Agrege(ikl) *roagr1(ikl) eta_SV(ikl,isn) = (1.-Agrege(ikl))*eta_SV(ikl,isn) . + Agrege(ikl) *etagr1(ikl) G1snSV(ikl,isn) = (1.-Agrege(ikl))*G1snSV(ikl,isn) . + Agrege(ikl) *G1agr1(ikl) G2snSV(ikl,isn) = (1.-Agrege(ikl))*G2snSV(ikl,isn) . + Agrege(ikl) *G2agr1(ikl) agsnSV(ikl,isn) = (1.-Agrege(ikl))*agsnSV(ikl,isn) . + Agrege(ikl) *agagr1(ikl) END DO C +--Above C + ^^^^^ DO ikl=1,knonv isn1(ikl)=max(i_thin(ikl),i_thin(ikl)+LIndsv(ikl)) END DO DO i= 1,nsno-1 DO ikl=1,knonv staggr = min(1,max(0,i +1 -isn1(ikl) )) istoSV(ikl,i) = (1.-staggr )*istoSV(ikl,i ) . + staggr*((1.-Agrege(ikl))*istoSV(ikl,i ) . + Agrege(ikl) *istoSV(ikl,i+1)) dzsnSV(ikl,i) = (1.-staggr )*dzsnSV(ikl,i ) . + staggr*((1.-Agrege(ikl))*dzsnSV(ikl,i ) . + Agrege(ikl) *dzsnSV(ikl,i+1)) TsisSV(ikl,i) = (1.-staggr )*TsisSV(ikl,i ) . + staggr*((1.-Agrege(ikl))*TsisSV(ikl,i ) . + Agrege(ikl) *TsisSV(ikl,i+1)) ro__SV(ikl,i) = (1.-staggr )*ro__SV(ikl,i ) . + staggr*((1.-Agrege(ikl))*ro__SV(ikl,i ) . + Agrege(ikl) *ro__SV(ikl,i+1)) eta_SV(ikl,i) = (1.-staggr )*eta_SV(ikl,i ) . + staggr*((1.-Agrege(ikl))*eta_SV(ikl,i ) . + Agrege(ikl) *eta_SV(ikl,i+1)) G1snSV(ikl,i) = (1.-staggr )*G1snSV(ikl,i ) . + staggr*((1.-Agrege(ikl))*G1snSV(ikl,i ) . + Agrege(ikl) *G1snSV(ikl,i+1)) G2snSV(ikl,i) = (1.-staggr )*G2snSV(ikl,i ) . + staggr*((1.-Agrege(ikl))*G2snSV(ikl,i ) . + Agrege(ikl) *G2snSV(ikl,i+1)) agsnSV(ikl,i) = (1.-staggr )*agsnSV(ikl,i ) . + staggr*((1.-Agrege(ikl))*agsnSV(ikl,i ) . + Agrege(ikl) *agsnSV(ikl,i+1)) END DO END DO DO ikl=1,knonv isn = min(isnoSV(ikl) +1,nsno) istoSV(ikl,isn) = (1.-Agrege(ikl))*istoSV(ikl,isn) dzsnSV(ikl,isn) = (1.-Agrege(ikl))*dzsnSV(ikl,isn) TsisSV(ikl,isn) = (1.-Agrege(ikl))*TsisSV(ikl,isn) ro__SV(ikl,isn) = (1.-Agrege(ikl))*ro__SV(ikl,isn) eta_SV(ikl,isn) = (1.-Agrege(ikl))*eta_SV(ikl,isn) G1snSV(ikl,isn) = (1.-Agrege(ikl))*G1snSV(ikl,isn) G2snSV(ikl,isn) = (1.-Agrege(ikl))*G2snSV(ikl,isn) agsnSV(ikl,isn) = (1.-Agrege(ikl))*agsnSV(ikl,isn) END DO C +--Constrains Splitting of too thick Layers C + ================================================= C +--Search the thickest non-zero Layer C + ---------------------------------- DO ikl=1,knonv dzthin(ikl) = 0. ! Arbitrary unrealistic END DO DO ikl=1,knonv DO isn=1,isnoSV(ikl) isno_n = isnoSV(ikl)-isn+1 ! Snow Normal.Profile iice_n = iiceSV(ikl)-isn ! Ice Normal.Profile iiceOK = min(1,max(0,iice_n +1)) ! Ice Switch dz_dif =( dzsnSV(ikl,isn) ! Actual Profile . - dz_max *((1-iiceOK)*isno_n*isno_n ! Theoretical Profile . + iiceOK *2. **iice_n) ) ! . /max(dzsnSV(ikl,isn),epsi) ! OKthin = max(zero, ! . sign(unun, ! . dz_dif-dzthin(ikl))) ! 1.=>New thickest Lay. . * max(0, ! 1 =>.le. isnoSV . min(1, ! . isnoSV(ikl)-isn +1 )) ! i_thin(ikl) = (1. - OKthin) * i_thin(ikl) ! Update thickest Lay. . + OKthin * isn ! Index dzthin(ikl) = (1. - OKthin) * dzthin(ikl) ! . + OKthin * dz_dif ! END DO isn=max(1,isnoSV(ikl)-3) if(dzsnSV(ikl,isn)>0.30) then ! surface layer > 30cm i_thin(ikl) = isn ! XF 04/07/2019 dzthin(ikl) = dzsnSV(ikl,isn) endif isn=max(1,isnoSV(ikl)-2) if(dzsnSV(ikl,isn)>0.10) then ! surface layer > 10cm i_thin(ikl) = isn ! XF 04/07/2019 dzthin(ikl) = dzsnSV(ikl,isn) endif isn=max(1,isnoSV(ikl)-1) if(dzsnSV(ikl,isn)>0.05) then ! surface layer > 5cm i_thin(ikl) = isn ! XF 04/07/2019 dzthin(ikl) = dzsnSV(ikl,isn) endif isn=max(1,isnoSV(ikl)) if(dzsnSV(ikl,isn)>0.02) then ! surface layer > 2cm i_thin(ikl) = isn ! XF 04/07/2019 dzthin(ikl) = dzsnSV(ikl,isn) endif END DO DO ikl=1,knonv ThickL = max(zero, ! 1. => a too thick . sign(unun,dzthin(ikl) ! Layer exists . -epsi )) ! . * max(0,1-max(0 , isnoSV(ikl) ! No spliting allowed . -nsno+1 )) ! if isno > nsno - 1 Agrege(ikl) = ThickL ! 1. => effective split . * max(0,1-max(0 , NLaysv(ikl) ! . +isnoSV(ikl) ! . -nsno+1 )) ! NLay_s(ikl) = ThickL ! Agregation . * max(0,1-max(0 , NLaysv(ikl) ! to allow Splitting . +isnoSV(ikl) ! at next Time Step . -nsno )) ! . -Agrege(ikl) ! NLay_s(ikl) = max(0 , NLay_s(ikl)) ! Agregation effective END DO C +--Rearranges the Layers C + --------------------- DO isn=nsno,2,-1 DO ikl=1,knonv IF (Agrege(ikl).gt.0..AND.i_thin(ikl).lt.isnoSV(ikl)) THEN staggr = min(1,max(0,isn-i_thin(ikl) -1)) . * min(1,max(0, isnoSV(ikl)-isn+2)) istoSV(ikl,isn) = staggr * istoSV(ikl ,isn-1) . + (1. - staggr) * istoSV(ikl ,isn ) dzsnSV(ikl,isn) = staggr * dzsnSV(ikl ,isn-1) . + (1. - staggr) * dzsnSV(ikl ,isn ) TsisSV(ikl,isn) = staggr * TsisSV(ikl ,isn-1) . + (1. - staggr) * TsisSV(ikl ,isn ) ro__SV(ikl,isn) = staggr * ro__SV(ikl ,isn-1) . + (1. - staggr) * ro__SV(ikl ,isn ) eta_SV(ikl,isn) = staggr * eta_SV(ikl ,isn-1) . + (1. - staggr) * eta_SV(ikl ,isn ) G1snSV(ikl,isn) = staggr * G1snSV(ikl ,isn-1) . + (1. - staggr) * G1snSV(ikl ,isn ) G2snSV(ikl,isn) = staggr * G2snSV(ikl ,isn-1) . + (1. - staggr) * G2snSV(ikl ,isn ) agsnSV(ikl,isn) = staggr * agsnSV(ikl ,isn-1) . + (1. - staggr) * agsnSV(ikl ,isn ) END IF END DO END DO DO ikl=1,knonv isn = i_thin(ikl) dzsnSV(ikl,isn) = 0.5*Agrege(ikl) *dzsnSV(ikl,isn) . + (1.-Agrege(ikl))*dzsnSV(ikl,isn) isn = min(i_thin(ikl) +1,nsno) istoSV(ikl,isn) = Agrege(ikl) *istoSV(ikl,isn-1) . + (1.-Agrege(ikl))*istoSV(ikl,isn) dzsnSV(ikl,isn) = Agrege(ikl) *dzsnSV(ikl,isn-1) . + (1.-Agrege(ikl))*dzsnSV(ikl,isn) TsisSV(ikl,isn) = Agrege(ikl) *TsisSV(ikl,isn-1) . + (1.-Agrege(ikl))*TsisSV(ikl,isn) ro__SV(ikl,isn) = Agrege(ikl) *ro__SV(ikl,isn-1) . + (1.-Agrege(ikl))*ro__SV(ikl,isn) eta_SV(ikl,isn) = Agrege(ikl) *eta_SV(ikl,isn-1) . + (1.-Agrege(ikl))*eta_SV(ikl,isn) G1snSV(ikl,isn) = Agrege(ikl) *G1snSV(ikl,isn-1) . + (1.-Agrege(ikl))*G1snSV(ikl,isn) G2snSV(ikl,isn) = Agrege(ikl) *G2snSV(ikl,isn-1) . + (1.-Agrege(ikl))*G2snSV(ikl,isn) agsnSV(ikl,isn) = Agrege(ikl) *agsnSV(ikl,isn-1) . + (1.-Agrege(ikl))*agsnSV(ikl,isn) isnoSV(ikl) = min(Agrege(ikl) +isnoSV(ikl),real(nsno)) iiceSV(ikl) = iiceSV(ikl) . + Agrege(ikl) *max(0,sign(1,iiceSV(ikl) . -isn +icemix)) . *max(0,sign(1,iiceSV(ikl) . -1 )) END DO C +--Constrains Agregation in case of too much Layers C + ================================================= C +--Search the thinest non-zero Layer C + ----------------------------------- DO ikl=1,knonv dzthin(ikl) = 0. ! Arbitrary unrealistic END DO ! Layer Thickness DO ikl=1,knonv DO isn=1,isnoSV(ikl)-3 ! no agregation of 3 first snowlayers ! XF 04/07/2019 isno_n = isnoSV(ikl)-isn+1 ! Snow Normal.Profile iice_n = iiceSV(ikl)-isn ! Ice Normal.Profile iiceOK = min(1,max(0,iice_n +1)) ! Ice Switch ! #vz dz_ref(isn) = ! ! #vz. dz_min *((1-iiceOK)*isno_n*isno_n ! Theoretical Profile ! #vz. + iiceOK * 2**iice_n) ! ! #vz. /max(1,isnoSV(ikl)) ! dz_dif = dz_min ! Actual Profile . - dzsnSV(ikl ,isn) ! . /max(epsi,((1-iiceOK)*isno_n*isno_n ! Theoretical Profile . + iiceOK *2. **iice_n)) ! ! #vz dzwdif(isn) = dz_dif ! OKthin = max(zero, ! . sign(unun, ! . dz_dif - dzthin(ikl)))! 1.=> New thinest Lay. . * max(0, ! 1 => .le. isnoSV . min(1, ! . isnoSV(ikl)-isn +1 )) ! i_thin(ikl) = (1. - OKthin) * i_thin(ikl) ! Update thinest Lay. . + OKthin * isn ! Index dzthin(ikl) = (1. - OKthin) * dzthin(ikl) ! . + OKthin * dz_dif ! END DO END DO C +--Index of the contiguous Layer to agregate C + ----------------------------------------- C + *************** call SISVAT_zCr C + *************** C +--Assign the 2 Layers to agregate C + ------------------------------- DO ikl=1,knonv isn = i_thin(ikl) if(LIndsv(ikl)>0) isn=min(isn, nsno-1) !cXF isagr1(ikl) = istoSV(ikl,isn) isagr2(ikl) = istoSV(ikl,isn+LIndsv(ikl)) dzagr1(ikl) = dzsnSV(ikl,isn) dzagr2(ikl) = dzsnSV(ikl,isn+LIndsv(ikl)) T_agr1(ikl) = TsisSV(ikl,isn) T_agr2(ikl) = TsisSV(ikl,isn+LIndsv(ikl)) roagr1(ikl) = ro__SV(ikl,isn) roagr2(ikl) = ro__SV(ikl,isn+LIndsv(ikl)) etagr1(ikl) = eta_SV(ikl,isn) etagr2(ikl) = eta_SV(ikl,isn+LIndsv(ikl)) G1agr1(ikl) = G1snSV(ikl,isn) G1agr2(ikl) = G1snSV(ikl,isn+LIndsv(ikl)) G2agr1(ikl) = G2snSV(ikl,isn) G2agr2(ikl) = G2snSV(ikl,isn+LIndsv(ikl)) agagr1(ikl) = agsnSV(ikl,isn) agagr2(ikl) = agsnSV(ikl,isn+LIndsv(ikl)) LstLay = min(1,max( 0, isnoSV(ikl)-1 )) Agrege(ikl) = min(1, . max(0, . NLaysv(ikl) +isnoSV(ikl)-nsno . +NLay_s(ikl) ) . *LstLay ) C + minimum uppermost layer thickness to guarantee a correct reproduction of the snow C + atmosphere coupling if(dzsnSV(ikl,max(1,isnoSV(ikl)-0))>0.02 .or. ! surface layers> 2-5-10 . dzsnSV(ikl,max(1,isnoSV(ikl)-1))>0.05 .or. ! XF 04/07/2019 . dzsnSV(ikl,max(1,isnoSV(ikl)-2))>0.10 .or. . dzsnSV(ikl,max(1,isnoSV(ikl)-3))>0.30 )then Agrege(ikl) = min(1, . max(0, . NLaysv(ikl) +isnoSV(ikl)+1-nsno ! nsno-1 layers ma . +NLay_s(ikl) ) . *LstLay ) endif isnoSV(ikl) = isnoSV(ikl) . -(1-LstLay)*max(zero, . sign(unun, eps_21 . -dzsnSV(ikl,1) )) isnoSV(ikl) =max( 0, isnoSV(ikl) ) WEagre(ikl) = 0. END DO DO isn=1,nsno DO ikl=1,knonv WEagre(ikl) = WEagre(ikl) + ro__SV(ikl,isn)*dzsnSV(ikl,isn) . *min(1,max(0,i_thin(ikl)+1-isn)) ENDDO ENDDO C +--Agregates C + --------- C + *************** call SISVAT_zAg . (isagr1,isagr2,WEagre . ,dzagr1,dzagr2,T_agr1,T_agr2 . ,roagr1,roagr2,etagr1,etagr2 . ,G1agr1,G1agr2,G2agr1,G2agr2 . ,agagr1,agagr2,Agrege . ) C + *************** C +--Rearranges the Layers C + --------------------- C +--New (agregated) Snow layer C + ^^^^^^^^^^^^^^^^^^^^^^^^^^ DO ikl=1,knonv isn = i_thin(ikl) isn = min(isn,isn+LIndsv(ikl)) isnoSV(ikl) = max(0.,isnoSV(ikl) -Agrege(ikl)) iiceSV(ikl) = iiceSV(ikl) . -max(0,sign(1,iiceSV(ikl) -isn +icemix)) . *Agrege(ikl) . *max(0,sign(1,iiceSV(ikl) -1 )) istoSV(ikl,isn) = (1.-Agrege(ikl))*istoSV(ikl,isn) . + Agrege(ikl) *isagr1(ikl) dzsnSV(ikl,isn) = (1.-Agrege(ikl))*dzsnSV(ikl,isn) . + Agrege(ikl) *dzagr1(ikl) TsisSV(ikl,isn) = (1.-Agrege(ikl))*TsisSV(ikl,isn) . + Agrege(ikl) *T_agr1(ikl) ro__SV(ikl,isn) = (1.-Agrege(ikl))*ro__SV(ikl,isn) . + Agrege(ikl) *roagr1(ikl) eta_SV(ikl,isn) = (1.-Agrege(ikl))*eta_SV(ikl,isn) . + Agrege(ikl) *etagr1(ikl) G1snSV(ikl,isn) = (1.-Agrege(ikl))*G1snSV(ikl,isn) . + Agrege(ikl) *G1agr1(ikl) G2snSV(ikl,isn) = (1.-Agrege(ikl))*G2snSV(ikl,isn) . + Agrege(ikl) *G2agr1(ikl) agsnSV(ikl,isn) = (1.-Agrege(ikl))*agsnSV(ikl,isn) . + Agrege(ikl) *agagr1(ikl) END DO C +--Above C + ^^^^^ DO ikl=1,knonv isn1(ikl)=max(i_thin(ikl),i_thin(ikl)+LIndsv(ikl)) END DO DO i= 1,nsno-1 DO ikl=1,knonv staggr = min(1,max(0,i +1 -isn1(ikl) )) istoSV(ikl,i) = (1.-staggr )*istoSV(ikl,i ) . + staggr*((1.-Agrege(ikl))*istoSV(ikl,i ) . + Agrege(ikl) *istoSV(ikl,i+1)) dzsnSV(ikl,i) = (1.-staggr )*dzsnSV(ikl,i ) . + staggr*((1.-Agrege(ikl))*dzsnSV(ikl,i ) . + Agrege(ikl) *dzsnSV(ikl,i+1)) TsisSV(ikl,i) = (1.-staggr )*TsisSV(ikl,i ) . + staggr*((1.-Agrege(ikl))*TsisSV(ikl,i ) . + Agrege(ikl) *TsisSV(ikl,i+1)) ro__SV(ikl,i) = (1.-staggr )*ro__SV(ikl,i ) . + staggr*((1.-Agrege(ikl))*ro__SV(ikl,i ) . + Agrege(ikl) *ro__SV(ikl,i+1)) eta_SV(ikl,i) = (1.-staggr )*eta_SV(ikl,i ) . + staggr*((1.-Agrege(ikl))*eta_SV(ikl,i ) . + Agrege(ikl) *eta_SV(ikl,i+1)) G1snSV(ikl,i) = (1.-staggr )*G1snSV(ikl,i ) . + staggr*((1.-Agrege(ikl))*G1snSV(ikl,i ) . + Agrege(ikl) *G1snSV(ikl,i+1)) G2snSV(ikl,i) = (1.-staggr )*G2snSV(ikl,i ) . + staggr*((1.-Agrege(ikl))*G2snSV(ikl,i ) . + Agrege(ikl) *G2snSV(ikl,i+1)) agsnSV(ikl,i) = (1.-staggr )*agsnSV(ikl,i ) . + staggr*((1.-Agrege(ikl))*agsnSV(ikl,i ) . + Agrege(ikl) *agsnSV(ikl,i+1)) END DO END DO DO ikl=1,knonv isn = min(isnoSV(ikl) +1,nsno) istoSV(ikl,isn) = (1.-Agrege(ikl))*istoSV(ikl,isn) dzsnSV(ikl,isn) = (1.-Agrege(ikl))*dzsnSV(ikl,isn) TsisSV(ikl,isn) = (1.-Agrege(ikl))*TsisSV(ikl,isn) ro__SV(ikl,isn) = (1.-Agrege(ikl))*ro__SV(ikl,isn) eta_SV(ikl,isn) = (1.-Agrege(ikl))*eta_SV(ikl,isn) G1snSV(ikl,isn) = (1.-Agrege(ikl))*G1snSV(ikl,isn) G2snSV(ikl,isn) = (1.-Agrege(ikl))*G2snSV(ikl,isn) agsnSV(ikl,isn) = (1.-Agrege(ikl))*agsnSV(ikl,isn) END DO return end