Index: /trunk/LMDZ.COMMON/libf/evolution/changelog.txt =================================================================== --- /trunk/LMDZ.COMMON/libf/evolution/changelog.txt (revision 3330) +++ /trunk/LMDZ.COMMON/libf/evolution/changelog.txt (revision 3331) @@ -315,2 +315,5 @@ == 16/05/2024 == JBC Improvement of a flag-like variable (more robust as a logical) to know where co2 ice was at the beginning. + +== 16/05/2024 == JBC +Update of the layering algorithm + corrections of wrong lines commited in r3330. Index: /trunk/LMDZ.COMMON/libf/evolution/layering_mod.F90 =================================================================== --- /trunk/LMDZ.COMMON/libf/evolution/layering_mod.F90 (revision 3330) +++ /trunk/LMDZ.COMMON/libf/evolution/layering_mod.F90 (revision 3331) @@ -474,5 +474,5 @@ htend_dust = tend_dust/rho_dust -if (htend_dust < 0.) then ! Dust lifting +if (htend_dust < 0.) then ! Dust lifting only if (abs(htend_co2ice) > eps .or. abs(htend_h2oice) > eps) error stop 'Situation not managed: dust lifting + CO2/H2O ice condensation/sublimation!' write(*,'(a)') ' Stratification -> Dust lifting' @@ -506,10 +506,10 @@ else -!------ Dust sedimentation +!------ Dust sedimentation only if (abs(htend_co2ice) < eps .and. abs(htend_h2oice) < eps) then write(*,'(a)') ' Stratification -> Dust sedimentation' ! New stratum at the layering top by sedimentation of dust thickness = htend_dust/(1. - dry_regolith_porosity) - if (thickness > 0.) then ! Only if the stratum is builiding up + if (thickness > 0.) then ! Only if the stratum is building up if (new_str) then call add_stratum(this,thickness,this%top%top_elevation + thickness,0.,0.,htend_dust/thickness,dry_regolith_porosity) @@ -526,5 +526,5 @@ ! New stratum at the layering top by condensation of CO2 and H2O ice thickness = htend_co2ice/(1. - co2ice_porosity) + htend_h2oice/(1. - h2oice_porosity) + htend_dust - if (thickness > 0.) then ! Only if the stratum is builiding up + if (thickness > 0.) then ! Only if the stratum is building up if (new_str) then call add_stratum(this,thickness,this%top%top_elevation + thickness,htend_co2ice/thickness,htend_h2oice/thickness,htend_dust/thickness,1. - (htend_co2ice/thickness + htend_h2oice/thickness + htend_dust/thickness)) @@ -550,7 +550,8 @@ h2subl_tot = 0. call subl_co2ice_layering(this,current1,h2subl,h_co2ice_old,h_h2oice_old,h_dust_old,new_lag1) - else if (h_co2ice_old < eps) then ! There is nothing to sublimate so we move to the underlying stratum - current1 => current1%down - new_lag1 = .true. + else if (h_co2ice_old < eps) then ! There is nothing to sublimate so we stop here + !current1 => current1%down ! Move to the underlying stratum + !new_lag1 = .true. + exit else ! Only a fraction can sublimate and so we move to the underlying stratum h2subl = h_co2ice_old @@ -567,5 +568,5 @@ ! New stratum at the layering top by condensation of H2O ice thickness = htend_h2oice/(1. - h2oice_porosity) + htend_dust - if (thickness > 0.) then ! Only if the stratum is builiding up + if (thickness > 0.) then ! Only if the stratum is building up if (new_str) then call add_stratum(this,thickness,this%top%top_elevation + thickness,0.,htend_h2oice/thickness,htend_dust/thickness,1. - (htend_h2oice/thickness + htend_dust/thickness)) @@ -591,7 +592,8 @@ h2subl_tot = 0. call subl_co2ice_layering(this,current1,h2subl,h_co2ice_old,h_h2oice_old,h_dust_old,new_lag1) - else if (h_co2ice_old < eps) then ! There is nothing to sublimate so we move to the underlying stratum - current1 => current1%down - new_lag1 = .true. + else if (h_co2ice_old < eps) then ! There is nothing to sublimate so we stop here + !current1 => current1%down ! Move to the underlying stratum + !new_lag1 = .true. + exit else ! Only a fraction can sublimate and so we move to the underlying stratum h2subl = h_co2ice_old @@ -617,7 +619,8 @@ h2subl_tot = 0. call subl_h2oice_layering(this,current2,h2subl,h_co2ice_old,h_h2oice_old,h_dust_old,new_lag2) - else if (h_h2oice_old < eps) then ! There is nothing to sublimate so we move to the underlying stratum - current2 => current2%down - new_lag2 = .true. + else if (h_h2oice_old < eps) then ! There is nothing to sublimate so we stop here + !current2 => current1%down ! Move to the underlying stratum + !new_lag2 = .true. + exit else ! Only a fraction can sublimate and so we move to the underlying stratum h2subl = h_h2oice_old @@ -682,5 +685,4 @@ endif - ! New stratum = dust lag h_lag = hsubl_dust/(1. - dry_lag_porosity) Index: /trunk/LMDZ.COMMON/libf/evolution/pem.F90 =================================================================== --- /trunk/LMDZ.COMMON/libf/evolution/pem.F90 (revision 3330) +++ /trunk/LMDZ.COMMON/libf/evolution/pem.F90 (revision 3331) @@ -792,5 +792,5 @@ do islope = 1,nslope if (is_co2ice_ini(ig,islope) .and. co2_ice(ig,islope) < 1.e-10 .and. .not. co2ice_disappeared(ig,islope)) then ! co2 ice disappeared, look for closest point without co2ice - co2ice_disappeared = .true. + co2ice_disappeared(ig,islope) = .true. if (latitude_deg(ig) > 0) then do ig_loop = ig,ngrid @@ -816,5 +816,4 @@ enddo endif - is_co2ice_ini(ig,islope) = .false. if ((co2_ice(ig,islope) < 1.e-10) .and. (h2o_ice(ig,islope) > frost_albedo_threshold)) then albedo(ig,1,islope) = albedo_h2o_frost Index: /trunk/LMDZ.MARS/deftank/field_def_physics_mars.xml =================================================================== --- /trunk/LMDZ.MARS/deftank/field_def_physics_mars.xml (revision 3330) +++ /trunk/LMDZ.MARS/deftank/field_def_physics_mars.xml (revision 3331) @@ -542,6 +542,5 @@ - + unit="kg.m-2.s-1" />