Changeset 3493 for trunk/LMDZ.COMMON/libf/evolution/pem.F90
- Timestamp:
- Nov 5, 2024, 5:19:11 PM (3 weeks ago)
- File:
-
- 1 edited
-
trunk/LMDZ.COMMON/libf/evolution/pem.F90 (modified) (11 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.COMMON/libf/evolution/pem.F90
r3492 r3493 54 54 use orbit_param_criterion_mod, only: orbit_param_criterion 55 55 use recomp_orb_param_mod, only: recomp_orb_param 56 use ice_table_mod, only: porefillingice_depth, porefillingice_thickness, end_ice_table_porefilling, &57 ini_ice_table _porefilling, icetable_equilibrium, icetable_dynamic, computeice_table_equilibrium,compute_massh2o_exchange_ssi56 use ice_table_mod, only: icetable_depth, icetable_thickness, end_ice_table, ice_porefilling, & 57 ini_ice_table, icetable_equilibrium, icetable_dynamic, computeice_table_equilibrium, compute_massh2o_exchange_ssi 58 58 use soil_thermalproperties_mod, only: update_soil_thermalproperties 59 59 use time_phylmdz_mod, only: daysec, dtphys … … 204 204 205 205 ! Variables for surface and soil 206 real, dimension(:,:), allocatable :: tsurf_avg ! Physic x SLOPE field: Averaged Surface Temperature [K] 207 real, dimension(:,:,:), allocatable :: tsoil_avg ! Physic x SOIL x SLOPE field: Averaged Soil Temperature [K] 208 real, dimension(:,:,:), allocatable :: tsurf_PCM_timeseries ! ngrid x SLOPE x TIMES field: Surface Temperature in timeseries [K] 209 real, dimension(:,:,:,:), allocatable :: tsoil_phys_PEM_timeseries ! IG x SOIL x SLOPE x TIMES field: Non averaged Soil Temperature [K] 210 real, dimension(:,:,:,:), allocatable :: tsoil_anom ! IG x SOIL x SLOPE x TIMES field: Amplitude between instataneous and yearly average soil temperature at the beginning [K] 211 real, dimension(:,:), allocatable :: tsurf_avg_yr1 ! Physic x SLOPE field: Averaged Surface Temperature of first call of the PCM [K] 212 real, dimension(:,:), allocatable :: Tsoil_locslope ! Physic x Soil: Intermediate when computing Tsoil [K] 213 real, dimension(:), allocatable :: Tsurf_locslope ! Physic x Soil: Intermediate surface temperature to compute Tsoil [K] 214 real, dimension(:,:,:,:), allocatable :: watersoil_density_timeseries ! Physic x Soil x Slope x Times water soil density, time series [kg /m^3] 215 real, dimension(:,:), allocatable :: watersurf_density_avg ! Physic x Slope, water surface density, yearly averaged [kg/m^3] 216 real, dimension(:,:,:,:), allocatable :: watersoil_density_PEM_timeseries ! Physic x Soil x Slope x Times, water soil density, time series [kg/m^3] 217 real, dimension(:,:,:), allocatable :: watersoil_density_PEM_avg ! Physic x Soil x Slopes, water soil density, yearly averaged [kg/m^3] 218 real, dimension(:,:), allocatable :: Tsurfavg_before_saved ! Surface temperature saved from previous time step [K] 219 real, dimension(:), allocatable :: delta_co2_adsorbded ! Physics: quantity of CO2 that is exchanged because of adsorption / desorption [kg/m^2] 220 real, dimension(:), allocatable :: delta_h2o_adsorbded ! Physics: quantity of H2O that is exchanged because of adsorption / desorption [kg/m^2] 221 real :: totmassco2_adsorbded ! Total mass of CO2 that is exchanged because of adsorption / desoprtion over the planets [kg] 222 real :: totmassh2o_adsorbded ! Total mass of H2O that is exchanged because of adsorption / desoprtion over the planets [kg] 223 logical :: bool_sublim ! logical to check if there is sublimation or not 224 logical, dimension(:,:), allocatable :: co2ice_disappeared ! logical to check if a co2 ice reservoir already disappeared at a previous timestep 225 real, dimension(:,:), allocatable :: porefillingice_thickness_prev_iter ! ngrid x nslope: Thickness of the ice table at the previous iteration [m] 226 real, dimension(:), allocatable :: delta_h2o_icetablesublim ! ngrid x Total mass of the H2O that has sublimated / condenses from the ice table [kg] 227 228 ! Dynamic ice table 229 real, dimension(:,:), allocatable :: ss_ice_pf ! the amout of porefilling in each layer in each grid [m^3/m^3] 230 real, dimension(:,:), allocatable :: ss_ice_pf_new ! the amout of porefilling in each layer in each grid after the ss module[m^3/m^3] 231 real, dimension(:,:), allocatable :: porefillingice_depth_new ! new pure ice table depth 206 real, dimension(:,:), allocatable :: tsurf_avg ! Physic x SLOPE field: Averaged Surface Temperature [K] 207 real, dimension(:,:,:), allocatable :: tsoil_avg ! Physic x SOIL x SLOPE field: Averaged Soil Temperature [K] 208 real, dimension(:,:,:), allocatable :: tsurf_PCM_timeseries ! ngrid x SLOPE x TIMES field: Surface Temperature in timeseries [K] 209 real, dimension(:,:,:,:), allocatable :: tsoil_phys_PEM_timeseries ! IG x SOIL x SLOPE x TIMES field: Non averaged Soil Temperature [K] 210 real, dimension(:,:,:,:), allocatable :: tsoil_anom ! IG x SOIL x SLOPE x TIMES field: Amplitude between instataneous and yearly average soil temperature at the beginning [K] 211 real, dimension(:,:), allocatable :: tsurf_avg_yr1 ! Physic x SLOPE field: Averaged Surface Temperature of first call of the PCM [K] 212 real, dimension(:,:), allocatable :: Tsoil_locslope ! Physic x Soil: Intermediate when computing Tsoil [K] 213 real, dimension(:), allocatable :: Tsurf_locslope ! Physic x Soil: Intermediate surface temperature to compute Tsoil [K] 214 real, dimension(:,:,:,:), allocatable :: watersoil_density_timeseries ! Physic x Soil x Slope x Times water soil density, time series [kg /m^3] 215 real, dimension(:,:), allocatable :: watersurf_density_avg ! Physic x Slope, water surface density, yearly averaged [kg/m^3] 216 real, dimension(:,:,:,:), allocatable :: watersoil_density_PEM_timeseries ! Physic x Soil x Slope x Times, water soil density, time series [kg/m^3] 217 real, dimension(:,:,:), allocatable :: watersoil_density_PEM_avg ! Physic x Soil x Slopes, water soil density, yearly averaged [kg/m^3] 218 real, dimension(:,:), allocatable :: Tsurfavg_before_saved ! Surface temperature saved from previous time step [K] 219 real, dimension(:), allocatable :: delta_co2_adsorbded ! Physics: quantity of CO2 that is exchanged because of adsorption / desorption [kg/m^2] 220 real, dimension(:), allocatable :: delta_h2o_adsorbded ! Physics: quantity of H2O that is exchanged because of adsorption / desorption [kg/m^2] 221 real :: totmassco2_adsorbded ! Total mass of CO2 that is exchanged because of adsorption / desoprtion over the planets [kg] 222 real :: totmassh2o_adsorbded ! Total mass of H2O that is exchanged because of adsorption / desoprtion over the planets [kg] 223 logical :: bool_sublim ! logical to check if there is sublimation or not 224 logical, dimension(:,:), allocatable :: co2ice_disappeared ! logical to check if a co2 ice reservoir already disappeared at a previous timestep 225 real, dimension(:,:), allocatable :: icetable_thickness_prev_iter ! ngrid x nslope: Thickness of the ice table at the previous iteration [m] 226 real, dimension(:), allocatable :: delta_h2o_icetablesublim ! ngrid x Total mass of the H2O that has sublimated / condenses from the ice table [kg] 227 real, dimension(:), allocatable :: porefill ! Pore filling (output) to compute the dynamic ice table 228 real :: ssi_depth ! Ice table depth (output) to compute the dynamic ice table 229 232 230 ! Some variables for the PEM run 233 231 real, parameter :: year_step = 1 ! Timestep for the pem … … 557 555 allocate(delta_co2_adsorbded(ngrid)) 558 556 allocate(co2ice_disappeared(ngrid,nslope)) 559 allocate( porefillingice_thickness_prev_iter(ngrid,nslope))557 allocate(icetable_thickness_prev_iter(ngrid,nslope)) 560 558 allocate(delta_h2o_icetablesublim(ngrid)) 561 559 allocate(delta_h2o_adsorbded(ngrid)) … … 583 581 call end_adsorption_h_PEM 584 582 call ini_adsorption_h_PEM(ngrid,nslope,nsoilmx_PEM) 585 call end_ice_table _porefilling586 call ini_ice_table _porefilling(ngrid,nslope)583 call end_ice_table 584 call ini_ice_table(ngrid,nslope,nsoilmx_PEM) 587 585 588 586 if (soil_pem) then … … 646 644 endif 647 645 648 call pemetat0("startpem.nc",ngrid,nsoilmx,nsoilmx_PEM,nslope,timelen,timestep,TI_PEM,tsoil_PEM, porefillingice_depth, &649 porefillingice_thickness,tsurf_avg_yr1,tsurf_avg,q_co2_PEM_phys,q_h2o_PEM_phys,ps_timeseries,&650 tsoil_phys_PEM_timeseries,tend_h2o_ice,tend_co2_ice,co2_ice,h2o_ice,global_avg_press_PCM,&651 watersurf_density_avg,watersoil_density_PEM_avg,co2_adsorbded_phys,delta_co2_adsorbded,&652 h2o_adsorbded_phys,delta_h2o_adsorbded,stratif)646 call pemetat0("startpem.nc",ngrid,nsoilmx,nsoilmx_PEM,nslope,timelen,timestep,TI_PEM,tsoil_PEM,icetable_depth, & 647 icetable_thickness,ice_porefilling,tsurf_avg_yr1,tsurf_avg,q_co2_PEM_phys,q_h2o_PEM_phys, & 648 ps_timeseries,tsoil_phys_PEM_timeseries,tend_h2o_ice,tend_co2_ice,co2_ice,h2o_ice, & 649 global_avg_press_PCM,watersurf_density_avg,watersoil_density_PEM_avg,co2_adsorbded_phys, & 650 delta_co2_adsorbded,h2o_adsorbded_phys,delta_h2o_adsorbded,stratif) 653 651 654 652 ! We save the places where h2o ice is sublimating … … 928 926 ! II_d.3 Update the ice table 929 927 if (icetable_equilibrium) then 930 write(*,*) "Compute ice table "931 porefillingice_thickness_prev_iter = porefillingice_thickness932 call computeice_table_equilibrium(ngrid,nslope,nsoilmx_PEM,watercaptag,watersurf_density_avg,watersoil_density_PEM_avg,TI_PEM(:,1,:), porefillingice_depth,porefillingice_thickness)933 call compute_massh2o_exchange_ssi(ngrid,nslope,nsoilmx_PEM, porefillingice_thickness_prev_iter,porefillingice_thickness,porefillingice_depth,tsurf_avg,tsoil_PEM,delta_h2o_icetablesublim) ! Mass of H2O exchange between the ssi and the atmosphere928 write(*,*) "Compute ice table (equilibrium method)" 929 icetable_thickness_prev_iter = icetable_thickness 930 call computeice_table_equilibrium(ngrid,nslope,nsoilmx_PEM,watercaptag,watersurf_density_avg,watersoil_density_PEM_avg,TI_PEM(:,1,:),icetable_depth,icetable_thickness) 931 call compute_massh2o_exchange_ssi(ngrid,nslope,nsoilmx_PEM,icetable_thickness_prev_iter,icetable_thickness,icetable_depth,tsurf_avg,tsoil_PEM,delta_h2o_icetablesublim) ! Mass of H2O exchange between the ssi and the atmosphere 934 932 else if (icetable_dynamic) then 935 write(*,*) "Compute ice table" 936 porefillingice_thickness_prev_iter = porefillingice_thickness 937 call dyn_ss_ice_m_wrapper(ngrid,nsoilmx,TI_PEM,ps,mmol(igcm_h2o_vap),tsoil_PEM,porefillingice_depth,ss_ice_pf,ss_ice_pf_new,porefillingice_depth_new) 938 939 call compute_massh2o_exchange_ssi(ngrid,nslope,nsoilmx_PEM,porefillingice_thickness_prev_iter,porefillingice_thickness,porefillingice_depth,tsurf_avg, tsoil_PEM,delta_h2o_icetablesublim) ! Mass of H2O exchange between the ssi and the atmosphere 933 write(*,*) "Compute ice table (dynamic method)" 934 allocate(porefill(nsoilmx_PEM)) 935 do ig = 1,ngrid 936 do islope = 1,nslope 937 call dyn_ss_ice_m(icetable_depth(ig,islope),tsurf_avg(ig,islope),tsoil_PEM(ig,:,islope),nsoilmx_PEM,TI_PEM(ig,1,nslope),ps(ig),sum(q_h2o_PEM_phys(ig,:))/size(q_h2o_PEM_phys,2),ice_porefilling(ig,:,islope),porefill,ssi_depth) 938 icetable_depth(ig,islope) = ssi_depth 939 ice_porefilling(ig,:,islope) = porefill 940 enddo 941 enddo 942 deallocate(porefill) 943 call compute_massh2o_exchange_ssi(ngrid,nslope,nsoilmx_PEM,icetable_thickness_prev_iter,icetable_thickness,icetable_depth,tsurf_avg, tsoil_PEM,delta_h2o_icetablesublim) ! Mass of H2O exchange between the ssi and the atmosphere 940 944 endif 941 945 942 946 ! II_d.4 Update the soil thermal properties 943 call update_soil_thermalproperties(ngrid,nslope,nsoilmx_PEM,tend_h2o_ice,h2o_ice,global_avg_press_new, porefillingice_depth,porefillingice_thickness,TI_PEM)947 call update_soil_thermalproperties(ngrid,nslope,nsoilmx_PEM,tend_h2o_ice,h2o_ice,global_avg_press_new,icetable_depth,icetable_thickness,TI_PEM) 944 948 945 949 ! II_d.5 Update the mass of the regolith adsorbed … … 952 956 totmassh2o_adsorbded = 0. 953 957 do ig = 1,ngrid 954 do islope = 1,nslope958 do islope = 1,nslope 955 959 do l = 1,nsoilmx_PEM 956 if (l ==1) then960 if (l == 1) then 957 961 totmassco2_adsorbded = totmassco2_adsorbded + co2_adsorbded_phys(ig,l,islope)*(layer_PEM(l))* & 958 962 subslope_dist(ig,islope)/cos(pi*def_slope_mean(islope)/180.)*cell_area(ig) … … 960 964 subslope_dist(ig,islope)/cos(pi*def_slope_mean(islope)/180.)*cell_area(ig) 961 965 else 962 totmassco2_adsorbded = totmassco2_adsorbded + co2_adsorbded_phys(ig,l,islope)*(layer_PEM(l) - layer_PEM(l -1))* &966 totmassco2_adsorbded = totmassco2_adsorbded + co2_adsorbded_phys(ig,l,islope)*(layer_PEM(l) - layer_PEM(l - 1))* & 963 967 subslope_dist(ig,islope)/cos(pi*def_slope_mean(islope)/180.)*cell_area(ig) 964 totmassh2o_adsorbded = totmassh2o_adsorbded + h2o_adsorbded_phys(ig,l,islope)*(layer_PEM(l) - layer_PEM(l -1))* &968 totmassh2o_adsorbded = totmassh2o_adsorbded + h2o_adsorbded_phys(ig,l,islope)*(layer_PEM(l) - layer_PEM(l - 1))* & 965 969 subslope_dist(ig,islope)/cos(pi*def_slope_mean(islope)/180.)*cell_area(ig) 966 970 endif … … 987 991 call writediagpem(ngrid,'tsurf_slope'//str2,'tsurf','K',2,tsurf(:,islope)) 988 992 if (icetable_equilibrium) then 989 call writediagpem(ngrid,'ssi_depth_slope'//str2,'ice table depth','m',2, porefillingice_depth(:,islope))990 call writediagpem(ngrid,'ssi_thick_slope'//str2,'ice table depth','m',2, porefillingice_thickness(:,islope))993 call writediagpem(ngrid,'ssi_depth_slope'//str2,'ice table depth','m',2,icetable_depth(:,islope)) 994 call writediagpem(ngrid,'ssi_thick_slope'//str2,'ice table depth','m',2,icetable_thickness(:,islope)) 991 995 else if (icetable_dynamic) then 992 call writediagpem(ngrid,'ssi_depth_slope'//str2,'ice table depth','m',2, porefillingice_depth(:,islope))993 call writediagpem(ngrid,'ssi_thick_slope'//str2,'ice table depth','m',2, porefillingice_thickness(:,islope))996 call writediagpem(ngrid,'ssi_depth_slope'//str2,'ice table depth','m',2,icetable_depth(:,islope)) 997 call writediagpem(ngrid,'ssi_thick_slope'//str2,'ice table depth','m',2,icetable_thickness(:,islope)) 994 998 endif 995 999 … … 1205 1209 call pemdem0("restartpem.nc",longitude,latitude,cell_area,ngrid,nslope,def_slope,subslope_dist) 1206 1210 call pemdem1("restartpem.nc",i_myear,nsoilmx_PEM,ngrid,nslope,tsoil_PEM, & 1207 TI_PEM, porefillingice_depth,porefillingice_thickness,&1211 TI_PEM,icetable_depth,icetable_thickness,ice_porefilling, & 1208 1212 co2_adsorbded_phys,h2o_adsorbded_phys,h2o_ice,stratif) 1209 1213 write(*,*) "restartpem.nc has been written" … … 1227 1231 deallocate(watersurf_density_avg,watersoil_density_timeseries,Tsurfavg_before_saved) 1228 1232 deallocate(tsoil_phys_PEM_timeseries,watersoil_density_PEM_timeseries,watersoil_density_PEM_avg) 1229 deallocate(delta_co2_adsorbded,delta_h2o_adsorbded,vmr_co2_PEM_phys,delta_h2o_icetablesublim, porefillingice_thickness_prev_iter)1233 deallocate(delta_co2_adsorbded,delta_h2o_adsorbded,vmr_co2_PEM_phys,delta_h2o_icetablesublim,icetable_thickness_prev_iter) 1230 1234 deallocate(is_co2ice_ini,co2ice_disappeared,ini_co2ice_sublim,ini_h2oice_sublim,stratif) 1231 1235 !----------------------------- END OUTPUT ------------------------------
Note: See TracChangeset
for help on using the changeset viewer.
