Changeset 5116 for LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis
- Timestamp:
- Jul 24, 2024, 2:54:37 PM (7 months ago)
- Location:
- LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis
- Files:
-
- 13 edited
-
inlandsis.f90 (modified) (18 diffs)
-
sisvat_bsn.f90 (modified) (4 diffs)
-
sisvat_gsn.f90 (modified) (7 diffs)
-
sisvat_qsn.f90 (modified) (9 diffs)
-
sisvat_qso.f90 (modified) (7 diffs)
-
sisvat_sno_albedo.f90 (modified) (4 diffs)
-
sisvat_ts2.f90 (modified) (8 diffs)
-
sisvat_tso.f90 (modified) (5 diffs)
-
sisvat_weq.f90 (modified) (4 diffs)
-
sisvat_zag.f90 (modified) (4 diffs)
-
sisvat_zcr.f90 (modified) (2 diffs)
-
sisvat_zsn.f90 (modified) (11 diffs)
-
surf_inlandsis_mod.F90 (modified) (16 diffs)
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis/inlandsis.f90
r5113 r5116 186 186 logical :: BloMod 187 187 logical :: debut 188 integer:: jjtime188 INTEGER :: jjtime 189 189 190 190 … … 195 195 ! --------- 196 196 197 real:: TBr_sv(klonv) ! Brightness Temperature198 real:: IRdwsv(klonv) ! DOWNward IR Flux199 real:: IRupsv(klonv) ! UPward IR Flux200 real:: d_Bufs,Bufs_N ! Buffer Snow Layer Increment201 real:: Buf_ro,Bros_N ! Buffer Snow Layer Density202 real:: BufPro ! Buffer Snow Layer Density203 real:: Buf_G1,BG1__N ! Buffer Snow Layer Dendr/Sphe[-]204 real:: Buf_G2,BG2__N ! Buffer Snow Layer Spher/Size[-]205 real:: Bdzssv(klonv) ! Buffer Snow Layer Thickness206 real:: z_snsv(klonv) ! Snow-Ice, current Thickness197 REAL :: TBr_sv(klonv) ! Brightness Temperature 198 REAL :: IRdwsv(klonv) ! DOWNward IR Flux 199 REAL :: IRupsv(klonv) ! UPward IR Flux 200 REAL :: d_Bufs,Bufs_N ! Buffer Snow Layer Increment 201 REAL :: Buf_ro,Bros_N ! Buffer Snow Layer Density 202 REAL :: BufPro ! Buffer Snow Layer Density 203 REAL :: Buf_G1,BG1__N ! Buffer Snow Layer Dendr/Sphe[-] 204 REAL :: Buf_G2,BG2__N ! Buffer Snow Layer Spher/Size[-] 205 REAL :: Bdzssv(klonv) ! Buffer Snow Layer Thickness 206 REAL :: z_snsv(klonv) ! Snow-Ice, current Thickness 207 207 208 208 … … 211 211 ! ----- 212 212 213 integer:: iwr214 integer:: ikl ,isn ,isl ,ist !215 integer:: ist__s,ist__w ! Soil/Water Body Identifier216 integer:: growth ! Seasonal Mask217 integer:: LISmsk ! Land+Ice / Open Sea Mask218 integer:: LSnMsk ! Snow-Ice / No Snow-Ice Mask219 integer:: IceMsk,IcIndx(klonv) ! Ice / No Ice Mask220 integer:: SnoMsk ! Snow / No Snow Mask221 real:: roSMin,roSMax,roSn_1,roSn_2,roSn_3 ! Fallen Snow Density (PAHAUT)222 real:: Dendr1,Dendr2,Dendr3 ! Fallen Snow Dendric.(GIRAUD)223 real:: Spher1,Spher2,Spher3,Spher4 ! Fallen Snow Spheric.(GIRAUD)224 real:: Polair ! Polar Snow Switch225 real:: PorSno,Salt_f,PorRef !213 INTEGER :: iwr 214 INTEGER :: ikl ,isn ,isl ,ist ! 215 INTEGER :: ist__s,ist__w ! Soil/Water Body Identifier 216 INTEGER :: growth ! Seasonal Mask 217 INTEGER :: LISmsk ! Land+Ice / Open Sea Mask 218 INTEGER :: LSnMsk ! Snow-Ice / No Snow-Ice Mask 219 INTEGER :: IceMsk,IcIndx(klonv) ! Ice / No Ice Mask 220 INTEGER :: SnoMsk ! Snow / No Snow Mask 221 REAL :: roSMin,roSMax,roSn_1,roSn_2,roSn_3 ! Fallen Snow Density (PAHAUT) 222 REAL :: Dendr1,Dendr2,Dendr3 ! Fallen Snow Dendric.(GIRAUD) 223 REAL :: Spher1,Spher2,Spher3,Spher4 ! Fallen Snow Spheric.(GIRAUD) 224 REAL :: Polair ! Polar Snow Switch 225 REAL :: PorSno,Salt_f,PorRef ! 226 226 ! #sw real PorVol,rWater ! 227 227 ! #sw real rusNEW,rdzNEW,etaNEW ! 228 real:: ro_new !229 real:: TaPole ! Maximum Polar Temperature230 real:: T__Min ! Minimum realistic Temperature231 real:: EmiSol ! Emissivity of Soil232 real:: EmiSno ! Emissivity of Snow233 real:: EmiWat ! Emissivity of a Water Area234 real:: vk2 ! Square of Von Karman Constant235 real:: u2star !(u*)**2236 real:: Z0mLnd ! Land Roughness Length228 REAL :: ro_new ! 229 REAL :: TaPole ! Maximum Polar Temperature 230 REAL :: T__Min ! Minimum realistic Temperature 231 REAL :: EmiSol ! Emissivity of Soil 232 REAL :: EmiSno ! Emissivity of Snow 233 REAL :: EmiWat ! Emissivity of a Water Area 234 REAL :: vk2 ! Square of Von Karman Constant 235 REAL :: u2star !(u*)**2 236 REAL :: Z0mLnd ! Land Roughness Length 237 237 ! #ZN real sqrrZ0 ! u*t/u* 238 real:: f_eff ! Marticorena & B. 1995 JGR (20)239 real:: A_Fact ! Fundamental * Roughness240 real:: Z0m_nu ! Smooth R Snow Roughness Length241 real:: Z0mBSn ! BSnow Roughness Length242 real:: Z0mBS0 ! Mimimum BSnow Roughness Length243 real:: Z0m_S0 ! Mimimum Snow Roughness Length244 real:: Z0m_S1 ! Maximum Snow Roughness Length238 REAL :: f_eff ! Marticorena & B. 1995 JGR (20) 239 REAL :: A_Fact ! Fundamental * Roughness 240 REAL :: Z0m_nu ! Smooth R Snow Roughness Length 241 REAL :: Z0mBSn ! BSnow Roughness Length 242 REAL :: Z0mBS0 ! Mimimum BSnow Roughness Length 243 REAL :: Z0m_S0 ! Mimimum Snow Roughness Length 244 REAL :: Z0m_S1 ! Maximum Snow Roughness Length 245 245 ! #SZ real Z0Sa_N ! Regime Snow Roughness Length 246 246 ! #SZ real Z0SaSi ! 1.IF Rgm Snow Roughness Length 247 247 ! #GL real Z0_GIM ! Mimimum GIMEX Roughness Length 248 real:: Z0_ICE ! Ice ISW Roughness Length249 real:: Z0m_Sn,Z0m_90 ! Snow Surface Roughness Length250 real:: SnoWat ! Snow Layer Switch251 real:: rstar,alors !252 real:: rstar0,rstar1,rstar2 !253 real:: SameOK ! 1. => Same Type of Grains254 real:: G1same ! Averaged G1, same Grains255 real:: G2same ! Averaged G2, same Grains256 real:: typ__1 ! 1. => Lay1 Type: Dendritic257 real:: zroNEW ! dz X ro, if fresh Snow258 real:: G1_NEW ! G1, if fresh Snow259 real:: G2_NEW ! G2, if fresh Snow260 real:: zroOLD ! dz X ro, if old Snow261 real:: G1_OLD ! G1, if old Snow262 real:: G2_OLD ! G2, if old Snow263 real:: SizNEW ! Size, if fresh Snow264 real:: SphNEW ! Spheric.,if fresh Snow265 real:: SizOLD ! Size, if old Snow266 real:: SphOLD ! Spheric.,if old Snow267 real:: Siz_av ! Averaged Grain Size268 real:: Sph_av ! Averaged Grain Spher.269 real:: Den_av ! Averaged Grain Dendr.270 real:: G1diff ! Averaged G1, diff. Grains271 real:: G2diff ! Averaged G2, diff. Grains272 real:: G1 ! Averaged G1273 real:: G2 ! Averaged G2274 real:: param ! Polynomial fit z0=f(T)275 real:: Z0_obs ! Fit Z0_obs=f(T) (m)276 real:: tamin ! min T of linear fit (K)277 real:: tamax ! max T of linear fit (K)278 real:: coefa,coefb,coefc,coefd ! Coefs for z0=f(T)279 real:: ta1,ta2,ta3 ! Air temperature thresholds280 real:: z01,z02,z03 ! z0 thresholds281 real:: tt_c,vv_c ! Critical param.282 real:: tt_tmp,vv_tmp,vv_virt ! Temporary variables283 real:: e_prad,e1pRad,A_Rad0,absg_V,absgnI,exdRad ! variables for SoSosv calculations284 real:: zm1, zm2, coefslope ! variables for surface temperature extrapolation248 REAL :: Z0_ICE ! Ice ISW Roughness Length 249 REAL :: Z0m_Sn,Z0m_90 ! Snow Surface Roughness Length 250 REAL :: SnoWat ! Snow Layer Switch 251 REAL :: rstar,alors ! 252 REAL :: rstar0,rstar1,rstar2 ! 253 REAL :: SameOK ! 1. => Same Type of Grains 254 REAL :: G1same ! Averaged G1, same Grains 255 REAL :: G2same ! Averaged G2, same Grains 256 REAL :: typ__1 ! 1. => Lay1 Type: Dendritic 257 REAL :: zroNEW ! dz X ro, if fresh Snow 258 REAL :: G1_NEW ! G1, if fresh Snow 259 REAL :: G2_NEW ! G2, if fresh Snow 260 REAL :: zroOLD ! dz X ro, if old Snow 261 REAL :: G1_OLD ! G1, if old Snow 262 REAL :: G2_OLD ! G2, if old Snow 263 REAL :: SizNEW ! Size, if fresh Snow 264 REAL :: SphNEW ! Spheric.,if fresh Snow 265 REAL :: SizOLD ! Size, if old Snow 266 REAL :: SphOLD ! Spheric.,if old Snow 267 REAL :: Siz_av ! Averaged Grain Size 268 REAL :: Sph_av ! Averaged Grain Spher. 269 REAL :: Den_av ! Averaged Grain Dendr. 270 REAL :: G1diff ! Averaged G1, diff. Grains 271 REAL :: G2diff ! Averaged G2, diff. Grains 272 REAL :: G1 ! Averaged G1 273 REAL :: G2 ! Averaged G2 274 REAL :: param ! Polynomial fit z0=f(T) 275 REAL :: Z0_obs ! Fit Z0_obs=f(T) (m) 276 REAL :: tamin ! min T of linear fit (K) 277 REAL :: tamax ! max T of linear fit (K) 278 REAL :: coefa,coefb,coefc,coefd ! Coefs for z0=f(T) 279 REAL :: ta1,ta2,ta3 ! Air temperature thresholds 280 REAL :: z01,z02,z03 ! z0 thresholds 281 REAL :: tt_c,vv_c ! Critical param. 282 REAL :: tt_tmp,vv_tmp,vv_virt ! Temporary variables 283 REAL :: e_prad,e1pRad,A_Rad0,absg_V,absgnI,exdRad ! variables for SoSosv calculations 284 REAL :: zm1, zm2, coefslope ! variables for surface temperature extrapolation 285 285 ! for Aeolian erosion and blowing snow 286 integer:: nit ,iit287 real:: Fac ! Correc. factor for drift ratio288 real:: dusuth,signus289 real:: sss__F,sss__N290 real:: sss__K,sss__G291 real:: us_127,us_227,us_327,us_427,us_527292 real:: VVa_OK, usuth0293 real:: ssstar294 real:: SblPom295 real:: rCd10n ! Square root of drag coefficient296 real:: DendOK ! Dendricity Switch297 real:: SaltOK ! Saltation Switch298 real:: MeltOK ! Saltation Switch (Melting Snow)299 real:: SnowOK ! Pack Top Switch300 real:: SaltM1,SaltM2,SaltMo,SaltMx ! Saltation Parameters301 real:: ShearX, ShearS ! Arg. Max Shear Stress302 real:: Por_BS ! Snow Porosity303 real:: Salt_us ! New thresh.friction velocity u*t304 real:: Fac_Mo,ArguSi,FacRho ! Numerical factors for u*t305 real:: SaltSI(klonv,0:nsno) ! Snow Drift Index !306 real:: MIN_Mo ! Minimum Mobility Fresh Fallen *307 character(len=3) :: qsalt_param ! Switch for saltation flux param.308 character(len=3) :: usth_param ! Switch for u*t param286 INTEGER :: nit ,iit 287 REAL :: Fac ! Correc. factor for drift ratio 288 REAL :: dusuth,signus 289 REAL :: sss__F,sss__N 290 REAL :: sss__K,sss__G 291 REAL :: us_127,us_227,us_327,us_427,us_527 292 REAL :: VVa_OK, usuth0 293 REAL :: ssstar 294 REAL :: SblPom 295 REAL :: rCd10n ! Square root of drag coefficient 296 REAL :: DendOK ! Dendricity Switch 297 REAL :: SaltOK ! Saltation Switch 298 REAL :: MeltOK ! Saltation Switch (Melting Snow) 299 REAL :: SnowOK ! Pack Top Switch 300 REAL :: SaltM1,SaltM2,SaltMo,SaltMx ! Saltation Parameters 301 REAL :: ShearX, ShearS ! Arg. Max Shear Stress 302 REAL :: Por_BS ! Snow Porosity 303 REAL :: Salt_us ! New thresh.friction velocity u*t 304 REAL :: Fac_Mo,ArguSi,FacRho ! Numerical factors for u*t 305 REAL :: SaltSI(klonv,0:nsno) ! Snow Drift Index ! 306 REAL :: MIN_Mo ! Minimum Mobility Fresh Fallen * 307 CHARACTER(LEN=3) :: qsalt_param ! Switch for saltation flux param. 308 CHARACTER(LEN=3) :: usth_param ! Switch for u*t param 309 309 310 310 … … 409 409 IF (BloMod) THEN 410 410 411 if (klonv==1) then412 if(isnoSV(1)>=2 .and. &411 if (klonv==1) THEN 412 IF(isnoSV(1)>=2 .and. & 413 413 TsisSV(1,max(1,isnoSV(1)))<273. .and. & 414 414 ro__SV(1,max(1,isnoSV(1)))<500. .and. & 415 eta_SV(1,max(1,isnoSV(1)))<epsi) then415 eta_SV(1,max(1,isnoSV(1)))<epsi) THEN 416 416 ! + ********** 417 417 CALL SISVAT_BSn … … 508 508 ! +--Threshold Friction Velocity 509 509 ! + ~~~~~~~~~~~~~~~~~~~~~~~~~~~ 510 if(ro__SV(ikl,isn)>300.) then510 IF(ro__SV(ikl,isn)>300.) THEN 511 511 Por_BS = 1.000 - ro__SV(ikl,isn) /ro_Ice 512 512 else … … 518 518 ! + Gallee et al., 2001 eq 5, p5 519 519 520 if (usth_param == "gal") then520 if (usth_param == "gal") THEN 521 521 Salt_us = (log(2.868) - log(1 + SaltMo)) * rCd10n/0.085 522 522 Salt_us = Salt_us * Fac_Mo … … 526 526 527 527 if (usth_param == "lis") then !Liston et al. 2007 528 if(ro__SV(ikl,isn)>300.) then528 IF(ro__SV(ikl,isn)>300.) THEN 529 529 Salt_us = 0.005*exp(0.013*ro__SV(ikl,isn)) 530 530 else … … 607 607 hSalSV(ikl) = 8.436e-2 * us__SV(ikl)**SblPom 608 608 609 if (qsalt_param == "pom") then609 if (qsalt_param == "pom") THEN 610 610 qSalSV(ikl) = (us__SV(ikl)**2 - usthSV(ikl)**2) *signus & 611 611 / (hSalSV(ikl) * gravit * us__SV(ikl) * 3.25) 612 612 endif 613 613 614 if (qsalt_param == "bin") then614 if (qsalt_param == "bin") THEN 615 615 qSalSV(ikl) = (us__SV(ikl) * us__SV(ikl) & 616 616 -usthSV(ikl) * usthSV(ikl))*signus & … … 662 662 ! #BS density_kotlyakov = .FALSE. !C.Amory BS 2018 663 663 ! + ... Fallen Snow Density, Adapted for Antarctica 664 if (is_ok_density_kotlyakov) then664 if (is_ok_density_kotlyakov) THEN 665 665 tt_tmp = TaT_SV(ikl)-TfSnow 666 666 !vv_tmp = VV10SV(ikl) … … 668 668 ! + ... [ A compromise between 669 669 ! + ... Kotlyakov (1961) and Lenaerts (2012, JGR, Part1) ] 670 if (tt_tmp>=-10) then670 if (tt_tmp>=-10) THEN 671 671 BufPro = max( rosMin, & 672 672 104. *sqrt( max( vv_tmp-6.0,0.0))) ! Kotlyakov (1961) … … 696 696 ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 697 697 698 if (BloMod) then698 if (BloMod) THEN 699 699 Bros_N = frsno 700 700 ro_new = ro__SV(ikl,max(1,isnoSV(ikl))) … … 894 894 895 895 896 if (discret_xf.AND.klonv==1) then 897 898 if(isnoSV(1)>=1.or.NLaysv(1)>=1) then 896 if (discret_xf.AND.klonv==1) THEN 897 IF(isnoSV(1)>=1.or.NLaysv(1)>=1) THEN 899 898 ! + ********** 900 899 CALL SISVAT_zSn … … 1081 1080 1082 1081 1083 if (iflag_temp_inlandsis == 0) then 1084 1082 if (iflag_temp_inlandsis == 0) THEN 1085 1083 CALL SISVAT_TSo 1086 1084 … … 1205 1203 ! Etienne: extrapolation from the two uppermost levels: 1206 1204 1207 if (isnoSV(ikl) >=2) then1205 if (isnoSV(ikl) >=2) THEN 1208 1206 zm1=-dzsnSV(ikl,isnoSV(ikl))/2. 1209 1207 zm2=-(dzsnSV(ikl,isnoSV(ikl)) + dzsnSV(ikl,isnoSV(ikl)-1)/2.) 1210 else if (isnoSV(ikl) == 1) then1208 else if (isnoSV(ikl) == 1) THEN 1211 1209 zm1=-dzsnSV(ikl,isnoSV(ikl))/2. 1212 1210 zm2=-(dzsnSV(ikl,isnoSV(ikl))+dz_dSV(0)/2.) … … 1236 1234 IF (SnoMod) THEN 1237 1235 1238 if (discret_xf .AND. klonv==1) then1239 if(isnoSV(1)>=1) then1236 if (discret_xf .AND. klonv==1) THEN 1237 IF(isnoSV(1)>=1) THEN 1240 1238 ! + ********** 1241 1239 CALL SISVAT_GSn … … 1336 1334 coefd = log(z03)-coefc*ta3 1337 1335 1338 if (TaT_SV(ikl) < ta1) then1336 if (TaT_SV(ikl) < ta1) THEN 1339 1337 Z0_obs = z01 1340 else if (TaT_SV(ikl)>=ta1 .and. TaT_SV(ikl)<ta2) then1338 else if (TaT_SV(ikl)>=ta1 .and. TaT_SV(ikl)<ta2) THEN 1341 1339 Z0_obs = exp(coefa*TaT_SV(ikl) + coefb) 1342 else if (TaT_SV(ikl)>=ta2 .and. TaT_SV(ikl)<ta3) then1340 else if (TaT_SV(ikl)>=ta2 .and. TaT_SV(ikl)<ta3) THEN 1343 1341 ! if st > 0, melting induce smooth surface 1344 1342 Z0_obs = exp(coefc*TaT_SV(ikl) + coefd) … … 1472 1470 !XF MAR is then too warm and not enough melt 1473 1471 1474 if(ro__SV(ikl,isnoSV(ikl))>50 & 1475 .and.ro__SV(ikl,isnoSV(ikl))<roSdSV)then 1476 1472 IF(ro__SV(ikl,isnoSV(ikl))>50 & 1473 .and.ro__SV(ikl,isnoSV(ikl))<roSdSV)THEN 1477 1474 Z0hnSV(ikl) = max(zero & 1478 1475 , sign(unun,zzsnsv(ikl,isnoSV(ikl))-epsi)) & … … 1497 1494 1498 1495 1499 end subroutineinlandsis1500 1501 1502 1503 1504 1505 1506 1507 1508 1509 1510 1496 END SUBROUTINE inlandsis 1497 1498 1499 1500 1501 1502 1503 1504 1505 1506 1507 -
LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis/sisvat_bsn.f90
r5113 r5116 33 33 34 34 35 integer:: ikl ,isn36 real:: h_mmWE ! Eroded Snow Layer Min Thickness37 real:: dbsaux(knonv) ! Drift Amount (Dummy Variable)38 real:: dzweqo,dzweqn,bsno_x ! Conversion variables for erosion39 real:: dz_new,rho_new40 real:: snofOK ! Threshd Snow Fall41 real:: Fac ! Correction factor for erosion42 real:: densif ! Densification rate if erosion35 INTEGER :: ikl ,isn 36 REAL :: h_mmWE ! Eroded Snow Layer Min Thickness 37 REAL :: dbsaux(knonv) ! Drift Amount (Dummy Variable) 38 REAL :: dzweqo,dzweqn,bsno_x ! Conversion variables for erosion 39 REAL :: dz_new,rho_new 40 REAL :: snofOK ! Threshd Snow Fall 41 REAL :: Fac ! Correction factor for erosion 42 REAL :: densif ! Densification rate if erosion 43 43 44 44 ! +--DATA … … 73 73 if((dzweqo-dzweqn)>0 .and. & 74 74 dzsnSV(ikl,isn)>0 .and. & 75 ro__SV(ikl,max(1,isnoSV(ikl)))<roBdSV) then 76 75 ro__SV(ikl,max(1,isnoSV(ikl)))<roBdSV) THEN 77 76 !characteristic time scale for drifting snow compaction set to 24h 78 77 !linear densification rate [kg/m3/s] over 24h … … 84 83 Fac = max(0.,min(1.,Fac)) 85 84 86 if (ro__SV(ikl,max(1,isnoSV(ikl)))>roBdSV) then85 if (ro__SV(ikl,max(1,isnoSV(ikl)))>roBdSV) THEN 87 86 densif=densif*Fac 88 87 endif … … 94 93 endif 95 94 96 if(dzsnSV(ikl,isn)>0 .and.dzsnSV(ikl,isn)<0.0001)then95 IF(dzsnSV(ikl,isn)>0 .and.dzsnSV(ikl,isn)<0.0001)THEN 97 96 dbs_SV(ikl) = dbs_SV(ikl)+ dzsnSV(ikl,isn)*ro__SV(ikl,isn) 98 97 dbs_Er(ikl) = dbs_Er(ikl)+ dzsnSV(ikl,isn)*ro__SV(ikl,isn) -
LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis/sisvat_gsn.f90
r5105 r5116 184 184 ! + ------ 185 185 186 integer:: dt__SV2186 INTEGER :: dt__SV2 187 187 188 188 … … 191 191 192 192 logical :: vector ! 193 integer:: ikl !194 integer:: isn ,isnp !195 integer:: istoOK !196 real:: G1_bak,G2_bak ! Old Values of G1, G2197 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 Gradient201 real:: sWater ! Water Content [%]202 real:: exp1Wa !203 real:: dDENDR ! Dendricity Increment204 real:: DENDRn ! Normalized Dendricity205 real:: SPHERn ! Normalized Sphericity206 real:: Wet_OK ! Wet Metamorphism Switch207 real:: OK__DE !208 real:: OK__wd ! New G*, from wet Dendritic209 real:: G1__wd ! New G1, from wet Dendritic210 real:: G2__wd ! New G2, from wet Dendritic211 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_new247 real:: zn_old,zn_new248 249 real :: epsi5 ! Alpha ev67 single precision250 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 Change254 real:: vsphe1 ! Max Sphericity255 real:: vsphe2 ! Low T Metamorphism Coeff.256 real:: vsphe3 ! Max.Sphericity (history=1)257 real:: vsphe4 ! Min.Sphericity=>history=1258 real:: vtang1,vtang2,vtang3,vtang4 ! Temperature Contribution259 real:: vtang5,vtang6,vtang7,vtang8 !260 real:: vtang9,vtanga,vtangb,vtangc !261 real:: vrang1,vrang2 ! Density Contribution262 real:: vgang1,vgang2,vgang3,vgang4 ! Grad(T) Contribution263 real:: vgang5,vgang6,vgang7,vgang8 !264 real:: vgang9,vganga,vgangb,vgangc !265 real:: vgran6 ! Max.Sphericity for Settling266 real:: vtelv1 ! Threshold | history = 2, 3267 real:: vvap1 ! Vapor Pressure Coefficient268 real:: vvap2 ! Vapor Pressure Exponent269 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 Coefficients273 real:: vvisc5,vvisc6,vvisc7 ! id., wet Snow274 real:: rovisc ! Wet Snow Density Influence275 real:: vdz3 ! Maximum Layer Densification276 real:: OK__ws ! New G2277 real:: G1__ws ! New G1, from wet Spheric278 real:: G2__ws ! New G2, from wet Spheric279 real:: husi_0,husi_1,husi_2,husi_3 ! Constants for New G2280 real:: vtail1,vtail2 ! Constants for New G2281 real:: frac_j ! Time Step [Day]282 283 real:: vdent1 ! Wet Snow Metamorphism284 integer:: nvdent1 ! (Coefficients for285 integer:: nvdent2 ! Dendricity)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 286 287 287 ! +--Snow Properties: IO … … 679 679 680 680 !XF 681 if(G1snSV(ikl,isn)<0.1) &681 IF(G1snSV(ikl,isn)<0.1) & 682 682 G2_hds = G2snSV(ikl,isn) + 1.d1 *AngSno*vfi *frac_j 683 683 !XF … … 739 739 ! + ~~~~~~~~~~~~~~~~~~~ 740 740 ! #vp IF (isn .le. isnoSV(ikl)) 741 ! #vp. write(47,471)isn ,isnoSV(ikl) ,741 ! #vp. WRITE(47,471)isn ,isnoSV(ikl) , 742 742 ! #vp. TsisSV(ikl,isn),ro__SV(ikl,isn),eta_SV(ikl,isn), 743 743 ! #vp. G1_bak ,G2_bak ,istoSV(ikl,isn), … … 860 860 ! #wp. .AND.istoSV(ikl,isn).eq. 0) 861 861 ! #wp. THEN 862 ! #wp write(6,*) ikl,isn,' G1,G2,hist,OK_Ang ',862 ! #wp WRITE(6,*) ikl,isn,' G1,G2,hist,OK_Ang ', 863 863 ! #wp. G1snSV(ikl,isn), G2snSV(ikl,isn),istoSV(ikl,isn),OK_Ang 864 864 ! #wp stop "Grains anguleux mal d?finis" … … 907 907 908 908 isn=1 909 if (dzsnSV(ikl,isn)>0.and.ro__SV(ikl,isn)>0) then909 if (dzsnSV(ikl,isn)>0.and.ro__SV(ikl,isn)>0) THEN 910 910 dzsnSV(ikl,isn) = dzsnSV(ikl,isn) +0.9999*(smb_old-smb_new) & 911 911 / ro__SV(ikl,isn) … … 923 923 924 924 925 end subroutinesisvat_gsn925 END SUBROUTINE sisvat_gsn -
LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis/sisvat_qsn.f90
r5113 r5116 61 61 use VARxSV 62 62 use VARySV 63 use surface_data, only: is_ok_slush,opt_runoff_ac63 use surface_data, ONLY: is_ok_slush,opt_runoff_ac 64 64 65 65 … … 86 86 ! + ================== 87 87 88 integer:: ikl ,isn !89 integer:: nh ! Non erodible Snow: up.lay.Index90 integer:: LayrOK ! 1 (0) if In(Above) Snow Pack91 integer:: k_face ! 1 (0) if Crystal(no) faceted92 integer:: LastOK ! 1 ==> 1! Snow Layer93 integer:: NOLayr ! 1 Layer Update94 integer:: noSnow(knonv) ! Nb of Layers Updater95 integer:: kSlush ! Slush Switch96 real:: dTSnow ! Temperature [C]97 real:: EExdum(knonv) ! Energy in Excess when no Snow98 real:: OKmelt ! 1 (0) if (no) Melting99 real:: EnMelt ! Energy in excess, for Melting100 real:: SnHLat ! Energy consumed in Melting101 real:: AdEnrg,B_Enrg ! Additional Energy from Vapor102 real:: dzVap0,dzVap1 ! Vaporized Thickness [m]103 real:: dzMelt(knonv) ! Melted Thickness [m]104 real:: rosDry ! Snow volumic Mass if no Water in105 real:: PorVol ! Pore volume106 real:: PClose ! Pore Hole Close OFF Switch107 real:: SGDiam ! Snow Grain Diameter108 real:: SGDmax ! Max. Snow Grain Diameter109 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 Freezing114 real:: WaFrez ! Water consumed in Melting115 real:: RapdOK ! 1. ==> Snow melts rapidly116 real:: ThinOK ! 1. ==> Snow Layer is thin117 real:: dzepsi ! Minim. Snow Layer Thickness (!)118 real:: dz_Min ! Minim. Snow Layer Thickness119 real:: z_Melt ! Last (thin) Layer Melting120 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:: Tsave88 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 127 128 128 ! +--OUTPUT of SISVAT Trace Statistics (see assignation in PHY_SISVAT) 129 129 ! + ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 130 integer:: isnnew,isinew,isnUpD,isnitr130 INTEGER :: isnnew,isinew,isnUpD,isnitr 131 131 132 132 ! OUTPUT in SISVAT at specified i,j,k,n (see assignation in PHY_SISVAT) … … 358 358 359 359 !XF 360 if(ro__SV(ikl,isn) >= roCdSV.and.ro__SV(ikl,1)<900) &360 IF(ro__SV(ikl,isn) >= roCdSV.and.ro__SV(ikl,1)<900) & 361 361 PClose = min(0.50,PClose * & 362 362 (1.-(ro_ice-ro__SV(ikl,isn))/(ro_ice-roCdSV))) … … 364 364 PClose = max(0.,min(1.,PClose)) 365 365 366 if(isn==1) then366 IF(isn==1) THEN 367 367 PClose = 1 368 368 ispiSV(ikl)= max(ispiSV(ikl),1) 369 369 endif 370 370 371 if(drr_SV(ikl) >0 .and.TsisSV(ikl,isn)>273.14) then371 IF(drr_SV(ikl) >0 .and.TsisSV(ikl,isn)>273.14) THEN 372 372 if((ro__SV(ikl,isn)>900.and.ro__SV(ikl,isn)<920).or. & 373 ro__SV(ikl,isn)>950) then373 ro__SV(ikl,isn)>950) THEN 374 374 dzsnSV(ikl,isn) = dzsnSV(ikl,isn)*ro__SV(ikl,isn)/ro_ice 375 375 ro__SV(ikl,isn) = ro_ice … … 384 384 ! . TsisSV(ikl,isn) >273.14 .and. 385 385 ! . TsisSV(ikl,isn+1)<273.15 .and. 386 ! . drr_SV(ikl) >0) then386 ! . drr_SV(ikl) >0) THEN 387 387 ! TsisSV(ikl,isn)=273.14 388 388 ! PClose = 1 … … 492 492 rusnew = rusnSV(ikl) * SWf_SV(ikl) 493 493 494 if(isnoSV(ikl)<=1 .OR. opt_runoff_ac) rusnew = 0.495 ! if(ivgtSV(ikl)>=1) rusnew = 0.494 IF(isnoSV(ikl)<=1 .OR. opt_runoff_ac) rusnew = 0. 495 !IF(ivgtSV(ikl)>=1) rusnew = 0. 496 496 497 497 ! #EU rusnew = 0. … … 525 525 ENDDO 526 526 527 if(zt<0.005+(TaT_SV(ikl)-TfSnow)/1000..and. &527 IF(zt<0.005+(TaT_SV(ikl)-TfSnow)/1000..and. & 528 528 isnoSV(ikl) >0 .and. & 529 529 TaT_SV(ikl) >=TfSnow .and. & 530 istoSV(ikl,isnoSV(ikl)) >1 ) then530 istoSV(ikl,isnoSV(ikl)) >1 ) THEN 531 531 DO isn=1,isnoSV(ikl) 532 532 drr_SV(ikl) = drr_SV(ikl) & … … 562 562 +zSlush ) & ! 563 563 / max(dzsnSV(ikl,isn) , epsi ) ! 564 if(ro_new<ro_Ice+20) then ! MAX 940kg/m3 !564 IF(ro_new<ro_Ice+20) then ! MAX 940kg/m3 ! 565 565 rusnSV(ikl) = rusnSV(ikl) - zSlush ! [mm] OR [kg/m2] 566 566 RuofSV(ikl,4)= max(0.,RuofSV(ikl,4) - zSlush/dt__SV) … … 715 715 716 716 717 end subroutinesisvat_qsn717 END SUBROUTINE sisvat_qsn -
LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis/sisvat_qso.f90
r5105 r5116 109 109 ! + ================== 110 110 111 integer:: isl ,jsl ,ist ,ikl !112 integer:: ikm ,ikp ,ik0 ,ik1 !113 integer:: ist__s,ist__w ! Soil/Water Body Identifier111 INTEGER :: isl ,jsl ,ist ,ikl ! 112 INTEGER :: ikm ,ikp ,ik0 ,ik1 ! 113 INTEGER :: ist__s,ist__w ! Soil/Water Body Identifier 114 114 ! #BP real z0soil ! Soil Surface Bumps Height [m] 115 115 ! #BP real z_Bump !(Partly)Bumpy Layers Height [m] … … 118 118 119 119 ! #BP real etBump(knonv) ! Bumps Layer Averaged Humidity 120 real:: etaMid ! Layer Interface's Humidity121 real:: Dhydif ! Hydraulic Diffusivity [m2/s]122 real:: eta__f ! Water Front Soil Water Content123 real:: Khyd_f ! Water Front Hydraulic Conduct.124 real:: Khydav ! Hydraulic Conductivity [m/s]125 real:: WgFlow ! Water gravitat. Flux [kg/m2/s]126 real:: Wg_MAX ! Water MAX.grav. Flux [kg/m2/s]127 real:: SatRat ! Saturation Flux [kg/m2/s]128 real:: WExces ! Saturat. Excess Flux [kg/m2/s]129 real:: SoRnOF(knonv) ! Soil Run OFF130 real:: Dhydtz(knonv,-nsol:0) ! Dhydif * dt / dz [m]131 real:: Elem_A,Elem_B,Elem_C ! Diagonal Coefficients132 real:: Diag_A(knonv,-nsol:0) ! A Diagonal133 real:: Diag_B(knonv,-nsol:0) ! B Diagonal134 real:: Diag_C(knonv,-nsol:0) ! C Diagonal135 real:: Term_D(knonv,-nsol:0) ! Independant Term136 real:: Aux__P(knonv,-nsol:0) ! P Auxiliary Variable137 real:: Aux__Q(knonv,-nsol:0) ! Q Auxiliary Variable138 real:: etaaux(knonv,-nsol:-nsol+1) ! Soil Water Content [m3/m3]139 real:: FreeDr ! Free Drainage Fraction (actual)140 real:: FreeD0 ! Free Drainage Fraction (1=Full)141 real:: aKdtSV3( 0:nsot, 0:nkhy) ! Khyd=a*eta+b: a * dt142 real:: bKdtSV3( 0:nsot, 0:nkhy) ! Khyd=a*eta+b: b * dt120 REAL :: etaMid ! Layer Interface's Humidity 121 REAL :: Dhydif ! Hydraulic Diffusivity [m2/s] 122 REAL :: eta__f ! Water Front Soil Water Content 123 REAL :: Khyd_f ! Water Front Hydraulic Conduct. 124 REAL :: Khydav ! Hydraulic Conductivity [m/s] 125 REAL :: WgFlow ! Water gravitat. Flux [kg/m2/s] 126 REAL :: Wg_MAX ! Water MAX.grav. Flux [kg/m2/s] 127 REAL :: SatRat ! Saturation Flux [kg/m2/s] 128 REAL :: WExces ! Saturat. Excess Flux [kg/m2/s] 129 REAL :: SoRnOF(knonv) ! Soil Run OFF 130 REAL :: Dhydtz(knonv,-nsol:0) ! Dhydif * dt / dz [m] 131 REAL :: Elem_A,Elem_B,Elem_C ! Diagonal Coefficients 132 REAL :: Diag_A(knonv,-nsol:0) ! A Diagonal 133 REAL :: Diag_B(knonv,-nsol:0) ! B Diagonal 134 REAL :: Diag_C(knonv,-nsol:0) ! C Diagonal 135 REAL :: Term_D(knonv,-nsol:0) ! Independant Term 136 REAL :: Aux__P(knonv,-nsol:0) ! P Auxiliary Variable 137 REAL :: Aux__Q(knonv,-nsol:0) ! Q Auxiliary Variable 138 REAL :: etaaux(knonv,-nsol:-nsol+1) ! Soil Water Content [m3/m3] 139 REAL :: FreeDr ! Free Drainage Fraction (actual) 140 REAL :: FreeD0 ! Free Drainage Fraction (1=Full) 141 REAL :: aKdtSV3( 0:nsot, 0:nkhy) ! Khyd=a*eta+b: a * dt 142 REAL :: bKdtSV3( 0:nsot, 0:nkhy) ! Khyd=a*eta+b: b * dt 143 143 144 144 ! Water (Mass) Budget … … 377 377 ikp = nkhy * eta_SV(ikl,isl+1) / etadSV(ist) 378 378 379 if(ikm<0.or.ik0<0.or.ikp<0)then379 IF(ikm<0.or.ik0<0.or.ikp<0)THEN 380 380 print *,"CRASH1 in sisvat_qso.f on pixel (i,j,n)", & 381 381 ii__SV(ikl),jj__SV(ikl),nn__SV(ikl) … … 423 423 ikp = nkhy * eta_SV(ikl,isl+1) / etadSV(ist) 424 424 425 if(ik0<0.or.ikp<0)then425 IF(ik0<0.or.ikp<0)THEN 426 426 print *,"CRASH2 in sisvat_qso.f on pixel (i,j,n)", & 427 427 ii__SV(ikl),jj__SV(ikl),nn__SV(ikl) … … 550 550 ! +--IO, for Verification 551 551 ! + ~~~~~~~~~~~~~~~~~~~~ 552 ! #WR write(6,6010)552 ! #WR WRITE(6,6010) 553 553 DO isl= 0,-nsol,-1 554 554 DO ikl= 1,knonv … … 557 557 Khydsv(ikl,isl) =(aKdtSV3(ist,ikp) *eta_SV(ikl,isl) & 558 558 +bKdtSV3(ist,ikp)) *2.0/dt__SV 559 ! #WR write(6,6011) ikl,isl,eta_SV(ikl,isl)*1.e3,559 ! #WR WRITE(6,6011) ikl,isl,eta_SV(ikl,isl)*1.e3, 560 560 ! #WR. ikp, aKdtSV3(ist,ikp),bKdtSV3(ist,ikp), 561 561 ! #WR. Khydsv(ikl,isl) … … 706 706 707 707 708 end subroutinesisvat_qso708 END SUBROUTINE sisvat_qso -
LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis/sisvat_sno_albedo.f90
r5113 r5116 67 67 use VARySV 68 68 use VARtSV 69 USE surface_data, only: iflag_albcalc,correc_alb69 USE surface_data, ONLY: iflag_albcalc,correc_alb 70 70 71 71 IMPLICIT NONE … … 73 73 74 74 ! + -- INPUT 75 integer:: jjtime75 INTEGER :: jjtime 76 76 77 77 ! +--Internal Variables 78 78 ! + ================== 79 79 80 real:: coalb1(knonv) ! weighted Coalbedo, Vis.81 real:: coalb2(knonv) ! weighted Coalbedo, nIR 182 real:: coalb3(knonv) ! weighted Coalbedo, nIR 283 real:: coalbm ! weighted Coalbedo, mean84 real:: sExt_1(knonv) ! Extinction Coeff., Vis.85 real:: sExt_2(knonv) ! Extinction Coeff., nIR 186 real:: sExt_3(knonv) ! Extinction Coeff., nIR 287 real:: SnOpSV(knonv, nsno) ! Snow Grain optical Size80 REAL :: coalb1(knonv) ! weighted Coalbedo, Vis. 81 REAL :: coalb2(knonv) ! weighted Coalbedo, nIR 1 82 REAL :: coalb3(knonv) ! weighted Coalbedo, nIR 2 83 REAL :: coalbm ! weighted Coalbedo, mean 84 REAL :: sExt_1(knonv) ! Extinction Coeff., Vis. 85 REAL :: sExt_2(knonv) ! Extinction Coeff., nIR 1 86 REAL :: sExt_3(knonv) ! Extinction Coeff., nIR 2 87 REAL :: SnOpSV(knonv, nsno) ! Snow Grain optical Size 88 88 ! #AG real agesno 89 89 90 integer:: isn ,ikl ,isn1, i91 real:: sbeta1,sbeta2,sbeta3,sbeta4,sbeta592 real:: AgeMax,AlbMin,HSnoSV,HIceSV,doptmx,SignG1,Sph_OK93 real:: dalbed,dalbeS,dalbeW94 real:: bsegal,czemax,csegal,csza95 real:: RoFrez,DiffRo,SignRo,SnowOK,OpSqrt96 real:: albSn1,albIc1,a_SnI1,a_SII197 real:: albSn2,albIc2,a_SnI2,a_SII298 real:: albSn3,albIc3,a_SnI3,a_SII399 real:: albSno,albIce,albSnI,albSII,albWIc,albmax100 real:: doptic,Snow_H,SIce_H,SnownH,SIcenH101 real:: exarg1,exarg2,exarg3,sign_0,sExt_0102 real:: albedo_old,albCor103 real:: ro_ave,dz_ave,minalb104 real:: l1min,l1max,l2min,l2max,l3min,l3max105 real:: l6min(6), l6max(6), albSn6(6), a_SII6(6)106 real:: lmintmp,lmaxtmp,albtmp90 INTEGER :: isn ,ikl ,isn1, i 91 REAL :: sbeta1,sbeta2,sbeta3,sbeta4,sbeta5 92 REAL :: AgeMax,AlbMin,HSnoSV,HIceSV,doptmx,SignG1,Sph_OK 93 REAL :: dalbed,dalbeS,dalbeW 94 REAL :: bsegal,czemax,csegal,csza 95 REAL :: RoFrez,DiffRo,SignRo,SnowOK,OpSqrt 96 REAL :: albSn1,albIc1,a_SnI1,a_SII1 97 REAL :: albSn2,albIc2,a_SnI2,a_SII2 98 REAL :: albSn3,albIc3,a_SnI3,a_SII3 99 REAL :: albSno,albIce,albSnI,albSII,albWIc,albmax 100 REAL :: doptic,Snow_H,SIce_H,SnownH,SIcenH 101 REAL :: exarg1,exarg2,exarg3,sign_0,sExt_0 102 REAL :: albedo_old,albCor 103 REAL :: ro_ave,dz_ave,minalb 104 REAL :: l1min,l1max,l2min,l2max,l3min,l3max 105 REAL :: l6min(6), l6max(6), albSn6(6), a_SII6(6) 106 REAL :: lmintmp,lmaxtmp,albtmp 107 107 108 108 ! +--Local DATA … … 464 464 ! prescription for each time step with NEMO values 465 465 466 ! #AO if (LSmask(ikl) .eq. 0 .and. coupling_ao .eq. .TRUE.) then467 ! #AO if (AOmask(ikl) .eq. 0) then466 ! #AO if (LSmask(ikl) .eq. 0 .and. coupling_ao .eq. .TRUE.) THEN 467 ! #AO if (AOmask(ikl) .eq. 0) THEN 468 468 ! #AO albisv(ikl) = (1.-AOmask(ikl))* albAOsisv(ikl) 469 469 ! #AO. +(AOmask(ikl)*albisv(ikl)) … … 569 569 570 570 571 end subroutinesnoptp571 END SUBROUTINE snoptp 572 572 573 573 -
LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis/sisvat_ts2.f90
r5113 r5116 95 95 ! + ================== 96 96 97 integer:: ig, jk, isl98 real:: mu99 real:: Tsrf(klonv) ! surface temperature as extrapolated from soil100 real:: mug(klonv) !hj coef top layers101 real:: ztherm_i(klonv), zdz2(klonv, -nsol:nsno), z1s102 real:: pfluxgrd(klonv), pcapcal(klonv), cal(klonv)103 real:: beta(klonv), dif_grnd(klonv)104 real:: C_coef(klonv, -nsol:nsno), D_coef(klonv, -nsol:nsno)97 INTEGER :: ig, jk, isl 98 REAL :: mu 99 REAL :: Tsrf(klonv) ! surface temperature as extrapolated from soil 100 REAL :: mug(klonv) !hj coef top layers 101 REAL :: ztherm_i(klonv), zdz2(klonv, -nsol:nsno), z1s 102 REAL :: pfluxgrd(klonv), pcapcal(klonv), cal(klonv) 103 REAL :: beta(klonv), dif_grnd(klonv) 104 REAL :: C_coef(klonv, -nsol:nsno), D_coef(klonv, -nsol:nsno) 105 105 106 106 REAL, DIMENSION(klonv) :: zx_mh, zx_nh, zx_oh … … 123 123 124 124 125 ! write(*,*)'T check'125 ! WRITE(*,*)'T check' 126 126 127 127 ! DO ig = 1,knonv … … 196 196 ! IF (knonv>0) THEN 197 197 ! DO ig=1,8 198 ! write(*,*)ig,'sisvat: Tsis ',TsisSV(ig,isnoSV(ig))199 ! write(*,*)'max-1 ',TsisSV(ig,isnoSV(ig)-1)200 ! write(*,*)'max-2 ',TsisSV(ig,isnoSV(ig)-2)201 ! write(*,*)'0 ',TsisSV(ig,0)202 !! write(*,*)min(max(isnoSV(ig),0),1),max(1-isnoSV(ig),0)198 ! WRITE(*,*)ig,'sisvat: Tsis ',TsisSV(ig,isnoSV(ig)) 199 ! WRITE(*,*)'max-1 ',TsisSV(ig,isnoSV(ig)-1) 200 ! WRITE(*,*)'max-2 ',TsisSV(ig,isnoSV(ig)-2) 201 ! WRITE(*,*)'0 ',TsisSV(ig,0) 202 !! WRITE(*,*)min(max(isnoSV(ig),0),1),max(1-isnoSV(ig),0) 203 203 ! ENDDO 204 204 ! END IF … … 216 216 217 217 IF (mug(ig) <= 0.05) THEN 218 write(*, *)'Attention mu low', mug(ig)218 WRITE(*, *)'Attention mu low', mug(ig) 219 219 ENDIF 220 220 IF (mug(ig) >= 0.98) THEN 221 write(*, *)'Attention mu high', mug(ig)221 WRITE(*, *)'Attention mu high', mug(ig) 222 222 ENDIF 223 223 … … 243 243 ENDDO 244 244 ENDDO 245 ! write(*,*)ig,'Tsis',TsisSV(ig,0)245 ! WRITE(*,*)ig,'Tsis',TsisSV(ig,0) 246 246 247 247 ! IF (indice == is_sic) THEN … … 329 329 DO ig = 1, knonv 330 330 IF (ps__SV(ig)<1.) THEN 331 ! write(*,*)'ig',ig,'ps',ps__SV(ig)331 ! WRITE(*,*)'ig',ig,'ps',ps__SV(ig) 332 332 ps__SV(ig) = max(ps__SV(ig), 1.e-8) 333 333 ENDIF 334 334 IF (p1l_SV(ig)<1.) THEN 335 ! write(*,*)'ig',ig,'p1l',p1l_SV(ig)335 ! WRITE(*,*)'ig',ig,'p1l',p1l_SV(ig) 336 336 p1l_SV(ig) = max(p1l_SV(ig), 1.e-8) 337 337 ENDIF 338 338 IF (TaT_SV(ig)<180.) THEN 339 ! write(*,*)'ig',ig,'TaT',TaT_SV(ig)339 ! WRITE(*,*)'ig',ig,'TaT',TaT_SV(ig) 340 340 TaT_SV(ig) = max(TaT_SV(ig), 180.) 341 341 ENDIF 342 342 IF (QaT_SV(ig)<1.e-8) THEN 343 ! write(*,*)'ig',ig,'QaT',QaT_SV(ig)343 ! WRITE(*,*)'ig',ig,'QaT',QaT_SV(ig) 344 344 QaT_SV(ig) = max(QaT_SV(ig), 1.e-8) 345 345 ENDIF 346 346 IF (Tsf_SV(ig)<100.) THEN 347 ! write(*,*)'ig',ig,'Tsf',Tsf_SV(ig)347 ! WRITE(*,*)'ig',ig,'Tsf',Tsf_SV(ig) 348 348 Tsf_SV(ig) = max(Tsf_SV(ig), 180.) 349 349 ENDIF 350 350 IF (Tsf_SV(ig)>500.) THEN 351 ! write(*,*)'ig',ig,'Tsf',Tsf_SV(ig)351 ! WRITE(*,*)'ig',ig,'Tsf',Tsf_SV(ig) 352 352 Tsf_SV(ig) = min(Tsf_SV(ig), 400.) 353 353 ENDIF 354 354 ! IF (Tsrf(ig).LT.1.) THEN 355 !! write(*,*)'ig',ig,'Tsrf',Tsrf(ig)355 !! WRITE(*,*)'ig',ig,'Tsrf',Tsrf(ig) 356 356 ! Tsrf(ig)=max(Tsrf(ig),TaT_SV(ig)-20.) 357 357 ! ENDIF 358 358 IF (cdH_SV(ig)<1.e-10) THEN 359 ! IF (ig.le.3) write(*,*)'ig',ig,'cdH',cdH_SV(ig)359 ! IF (ig.le.3) WRITE(*,*)'ig',ig,'cdH',cdH_SV(ig) 360 360 cdH_SV(ig) = .5 361 361 ENDIF … … 370 370 zx_qs = r2es * FOEEW(Tsf_SV(ig), zdelta) / ps__SV(ig) 371 371 zx_qs = MIN(0.5, zx_qs) 372 ! write(*,*)'zcor',retv*zx_qs372 !WRITE(*,*)'zcor',retv*zx_qs 373 373 zcor = 1. / (1. - retv * zx_qs) 374 374 zx_qs = zx_qs * zcor … … 470 470 ENDDO 471 471 472 end subroutinesisvat_ts2472 END SUBROUTINE sisvat_ts2 -
LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis/sisvat_tso.f90
r5113 r5116 94 94 ! + ================== 95 95 96 integer:: ikl ,isl ,jsl ,ist !97 integer:: ist__s,ist__w ! Soil/Water Body Identifier98 integer:: islsgn ! Soil/Snow Surfac.Identifier99 real:: eps__3 ! Arbitrary Low Number100 real:: etaMid,psiMid ! Layer Interface's Humidity101 real:: mu_eta ! Soil thermal Conductivity102 real:: mu_exp ! arg Soil thermal Conductivity103 real:: mu_min ! Min Soil thermal Conductivity104 real:: mu_max ! Max Soil thermal Conductivity105 real:: mu_sno(knonv),mu_aux ! Snow thermal Conductivity106 real:: mu__dz(knonv,-nsol:nsno+1) ! mu_(eta,sno) / dz107 real:: dtC_sv(knonv,-nsol:nsno) ! dt / C108 real:: IRs__D(knonv) ! UpwardIR Previous Iter.Contr.109 real:: dIRsdT(knonv) ! UpwardIR T Derivat.110 real:: f_HSHL(knonv) ! Factor common to HS and HL111 real:: dRidTs(knonv) ! d(Rib)/d(Ts)112 real:: HS___D(knonv) ! Sensible Heat Flux Atm.Contr.113 real:: f___HL(knonv) !114 real:: HL___D(knonv) ! Latent Heat Flux Atm.Contr.96 INTEGER :: ikl ,isl ,jsl ,ist ! 97 INTEGER :: ist__s,ist__w ! Soil/Water Body Identifier 98 INTEGER :: islsgn ! Soil/Snow Surfac.Identifier 99 REAL :: eps__3 ! Arbitrary Low Number 100 REAL :: etaMid,psiMid ! Layer Interface's Humidity 101 REAL :: mu_eta ! Soil thermal Conductivity 102 REAL :: mu_exp ! arg Soil thermal Conductivity 103 REAL :: mu_min ! Min Soil thermal Conductivity 104 REAL :: mu_max ! Max Soil thermal Conductivity 105 REAL :: mu_sno(knonv),mu_aux ! Snow thermal Conductivity 106 REAL :: mu__dz(knonv,-nsol:nsno+1) ! mu_(eta,sno) / dz 107 REAL :: dtC_sv(knonv,-nsol:nsno) ! dt / C 108 REAL :: IRs__D(knonv) ! UpwardIR Previous Iter.Contr. 109 REAL :: dIRsdT(knonv) ! UpwardIR T Derivat. 110 REAL :: f_HSHL(knonv) ! Factor common to HS and HL 111 REAL :: dRidTs(knonv) ! d(Rib)/d(Ts) 112 REAL :: HS___D(knonv) ! Sensible Heat Flux Atm.Contr. 113 REAL :: f___HL(knonv) ! 114 REAL :: HL___D(knonv) ! Latent Heat Flux Atm.Contr. 115 115 REAL :: TSurf0(knonv),dTSurf ! Previous Surface Temperature 116 real:: qsatsg(knonv) !,den_qs,arg_qs ! Soil Saturat. Spec. Humidity117 real:: dqs_dT(knonv) ! d(qsatsg)/dTv118 real:: Psi( knonv) ! 1st Soil Layer Water Potential119 real:: RHuSol(knonv) ! Soil Surface Relative Humidity120 real:: etaSol ! Soil Surface Humidity121 real:: d__eta ! Soil Surface Humidity Increm.122 real:: Elem_A,Elem_C ! Diagonal Coefficients123 real:: Diag_A(knonv,-nsol:nsno) ! A Diagonal124 real:: Diag_B(knonv,-nsol:nsno) ! B Diagonal125 real:: Diag_C(knonv,-nsol:nsno) ! C Diagonal126 real:: Term_D(knonv,-nsol:nsno) ! Independant Term127 real:: Aux__P(knonv,-nsol:nsno) ! P Auxiliary Variable128 real:: Aux__Q(knonv,-nsol:nsno) ! Q Auxiliary Variable129 real:: Ts_Min,Ts_Max ! Temperature Limits116 REAL :: qsatsg(knonv) !,den_qs,arg_qs ! Soil Saturat. Spec. Humidity 117 REAL :: dqs_dT(knonv) ! d(qsatsg)/dTv 118 REAL :: Psi( knonv) ! 1st Soil Layer Water Potential 119 REAL :: RHuSol(knonv) ! Soil Surface Relative Humidity 120 REAL :: etaSol ! Soil Surface Humidity 121 REAL :: d__eta ! Soil Surface Humidity Increm. 122 REAL :: Elem_A,Elem_C ! Diagonal Coefficients 123 REAL :: Diag_A(knonv,-nsol:nsno) ! A Diagonal 124 REAL :: Diag_B(knonv,-nsol:nsno) ! B Diagonal 125 REAL :: Diag_C(knonv,-nsol:nsno) ! C Diagonal 126 REAL :: Term_D(knonv,-nsol:nsno) ! Independant Term 127 REAL :: Aux__P(knonv,-nsol:nsno) ! P Auxiliary Variable 128 REAL :: Aux__Q(knonv,-nsol:nsno) ! Q Auxiliary Variable 129 REAL :: Ts_Min,Ts_Max ! Temperature Limits 130 130 ! #e1 real Exist0 ! Existing Layer Switch 131 real:: psat_wat, psat_ice, sp ! computation of qsat132 133 integer:: nt_srf,it_srf,itEuBk ! HL: Surface Scheme131 REAL :: psat_wat, psat_ice, sp ! computation of qsat 132 133 INTEGER :: nt_srf,it_srf,itEuBk ! HL: Surface Scheme 134 134 parameter(nt_srf=10) ! 10 before 135 real:: agpsrf,xgpsrf,dt_srf,dt_ver !136 real:: etaBAK(knonv) !137 real:: etaNEW(knonv) !138 real:: etEuBk(knonv) !139 real:: fac_dt(knonv),faceta(knonv) !140 real:: PsiArg(knonv),SHuSol(knonv) !135 REAL :: agpsrf,xgpsrf,dt_srf,dt_ver ! 136 REAL :: etaBAK(knonv) ! 137 REAL :: etaNEW(knonv) ! 138 REAL :: etEuBk(knonv) ! 139 REAL :: fac_dt(knonv),faceta(knonv) ! 140 REAL :: PsiArg(knonv),SHuSol(knonv) ! 141 141 142 142 … … 436 436 * exp (6827.*(1. /273.16-1./TsisSV(ikl,isl))) 437 437 438 if(TsisSV(ikl,isl)<=273.16) then438 IF(TsisSV(ikl,isl)<=273.16) THEN 439 439 qsatsg(ikl) = 0.622 * psat_ice / (sp - 0.378 * psat_ice) 440 440 else … … 466 466 END DO 467 467 468 if(ist__s==1) then ! to reduce computer time468 IF(ist__s==1) then ! to reduce computer time 469 469 470 470 DO it_srf=1,nt_srf ! … … 615 615 TsisSV(ikl,isl) = Aux__Q(ikl,isl) *TsisSV(ikl,isl+1) & 616 616 +TsisSV(ikl,isl) 617 if(isl==0.and.isnoSV(ikl)==0) then 618 617 IF(isl==0.and.isnoSV(ikl)==0) THEN 619 618 TsisSV(ikl,isl) = min(TaT_SV(ikl)+30,TsisSV(ikl,isl)) 620 619 TsisSV(ikl,isl) = max(TaT_SV(ikl)-30,TsisSV(ikl,isl)) … … 673 672 674 673 675 end subroutinesisvat_tso674 END SUBROUTINE sisvat_tso -
LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis/sisvat_weq.f90
r5105 r5116 30 30 31 31 32 character(len=6) :: labWEq33 integer:: istart32 CHARACTER(LEN=6) :: labWEq 33 INTEGER :: istart 34 34 35 35 logical :: logWEq … … 40 40 ! + ================ 41 41 42 integer:: ikl ,isn43 real:: SnoWEQ,IceWEQ42 INTEGER :: ikl ,isn 43 REAL :: SnoWEQ,IceWEQ 44 44 45 45 … … 85 85 86 86 !! IF (istart.eq.1) THEN 87 !! write(45,45)dahost,i___SV(lwriSV(1)),j___SV(lwriSV(1)),87 !! WRITE(45,45)dahost,i___SV(lwriSV(1)),j___SV(lwriSV(1)), 88 88 !! . n___SV(lwriSV(1)) 89 89 !! 45 format(a18,10('-'),'Pt.',3i4,60('-')) 90 90 !! END IF 91 91 92 !! write(45,450) labWEq,IceWEQ,iiceSV(ikl),SnoWEQ92 !! WRITE(45,450) labWEq,IceWEQ,iiceSV(ikl),SnoWEQ 93 93 !! . ,IceWEQ+SnoWEQ,isnoSV(ikl) 94 94 !! . ,drr_SV(ikl)*dt__SV … … 102 102 103 103 104 end subroutinesisvat_weq104 END SUBROUTINE sisvat_weq -
LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis/sisvat_zag.f90
r5105 r5116 63 63 ! + ----- 64 64 65 integer:: isagrb(knonv) ! 2nd Layer History66 real:: dzagrb(knonv) ! 2nd Layer Thickness67 real:: T_agrb(knonv) ! 2nd Layer Temperature68 real:: roagrb(knonv) ! 2nd Layer Density69 real:: etagrb(knonv) ! 2nd Layer Water Content70 real:: G1agrb(knonv) ! 2nd Layer Dendricity/Spher.71 real:: G2agrb(knonv) ! 2nd Layer Sphericity/Size72 real:: agagrb(knonv) ! 2nd Layer Age65 INTEGER :: isagrb(knonv) ! 2nd Layer History 66 REAL :: dzagrb(knonv) ! 2nd Layer Thickness 67 REAL :: T_agrb(knonv) ! 2nd Layer Temperature 68 REAL :: roagrb(knonv) ! 2nd Layer Density 69 REAL :: etagrb(knonv) ! 2nd Layer Water Content 70 REAL :: G1agrb(knonv) ! 2nd Layer Dendricity/Spher. 71 REAL :: G2agrb(knonv) ! 2nd Layer Sphericity/Size 72 REAL :: agagrb(knonv) ! 2nd Layer Age 73 73 74 74 … … 76 76 ! + ------------ 77 77 78 integer:: isagra(knonv) ! 1st Layer History79 real:: WEagra(knonv) ! 1st Layer Height [mm w.e.]80 real:: Agreg1(knonv) ! 1. ===> Agregates81 real:: dzagra(knonv) ! 1st Layer Thickness82 real:: T_agra(knonv) ! 1st Layer Temperature83 real:: roagra(knonv) ! 1st Layer Density84 real:: etagra(knonv) ! 1st Layer Water Content85 real:: G1agra(knonv) ! 1st Layer Dendricity/Spher.86 real:: G2agra(knonv) ! 1st Layer Sphericity/Size87 real:: agagra(knonv) ! 1st Layer Age78 INTEGER :: isagra(knonv) ! 1st Layer History 79 REAL :: WEagra(knonv) ! 1st Layer Height [mm w.e.] 80 REAL :: Agreg1(knonv) ! 1. ===> Agregates 81 REAL :: dzagra(knonv) ! 1st Layer Thickness 82 REAL :: T_agra(knonv) ! 1st Layer Temperature 83 REAL :: roagra(knonv) ! 1st Layer Density 84 REAL :: etagra(knonv) ! 1st Layer Water Content 85 REAL :: G1agra(knonv) ! 1st Layer Dendricity/Spher. 86 REAL :: G2agra(knonv) ! 1st Layer Sphericity/Size 87 REAL :: agagra(knonv) ! 1st Layer Age 88 88 89 89 … … 91 91 ! + ================== 92 92 93 integer:: ikl94 integer:: nh ! Averaged Snow History95 integer:: nh__OK ! 1=>Conserve Snow History96 real:: rh !97 real:: dz ! Thickness98 real:: dzro_1 ! Thickness X Density, Lay.199 real:: dzro_2 ! Thickness X Density, Lay.2100 real:: dzro ! Thickness X Density, Aver.101 real:: ro ! Averaged Density102 real:: wn ! Averaged Water Content103 real:: tn ! Averaged Temperature104 real:: ag ! Averaged Snow Age105 real:: SameOK ! 1. => Same Type of Grains106 real:: G1same ! Averaged G1, same Grains107 real:: G2same ! Averaged G2, same Grains108 real:: typ__1 ! 1. => Lay1 Type: Dendritic109 real:: zroNEW ! dz X ro, if fresh Snow110 real:: G1_NEW ! G1, if fresh Snow111 real:: G2_NEW ! G2, if fresh Snow112 real:: zroOLD ! dz X ro, if old Snow113 real:: G1_OLD ! G1, if old Snow114 real:: G2_OLD ! G2, if old Snow115 real:: SizNEW ! Size, if fresh Snow116 real:: SphNEW ! Spheric.,if fresh Snow117 real:: SizOLD ! Size, if old Snow118 real:: SphOLD ! Spheric.,if old Snow119 real:: Siz_av ! Averaged Grain Size120 real:: Sph_av ! Averaged Grain Spher.121 real:: Den_av ! Averaged Grain Dendr.122 real:: DendOK ! 1. => Average is Dendr.123 real:: G1diff ! Averaged G1, diff. Grains124 real:: G2diff ! Averaged G2, diff. Grains125 real:: G1 ! Averaged G1126 real:: G2 ! Averaged G293 INTEGER :: ikl 94 INTEGER :: nh ! Averaged Snow History 95 INTEGER :: nh__OK ! 1=>Conserve Snow History 96 REAL :: rh ! 97 REAL :: dz ! Thickness 98 REAL :: dzro_1 ! Thickness X Density, Lay.1 99 REAL :: dzro_2 ! Thickness X Density, Lay.2 100 REAL :: dzro ! Thickness X Density, Aver. 101 REAL :: ro ! Averaged Density 102 REAL :: wn ! Averaged Water Content 103 REAL :: tn ! Averaged Temperature 104 REAL :: ag ! Averaged Snow Age 105 REAL :: SameOK ! 1. => Same Type of Grains 106 REAL :: G1same ! Averaged G1, same Grains 107 REAL :: G2same ! Averaged G2, same Grains 108 REAL :: typ__1 ! 1. => Lay1 Type: Dendritic 109 REAL :: zroNEW ! dz X ro, if fresh Snow 110 REAL :: G1_NEW ! G1, if fresh Snow 111 REAL :: G2_NEW ! G2, if fresh Snow 112 REAL :: zroOLD ! dz X ro, if old Snow 113 REAL :: G1_OLD ! G1, if old Snow 114 REAL :: G2_OLD ! G2, if old Snow 115 REAL :: SizNEW ! Size, if fresh Snow 116 REAL :: SphNEW ! Spheric.,if fresh Snow 117 REAL :: SizOLD ! Size, if old Snow 118 REAL :: SphOLD ! Spheric.,if old Snow 119 REAL :: Siz_av ! Averaged Grain Size 120 REAL :: Sph_av ! Averaged Grain Spher. 121 REAL :: Den_av ! Averaged Grain Dendr. 122 REAL :: DendOK ! 1. => Average is Dendr. 123 REAL :: G1diff ! Averaged G1, diff. Grains 124 REAL :: G2diff ! Averaged G2, diff. Grains 125 REAL :: G1 ! Averaged G1 126 REAL :: G2 ! Averaged G2 127 127 128 128 … … 236 236 237 237 238 end subroutinesisvat_zag238 END SUBROUTINE sisvat_zag -
LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis/sisvat_zcr.f90
r5105 r5116 49 49 ! + ================== 50 50 ! + 51 integer:: ikl ,isn ,is0 ,is152 integer:: isno_1 ! Switch: ! Snow Layer over Ice53 real:: Dtyp_0,Dtyp_1 ! Snow Grains Difference Measure54 real:: DenSph ! 1. when contiguous spheric51 INTEGER :: ikl ,isn ,is0 ,is1 52 INTEGER :: isno_1 ! Switch: ! Snow Layer over Ice 53 REAL :: Dtyp_0,Dtyp_1 ! Snow Grains Difference Measure 54 REAL :: DenSph ! 1. when contiguous spheric 55 55 ! + ! and dendritic Grains 56 real:: DendOK ! 1. when dendritic Grains57 real:: dTypMx ! Grain Type Differ.58 real:: dTypSp ! Sphericity Weight59 real:: dTypRo ! Density Weight60 real:: dTypDi ! Grain Diam.Weight61 real:: dTypHi ! History Weight56 REAL :: DendOK ! 1. when dendritic Grains 57 REAL :: dTypMx ! Grain Type Differ. 58 REAL :: dTypSp ! Sphericity Weight 59 REAL :: dTypRo ! Density Weight 60 REAL :: dTypDi ! Grain Diam.Weight 61 REAL :: dTypHi ! History Weight 62 62 63 63 … … 173 173 ! + 174 174 175 end subroutinesisvat_zcr175 END SUBROUTINE sisvat_zcr -
LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis/sisvat_zsn.f90
r5113 r5116 52 52 use VARxSV 53 53 use VARySV 54 use surface_data, only: ok_zsn_ii54 use surface_data, ONLY: ok_zsn_ii 55 55 56 56 IMPLICIT NONE … … 60 60 ! + ================== 61 61 62 integer:: ikl ,isn ,i !63 64 65 66 67 68 69 integer:: NLay_s(knonv) ! Split Snow Layer Switch70 integer:: isagr1(knonv) ! 1st Layer History71 integer:: isagr2(knonv) ! 2nd Layer History72 integer:: LstLay ! 0 ====> isnoSV = 173 integer:: isno_n ! Snow Normal.Profile74 integer:: iice_n ! Ice Normal.Profile75 integer:: iiceOK ! Ice Switch76 integer:: icemix ! 0 ====> Agregated Snow+Ice=Snow62 INTEGER :: ikl ,isn ,i ! 63 64 65 66 67 68 69 INTEGER :: NLay_s(knonv) ! Split Snow Layer Switch 70 INTEGER :: isagr1(knonv) ! 1st Layer History 71 INTEGER :: isagr2(knonv) ! 2nd Layer History 72 INTEGER :: LstLay ! 0 ====> isnoSV = 1 73 INTEGER :: isno_n ! Snow Normal.Profile 74 INTEGER :: iice_n ! Ice Normal.Profile 75 INTEGER :: iiceOK ! Ice Switch 76 INTEGER :: icemix ! 0 ====> Agregated Snow+Ice=Snow 77 77 ! + ! 1 Ice 78 integer:: isn1 (knonv) ! 1st layer to stagger79 real:: staggr ! stagger Switch80 81 real:: WEagre(knonv) ! Snow Water Equivalent Thickness82 real:: dzthin(knonv) ! Thickness of the thinest layer83 real:: OKthin ! Swich ON a new thinest layer84 real:: dz_dif ! difference from ideal discret.85 real:: thickL ! Thick Layer Indicator86 real:: OK_ICE ! Swich ON uppermost Ice Layer87 88 real:: Agrege(knonv) ! 1. when Agregation constrained89 real:: dzepsi ! Min Single Snw Layer Thickness90 real:: dzxmin ! Min Acceptable Layer Thickness91 real:: dz_min ! Min Layer Thickness92 real:: dz_max ! Max Layer Thickness93 real:: dzagr1(knonv) ! 1st Layer Thickness94 real:: dzagr2(knonv) ! 2nd Layer Thickness95 real:: T_agr1(knonv) ! 1st Layer Temperature96 real:: T_agr2(knonv) ! 2nd Layer Temperature97 real:: roagr1(knonv) ! 1st Layer Density98 real:: roagr2(knonv) ! 2nd Layer Density99 real:: etagr1(knonv) ! 1st Layer Water Content100 real:: etagr2(knonv) ! 2nd Layer Water Content101 real:: G1agr1(knonv) ! 1st Layer Dendricity/Spher.102 real:: G1agr2(knonv) ! 2nd Layer Dendricity/Spher.103 real:: G2agr1(knonv) ! 1st Layer Sphericity/Size104 real:: G2agr2(knonv) ! 2nd Layer Sphericity/Size105 real:: agagr1(knonv) ! 1st Layer Age106 real:: agagr2(knonv) ! 2nd Layer Age78 INTEGER :: isn1 (knonv) ! 1st layer to stagger 79 REAL :: staggr ! stagger Switch 80 81 REAL :: WEagre(knonv) ! Snow Water Equivalent Thickness 82 REAL :: dzthin(knonv) ! Thickness of the thinest layer 83 REAL :: OKthin ! Swich ON a new thinest layer 84 REAL :: dz_dif ! difference from ideal discret. 85 REAL :: thickL ! Thick Layer Indicator 86 REAL :: OK_ICE ! Swich ON uppermost Ice Layer 87 88 REAL :: Agrege(knonv) ! 1. when Agregation constrained 89 REAL :: dzepsi ! Min Single Snw Layer Thickness 90 REAL :: dzxmin ! Min Acceptable Layer Thickness 91 REAL :: dz_min ! Min Layer Thickness 92 REAL :: dz_max ! Max Layer Thickness 93 REAL :: dzagr1(knonv) ! 1st Layer Thickness 94 REAL :: dzagr2(knonv) ! 2nd Layer Thickness 95 REAL :: T_agr1(knonv) ! 1st Layer Temperature 96 REAL :: T_agr2(knonv) ! 2nd Layer Temperature 97 REAL :: roagr1(knonv) ! 1st Layer Density 98 REAL :: roagr2(knonv) ! 2nd Layer Density 99 REAL :: etagr1(knonv) ! 1st Layer Water Content 100 REAL :: etagr2(knonv) ! 2nd Layer Water Content 101 REAL :: G1agr1(knonv) ! 1st Layer Dendricity/Spher. 102 REAL :: G1agr2(knonv) ! 2nd Layer Dendricity/Spher. 103 REAL :: G2agr1(knonv) ! 1st Layer Sphericity/Size 104 REAL :: G2agr2(knonv) ! 2nd Layer Sphericity/Size 105 REAL :: agagr1(knonv) ! 1st Layer Age 106 REAL :: agagr2(knonv) ! 2nd Layer Age 107 107 108 108 … … 130 130 131 131 DO ikl=1,knonv 132 if(isnoSV(ikl)<=2) dz_min=max(0.0050,dz_min)132 IF(isnoSV(ikl)<=2) dz_min=max(0.0050,dz_min) 133 133 134 134 dzepsi=0.0015 135 if(ro__SV(ikl,isnoSV(ikl))>920) dzepsi=0.0020135 IF(ro__SV(ikl,isnoSV(ikl))>920) dzepsi=0.0020 136 136 137 137 dzthin(ikl) = 0. ! Arbitrary unrealistic … … 240 240 DO ikl=1,knonv 241 241 isn = i_thin(ikl) 242 if(LIndsv(ikl)>0) isn=min(nsno-1,isn) ! cXF242 IF(LIndsv(ikl)>0) isn=min(nsno-1,isn) ! cXF 243 243 isagr1(ikl) = istoSV(ikl,isn) 244 244 isagr2(ikl) = istoSV(ikl,isn+LIndsv(ikl)) … … 410 410 411 411 isn=max(1,isnoSV(ikl)-3) 412 if(dzsnSV(ikl,isn)>0.30) then ! surface layer > 30cm412 IF(dzsnSV(ikl,isn)>0.30) then ! surface layer > 30cm 413 413 i_thin(ikl) = isn ! XF 04/07/2019 414 414 dzthin(ikl) = dzsnSV(ikl,isn) … … 416 416 417 417 isn=max(1,isnoSV(ikl)-2) 418 if(dzsnSV(ikl,isn)>0.10) then ! surface layer > 10cm418 IF(dzsnSV(ikl,isn)>0.10) then ! surface layer > 10cm 419 419 i_thin(ikl) = isn ! XF 04/07/2019 420 420 dzthin(ikl) = dzsnSV(ikl,isn) … … 422 422 423 423 isn=max(1,isnoSV(ikl)-1) 424 if(dzsnSV(ikl,isn)>0.05) then ! surface layer > 5cm424 IF(dzsnSV(ikl,isn)>0.05) then ! surface layer > 5cm 425 425 i_thin(ikl) = isn ! XF 04/07/2019 426 426 dzthin(ikl) = dzsnSV(ikl,isn) … … 428 428 429 429 isn=max(1,isnoSV(ikl)) 430 if(dzsnSV(ikl,isn)>0.02) then ! surface layer > 2cm430 IF(dzsnSV(ikl,isn)>0.02) then ! surface layer > 2cm 431 431 i_thin(ikl) = isn ! XF 04/07/2019 432 432 dzthin(ikl) = dzsnSV(ikl,isn) … … 572 572 DO ikl=1,knonv 573 573 isn = i_thin(ikl) 574 if(LIndsv(ikl)>0) isn=min(isn, nsno-1) !cXF574 IF(LIndsv(ikl)>0) isn=min(isn, nsno-1) !cXF 575 575 isagr1(ikl) = istoSV(ikl,isn) 576 576 isagr2(ikl) = istoSV(ikl,isn+LIndsv(ikl)) … … 598 598 ! + minimum uppermost layer thickness to guarantee a correct reproduction of the snow 599 599 ! + atmosphere coupling 600 if(dzsnSV(ikl,max(1,isnoSV(ikl)-0))>0.02 .or. & ! surface layers> 2-5-10600 IF(dzsnSV(ikl,max(1,isnoSV(ikl)-0))>0.02 .or. & ! surface layers> 2-5-10 601 601 dzsnSV(ikl,max(1,isnoSV(ikl)-1))>0.05 .or. & ! XF 04/07/2019 602 602 dzsnSV(ikl,max(1,isnoSV(ikl)-2))>0.10 .or. & 603 dzsnSV(ikl,max(1,isnoSV(ikl)-3))>0.30 ) then603 dzsnSV(ikl,max(1,isnoSV(ikl)-3))>0.30 )THEN 604 604 Agrege(ikl) = min(1, & 605 605 max(0, & … … 739 739 740 740 741 end subroutinesisvat_zsn741 END SUBROUTINE sisvat_zsn -
LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis/surf_inlandsis_mod.F90
r5113 r5116 165 165 real, dimension(klon) :: mean_dens 166 166 ! lat_scale : temperature lapse rate against latitude [K degree-1] 167 real:: lat_scale167 REAL :: lat_scale 168 168 ! sh_scale : temperature lapse rate against altitude [K km-1] 169 real:: sh_scale169 REAL :: sh_scale 170 170 ! variables for density profile 171 171 ! E0, E1 : exponent 172 real:: E0, E1172 REAL :: E0, E1 173 173 ! depth at which 550 kg m-3 is reached [m] 174 real:: z550174 REAL :: z550 175 175 ! depths of snow layers 176 real:: depth, snow_depth, distup176 REAL :: depth, snow_depth, distup 177 177 ! number of initial snow layers 178 integer:: nb_snow_layer178 INTEGER :: nb_snow_layer 179 179 ! For density calc. 180 real:: alpha0, alpha1, ln_smb180 REAL :: alpha0, alpha1, ln_smb 181 181 ! theoritical densities [kg m-3] 182 real:: rho0, rho1, rho1_550182 REAL :: rho0, rho1, rho1_550 183 183 ! constants for density profile 184 184 ! C0, C1 : constant, 0.07 for z <= 550 kg m-3 … … 228 228 klonv = klon 229 229 knonv = knon 230 write(*, *) 'ikl, lon and lat in INLANDSIS'230 WRITE(*, *) 'ikl, lon and lat in INLANDSIS' 231 231 232 232 DO ikl = 1, knon 233 233 i=ikl2i(ikl) 234 write(*, *) 'ikl=', ikl, 'rlon=', rlon(i), 'rlat=', rlat(i)234 WRITE(*, *) 'ikl=', ikl, 'rlon=', rlon(i), 'rlat=', rlat(i) 235 235 END DO 236 236 … … 260 260 ! +--Soil layer thickness . Compute soil discretization (as for LMDZ) 261 261 ! + ---------------------------------------------------------------- 262 ! write(*,'(/a)') 'Start SISVAT init: soil discretization ', nsoilmx262 ! WRITE(*,'(/a)') 'Start SISVAT init: soil discretization ', nsoilmx 263 263 CALL get_soil_levels(dz1, dz2, lambda) 264 264 … … 368 368 ! with a moist-adiabatic lapse rate of 5 °C km-1 everywhere except for Antarctica, 369 369 ! for Antarctica, a dry-adiabatic lapse rate of 9.8 °C km-1 is assumed. 370 if (lati > 60.) then370 if (lati > 60.) THEN 371 371 ! CA todo : add longitude bounds 372 372 ! Greenland mean temperature : function of altitude and latitude … … 379 379 ! surface density: Fausto et al. 2018, https://doi.org/10.3389/feart.2018.00051 380 380 mean_dens(ikl) = 315. 381 else if (lati < -60.) then381 else if (lati < -60.) THEN 382 382 ! Antarctica mean temperature : function of altitude and latitude 383 383 ! for altitudes 0. to 500. m, lat_scale varies from 1.3 to 0.6 °C °lat-1 … … 399 399 else 400 400 401 ! write(*, *) 'Attention: temperature initialization is only defined for Greenland and Antarctica'401 ! WRITE(*, *) 'Attention: temperature initialization is only defined for Greenland and Antarctica' 402 402 403 403 mean_dens(ikl) =350. … … 463 463 rho0 = exp(E0 * depth) / (rho_ice / mean_dens(ikl) - 1 + exp(E0 * depth)) * rho_ice 464 464 rho1 = exp(E1 * depth) / (rho_ice / mean_dens(ikl) - 1 + exp(E1 * depth)) * rho_ice 465 if (depth <= z550) then465 if (depth <= z550) THEN 466 466 ro__SV(ikl, isl) = exp(E0 * depth) / (rho_ice / mean_dens(ikl) - 1 + exp(E0 * depth)) * rho_ice 467 467 else … … 499 499 open(unit = un_outfor, status = 'replace', file = fn_outfor) 500 500 ikl = gp_outfor ! index sur la grille land ice 501 write(un_outfor, *) fn_outfor, ikl, dt__SV, rlon(ikl2i(ikl)), rlat(ikl2i(ikl))502 write(un_outfor, *) 'nsnow - albedo - z0m - z0h , dz [m,30], temp [K,41], rho [kg/m3,41], eta [kg/kg,41] &501 WRITE(un_outfor, *) fn_outfor, ikl, dt__SV, rlon(ikl2i(ikl)), rlat(ikl2i(ikl)) 502 WRITE(un_outfor, *) 'nsnow - albedo - z0m - z0h , dz [m,30], temp [K,41], rho [kg/m3,41], eta [kg/kg,41] & 503 503 G1 [-,30], G2 [-,30], agesnow [d,30], history [-,30], DOP [m,30]' 504 504 END IF … … 553 553 ! => Upper bound for eroded snow mass 554 554 ! uss_SV(ikl) = SLussl(i,j,n) ! u*qs* (only for Tv in sisvatesbl.f) 555 ! #BS if(dsn_SV(ikl)>eps12.and.erprev(i,j,n).gt.eps9) then555 ! #BS IF(dsn_SV(ikl)>eps12.and.erprev(i,j,n).gt.eps9) THEN 556 556 ! #BS dsnbSV(ikl) =1.0-min(qsHY(i,j,kB) !BS neglib. at kb ~100 magl) 557 557 ! #BS. /max(qshy(i,j,mz),eps9),unun) … … 564 564 ! will be used for characterizing the Buffer Layer 565 565 ! (see update of Bros_N, G1same, G2same, zroOLD, zroNEW) 566 ! #BS if(n==1) qbs_HY(i,j) = dsnbSV(ikl)566 ! #BS IF(n==1) qbs_HY(i,j) = dsnbSV(ikl) 567 567 qsnoSV(ikl) = snow_cont_air(ikl) 568 568 … … 642 642 toicSV(ikl) = toicSV(ikl) - sn_add 643 643 ELSE 644 write(*, *) 'Attention, bare ice... point ', ikl644 WRITE(*, *) 'Attention, bare ice... point ', ikl 645 645 isnoSV(ikl) = 1 646 646 istoSV(ikl, 1) = 0 … … 713 713 IF (ok_outfor) THEN 714 714 ikl= gp_outfor 715 write(un_outfor, *) '+++++++++++', rlon(ikl2i(ikl)), rlat(ikl2i(ikl)),alt(ikl),'+++++++++++'716 write(un_outfor, *) isnoSV(ikl), alb_SV(ikl), Z0m_SV(ikl), Z0h_SV(ikl),HSs_sv(ikl),HLs_sv(ikl),alb1(ikl),alb2(ikl)717 write(un_outfor, *) dzsnSV(ikl, :)718 write(un_outfor, *) TsisSV(ikl, :)719 write(un_outfor, *) ro__SV(ikl, :)720 write(un_outfor, *) eta_SV(ikl, :)721 write(un_outfor, *) G1snSV(ikl, :)722 write(un_outfor, *) G2snSV(ikl, :)723 write(un_outfor, *) agsnSV(ikl, :)724 write(un_outfor, *) istoSV(ikl, :)725 write(un_outfor, *) DOPsnSV(ikl, :)715 WRITE(un_outfor, *) '+++++++++++', rlon(ikl2i(ikl)), rlat(ikl2i(ikl)),alt(ikl),'+++++++++++' 716 WRITE(un_outfor, *) isnoSV(ikl), alb_SV(ikl), Z0m_SV(ikl), Z0h_SV(ikl),HSs_sv(ikl),HLs_sv(ikl),alb1(ikl),alb2(ikl) 717 WRITE(un_outfor, *) dzsnSV(ikl, :) 718 WRITE(un_outfor, *) TsisSV(ikl, :) 719 WRITE(un_outfor, *) ro__SV(ikl, :) 720 WRITE(un_outfor, *) eta_SV(ikl, :) 721 WRITE(un_outfor, *) G1snSV(ikl, :) 722 WRITE(un_outfor, *) G2snSV(ikl, :) 723 WRITE(un_outfor, *) agsnSV(ikl, :) 724 WRITE(un_outfor, *) istoSV(ikl, :) 725 WRITE(un_outfor, *) DOPsnSV(ikl, :) 726 726 ENDIF 727 727 … … 771 771 fz(rk) = fz1 * (dalph_soil**rk - 1.) / (dalph_soil - 1.) 772 772 773 ! write(*,*)'Start soil level computation'773 ! WRITE(*,*)'Start soil level computation' 774 774 !----------------------------------------------------------------------- 775 775 ! Calculation of some constants … … 1190 1190 ENDIF 1191 1191 ENDDO 1192 write(*, *)'Read ', fichnom, ' finished!!'1192 WRITE(*, *)'Read ', fichnom, ' finished!!' 1193 1193 1194 1194 !********************************************************************************* … … 1216 1216 1217 1217 END DO 1218 write(*, *)'Copy histo', ikl1218 WRITE(*, *)'Copy histo', ikl 1219 1219 1220 1220 DO isn = 1, isnoSV(ikl) !nsno 1221 1221 snopts = snopts + 1 1222 1222 IF (isto(i, isn) > 10.) THEN !hj check 1223 write(*, *)'Irregular isto', ikl, i, isn, isto(i, isn)1223 WRITE(*, *)'Irregular isto', ikl, i, isn, isto(i, isn) 1224 1224 isto(i, isn) = 1. 1225 1225 ENDIF … … 1245 1245 errro = errro + 1 1246 1246 ENDIF 1247 write(*, *)ikl, i, isn, Tsis(i, isn), G1sn(i, isn)1247 WRITE(*, *)ikl, i, isn, Tsis(i, isn), G1sn(i, isn) 1248 1248 G1snSV(ikl, isn) = G1sn(i, isn) ! [-] [-] 1249 1249 G2snSV(ikl, isn) = G2sn(i, isn) ! [-] [0.0001 m]
Note: See TracChangeset
for help on using the changeset viewer.
