Changeset 2260
- Timestamp:
- Mar 13, 2020, 7:32:42 PM (5 years ago)
- Location:
- trunk/LMDZ.MARS
- Files:
-
- 9 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.MARS/README
r2258 r2260 2888 2888 values. 2889 2889 2890 == 13/03/2020 == JN 2891 Adding a new 2D variable "watercap", perennial ground water ice. 2892 qsurf(ig,igcm_h2o_ice) can no longer be negative. 2893 When watercaptag=true, ie when there is an infinite amount of water avalaible for sublimation, 2894 and there is no more h2o frost, sublimation digs into watercap reservoir. 2895 For now watercap can't grow, even on watercaptag, and has same albedo than qsurf(h2o_ice). 2896 When using old start file, watercap is automatically initiated as 0. 2897 Added watercap in newstart.F aswell. 2898 2899 2900 -
trunk/LMDZ.MARS/libf/dynphy_lonlat/phymars/newstart.F
r2258 r2260 28 28 use surfdat_h, only: phisfi, z0, zmea, zstd, zsig, zgam, zthe, 29 29 & albedodat, z0_default, qsurf, tsurf, 30 & co2ice, emis, hmons, summit, base 30 & co2ice, emis, hmons, summit, base, watercap 31 31 use comsoil_h, only: inertiedat, layer, mlayer, nsoilmx, tsoil 32 32 use control_mod, only: day_step, iphysiq, anneeref, planet_type … … 418 418 & date,tsurf,tsoil,emis,q2, 419 419 & t,ucov,vcov,ps,co2ice,teta,phisold_newgrid,q,qsurf, 420 & tauscaling,totcloudfrac,surfith,nid )420 & tauscaling,totcloudfrac,surfith,nid,watercap) 421 421 write(*,*) "OK, read start_archive file" 422 422 ! copy soil thermal inertia … … 440 440 & day_ini,time,tsurf,tsoil,albedo,emis, 441 441 & q2,qsurf,co2ice,tauscaling,totcloudfrac, 442 & wstar,mem_Mccn_co2,mem_Nccn_co2,mem_Mh2o_co2 )442 & wstar,mem_Mccn_co2,mem_Nccn_co2,mem_Mh2o_co2,watercap) 443 443 444 444 ! copy albedo and soil thermal inertia … … 1752 1752 & tsurf,tsoil,co2ice,albedo,emis,q2,qsurf,tauscaling, 1753 1753 & totcloudfrac,wstar, 1754 & mem_Mccn_co2,mem_Nccn_co2,mem_Mh2o_co2 )1754 & mem_Mccn_co2,mem_Nccn_co2,mem_Mh2o_co2,watercap) 1755 1755 1756 1756 c======================================================================= -
trunk/LMDZ.MARS/libf/phymars/co2cloud.F
r2257 r2260 840 840 No=Nccnco2*tauscaling(ig) 841 841 Rn=-dlog(riceco2(ig,l)) 842 n_derf = erf( (rb_cldco2(1)+Rn) *dev2)842 n_derf = derf( (rb_cldco2(1)+Rn) *dev2) 843 843 Qext1bins2(ig,l)=0. 844 844 do i = 1, nbinco2_cld 845 845 n_aer(i) = -0.5 * No * n_derf !! this ith previously computed 846 n_derf = erf((rb_cldco2(i+1)+Rn) *dev2)846 n_derf = derf((rb_cldco2(i+1)+Rn) *dev2) 847 847 n_aer(i) = n_aer(i) + 0.5 * No * n_derf 848 848 Qext1bins2(ig,l)=Qext1bins2(ig,l)+Qext1bins(i)*n_aer(i) -
trunk/LMDZ.MARS/libf/phymars/improvedco2clouds_mod.F
r2257 r2260 552 552 c Ic_rice Mass transfer rate (kg/s) for a rice particle >0 si croissance ! 553 553 554 if ( (Ic_rice.ne.Ic_rice) ! will be true if it is Nan 555 & .or. Ic_rice .eq. 0.) then 554 if (isnan(Ic_rice) .or. Ic_rice .eq. 0.) then 556 555 Ic_rice=0. 557 556 subpdtcloudco2(ig,l)=-sum_subpdt(ig,l) -
trunk/LMDZ.MARS/libf/phymars/phyetat0_mod.F90
r2220 r2260 8 8 day_ini,time0,tsurf,tsoil,albedo,emis,q2,qsurf,co2ice, & 9 9 tauscaling,totcloudfrac,wstar,mem_Mccn_co2,mem_Nccn_co2,& 10 mem_Mh2o_co2 )10 mem_Mh2o_co2,watercap) 11 11 12 12 use tracer_mod, only: noms ! tracer names … … 61 61 real,intent(out) :: mem_Nccn_co2(ngrid,nlay) ! Memory of CCN number of H2O and dust used by CO2 62 62 real,intent(out) :: mem_Mh2o_co2(ngrid,nlay) ! Memory of H2O mass integred into CO2 crystal 63 real,intent(out) :: watercap(ngrid) ! h2o_ice_cover 63 64 !====================================================================== 64 65 ! Local variables: … … 446 447 write(*,*) 'phyetat0: loading surface tracer', & 447 448 ' h2o_ice instead of h2o_vap' 449 write(*,*) 'iq = ', iq 448 450 endif 449 451 call get_field(txt,qsurf(:,iq),found,indextime) … … 458 460 endif ! of if (nq.ge.1) 459 461 462 463 ! Surface water ice : 464 write(*,*) "qsurf(vap) : ", qsurf(:,2) 465 write(*,*) "qsurf(ice) : ", qsurf(:,3) 466 467 468 call get_field("watercap",watercap,found,indextime) 469 if (.not.found) then 470 write(*,*) "phyetat0: Failed loading <watercap> : ", & 471 "<watercap> is set to zero" 472 watercap(:)=0 473 else 474 write(*,*) "phyetat0: Surface water ice <watercap> range:", & 475 minval(watercap), maxval(watercap) 476 endif 477 478 write(*,*) 'Surface water ice cannot be allowed', & 479 ' to be negative if not on watercap !' 480 ! tracer on surface 481 if (nq.ge.1) then 482 do iq=1,nq 483 txt=noms(iq) 484 ! if ((txt.eq."h2o_vap").or.(txt.eq."h2o_ice")) then 485 if (txt.eq."h2o_ice") then 486 do ig=1,ngrid 487 if (qsurf(ig,iq).lt.0.0) then 488 watercap(ig) = qsurf(ig,iq) 489 qsurf(ig,iq) = 0.0 490 ! else if (qsurf(ig,iq)>0.5) then 491 ! watercap(ig) = watercap(ig) + qsurf(ig,iq) - 0.5 492 ! qsurf(ig,iq)=0.5 493 end if 494 end do 495 endif 496 ! if (txt.eq."h2o_vap") then 497 ! do ig=1,ngrid 498 ! qsurf(ig,iq)=0.0 499 ! else if (qsurf(ig,iq)>0.5) then 500 ! watercap(ig) = watercap(ig) + qsurf(ig,iq) - 0.5 501 ! qsurf(ig,iq)=0.5 502 ! end do 503 ! end if 504 enddo 505 endif 506 ! tracer on surface 507 !if (nq.ge.1) then 508 ! do iq=1,nq 509 ! txt=noms(iq) 510 ! if (txt.eq."h2o_vap") then 511 ! There is no surface tracer for h2o_vap; 512 ! "h2o_ice" should be loaded instead 513 ! txt="h2o_ice" 514 ! write(*,*) 'phyetat0: loading surface tracer', & 515 ! ' h2o_ice instead of h2o_vap' 516 !do ig=1,ngrid 517 ! if (qsurf(ig,igcm_h2o_vap).lt.0.0) then 518 ! watercap(ig) = watercap(ig) + qsurf(ig,igcm_h2o_vap) 519 ! qsurf(ig,igcm_h2o_vap)=0.0 520 ! else if (qsurf(ig,iq)>0.5) then 521 ! watercap(ig) = watercap(ig) + qsurf(ig,iq) - 0.5 522 ! qsurf(ig,iq)=0.5 523 ! end if 524 !end do 525 ! endif 526 ! enddo 527 !endif 528 529 530 531 !do ig=1,ngrid 532 ! if (qsurf(ig,'h2o_ice').lt.0.0) then 533 ! watercap(ig) = watercap(ig) + qsurf(ig,'h2o_ice') 534 ! qsurf(ig,'h2o_ice')=0.0 535 ! else if (qsurf(ig,iq)>0.5) then 536 ! watercap(ig) = watercap(ig) + qsurf(ig,iq) - 0.5 537 ! qsurf(ig,iq)=0.5 538 ! end if 539 !end do 540 460 541 ! Call to soil_settings, in order to read soil temperatures, 461 542 ! as well as thermal inertia and volumetric heat capacity -
trunk/LMDZ.MARS/libf/phymars/phyredem.F90
r2220 r2260 150 150 phystep,time,tsurf,tsoil,co2ice,albedo,emis,q2,qsurf,& 151 151 tauscaling,totcloudfrac,wstar, & 152 mem_Mccn_co2,mem_Nccn_co2,mem_Mh2o_co2 )152 mem_Mccn_co2,mem_Nccn_co2,mem_Mh2o_co2, watercap) 153 153 ! write time-dependent variable to restart file 154 154 use iostart, only : open_restartphy, close_restartphy, & … … 181 181 real,intent(in) :: mem_Nccn_co2(ngrid,nlay) ! CCN number of H2O and dust used by CO2 182 182 real,intent(in) :: mem_Mh2o_co2(ngrid,nlay) ! H2O mass integred into CO2 crystal 183 real,intent(in) :: watercap(ngrid) 183 184 184 185 integer :: iq … … 198 199 ! CO2 ice layer 199 200 call put_field("co2ice","CO2 ice cover",co2ice,time) 201 202 ! Water ice layer 203 call put_field("watercap","H2O ice cover",watercap,time) 200 204 201 205 ! Surface temperature -
trunk/LMDZ.MARS/libf/phymars/physiq_mod.F
r2252 r2260 45 45 & tsurf, co2ice, emis, 46 46 & capcal, fluxgrd, qsurf, 47 & hmons,summit,base 47 & hmons,summit,base,watercap 48 48 use comsaison_h, only: dist_sol, declin, mu0, fract, local_time 49 49 use slope_mod, only: theta_sl, psi_sl … … 331 331 real zdqmoldiff(ngrid,nlayer,nq) 332 332 333 REAL dwatercap(ngrid), dwatercap_dif(ngrid) ! (kg/m-2) 334 333 335 c Local variable for local intermediate calcul: 334 336 REAL zflubid(ngrid) … … 512 514 & q2,qsurf,co2ice,tauscaling,totcloudfrac,wstar, 513 515 & mem_Mccn_co2,mem_Nccn_co2, 514 & mem_Mh2o_co2 )516 & mem_Mh2o_co2,watercap) 515 517 516 518 if (pday.ne.day_ini) then … … 539 541 PRINT*,'check: qsurf ',qsurf(1,1),qsurf(ngrid,nq) 540 542 PRINT*,'check: co2 ',co2ice(1),co2ice(ngrid) 543 PRINT*,'check: watercap ',watercap(1),watercap(ngrid) 541 544 !!! 542 545 day_ini = pday … … 675 678 dsords(:,:)=0. 676 679 dsotop(:,:)=0. 680 dwatercap(:)=0 677 681 678 682 #ifdef DUSTSTORM … … 1241 1245 & zdqdif,zdqsdif,wstar,zcdv,zcdh,hfmax_th, 1242 1246 & zcondicea_co2microp,sensibFlux, 1243 & dustliftday,local_time )1247 & dustliftday,local_time,watercap,dwatercap_dif) 1244 1248 1245 1249 DO ig=1,ngrid 1246 1250 zdtsurf(ig)=zdtsurf(ig)+zdtsdif(ig) 1251 dwatercap(ig)=dwatercap(ig)+dwatercap_dif(ig) 1247 1252 ENDDO 1248 1253 … … 2008 2013 2009 2014 c----------------------------------------------------------------------- 2015 c J. Naar : Surface and sub-surface water ice 2016 c----------------------------------------------------------------------- 2017 c 2018 c 2019 c Increment Watercap (surface h2o reservoirs): 2020 c ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 2021 2022 DO ig=1,ngrid 2023 watercap(ig)=watercap(ig)+ptimestep*dwatercap(ig) 2024 ENDDO 2025 2026 c Test watercap 2027 c DO ig=1,ngrid 2028 c qsurf(iq,igcm_h2o_ice)=qsurf(iq,igcm_h2o_ice)+watercap(ig) 2029 c ENDDO 2030 2031 c----------------------------------------------------------------------- 2010 2032 c 13. Write output files 2011 2033 c ---------------------- … … 2185 2207 . tsurf,tsoil,co2ice,albedo,emis, 2186 2208 . q2,qsurf,tauscaling,totcloudfrac,wstar, 2187 . mem_Mccn_co2,mem_Nccn_co2,mem_Mh2o_co2 )2209 . mem_Mccn_co2,mem_Nccn_co2,mem_Mh2o_co2,watercap) 2188 2210 2189 2211 ENDIF … … 2368 2390 call wstats(ngrid,"co2ice","CO2 ice cover", 2369 2391 & "kg.m-2",2,co2ice) 2392 call wstats(ngrid,"watercap","H2O ice cover", 2393 & "kg.m-2",2,watercap) 2370 2394 call wstats(ngrid,"tauref","reference dod at 610 Pa","NU", 2371 2395 & 2,tauref) … … 2688 2712 call WRITEDIAGFI(ngrid,"co2ice","co2 ice thickness" 2689 2713 & ,"kg.m-2",2,co2ice) 2714 call WRITEDIAGFI(ngrid,"watercap","Water ice thickness" 2715 & ,"kg.m-2",2,watercap) 2690 2716 2691 2717 call WRITEDIAGFI(ngrid,"temp7","temperature in layer 7", -
trunk/LMDZ.MARS/libf/phymars/surfdat_h.F90
r2079 r2260 31 31 REAL,SAVE,ALLOCATABLE :: fluxgrd(:) ! surface conduction flux (W.m-2) 32 32 REAL,ALLOCATABLE,SAVE :: qsurf(:,:) ! tracer on surface (e.g. kg.m-2) 33 REAL,SAVE,ALLOCATABLE :: watercap(:) ! Surface water ice (kg.m-2) 33 34 34 35 contains … … 53 54 allocate(tsurf(ngrid)) 54 55 allocate(co2ice(ngrid)) 56 allocate(watercap(ngrid)) 55 57 allocate(emis(ngrid)) 56 58 allocate(capcal(ngrid)) … … 80 82 if (allocated(tsurf)) deallocate(tsurf) 81 83 if (allocated(co2ice)) deallocate(co2ice) 84 if (allocated(watercap)) deallocate(watercap) 82 85 if (allocated(emis)) deallocate(emis) 83 86 if (allocated(capcal)) deallocate(capcal) -
trunk/LMDZ.MARS/libf/phymars/vdifc_mod.F
r2218 r2260 13 13 $ pdqdif,pdqsdif,wstar,zcdv_true,zcdh_true, 14 14 $ hfmax,pcondicea_co2microp,sensibFlux, 15 $ dustliftday,local_time )15 $ dustliftday,local_time,watercap, dwatercap_dif) 16 16 17 17 use tracer_mod, only: noms, igcm_dust_mass, igcm_dust_number, … … 142 142 REAL kmixmin 143 143 144 c Argument added for surface water ice budget: 145 REAL,INTENT(IN) :: watercap(ngrid) 146 REAL,INTENT(OUT) :: dwatercap_dif(ngrid) 147 144 148 c Mass-variation scheme : 145 149 c ~~~~~~~ … … 213 217 pdqdif(1:ngrid,1:nlay,1:nq)=0 214 218 pdqsdif(1:ngrid,1:nq)=0 219 dwatercap_dif(1:ngrid)=0 215 220 216 221 c ** calcul de rho*dz et dt*rho/dz=dt*rho**2 g/dp … … 838 843 839 844 DO ig=1,ngrid 840 if(.not.watercaptag(ig)) then 841 if ((-pdqsdif(ig,igcm_h2o_ice)*ptimestep) 842 & .gt.pqsurf(ig,igcm_h2o_ice)) then 843 c write(*,*)'on sublime plus que qsurf!' 844 pdqsdif(ig,igcm_h2o_ice)= 845 if ((-pdqsdif(ig,igcm_h2o_ice)*ptimestep) 846 & .gt.(pqsurf(ig,igcm_h2o_ice))) then 847 c write(*,*)'subliming more than available frost: qsurf!' 848 if(watercaptag(ig)) then 849 c dwatercap_dif same sign as pdqsdif (pos. toward ground) 850 dwatercap_dif(ig) = -pdqsdif(ig,igcm_h2o_ice) 851 & - pqsurf(ig,igcm_h2o_ice)/ptimestep 852 pdqsdif(ig,igcm_h2o_ice)= 845 853 & -pqsurf(ig,igcm_h2o_ice)/ptimestep 854 c write(*,*) "subliming more than available frost: qsurf!" 855 else ! (case not.watercaptag(ig)) 856 pdqsdif(ig,igcm_h2o_ice)= 857 & -pqsurf(ig,igcm_h2o_ice)/ptimestep 858 dwatercap_dif(ig) = 0.0 859 860 endif 846 861 c write(*,*)'flux vers le sol=',pdqsdif(ig,nq) 847 862 z1(ig)=1./(za(ig,1)+ zb(ig,2)*(1.-zd(ig,2))) 848 863 zc(ig,1)=(za(ig,1)*zq(ig,1,igcm_h2o_vap)+ 849 $ zb(ig,2)*zc(ig,2) + 850 $ (-pdqsdif(ig,igcm_h2o_ice)) *ptimestep) *z1(ig) 864 $ zb(ig,2)*zc(ig,2) + (dwatercap_dif(ig)- 865 c $ zb(ig,2)*zc(ig,2) + (- 866 $ pdqsdif(ig,igcm_h2o_ice)) *ptimestep) *z1(ig) 851 867 zq1temp(ig)=zc(ig,1) 852 endif853 868 endif ! if (.not.watercaptag(ig)) 854 869 c Starting upward calculations for water : … … 871 886 lh=(2834.3-0.28*(tsrf_lh(ig)-To)- 872 887 & 0.004*(tsrf_lh(ig)-To)*(tsrf_lh(ig)-To))*1.e+3 873 pdtsrf(ig)= pdtsrf(ig) + pdqsdif(ig,igcm_h2o_ice)*lh874 & 888 pdtsrf(ig)= pdtsrf(ig) + (dwatercap_dif(ig)+ 889 & pdqsdif(ig,igcm_h2o_ice))*lh /pcapcal(ig) 875 890 endif 876 891 877 892 if(pqsurf(ig,igcm_h2o_ice) 878 & + pdqsdif(ig,igcm_h2o_ice)*ptimestep893 & +(dwatercap_dif(ig)+pdqsdif(ig,igcm_h2o_ice))*ptimestep 879 894 & .gt.frost_albedo_threshold) ! if there is still ice, T cannot exceed To 880 895 & pdtsrf(ig) = min(pdtsrf(ig),(To-ptsrf(ig))/ptimestep) ! ice melt case
Note: See TracChangeset
for help on using the changeset viewer.