module phyredem implicit none !======================================================================= contains !======================================================================= subroutine physdem0(filename,lonfi,latfi,nsoil,ngrid,nlay,nq, & phystep,day_ini,time,airefi, & alb,ith,def_slope, & subslope_dist) ! create physics restart file and write time-independent variables use comsoil_h, only: inertiedat, volcapa, mlayer use geometry_mod, only: cell_area use surfdat_h, only: zmea, zstd, zsig, zgam, zthe, & z0_default, albedice, emisice, emissiv, & iceradius, dtemisice, phisfi, z0, & hmons,summit,base,watercaptag use dimradmars_mod, only: tauvis use iostart, only : open_restartphy, close_restartphy, & put_var, put_field, length use mod_grid_phy_lmdz, only : klon_glo use planete_h, only: aphelie, emin_turb, lmixmin, obliquit, & peri_day, periheli, year_day use comcstfi_h, only: g, mugaz, omeg, rad, rcp use time_phylmdz_mod, only: daysec use comslope_mod, ONLY: nslope USE paleoclimate_mod, ONLY: paleoclimate, h2o_ice_depth, lag_co2_ice, d_coef implicit none character(len=*), intent(in) :: filename real,intent(in) :: lonfi(ngrid) real,intent(in) :: latfi(ngrid) integer,intent(in) :: nsoil integer,intent(in) :: ngrid integer,intent(in) :: nlay integer,intent(in) :: nq real,intent(in) :: phystep real,intent(in) :: day_ini real,intent(in) :: time real,intent(in) :: airefi(ngrid) real,intent(in) :: alb(ngrid) real,intent(in) :: ith(ngrid,nsoil) ! thermal inertia for present day real, intent(in) :: def_slope(nslope+1) !boundaries for bining of the slopes real, intent(in) :: subslope_dist(ngrid,nslope) !undermesh statistics real :: tab_cntrl(length) ! nb "length=100" defined in iostart module integer :: ig real :: watercaptag_tmp(ngrid) ! Create physics start file call open_restartphy(filename) ! Build tab_cntrl(:) array tab_cntrl(:)=0.0 !cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc ! Fill control array tab_cntrl(:) with parameters for this run !cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc ! Informations on the physics grid tab_cntrl(1) = float(klon_glo) ! total number of nodes on physics grid tab_cntrl(2) = float(nlay) ! number of atmospheric layers tab_cntrl(3) = day_ini + int(time) ! initial day tab_cntrl(4) = time -int(time) ! initial time of day ! Informations about Mars, used by dynamics and physics tab_cntrl(5) = rad ! radius of Mars (m) ~3397200 tab_cntrl(6) = omeg ! rotation rate (rad.s-1) tab_cntrl(7) = g ! gravity (m.s-2) ~3.72 tab_cntrl(8) = mugaz ! Molar mass of the atmosphere (g.mol-1) ~43.49 tab_cntrl(9) = rcp ! = r/cp ~0.256793 (=kappa dans dynamique) tab_cntrl(10) = daysec ! length of a sol (s) ~88775 tab_cntrl(11) = phystep ! time step in the physics tab_cntrl(12) = 0. tab_cntrl(13) = 0. ! Informations about Mars, only for physics tab_cntrl(14) = year_day ! length of year (sols) ~668.6 tab_cntrl(15) = periheli ! min. Sun-Mars distance (Mkm) ~206.66 tab_cntrl(16) = aphelie ! max. SUn-Mars distance (Mkm) ~249.22 tab_cntrl(17) = peri_day ! date of perihelion (sols since N. spring) tab_cntrl(18) = obliquit ! Obliquity of the planet (deg) ~23.98 ! Boundary layer and turbulence tab_cntrl(19) = z0_default ! default surface roughness (m) ~0.01 tab_cntrl(20) = lmixmin ! mixing length ~100 tab_cntrl(21) = emin_turb ! minimal energy ~1.e-8 ! Optical properties of polar caps and ground emissivity tab_cntrl(22) = albedice(1) ! Albedo of northern cap ~0.5 tab_cntrl(23) = albedice(2) ! Albedo of southern cap ~0.5 tab_cntrl(24) = emisice(1) ! Emissivity of northern cap ~0.95 tab_cntrl(25) = emisice(2) ! Emissivity of southern cap ~0.95 tab_cntrl(26) = emissiv ! Emissivity of martian soil ~.95 tab_cntrl(31) = iceradius(1) ! mean scat radius of CO2 snow (north) tab_cntrl(32) = iceradius(2) ! mean scat radius of CO2 snow (south) tab_cntrl(33) = dtemisice(1) ! time scale for snow metamorphism (north) tab_cntrl(34) = dtemisice(2) ! time scale for snow metamorphism (south) ! dust aerosol properties tab_cntrl(27) = tauvis ! mean visible optical depth ! Soil properties: tab_cntrl(35) = volcapa ! soil volumetric heat capacity ! Write the controle array call put_var("controle","Control parameters",tab_cntrl) ! Write the mid-layer depths call put_var("soildepth","Soil mid-layer depth",mlayer) ! Write longitudes call put_field("longitude","Longitudes of physics grid",lonfi) ! Write latitudes call put_field("latitude","Latitudes of physics grid",latfi) ! Write mesh areas call put_field("area","Mesh area",cell_area) ! Write surface geopotential call put_field("phisfi","Geopotential at the surface",phisfi) ! Write surface albedo call put_field("albedodat","Albedo of bare ground",alb) ! Subgrid topogaphy variables call put_field("ZMEA","Relief: mean relief",zmea) call put_field("ZSTD","Relief: standard deviation",zstd) call put_field("ZSIG","Relief: sigma parameter",zsig) call put_field("ZGAM","Relief: gamma parameter",zgam) call put_field("ZTHE","Relief: theta parameter",zthe) call put_field("hmons","Relief: hmons parameter (summit - base)",hmons) call put_field("summit","Relief: altitude from the aeroid of the summit of the highest subgrid topography",summit) call put_field("base","Relief: altitude from the aeroid of the base of the highest subgrid topography",base) ! Soil thermal inertia call put_field("inertiedat","Soil thermal inertia - present day",ith) ! Surface roughness length call put_field("z0","Surface roughness length",z0) ! Water cap watercaptag_tmp = 0 where (watercaptag) watercaptag_tmp = 1 call put_field("watercaptag","Infinite water reservoir",watercaptag_tmp) ! Sub grid scale slope parametrization call put_var("def_slope","slope criterium stages",def_slope) call put_field("subslope_dist","under mesh slope distribution",subslope_dist) ! Paleoclimate outputs if (paleoclimate) then call put_field("h2o_ice_depth","Depth to the fisrt H2O ice",h2o_ice_depth) call put_field("lag_co2_ice","Depth of the CO2 lags",lag_co2_ice) call put_field("d_coef","Diffusion coefficent",d_coef) endif ! Close file call close_restartphy end subroutine physdem0 !======================================================================= subroutine physdem1(filename,nsoil,ngrid,nlay,nq,nqsoil, & phystep,time,tsurf,tsoil,inertiesoil, & albedo,emis,q2,qsurf,qsoil,& tauscaling,totcloudfrac,wstar, & watercap,perennial_co2ice) ! write time-dependent variable to restart file use iostart, only : open_restartphy, close_restartphy, & put_var, put_field use tracer_mod, only: noms ! tracer names use nonoro_gwd_ran_mod, only: du_nonoro_gwd, dv_nonoro_gwd use compute_dtau_mod, only: dtau use dust_rad_adjust_mod, only: dust_rad_adjust_prev,dust_rad_adjust_next use dust_param_mod, only: dustscaling_mode use comsoil_h,only: flux_geo,adsorption_soil,igcm_h2o_vap_soil, & igcm_h2o_ice_soil,igcm_h2o_vap_ads use comslope_mod, only: nslope use paleoclimate_mod, only: paleoclimate implicit none include "callkeys.h" character(len=*),intent(in) :: filename integer,intent(in) :: nsoil integer,intent(in) :: ngrid integer,intent(in) :: nlay integer,intent(in) :: nq integer,intent(in) :: nqsoil real,intent(in) :: phystep real,intent(in) :: time real,intent(in) :: tsurf(ngrid,nslope) real,intent(in) :: tsoil(ngrid,nsoil,nslope) real,intent(in) :: inertiesoil(ngrid,nsoil,nslope) real,intent(in) :: albedo(ngrid,2,nslope) real,intent(in) :: emis(ngrid,nslope) real,intent(in) :: q2(ngrid,nlay+1) real,intent(in) :: qsurf(ngrid,nq,nslope) real,intent(in) :: qsoil(ngrid,nsoil,nqsoil,nslope) real,intent(in) :: tauscaling(ngrid) real,intent(in) :: totcloudfrac(ngrid) real,intent(in) :: wstar(ngrid) real,intent(in) :: watercap(ngrid,nslope) real,intent(in) :: perennial_co2ice(ngrid,nslope) integer :: iq character(len=30) :: txt ! to store some text ! indexes of water vapour & water ice tracers (if any): integer :: i_h2o_vap=0 integer :: i_h2o_ice=0 integer :: i_hdo_vap=0 integer :: i_hdo_ice=0 ! Open file call open_restartphy(filename) ! First variable to write must be "Time", in order to correctly ! set time counter in file call put_var("Time","Temps de simulation",time) ! Water ice layer call put_field("watercap","H2O ice cover",watercap,time) ! Perennial CO2 ice layer if (paleoclimate) call put_field("perennial_co2ice","CO2 ice cover",perennial_co2ice,time) ! Surface temperature call put_field("tsurf","Surface temperature",tsurf,time) ! Soil temperature call put_field("inertiesoil","Soil thermal inertia",inertiesoil,time) ! Soil temperature call put_field("tsoil","Soil temperature",tsoil,time) ! Albedo of the surface call put_field("albedo","Surface albedo",albedo(:,1,:),time) ! Emissivity of the surface call put_field("emis","Surface emissivity",emis,time) ! Planetary Boundary Layer call put_field("q2","pbl wind variance",q2,time) ! Sub-grid cloud fraction call put_field("totcloudfrac","Total cloud fraction",totcloudfrac,time) ! Dust conversion factor ! Only to be read by newstart to convert to actual dust values ! Or by any user who wants to reconstruct dust, opacity from the start files. call put_field("tauscaling","dust conversion factor",tauscaling,time) ! Radiative scaling coefficients if (dustscaling_mode==2) then call put_field("dust_rad_adjust_prev", & "radiative dust adjustement factor prev. sol", & dust_rad_adjust_prev,time) call put_field("dust_rad_adjust_next", & "radiative dust adjustement factor next sol", & dust_rad_adjust_next,time) endif if (dustinjection.gt.0) then call put_field("dtau","dust opacity difference between GCM and scenario",& dtau,time) endif if (calltherm) then call put_field("wstar","Max vertical velocity in thermals",wstar,time) endif ! Tracers on the surface ! preliminary stuff: look for water vapour & water ice tracers (if any) do iq=1,nq if (noms(iq).eq."h2o_vap") then i_h2o_vap=iq endif if (noms(iq).eq."h2o_ice") then i_h2o_ice=iq endif ! look for HDO vapour & HDO ice tracers (if any) if (noms(iq).eq."hdo_vap") then i_hdo_vap=iq endif if (noms(iq).eq."hdo_ice") then i_hdo_ice=iq endif enddo if (nq.gt.0) then do iq=1,nq txt=noms(iq) ! Exception: there is no water vapour surface tracer if (txt.eq."h2o_vap") then write(*,*)"physdem1: skipping water vapour tracer" if (i_h2o_ice.eq.i_h2o_vap) then ! then there is no "water ice" tracer; but still ! there is some water ice on the surface write(*,*)" writing water ice instead" txt="h2o_ice" else ! there is a "water ice" tracer which has been / will be ! delt with in due time cycle endif ! of if (igcm_h2o_ice.eq.igcm_h2o_vap) endif ! of if (txt.eq."h2o_vap") if (txt.eq."hdo_vap") then write(*,*)"physdem1: skipping HDO vapour tracer" if (i_hdo_ice.eq.i_hdo_vap) then ! then there is no "water ice" tracer; but still ! there is some water ice on the surface write(*,*)" writing HDO ice instead" txt="hdo_ice" else ! there is a "water ice" tracer which has been / will be ! delt with in due time cycle endif ! of if (igcm_hdo_ice.eq.igcm_hdo_vap) endif ! of if (txt.eq."hdo_vap") ! co2_ice has been added to qsurf(:,igcm_co2) in co2condens4micro if (txt.eq."co2_ice") then write(*,*)"physdem1: skipping co2_ice tracer" cycle end if call put_field(trim(txt),"tracer on surface",qsurf(:,iq,:),time) enddo endif ! Non-orographic gavity waves if (calllott_nonoro) then call put_field("du_nonoro_gwd","Zonal wind tendency due to GW",du_nonoro_gwd,time) call put_field("dv_nonoro_gwd","Meridional wind tendency due to GW",dv_nonoro_gwd,time) endif ! Geothermal Flux call put_field('flux_geo','Geothermal flux',flux_geo,time) ! Adsorption if (adsorption_soil) then call put_field("h2o_vap_soil","subsurface water vapour", & qsoil(:,:,igcm_h2o_vap_soil,:), time) call put_field("h2o_ice_soil","subsurface water ice", & qsoil(:,:,igcm_h2o_ice_soil,:), time) call put_field("h2o_vap_ads", "adsorbed water", & qsoil(:,:,igcm_h2o_vap_ads,:), time) endif ! Close file call close_restartphy end subroutine physdem1 end module phyredem