Index: trunk/LMDZ.COMMON/libf/evolution/allocation.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/allocation.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/allocation.F90 (revision 4170) @@ -28,4 +28,5 @@ use orbit, only: ini_orbit, end_orbit use output, only: ini_output, end_output +use ice_table, only: ini_ice_table, end_ice_table use planet, only: ini_planet, end_planet @@ -80,4 +81,5 @@ call ini_orbit() call ini_output() +call ini_ice_table() END SUBROUTINE ini_allocation @@ -106,4 +108,5 @@ ! CODE ! ---- +call end_ice_table() call end_output() call end_orbit() Index: trunk/LMDZ.COMMON/libf/evolution/atmosphere.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/atmosphere.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/atmosphere.F90 (revision 4170) @@ -118,5 +118,5 @@ ! ! DESCRIPTION -! Setter for 'atmosphere' configuration parameters. +! Setter for "atmosphere" configuration parameters. ! ! AUTHORS & DATE @@ -155,5 +155,5 @@ ! ! DESCRIPTION -! Setter for 'set_CO2cond_ps_PCM' parameter. +! Setter for 'CO2cond_ps_PCM' parameter. ! ! AUTHORS & DATE Index: trunk/LMDZ.COMMON/libf/evolution/backup.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/backup.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/backup.F90 (revision 4170) @@ -36,5 +36,5 @@ ! ! DESCRIPTION -! Setter for backup configuration parameters. +! Setter for "backup" configuration parameters. ! ! AUTHORS & DATE @@ -70,5 +70,5 @@ !======================================================================= SUBROUTINE save_clim_state(h2o_ice,co2_ice,tsurf_avg,tsurf_dev,tsoil_avg,tsoil_dev,ps_avg,ps_dev,ps_avg_glob,ps_avg_glob_ini, & - icetable_depth,icetable_thickness,ice_porefilling,h2o_ads_reg,co2_ads_reg,layerings_map,backup_idt) + icetable_depth,icetable_thickness,ice_porefilling,h2oice_depth,flux_ssice_avg,h2o_ads_reg,co2_ads_reg,layerings_map,backup_idt) !----------------------------------------------------------------------- ! NAME @@ -95,4 +95,5 @@ use atmosphere, only: build4PCM_atmosphere use tracers, only: build4PCM_tracers, nq +use ice_table, only: build4PCM_ssice use clim_state_rec, only: write_restart, write_restartfi, write_restartevo use layered_deposits, only: layering @@ -104,15 +105,15 @@ ! ARGUMENTS ! --------- -real(dp), dimension(:), intent(in) :: ps_avg, ps_dev -real(dp), intent(in) :: ps_avg_glob, ps_avg_glob_ini -real(dp), dimension(:,:), intent(in) :: tsurf_avg, tsurf_dev, icetable_depth, icetable_thickness -real(dp), dimension(:,:,:), intent(in) :: tsoil_avg, tsoil_dev, ice_porefilling, h2o_ads_reg, co2_ads_reg -type(layering), dimension(:,:), intent(in) :: layerings_map -integer(di), optional, intent(in) :: backup_idt +real(dp), dimension(:), intent(in) :: ps_avg, ps_dev +real(dp), intent(in) :: ps_avg_glob, ps_avg_glob_ini +real(dp), dimension(:,:), intent(in) :: tsurf_avg, tsurf_dev, icetable_depth, icetable_thickness, h2oice_depth, flux_ssice_avg +real(dp), dimension(:,:,:), intent(in) :: tsoil_avg, tsoil_dev, ice_porefilling, h2o_ads_reg, co2_ads_reg +type(layering), dimension(:,:), intent(in) :: layerings_map +integer(di), intent(in), optional :: backup_idt real(dp), dimension(:,:), intent(inout) :: h2o_ice, co2_ice ! LOCAL VARIABLES ! --------------- -real(dp), dimension(ngrid,nslope) :: h2o_ice4PCM, co2_ice4PCM, tsurf4PCM, flux_geo4PCM, albedo4PCM, emissivity4PCM +real(dp), dimension(ngrid,nslope) :: h2o_ice4PCM, co2_ice4PCM, tsurf4PCM, flux_geo4PCM, albedo4PCM, emissivity4PCM, flux_ssice4PCM, h2oice_depth4PCM real(dp), dimension(ngrid,nlayer) :: teta4PCM, air_mass4PCM real(dp), dimension(ngrid,nsoil_PCM,nslope) :: tsoil4PCM, inertiesoil4PCM @@ -130,6 +131,11 @@ call build4PCM_tsurf(tsurf_avg,tsurf_dev,tsurf4PCM) -! Build soil for the PCM -if (do_soil) call build4PCM_soil(tsoil_avg,tsoil_dev,inertiesoil4PCM,tsoil4PCM,flux_geo4PCM) +if (do_soil) then + ! Build soil for the PCM + call build4PCM_soil(tsoil_avg,tsoil_dev,inertiesoil4PCM,tsoil4PCM,flux_geo4PCM) + + ! Build subsurface ice for the PCM + call build4PCM_ssice(icetable_depth,h2oice_depth,h2oice_depth4PCM) +end if ! Build atmosphere for the PCM @@ -144,5 +150,5 @@ ! Write restart files call write_restartevo(h2o_ice,co2_ice,tsoil_avg,TI,icetable_depth,icetable_thickness,ice_porefilling,h2o_ads_reg,co2_ads_reg,layerings_map) -call write_restartfi(is_h2o_perice,h2o_ice4PCM,co2_ice4PCM,tsurf4PCM,tsoil4PCM,inertiesoil4PCM,albedo4PCM,emissivity4PCM,flux_geo4PCM) +call write_restartfi(is_h2o_perice,h2o_ice4PCM,co2_ice4PCM,tsurf4PCM,tsoil4PCM,inertiesoil4PCM,albedo4PCM,emissivity4PCM,flux_geo4PCM,h2oice_depth4PCM) call write_restart(ps4PCM,pa4PCM,preff4PCM,q4PCM,teta4PCM,air_mass4PCM) @@ -191,6 +197,5 @@ ! CODE ! ---- -suffix = '' -suffix = suffix//'_ts'//int2str(backup_idt) +suffix = '_ts'//int2str(backup_idt) call print_msg('> Backup of "restart" files at dt = '//int2str(backup_idt),LVL_NFO) Index: trunk/LMDZ.COMMON/libf/evolution/changelog.txt =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/changelog.txt (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/changelog.txt (revision 4170) @@ -950,2 +950,8 @@ - Add lifecycle helpers for clear allocation/deallocation logic. - Simplify string suffix for slopes variables. + +== 03/04/2026 == JBC +- Deletion of 'flux_ssice' ('zqdsdif_tot') from the "startfi.nc" as it is not needed. +- Using the yearly average flux for the sublimating subsurface ice instead of the value at last timestep. +- Keeping a clear separation between the subsurface ice flux and the surface ice tendency. +- Making sure that subsurface ice depth is well given to the PCM at the end of the PEM (ice table VS layering). Index: trunk/LMDZ.COMMON/libf/evolution/clim_state_init.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/clim_state_init.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/clim_state_init.F90 (revision 4170) @@ -205,4 +205,5 @@ use surf_ice, only: set_is_h2o_perice_PCM, set_co2_perice_PCM use soil, only: set_TI_PCM, set_inertiedat_PCM +use ice_table, only: set_coef_ssdif_PCM use display, only: print_msg, LVL_NFO @@ -276,4 +277,7 @@ call set_tsoil_PCM(tmp3d) +call get_var_nc('coef_ssdif',tmp2d) +call set_coef_ssdif_PCM(tmp2d) + deallocate(tmp2d); allocate(tmp2d(ngrid,nsoil_PCM)) call get_var_nc('inertiedat',tmp2d) @@ -282,9 +286,4 @@ call get_var_nc('inertiesoil',tmp3d) call set_TI_PCM(tmp3d) - -! To do? -! h2oice_depth -! d_coef -! lag_co2_ice ! Close/Deallocate @@ -398,7 +397,7 @@ !!! transition delta = depth_breccia - TI(i,index_breccia + 1,islope) = sqrt((layer(index_breccia + 1) - layer(index_breccia))/ & - (((delta - layer(index_breccia))/(TI(i,index_breccia,islope)**2)) + & - ((layer(index_breccia + 1) - delta)/(TI_breccia**2)))) + TI(i,index_breccia + 1,islope) = sqrt((layer(index_breccia + 1) - layer(index_breccia))/ & + (((delta - layer(index_breccia))/(TI(i,index_breccia,islope)**2)) + & + ((layer(index_breccia + 1) - delta)/(TI_breccia**2)))) do k = index_breccia + 2,index_bedrock TI(i,k,islope) = TI_breccia @@ -411,7 +410,7 @@ !! transition delta = depth_bedrock - TI(i,index_bedrock + 1,islope) = sqrt((layer(index_bedrock + 1) - layer(index_bedrock))/ & - (((delta - layer(index_bedrock))/(TI(i,index_bedrock,islope)**2)) + & - ((layer(index_bedrock + 1) - delta)/(TI_breccia**2)))) + TI(i,index_bedrock + 1,islope) = sqrt((layer(index_bedrock + 1) - layer(index_bedrock))/ & + (((delta - layer(index_bedrock))/(TI(i,index_bedrock,islope)**2)) + & + ((layer(index_bedrock + 1) - delta)/(TI_breccia**2)))) do k = index_bedrock + 2,nsoil TI(i,k,islope) = TI_bedrock @@ -427,7 +426,7 @@ do i = 1,ngrid if (inertiedat(i,index_breccia) < TI_breccia) then - inertiedat(i,index_breccia + 1) = sqrt((layer(index_breccia + 1) - layer(index_breccia))/ & - (((delta - layer(index_breccia))/(inertiedat(i,index_breccia)**2)) + & - ((layer(index_breccia + 1) - delta)/(TI_breccia**2)))) + inertiedat(i,index_breccia + 1) = sqrt((layer(index_breccia + 1) - layer(index_breccia))/ & + (((delta - layer(index_breccia))/(inertiedat(i,index_breccia)**2)) + & + ((layer(index_breccia + 1) - delta)/(TI_breccia**2)))) do k = index_breccia + 2,index_bedrock inertiedat(i,k) = TI_breccia @@ -442,7 +441,7 @@ delta = depth_bedrock do i = 1,ngrid - inertiedat(i,index_bedrock + 1) = sqrt((layer(index_bedrock + 1) - layer(index_bedrock))/ & - (((delta - layer(index_bedrock))/(inertiedat(i,index_bedrock)**2)) + & - ((layer(index_bedrock + 1) - delta)/(TI_bedrock**2)))) + inertiedat(i,index_bedrock + 1) = sqrt((layer(index_bedrock + 1) - layer(index_bedrock))/ & + (((delta - layer(index_bedrock))/(inertiedat(i,index_bedrock)**2)) + & + ((layer(index_bedrock + 1) - delta)/(TI_bedrock**2)))) end do Index: trunk/LMDZ.COMMON/libf/evolution/clim_state_rec.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/clim_state_rec.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/clim_state_rec.F90 (revision 4170) @@ -180,5 +180,5 @@ !======================================================================= -SUBROUTINE write_restartfi(is_h2o_perice,h2o_ice4PCM,co2_ice4PCM,tsurf4PCM,tsoil4PCM,inertiesoil4PCM,albedo4PCM,emissivity4PCM,flux_geo4PCM) +SUBROUTINE write_restartfi(is_h2o_perice,h2o_ice4PCM,co2_ice4PCM,tsurf4PCM,tsoil4PCM,inertiesoil4PCM,albedo4PCM,emissivity4PCM,flux_geo4PCM,h2oice_depth4PCM) !----------------------------------------------------------------------- ! NAME @@ -208,5 +208,5 @@ ! ARGUMENTS ! --------- -real(dp), dimension(:,:), intent(in) :: h2o_ice4PCM, co2_ice4PCM, albedo4PCM, emissivity4PCM, tsurf4PCM, flux_geo4PCM +real(dp), dimension(:,:), intent(in) :: h2o_ice4PCM, co2_ice4PCM, albedo4PCM, emissivity4PCM, tsurf4PCM, flux_geo4PCM, h2oice_depth4PCM real(dp), dimension(:,:,:), intent(in) :: tsoil4PCM, inertiesoil4PCM logical(k4), dimension(:), intent(in) :: is_h2o_perice @@ -240,8 +240,8 @@ ! Orbital parameters (controle) -controle(18) = obliquity ! Obliquity +controle(18) = obliquity ! Obliquity controle(15) = perihelion ! Perihelion -controle(16) = aphelion ! Aphelion -controle(17) = date_peri ! Date of perihelion +controle(16) = aphelion ! Aphelion +controle(17) = date_peri ! Date of perihelion call put_var_nc('controle',controle) deallocate(controle) @@ -259,10 +259,7 @@ call put_var_nc('emis',emissivity4PCM,1) call put_var_nc('flux_geo',flux_geo4PCM,1) - -! h2oice_depth -! lag_co2_ice -! zdqsdif_ssi_tot -! d_coef - +call put_var_nc('h2oice_depth',h2oice_depth4PCM) + +! Close call close_nc('re'//startfi_name) Index: trunk/LMDZ.COMMON/libf/evolution/config.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/config.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/config.F90 (revision 4170) @@ -172,5 +172,5 @@ call getin('h2oice_huge_ini',h2oice_huge_ini_l) -threshold_h2oice_cap_l = 460._dp ! kg.m-2 (= 0.5 m) +threshold_h2oice_cap_l = 460._dp ! kg/m2 (= 0.5 m) call getin('threshold_h2oice_cap',threshold_h2oice_cap_l) @@ -186,5 +186,5 @@ call getin('do_layering',do_layering_l) -d_dust_l = 5.78e-2_dp ! kg.m-2.y-1 (= 1.e-9 kg.m-2.s-1) +d_dust_l = 5.78e-2_dp ! kg/m2/y (= 1.e-9 kg/m2/s) call getin('d_dust',d_dust_l) Index: trunk/LMDZ.COMMON/libf/evolution/display.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/display.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/display.F90 (revision 4170) @@ -50,5 +50,5 @@ ! ! DESCRIPTION -! Setter for 'display' configuration parameters. +! Setter for "display" configuration parameters. ! ! AUTHORS & DATE Index: trunk/LMDZ.COMMON/libf/evolution/evolution.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/evolution.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/evolution.F90 (revision 4170) @@ -77,5 +77,5 @@ ! ! DESCRIPTION -! Setter for 'evolution' configuration parameters. +! Setter for "evolution" configuration parameters. ! ! AUTHORS & DATE Index: trunk/LMDZ.COMMON/libf/evolution/glaciers.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/glaciers.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/glaciers.F90 (revision 4170) @@ -39,5 +39,5 @@ ! ! DESCRIPTION -! Setter for 'glaciers' configuration parameters. +! Setter for "glaciers" configuration parameters. ! ! AUTHORS & DATE Index: trunk/LMDZ.COMMON/libf/evolution/ice_table.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/ice_table.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/ice_table.F90 (revision 4170) @@ -25,6 +25,7 @@ ! PARAMETERS ! ---------- -logical(k4), protected :: icetable_equilibrium ! Flag to compute the icetable depth according to equilibrium -logical(k4), protected :: icetable_dynamic ! Flag to compute the icetable depth according to the dynamic method +logical(k4), protected :: icetable_equilibrium ! Flag to compute the icetable depth according to equilibrium +logical(k4), protected :: icetable_dynamic ! Flag to compute the icetable depth according to the dynamic method +real(dp), dimension(:,:), allocatable, protected :: coef_ssdif_PCM ! Diffusion coefficient for subsurface water contains @@ -32,4 +33,63 @@ !======================================================================= +SUBROUTINE ini_ice_table() +!----------------------------------------------------------------------- +! NAME +! ini_ice_table +! +! DESCRIPTION +! Initialize the parameters of module 'ice_table'. +! +! AUTHORS & DATE +! JB Clement, 03/2026 +! +! NOTES +! +!----------------------------------------------------------------------- + +! DEPENDENCIES +! ------------ +use geometry, only: ngrid, nslope + +! DECLARATION +! ----------- +implicit none + +! CODE +! ---- +allocate(coef_ssdif_PCM(ngrid,nslope)) +coef_ssdif_PCM(:,:) = 0._dp + +END SUBROUTINE ini_ice_table +!======================================================================= + +!======================================================================= +SUBROUTINE end_ice_table() +!----------------------------------------------------------------------- +! NAME +! end_ice_table +! +! DESCRIPTION +! Deallocate ice_table arrays. +! +! AUTHORS & DATE +! JB Clement, 03/2026 +! +! NOTES +! +!----------------------------------------------------------------------- + +! DECLARATION +! ----------- +implicit none + +! CODE +! ---- +if (allocated(coef_ssdif_PCM)) deallocate(coef_ssdif_PCM) + +END SUBROUTINE end_ice_table +!======================================================================= + +!======================================================================= SUBROUTINE set_ice_table_config(icetable_equilibrium_in,icetable_dynamic_in) !----------------------------------------------------------------------- @@ -38,5 +98,5 @@ ! ! DESCRIPTION -! Setter for 'ice_table' configuration parameters. +! Setter for "ice_table" configuration parameters. ! ! AUTHORS & DATE @@ -72,4 +132,35 @@ END SUBROUTINE set_ice_table_config +!======================================================================= + +!======================================================================= +SUBROUTINE set_coef_ssdif_PCM(coef_ssdif_PCM_in) +!----------------------------------------------------------------------- +! NAME +! set_coef_ssdif_PCM +! +! DESCRIPTION +! Setter for 'coef_ssdif_PCM' parameter. +! +! AUTHORS & DATE +! JB Clement, 03/2026 +! +! NOTES +! +!----------------------------------------------------------------------- + +! DECLARATION +! ----------- +implicit none + +! ARGUMENTS +! --------- +real(dp), dimension(:,:), intent(in) :: coef_ssdif_PCM_in + +! CODE +! ---- +coef_ssdif_PCM(:,:) = coef_ssdif_PCM_in(:,:) + +END SUBROUTINE set_coef_ssdif_PCM !======================================================================= @@ -609,3 +700,46 @@ !======================================================================= +!======================================================================= +SUBROUTINE build4PCM_ssice(icetable_depth,h2oice_depth,h2oice_depth4PCM) +!----------------------------------------------------------------------- +! NAME +! build4PCM_ssice +! +! DESCRIPTION +! Reconstructs subsurface ice for the PCM. +! +! AUTHORS & DATE +! JB Clement, 03/2026 +! +! NOTES +! +!----------------------------------------------------------------------- + +! DEPENDENCIES +! ------------ +use layered_deposits, only: do_layering +use display, only: print_msg, LVL_NFO + +! DECLARATION +! ----------- +implicit none + +! ARGUMENTS +! --------- +real(dp), dimension(:,:), intent(in) :: icetable_depth, h2oice_depth +real(dp), dimension(:,:), intent(out) :: h2oice_depth4PCM + +! CODE +! ---- +call print_msg('> Building subsurface ice for the PCM',LVL_NFO) +h2oice_depth4PCM(:,:) = -999. ! Default initialization (no subsurface ice) +if (do_layering) then + h2oice_depth4PCM(:,:) = h2oice_depth(:,:) +else if (icetable_equilibrium .or. icetable_dynamic) then + h2oice_depth4PCM(:,:) = icetable_depth(:,:) +end if + +END SUBROUTINE build4PCM_ssice +!======================================================================= + END MODULE ice_table Index: trunk/LMDZ.COMMON/libf/evolution/layered_deposits.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/layered_deposits.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/layered_deposits.F90 (revision 4170) @@ -33,5 +33,5 @@ logical(k4), protected, private :: impose_dust_ratio ! Impose dust-to-ice ratio instead of dust tendency (see 'dust2ice_ratio') real(dp), protected, private :: dust2ice_ratio ! Dust-to-ice ratio when ice condenses (10% by default) -real(dp), protected, private :: d_dust ! Tendency of dust [kg.m-2.y-1] +real(dp), protected, private :: d_dust ! Tendency of dust [kg/m2/y] real(dp), parameter, private :: dry_lag_porosity = 0.2_dp ! Porosity of dust lag real(dp), parameter, private :: wet_lag_porosity = 0.15_dp ! Porosity of water ice lag @@ -85,5 +85,5 @@ ! ! DESCRIPTION -! Setter for 'layered_deposits' configuration parameters. +! Setter for "layered_deposits" configuration parameters. ! ! AUTHORS & DATE @@ -749,5 +749,5 @@ h2o_ice(ig,islope) = layerings_map(ig,islope)%top%h_h2oice*rho_h2oice else - call subsurface_ice_layering(layerings_map(ig,islope),h2oice_depth(ig,islope),h2o_ice(ig,islope),co2_ice(ig,islope)) + call get_ssice_layering(layerings_map(ig,islope),h2oice_depth(ig,islope),h2o_ice(ig,islope),co2_ice(ig,islope)) end if end do @@ -1034,8 +1034,8 @@ !======================================================================= -SUBROUTINE subsurface_ice_layering(this,h_toplag,h2o_ice,co2_ice) -!----------------------------------------------------------------------- -! NAME -! subsurface_ice_layering +SUBROUTINE get_ssice_layering(this,h_toplag,h2o_ice,co2_ice) +!----------------------------------------------------------------------- +! NAME +! get_ssice_layering ! ! DESCRIPTION @@ -1095,5 +1095,5 @@ end if -END SUBROUTINE subsurface_ice_layering +END SUBROUTINE get_ssice_layering !======================================================================= Index: trunk/LMDZ.COMMON/libf/evolution/orbit.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/orbit.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/orbit.F90 (revision 4170) @@ -59,5 +59,5 @@ ! ! DESCRIPTION -! Setter for 'orbit' configuration parameters. +! Setter for "orbit" configuration parameters. ! ! AUTHORS & DATE Index: trunk/LMDZ.COMMON/libf/evolution/output.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/output.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/output.F90 (revision 4170) @@ -100,5 +100,5 @@ ! ! DESCRIPTION -! Setter for 'output' configuration parameters. +! Setter for "output" configuration parameters. ! ! AUTHORS & DATE Index: trunk/LMDZ.COMMON/libf/evolution/pem.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/pem.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/pem.F90 (revision 4170) @@ -52,5 +52,5 @@ use surf_ice, only: evolve_co2ice, evolve_h2oice, balance_h2oice_reservoirs use surf_temp, only: tsurf_PCM, adapt_tsurf2disappearedice -use tendencies, only: compute_tendice, evolve_tend_co2ice, evolve_tend_h2oice +use tendencies, only: compute_tendice, evolve_tend_co2ice, evolve_flux_ssice use tracers, only: adapt_tracers2pressure use utility, only: real2str @@ -72,8 +72,8 @@ logical(k4), dimension(:,:), allocatable :: is_co2ice_flow ! Flag for location of CO2 glacier flow logical(k4), dimension(:,:), allocatable :: is_h2oice_flow ! Flag for location of H2O glacier flow -real(dp), dimension(:,:,:), allocatable :: minPCM_h2operice ! Minimum of H2O perennial ice over the last PCM year [kg.m-2] -real(dp), dimension(:,:,:), allocatable :: minPCM_co2perice ! Minimum of CO2 perennial ice over the last PCM year [kg.m-2] -real(dp), dimension(:,:,:), allocatable :: minPCM_h2ofrost ! Minimum of H2O frost over the last PCM year [kg.m-2] -real(dp), dimension(:,:,:), allocatable :: minPCM_co2frost ! Minimum of CO2 frost over the last PCM year [kg.m-2] +real(dp), dimension(:,:,:), allocatable :: minPCM_h2operice ! Minimum of H2O perennial ice over the last PCM year [kg/m2] +real(dp), dimension(:,:,:), allocatable :: minPCM_co2perice ! Minimum of CO2 perennial ice over the last PCM year [kg/m2] +real(dp), dimension(:,:,:), allocatable :: minPCM_h2ofrost ! Minimum of H2O frost over the last PCM year [kg/m2] +real(dp), dimension(:,:,:), allocatable :: minPCM_co2frost ! Minimum of CO2 frost over the last PCM year [kg/m2] ! Surface-related: real(dp), dimension(:,:), allocatable :: tsurf_avg_yr1 ! Average surface temperature of the second to last PCM run [K] @@ -141,5 +141,5 @@ call load_xios_data(ps_avg,ps_ts,tsurf_avg,tsurf_avg_yr1,tsoil_avg,tsoil_ts,h2o_surfdensity_avg,h2o_soildensity_avg, & - q_h2o_ts,q_co2_ts,minPCM_h2operice,minPCM_co2perice,minPCM_h2ofrost,minPCM_co2frost) + q_h2o_ts,q_co2_ts,minPCM_h2operice,minPCM_co2perice,minPCM_h2ofrost,minPCM_co2frost,flux_ssice_avg) ! Initiate soil settings and TI @@ -169,7 +169,7 @@ call print_msg('> Computing surface ice tendencies',LVL_NFO) call compute_tendice(minPCM_h2operice,minPCM_h2ofrost,h2o_ice,d_h2oice) -call print_msg('H2O ice tendencies [kg/m2/yr] (min|max): '//real2str(minval(d_h2oice))//' | '//real2str(maxval(d_h2oice)),LVL_NFO) +call print_msg('H2O ice tendencies [kg/m2/y] (min|max): '//real2str(minval(d_h2oice))//' | '//real2str(maxval(d_h2oice)),LVL_NFO) call compute_tendice(minPCM_co2perice,minPCM_co2frost,co2_ice,d_co2ice) -call print_msg('CO2 ice tendencies [kg/m2/yr] (min|max): '//real2str(minval(d_co2ice))//' | '//real2str(maxval(d_co2ice)),LVL_NFO) +call print_msg('CO2 ice tendencies [kg/m2/y] (min|max): '//real2str(minval(d_co2ice))//' | '//real2str(maxval(d_co2ice)),LVL_NFO) deallocate(minPCM_h2operice,minPCM_co2perice,minPCM_h2ofrost,minPCM_co2frost) @@ -220,6 +220,4 @@ ! Conversion to surface ice call layering2surfice(layerings_map,h2o_ice,co2_ice,h2oice_depth) - ! We put the sublimating tendency coming from subsurface ice into the overall tendency - !where (h2oice_depth > 0. .and. zdqsdif_ssi_tot < -minieps) d_h2oice = zdqsdif_ssi_tot end if call allocate_loop_state() @@ -227,5 +225,5 @@ idt = 0 do while (n_yr_run < min(nmax_yr_runorb,nmax_yr_run) .and. n_yr_sim < ntot_yr_sim) - call print_msg('**** Iteration of the PEM run [plnt yr]: '//real2str(n_yr_run + dt),LVL_NFO) + call print_msg('**** Iteration of the PEM run [plnt y]: '//real2str(n_yr_run + dt),LVL_NFO) ! Evolve global surface pressure call evolve_pressure(d_co2ice,delta_co2_ads,do_sorption,ps_avg_glob_old,ps_avg_glob,ps_avg) @@ -237,5 +235,8 @@ allocate(zshift_surf(ngrid,nslope),zlag(ngrid,nslope)) if (do_layering) then + allocate(d_h2oice_new(ngrid,nslope)) h2oice_depth_old(:,:) = h2oice_depth(:,:) + ! The sublimating tendency coming from subsurface ice is given through the surface ice tendency + where (h2oice_depth(:,:) > 0. .and. flux_ssice_avg(:,:) < 0._dp) d_h2oice(:,:) = flux_ssice_avg(:,:) do islope = 1,nslope do i = 1,ngrid @@ -247,5 +248,4 @@ call layering2surfice(layerings_map,h2o_ice,co2_ice,h2oice_depth) ! Balance H2O ice reservoirs - allocate(d_h2oice_new(ngrid,nslope)) call stopping_crit_h2o(delta_h2o_ads,delta_icetable,h2o_ice,d_h2oice,S_atm_2_h2o,S_h2o_2_atm,S_atm_2_h2oice,S_h2oice_2_atm,stopcrit) call balance_h2oice_reservoirs(S_atm_2_h2o,S_h2o_2_atm,S_atm_2_h2oice,S_h2oice_2_atm,h2o_ice,d_h2oice,d_h2oice_new) @@ -253,4 +253,5 @@ else zlag(:,:) = 0._dp + zshift_surf(:,:) = 0._dp call evolve_h2oice(delta_h2o_ads,delta_icetable,h2o_ice,d_h2oice,zshift_surf,stopcrit) call evolve_co2ice(co2_ice,d_co2ice,zshift_surf) @@ -306,8 +307,8 @@ write(num,'("_slope",i2.2)') islope end if - call write_diagevo('h2oice'//trim(num),'H2O ice','kg.m-2',h2o_ice(:,islope),(/dim_ngrid/)) - call write_diagevo('co2ice'//trim(num),'CO2 ice','kg.m-2',co2_ice(:,islope),(/dim_ngrid/)) - call write_diagevo('d_h2oice'//trim(num),'H2O ice tendency','kg.m-2.yr-1',d_h2oice(:,islope),(/dim_ngrid/)) - call write_diagevo('d_co2ice'//trim(num),'CO2 ice tendency','kg.m-2.yr-1',d_co2ice(:,islope),(/dim_ngrid/)) + call write_diagevo('h2oice'//trim(num),'H2O ice','kg/m2',h2o_ice(:,islope),(/dim_ngrid/)) + call write_diagevo('co2ice'//trim(num),'CO2 ice','kg/m2',co2_ice(:,islope),(/dim_ngrid/)) + call write_diagevo('d_h2oice'//trim(num),'H2O ice tendency','kg/m2/y',d_h2oice(:,islope),(/dim_ngrid/)) + call write_diagevo('d_co2ice'//trim(num),'CO2 ice tendency','kg/m2/y',d_co2ice(:,islope),(/dim_ngrid/)) if (co2ice_flow) then call write_diagevo('Flow_co2ice'//trim(num),'CO2 ice flow location','T/F',merge(1._dp,0._dp,is_co2ice_flow(:,islope)),(/dim_ngrid/)) @@ -328,6 +329,6 @@ call write_diagevo('inertiesoil'//trim(num),'Thermal inertia','SI',TI(:,:,islope),(/dim_ngrid,dim_nsoil/)) if (do_sorption) then - call write_diagevo('co2_ads_reg'//trim(num),'CO2 adsorbed in regolith','kg.m-2',co2_ads_reg(:,:,islope),(/dim_ngrid,dim_nsoil/)) - call write_diagevo('h2o_ads_reg'//trim(num),'H2O adsorbed in regolith','kg.m-2',h2o_ads_reg(:,:,islope),(/dim_ngrid,dim_nsoil/)) + call write_diagevo('co2_ads_reg'//trim(num),'CO2 adsorbed in regolith','kg/m2',co2_ads_reg(:,:,islope),(/dim_ngrid,dim_nsoil/)) + call write_diagevo('h2o_ads_reg'//trim(num),'H2O adsorbed in regolith','kg/m2',h2o_ads_reg(:,:,islope),(/dim_ngrid,dim_nsoil/)) end if end if @@ -361,17 +362,9 @@ ! Evolve the tendencies call evolve_tend_co2ice(d_co2ice_ini,co2_ice,emissivity_PCM,q_co2_ts_ini,q_co2_ts,ps_ts,ps_avg_glob_ini,ps_avg_glob,d_co2ice) -!~ if (do_layering) then -!~ do i = 1,ngrid -!~ do islope = 1,nslope -!~ if (is_h2oice_sublim_ini(i,islope) .and. h2oice_depth(i,islope) > 0._dp) call evolve_tend_h2oice(h2oice_depth_old(i,islope),h2oice_depth(i,islope),tsurf_avg(i,islope),tsoil_ts_old(i,:,islope,:),tsoil_ts(i,:,islope,:),d_h2oice(i,islope)) -!~ end do -!~ end do -! else -! do i = 1,ngrid -! do islope = 1,nslope -! call evolve_tend_h2oice(icetable_depth_old(i,islope),icetable_depth(i,islope),tsurf_avg(i,islope),tsoil_ts_old(i,:,islope,:),tsoil_ts(i,:,islope,:),d_h2oice(i,islope)) -! end do -! end do -!~ end if + if (do_layering) then + call evolve_flux_ssice(h2oice_depth_old,h2oice_depth,tsurf_avg,tsoil_ts_old,tsoil_ts,flux_ssice_avg) + else if (icetable_equilibrium .or. icetable_dynamic) then + call evolve_flux_ssice(icetable_depth_old,icetable_depth,tsurf_avg,tsoil_ts_old,tsoil_ts,flux_ssice_avg) + end if ! Increment time @@ -383,5 +376,5 @@ if (backup_rate > 0) then if (mod(idt,backup_rate) == 0) call save_clim_state(h2o_ice,co2_ice,tsurf_avg,tsurf_dev,tsoil_avg,tsoil_dev,ps_avg,ps_dev,ps_avg_glob,ps_avg_glob_ini, & - icetable_depth,icetable_thickness,ice_porefilling,h2o_ads_reg,co2_ads_reg,layerings_map,idt) + icetable_depth,icetable_thickness,ice_porefilling,h2oice_depth,flux_ssice_avg,h2o_ads_reg,co2_ads_reg,layerings_map,idt) end if @@ -415,5 +408,5 @@ call save_clim_state(h2o_ice,co2_ice,tsurf_avg,tsurf_dev,tsoil_avg,tsoil_dev,ps_avg,ps_dev,ps_avg_glob,ps_avg_glob_ini, & - icetable_depth,icetable_thickness,ice_porefilling,h2o_ads_reg,co2_ads_reg,layerings_map) + icetable_depth,icetable_thickness,ice_porefilling,h2oice_depth,flux_ssice_avg,h2o_ads_reg,co2_ads_reg,layerings_map) ! Update the duration information of the workflow Index: trunk/LMDZ.COMMON/libf/evolution/planet.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/planet.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/planet.F90 (revision 4170) @@ -38,6 +38,6 @@ ! Ice-related: -real(dp), dimension(:,:), allocatable :: h2o_ice ! H2O ice [kg.m-2] -real(dp), dimension(:,:), allocatable :: co2_ice ! CO2 ice [kg.m-2] +real(dp), dimension(:,:), allocatable :: h2o_ice ! H2O ice [kg/m2] +real(dp), dimension(:,:), allocatable :: co2_ice ! CO2 ice [kg/m2] real(dp) :: h2oice_sublim_coverage_ini ! Initial surface area of sublimating H2O ice [m2] real(dp) :: co2ice_sublim_coverage_ini ! Initial surface area of sublimating CO2 ice [m2] @@ -73,4 +73,5 @@ real(dp), dimension(:,:), allocatable :: icetable_depth_old ! Old depth of the ice table [m] real(dp), dimension(:), allocatable :: delta_icetable ! Total mass of the H2O exchanged with the ice table [kg] +real(dp), dimension(:,:), allocatable :: flux_ssice_avg ! Average of total flux exchanged with subsurface ice [kg/m2/y] ! Tracer-related: @@ -153,4 +154,5 @@ allocate(q_h2o_ts(ngrid,nday)) allocate(q_co2_ts(ngrid,nday)) +allocate(flux_ssice_avg(ngrid,nslope)) ps_avg(:) = 0._dp @@ -163,4 +165,5 @@ q_h2o_ts(:,:) = 0._dp q_co2_ts(:,:) = 0._dp +flux_ssice_avg(:,:) = 0._dp END SUBROUTINE allocate_xios_state @@ -422,4 +425,5 @@ if (allocated(d_h2oice)) deallocate(d_h2oice) if (allocated(layerings_map)) deallocate(layerings_map) +if (allocated(flux_ssice_avg)) deallocate(flux_ssice_avg) END SUBROUTINE end_planet Index: trunk/LMDZ.COMMON/libf/evolution/soil.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/soil.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/soil.F90 (revision 4170) @@ -148,5 +148,5 @@ ! ! DESCRIPTION -! Setter for 'soil' configuration parameters. +! Setter for "soil" configuration parameters. ! ! AUTHORS & DATE Index: trunk/LMDZ.COMMON/libf/evolution/sorption.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/sorption.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/sorption.F90 (revision 4170) @@ -38,5 +38,5 @@ ! ! DESCRIPTION -! Setter for 'sorption' configuration parameters. +! Setter for "sorption" configuration parameters. ! ! AUTHORS & DATE Index: trunk/LMDZ.COMMON/libf/evolution/stopping_crit.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/stopping_crit.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/stopping_crit.F90 (revision 4170) @@ -55,5 +55,5 @@ ! ! DESCRIPTION -! Setter for 'stopping_crit' configuration parameters. +! Setter for "stopping_crit" configuration parameters. ! ! AUTHORS & DATE Index: trunk/LMDZ.COMMON/libf/evolution/surf_ice.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/surf_ice.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/surf_ice.F90 (revision 4170) @@ -26,8 +26,8 @@ real(dp), parameter :: rho_co2ice = 1650._dp ! Density of CO2 ice [kg.m-3] real(dp), parameter :: rho_h2oice = 920._dp ! Density of H2O ice [kg.m-3] -real(dp), protected :: threshold_h2oice_cap ! Threshold to consider H2O ice as infinite reservoir [kg.m-2] +real(dp), protected :: threshold_h2oice_cap ! Threshold to consider H2O ice as infinite reservoir [kg/m2] real(dp), protected :: h2oice_huge_ini ! Initial value for huge reservoirs of H2O ice [kg/m^2] -logical(k4), dimension(:), allocatable, protected :: is_h2o_perice_PCM ! Flag for the presence of perennial H2O ice in the PCM at the beginning [kg.m-2] -real(dp), dimension(:,:), allocatable, protected :: co2_perice_PCM ! Perennial CO2 ice in the PCM at the beginning [kg.m-2] +logical(k4), dimension(:), allocatable, protected :: is_h2o_perice_PCM ! Flag for the presence of perennial H2O ice in the PCM at the beginning [kg/m2] +real(dp), dimension(:,:), allocatable, protected :: co2_perice_PCM ! Perennial CO2 ice in the PCM at the beginning [kg/m2] contains @@ -103,5 +103,5 @@ ! ! DESCRIPTION -! Setter for 'surf_ice' configuration parameters. +! Setter for "surf_ice" configuration parameters. ! ! AUTHORS & DATE @@ -306,12 +306,12 @@ call print_msg('> Evolution of CO2 ice',LVL_NFO) -co2_ice_old = co2_ice -co2_ice = co2_ice + d_co2ice*dt +co2_ice_old(:,:) = co2_ice(:,:) +co2_ice(:,:) = co2_ice(:,:) + d_co2ice(:,:)*dt where (co2_ice < 0._dp) - co2_ice = 0._dp - d_co2ice = -co2_ice_old/dt + co2_ice(:,:) = 0._dp + d_co2ice(:,:) = -co2_ice_old(:,:)/dt end where -zshift_surf = d_co2ice*dt/rho_co2ice +zshift_surf(:,:) = zshift_surf(:,:) + d_co2ice(:,:)*dt/rho_co2ice END SUBROUTINE evolve_co2ice @@ -379,5 +379,5 @@ h2o_ice(:,:) = h2o_ice(:,:) + d_h2oice_new(:,:)*dt -zshift_surf(:,:) = d_h2oice_new(:,:)*dt/rho_h2oice +zshift_surf(:,:) = zshift_surf(:,:) + d_h2oice_new(:,:)*dt/rho_h2oice END SUBROUTINE evolve_h2oice Index: trunk/LMDZ.COMMON/libf/evolution/surface.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/surface.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/surface.F90 (revision 4170) @@ -91,5 +91,5 @@ if (allocated(albedo_PCM)) deallocate(albedo_PCM) if (allocated(albedodat_PCM)) deallocate(albedodat_PCM) -if (allocated(emissivity_PCM)) deallocate(emissivity_PCM ) +if (allocated(emissivity_PCM)) deallocate(emissivity_PCM) END SUBROUTINE end_surface Index: trunk/LMDZ.COMMON/libf/evolution/tendencies.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/tendencies.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/tendencies.F90 (revision 4170) @@ -114,5 +114,5 @@ ! ---- call print_msg("> Evolving the CO2 ice tendency according to the new pressure",LVL_NFO) -call print_msg('Old CO2 ice tendencies [kg/m2/yr] (min|max): '//real2str(minval(d_co2ice))//' | '//real2str(maxval(d_co2ice)),LVL_NFO) +call print_msg('Old CO2 ice tendencies [kg/m2/y] (min|max): '//real2str(minval(d_co2ice))//' | '//real2str(maxval(d_co2ice)),LVL_NFO) ! Compute volume mixing ratio @@ -135,5 +135,5 @@ end do end do -call print_msg('New CO2 ice tendencies [kg/m2/yr] (min|max): '//real2str( minval(d_co2ice))//' | '//real2str(maxval(d_co2ice)),LVL_NFO) +call print_msg('New CO2 ice tendencies [kg/m2/y] (min|max): '//real2str( minval(d_co2ice))//' | '//real2str(maxval(d_co2ice)),LVL_NFO) END SUBROUTINE evolve_tend_co2ice @@ -141,8 +141,8 @@ !======================================================================= -SUBROUTINE evolve_tend_h2oice(h2oice_depth_old,h2oice_depth_new,tsurf,tsoil_ts_old,tsoil_ts_new,d_h2oice) -!----------------------------------------------------------------------- -! NAME -! evolve_tend_h2oice +SUBROUTINE evolve_flux_ssice(h2oice_depth_old,h2oice_depth_new,tsurf,tsoil_ts_old,tsoil_ts_new,flux_ssice_avg) +!----------------------------------------------------------------------- +! NAME +! evolve_flux_ssice ! ! DESCRIPTION @@ -158,6 +158,8 @@ ! DEPENDENCIES ! ------------ +use geometry, only: ngrid, nslope use soil_temp, only: itp_tsoil use subsurface_ice, only: psv +use ice_table, only: coef_ssdif_PCM use display, only: print_msg, LVL_NFO @@ -168,46 +170,49 @@ ! ARGUMENTS ! --------- -real(dp), intent(in) :: h2oice_depth_old ! Old H2O ice depth -real(dp), intent(in) :: h2oice_depth_new ! New H2O ice depth -real(dp), intent(in) :: tsurf ! Surface temperature -real(dp), dimension(:,:), intent(in) :: tsoil_ts_old ! Old soil temperature time series -real(dp), dimension(:,:), intent(in) :: tsoil_ts_new ! New soil temperature time series -real(dp), intent(inout) :: d_h2oice ! Evolution of perennial ice +real(dp), dimension(:,:), intent(in) :: h2oice_depth_old ! Old H2O ice depth +real(dp), dimension(:,:), intent(in) :: h2oice_depth_new ! New H2O ice depth +real(dp), dimension(:,:), intent(in) :: tsurf ! Surface temperature +real(dp), dimension(:,:,:,:), intent(in) :: tsoil_ts_old ! Old soil temperature time series +real(dp), dimension(:,:,:,:), intent(in) :: tsoil_ts_new ! New soil temperature time series +real(dp), dimension(:,:), intent(inout) :: flux_ssice_avg ! Tendency of sub-surface ice ! LOCAL VARIABLES ! --------------- +real(dp), parameter :: zcdv = 0.0325_dp ! Drag coefficient real(dp) :: Rz_old, Rz_new, R_dec, hum_dec, psv_max_old, psv_max_new -integer(di) :: iday -real(dp), parameter :: coef_diff = 4.e-4_dp ! Diffusion coefficient -real(dp), parameter :: zcdv = 0.0325_dp ! Drag coefficient +integer(di) :: i, islope, iday ! CODE ! ---- -call print_msg("> Updating the H2O sub-surface ice tendency due to lag layer",LVL_NFO) - -! Higher resistance due to growing lag layer (higher depth) -Rz_old = h2oice_depth_old*zcdv/coef_diff ! Old resistance from PCM -Rz_new = h2oice_depth_new*zcdv/coef_diff ! New resistance based on new depth -R_dec = Rz_old/Rz_new ! Decrease because of resistance - -! The maxmimum of the daily averages over one year for the saturation vapor pressure at the ice table location -psv_max_old = 0._dp -psv_max_new = 0._dp -do iday = 1,size(tsoil_ts_old,2) - psv_max_old = max(psv_max_old,psv(itp_tsoil(tsoil_ts_old(:,iday),tsurf,h2oice_depth_old))) - psv_max_new = max(psv_max_new,psv(itp_tsoil(tsoil_ts_new(:,iday),tsurf,h2oice_depth_new))) -end do - -! Lower humidity due to growing lag layer (higher depth) -if (abs(psv_max_old) < tol) then - hum_dec = 1._dp -else - hum_dec = psv_max_new/psv_max_old ! Decrease because of lower water vapor pressure at the new depth -end if - -! Flux correction (decrease) -d_h2oice = d_h2oice*R_dec*hum_dec - -END SUBROUTINE evolve_tend_h2oice +call print_msg("> Updating the H2O sub-surface ice tendency due to dust lag layer thickness",LVL_NFO) + +do i = 1,ngrid + do islope = 1,nslope + ! Higher resistance due to growing lag layer (higher depth) + Rz_old = h2oice_depth_old(i,islope)*zcdv/coef_ssdif_PCM(i,islope) ! Old resistance from PCM + Rz_new = h2oice_depth_new(i,islope)*zcdv/coef_ssdif_PCM(i,islope) ! New resistance based on new depth + R_dec = Rz_old/Rz_new ! Decrease because of resistance + + ! The maxmimum of the daily averages over one year for the saturation vapor pressure at the ice table location + psv_max_old = 0._dp + psv_max_new = 0._dp + do iday = 1,size(tsoil_ts_old,2) + psv_max_old = max(psv_max_old,psv(itp_tsoil(tsoil_ts_old(i,:,islope,iday),tsurf(i,islope),h2oice_depth_old(i,islope)))) + psv_max_new = max(psv_max_new,psv(itp_tsoil(tsoil_ts_new(i,:,islope,iday),tsurf(i,islope),h2oice_depth_new(i,islope)))) + end do + + ! Lower humidity due to growing lag layer (higher depth) + if (abs(psv_max_old) < tol) then + hum_dec = 1._dp + else + hum_dec = psv_max_new/psv_max_old ! Decrease because of lower water vapor pressure at the new depth + end if + + ! Flux correction (decrease) + flux_ssice_avg(i,islope) = flux_ssice_avg(i,islope)*R_dec*hum_dec + enddo +enddo + +END SUBROUTINE evolve_flux_ssice !======================================================================= Index: trunk/LMDZ.COMMON/libf/evolution/workflow_status.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/workflow_status.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/workflow_status.F90 (revision 4170) @@ -152,6 +152,5 @@ n_id_PCM = n_pcm_runs end if -headline = '' -headline = headline//'cycle=PCM('//int2str(id_PCM_1) +headline = 'cycle=PCM('//int2str(id_PCM_1) do i = 1,n_id_PCM - 1 headline = headline//'+'//int2str(id_PCM_1 + i) Index: trunk/LMDZ.COMMON/libf/evolution/xios_data.F90 =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/xios_data.F90 (revision 4169) +++ trunk/LMDZ.COMMON/libf/evolution/xios_data.F90 (revision 4170) @@ -27,5 +27,5 @@ !======================================================================= SUBROUTINE load_xios_data(ps_avg,ps_ts,tsurf_avg,tsurf_avg_yr1,tsoil_avg,tsoil_ts,h2o_surfdensity_avg,h2o_soildensity_avg, & - q_h2o_ts,q_co2_ts,minPCM_h2operice,minPCM_co2perice,minPCM_h2ofrost,minPCM_co2frost) + q_h2o_ts,q_co2_ts,minPCM_h2operice,minPCM_co2perice,minPCM_h2ofrost,minPCM_co2frost,flux_ssice_avg) !----------------------------------------------------------------------- ! NAME @@ -59,5 +59,5 @@ ! --------- real(dp), dimension(:), intent(out) :: ps_avg -real(dp), dimension(:,:), intent(out) :: tsurf_avg, tsurf_avg_yr1, h2o_surfdensity_avg, ps_ts, q_h2o_ts, q_co2_ts +real(dp), dimension(:,:), intent(out) :: tsurf_avg, tsurf_avg_yr1, h2o_surfdensity_avg, ps_ts, q_h2o_ts, q_co2_ts, flux_ssice_avg real(dp), dimension(:,:,:), intent(out) :: tsoil_avg, h2o_soildensity_avg real(dp), dimension(:,:,:,:), intent(out) :: tsoil_ts @@ -130,4 +130,5 @@ call get_var_nc('perennial_co2ice'//trim(num),var_read_2d); call lonlat2vect(var_read_2d,minPCM_co2perice(:,islope,2)) if (do_soil) then + call get_var_nc('flux_ssice'//trim(num),var_read_2d); call lonlat2vect(var_read_2d,flux_ssice_avg(:,islope)) call get_var_nc('soiltemp'//trim(num),var_read_3d) do isoil = 1,nsoil_PCM Index: trunk/LMDZ.MARS/libf/phymars/phyetat0_mod.F90 =================================================================== --- trunk/LMDZ.MARS/libf/phymars/phyetat0_mod.F90 (revision 4169) +++ trunk/LMDZ.MARS/libf/phymars/phyetat0_mod.F90 (revision 4170) @@ -37,5 +37,5 @@ use comsoil_h, only: flux_geo use comslope_mod, only: nslope, major_slope -use paleoclimate_mod, only: paleoclimate, h2oice_depth, co2ice_depth, coef_ssdif, zdqsdif_ssi_tot +use paleoclimate_mod, only: paleoclimate, h2oice_depth, co2ice_depth, coef_ssdif use comcstfi_h, only: pi use geometry_mod, only: latitude @@ -914,12 +914,4 @@ h2oice_depth(:,:) = -999. endif - -! Total flux with SSI - call get_field("flux_ssice",zdqsdif_ssi_tot,found,indextime) - if (.not.found) then - write(*,*) "phyetat0: Failed loading : ", & - " is set as 0 (no subsurface ice)" - zdqsdif_ssi_tot(:,:) = 0. - endif ! Diffusion coeficent @@ -950,5 +942,4 @@ h2oice_depth(:,:) = -999. co2ice_depth(:,:) = -999. - zdqsdif_ssi_tot(:,:) = 0. coef_ssdif(:,:) = 4.e-4 perennial_co2ice(:,:) = 0. @@ -958,5 +949,4 @@ h2oice_depth(:,:) = -999. co2ice_depth(:,:) = -999. - zdqsdif_ssi_tot(:,:) = 0. coef_ssdif(:,:) = 4.e-4 perennial_co2ice(:,:) = 0. Index: trunk/LMDZ.MARS/libf/phymars/phyredem.F90 =================================================================== --- trunk/LMDZ.MARS/libf/phymars/phyredem.F90 (revision 4169) +++ trunk/LMDZ.MARS/libf/phymars/phyredem.F90 (revision 4170) @@ -210,5 +210,5 @@ use callkeys_mod, only: calltherm, dustinjection, calllott_nonoro use callkeys_mod, only: CLFvarying - use paleoclimate_mod, only: paleoclimate, h2oice_depth, co2ice_depth, coef_ssdif, zdqsdif_ssi_tot + use paleoclimate_mod, only: paleoclimate, h2oice_depth, co2ice_depth, coef_ssdif implicit none @@ -390,5 +390,4 @@ call put_field("h2oice_depth","Depth of the shallowest H2O ice layer",h2oice_depth) call put_field("co2ice_depth","Depth of the shallowest CO2 ice layer",co2ice_depth) - call put_field("flux_ssice","Total flux exchanged with subsurface water ice",zdqsdif_ssi_tot) call put_field("coef_ssdif","Diffusion coefficent for subsurface water",coef_ssdif) endif