Index: trunk/LMDZ.MARS/libf/phymars/soilwater.F90 =================================================================== --- trunk/LMDZ.MARS/libf/phymars/soilwater.F90 (revision 3120) +++ trunk/LMDZ.MARS/libf/phymars/soilwater.F90 (revision 3121) @@ -80,5 +80,5 @@ real, intent(in) :: pc(ngrid, nlayer) ! Coefficients zc real, intent(in) :: pd(ngrid, nlayer) ! Coefficients zd -real, intent(in) :: pdqsdifpot(ngrid, nq) ! Potential pdqsdif (without subsurface - atmosphere exchange) +real, intent(in) :: pdqsdifpot(ngrid) ! Potential pdqsdif (without subsurface - atmosphere exchange) real, intent(in) :: pplev(ngrid, nlayer+1) ! Atmospheric pressure levels @@ -90,5 +90,5 @@ ! ---------------------------------------------------------------------- real, intent(inout) :: pqsoil(ngrid, nsoil, nqsoil) ! Subsurface tracers (water vapor and ice) -real, intent(in) :: pqsurf(ngrid, nq) ! Water ice on the surface +real, intent(in) :: pqsurf(ngrid) ! Water ice on the surface ! (in kg.m - 3 : m - 3 of pore air for water vapor and m - 3 of regolith for water ice and adsorbed water) ! Outputs : @@ -128,5 +128,4 @@ real*8 :: alpha0(ngrid) real*8 :: beta0(ngrid) -real :: co2ice(ngrid) ! CO2 ice on the surface ! Adsorbtion/Desorption variables and parameters @@ -247,4 +246,5 @@ integer :: sublim_n ! number of sublimation loop runs integer :: satur_n ! number of saturation loop runs + @@ -286,5 +286,4 @@ ! ================= -co2ice = pqsurf(:,igcm_co2) if (firstcall_soil) then @@ -648,5 +647,5 @@ alpha0(ig) = porosity_ice_free(ig, 1) * pb(ig, 1) / delta0(ig) - beta0(ig) = ( pb(ig, 2) * pc(ig, 2) + pq(ig, 1, igcm_h2o_vap) * pa(ig, 1) + pqsurf(ig, igcm_h2o_ice) ) & + beta0(ig) = ( pb(ig, 2) * pc(ig, 2) + pq(ig, 1, igcm_h2o_vap) * pa(ig, 1) + pqsurf(ig) ) & / delta0(ig) @@ -818,5 +817,5 @@ flux(ig, 0) = porosity_ice_free(ig, 1) * pb(ig, 1) / ptimestep & * ( znsoil(ig, 1) / rho(ig) - beta0(ig) - alpha0(ig) * znsoil(ig, 1) / rho(ig) ) - elseif (pqsurf(ig, igcm_h2o_ice).gt.0.) then + elseif (pqsurf(ig).gt.0.) then ! assume that the surface is covered by waterice (if it is co2ice it should not call this subroutine at all) flux(ig, 0) = D_mid(ig, 1) / midlayer_dz(ig, 0) * ( znsoil(ig, 1) - qsat_surf(ig) * rhosurf(ig) ) @@ -824,5 +823,5 @@ ! check if the ground would take up water from the surface but there is none - if ((.not.exchange(ig)).and.(pqsurf(ig, igcm_h2o_ice).eq.0.).and.(flux(ig, 0).lt.0.)) then + if ((.not.exchange(ig)).and.(pqsurf(ig).eq.0.).and.(flux(ig, 0).lt.0.)) then flux(ig, 0) = 0. endif @@ -926,5 +925,5 @@ ! calculate the temporty mixing ratio above the surface zq1temp2(ig) = beta0(ig) + alpha0(ig) * znsoil(ig, 1) / rho(ig) - + ! calculate the flux from the subsurface zdqsdifrego(ig) = porosity_ice_free(ig, 1) * pb(ig, 1) / ptimestep * (zq1temp2(ig) - znsoil(ig, 1) / rho(ig) ) @@ -933,5 +932,5 @@ ! calculate the flux from the subsurface zdqsdifrego(ig) = -D_mid(ig, 1) / midlayer_dz(ig, 0) * (qsat_surf(ig) * rhosurf(ig) - znsoil(ig, 1)) ! faut - il faire intervenir la porosite ? - if ((pqsurf(ig, igcm_h2o_ice).eq.0.).and.(zdqsdifrego(ig).lt.0.)) then + if ((pqsurf(ig).eq.0.).and.(zdqsdifrego(ig).lt.0.)) then zdqsdifrego(ig) = 0. endif @@ -946,9 +945,9 @@ if ( (.not.exchange(ig)) & .and. ( -(zdqsdifrego(ig) * ptimestep) & - .gt.( pqsurf(ig, igcm_h2o_ice) + pdqsdifpot(ig, igcm_h2o_ice) * ptimestep) ) & - .and.( (pqsurf(ig, igcm_h2o_ice) + pdqsdifpot(ig, igcm_h2o_ice) * ptimestep).gt.0. ) ) then + .gt.( pqsurf(ig) + pdqsdifpot(ig) * ptimestep) ) & + .and.( (pqsurf(ig) + pdqsdifpot(ig) * ptimestep).gt.0. ) ) then ! calculate a new flux from the subsurface - zdqsdifrego(ig) = -( pqsurf(ig, igcm_h2o_ice) + pdqsdifpot(ig, igcm_h2o_ice) * ptimestep ) / ptimestep + zdqsdifrego(ig) = -( pqsurf(ig) + pdqsdifpot(ig) * ptimestep ) / ptimestep ! ! check case where there is CO2 ice on the surface but qsurf = 0 Index: trunk/LMDZ.MARS/libf/phymars/vdifc_mod.F =================================================================== --- trunk/LMDZ.MARS/libf/phymars/vdifc_mod.F (revision 3120) +++ trunk/LMDZ.MARS/libf/phymars/vdifc_mod.F (revision 3121) @@ -955,11 +955,7 @@ c make_tsub : sous pas de temps adaptatif c la subroutine est a la fin du fichier - if (adsorption_soil) then - nsubtimestep(:) = 1 - else call make_tsub(ngrid,pdtsrf(:,islope),zqsurf, & ptimestep,dtmax,watercaptag, & nsubtimestep) - endif c Calculation for turbulent exchange with the surface (for ice) c initialization of ztsrf, which is surface temperature in @@ -1030,5 +1026,4 @@ if (adsorption_soil) then - call soilwater(1,nlay,nq,nsoil, nqsoil, & ztsrf(ig),ptsoil(ig,:,islope),subtimestep, @@ -1047,5 +1042,4 @@ zdqsdif_tot(ig)= & -zqsurf(ig)/subtimestep - zdqsdif_tot(ig) = zdqsdif_surf(ig) else zdqsdif_tot(ig) = zdqsdif_surf(ig) + @@ -1055,5 +1049,5 @@ endif ! adsorption - if(.not.watercaptag(ig)) then + if(.not.watercaptag(ig).and.(.not.adsorption_soil)) then if ((-zdqsdif_tot(ig)*subtimestep) & .gt.(zqsurf(ig))) then @@ -1185,5 +1179,5 @@ c Actualisation de h2o_vap dans le premier niveau zq_tmp_vap(ig,1,igcm_h2o_vap)=zq1temp(ig) -c Take into account the H2O latent heat impact on the surface temperature +c Take into account the H2O latent heat impact on the surface temperature if (latentheat_surfwater) then lh=(2834.3-0.28*(ztsrf(ig)-To)- @@ -1198,5 +1192,4 @@ ENDDO c Subtimestep water budget : - ztsrf(ig) = ztsrf(ig)+(pdtsrf(ig,islope) & + zdtsrf(ig,islope))*subtimestep