Changeset 2855


Ignore:
Timestamp:
Dec 22, 2022, 5:43:14 PM (2 years ago)
Author:
llange
Message:

PEM
Documentation of the main subroutines, and variables.
Unused programs have been removed.
LL

Location:
trunk/LMDZ.COMMON/libf/evolution
Files:
3 deleted
15 edited

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.COMMON/libf/evolution/computeice_table.F90

    r2849 r2855  
    11   SUBROUTINE computeice_table(ngrid,nslope,nsoil_PEM,rhowatersurf_ave,rhowatersoil_ave,ice_table)
     2!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3!!!
     4!!! Purpose: Compute the ice table depth knowing the yearly average water
     5!!! density at the surface and at depth.
     6!!! Computations are made following the methods in Schorgofer et al., 2005
     7!!!
     8!!! Author: LL
     9!!!
     10!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     11
    212#ifndef CPP_STD
    3     USE comsoil_h, only:  inertiedat, volcapa
    4     USE vertical_layers_mod, ONLY: ap,bp
    5     USE comsoil_h_PEM, only: layer_PEM
     13    USE comsoil_h_PEM, only: layer_PEM                             ! Depth of the vertical grid
    614    implicit none
    715
    8     integer,intent(in) :: ngrid,nslope,nsoil_PEM
    9     real,intent(in) :: rhowatersurf_ave(ngrid,nslope)    ! surface temperature, timeseries [K]
    10     real,intent(in) :: rhowatersoil_ave(ngrid,nsoil_PEM,nslope)    ! soil water density, yearly average [kg/m^3]
    11     real,intent(inout) :: ice_table(ngrid,nslope)                 ! ice table [m]
    12     real :: z1,z2
    13     integer ig, islope,isoil
    14     real :: diff_rho(nsoil_PEM)                    ! difference of vapor content
     16    integer,intent(in) :: ngrid,nslope,nsoil_PEM                   ! Size of the physical grid, number of subslope, number of soil layer in the PEM
     17    real,intent(in) :: rhowatersurf_ave(ngrid,nslope)              ! Water density at the surface, yearly averaged [kg/m^3]
     18    real,intent(in) :: rhowatersoil_ave(ngrid,nsoil_PEM,nslope)    ! Water density at depth, computed from clapeyron law's (Murchy and Koop 2005), yearly averaged [kg/m^3]
     19    real,intent(inout) :: ice_table(ngrid,nslope)                  ! ice table depth [m]
     20    real :: z1,z2                                                  ! intermediate variables used when doing a linear interpolation between two depths to find the root
     21    integer ig, islope,isoil                                       ! loop variables
     22    real :: diff_rho(nsoil_PEM)                                    ! difference of water vapor density between the surface and at depth [kg/m^3]
    1523
    1624
     
    2432         enddo
    2533
    26          if(diff_rho(1) > 0) then
     34         if(diff_rho(1) > 0) then                                   ! ice is at the surface
    2735           ice_table(ig,islope) = 0.
    2836         else
    29            do isoil = 1,nsoil_PEM -1
     37           do isoil = 1,nsoil_PEM -1                                ! general case, we find the ice table depth by doing a linear approximation between the two depth, and then solve the first degree equation to find the root
    3038             if((diff_rho(isoil).lt.0).and.(diff_rho(isoil+1).gt.0.)) then
    3139               z1 = (diff_rho(isoil) - diff_rho(isoil+1))/(layer_PEM(isoil) - layer_PEM(isoil+1))
     
    3846        enddo
    3947      enddo
    40 
    41    
    42 
    4348!=======================================================================
    4449      RETURN
    45 
    4650#endif
    4751      END
  • trunk/LMDZ.COMMON/libf/evolution/comsoil_h_PEM.F90

    r2794 r2855  
    22
    33implicit none
    4 ! nsoilmx : number of subterranean layers
    5 !EM: old soil routine:      integer, parameter :: nsoilmx = 10
    6   integer, parameter :: nsoilmx_PEM = 27
    7   integer, parameter :: n_1km = 23
     4  integer, parameter :: nsoilmx_PEM = 27               ! number of layers in the PEM
     5  integer, parameter :: n_1km = 23                     ! index at which we have overcome 1km
    86
    9   real,save,allocatable,dimension(:) :: layer_PEM      ! soil layer depths
    10   real,save,allocatable,dimension(:) :: mlayer_PEM     ! soil mid-layer depths
    11   real,save,allocatable,dimension(:,:,:) :: TI_PEM ! soil thermal inertia
    12     real,save,allocatable,dimension(:,:) :: inertiedat_PEM ! soil thermal inertia
     7  real,save,allocatable,dimension(:) :: layer_PEM      ! soil layer depths   [m]
     8  real,save,allocatable,dimension(:) :: mlayer_PEM     ! soil mid-layer depths [m]
     9  real,save,allocatable,dimension(:,:,:) :: TI_PEM     ! soil thermal inertia [SI]
     10    real,save,allocatable,dimension(:,:) :: inertiedat_PEM ! soil thermal inertia saved as reference for current climate [SI]
    1311  ! variables (FC: built in firstcall in soil.F)
    14   REAL,SAVE,ALLOCATABLE :: tsoil_PEM(:,:,:)       ! sub-surface temperatures (K)
    15   real,save,allocatable :: mthermdiff_PEM(:,:)  ! (FC) mid-layer thermal diffusivity
    16   real,save,allocatable :: thermdiff_PEM(:,:)   ! (FC) inter-layer thermal diffusivity
    17   real,save,allocatable :: coefq_PEM(:)         ! (FC) q_{k+1/2} coefficients
    18   real,save,allocatable :: coefd_PEM(:,:)       ! (FC) d_k coefficients
    19   real,save,allocatable :: alph_PEM(:,:,:)        ! (FC) alpha_k coefficients
    20   real,save,allocatable :: beta_PEM(:,:,:)        ! beta_k coefficients
    21   real,save :: mu_PEM
    22   real,parameter :: fluxgeo = 30e-3 !W/m^2
    23   real, save, allocatable :: co2_adsorbded_phys(:,:,:)  ! co2 that is in the regolith (kg/m^2)
     12  REAL,SAVE,ALLOCATABLE :: tsoil_PEM(:,:,:)       ! sub-surface temperatures [K]
     13  real,save,allocatable :: mthermdiff_PEM(:,:)  ! (FC) mid-layer thermal diffusivity [SI]
     14  real,save,allocatable :: thermdiff_PEM(:,:)   ! (FC) inter-layer thermal diffusivity [SI]
     15  real,save,allocatable :: coefq_PEM(:)         ! (FC) q_{k+1/2} coefficients [SI]
     16  real,save,allocatable :: coefd_PEM(:,:)       ! (FC) d_k coefficients [SI]
     17  real,save,allocatable :: alph_PEM(:,:,:)      ! (FC) alpha_k coefficients [SI]
     18  real,save,allocatable :: beta_PEM(:,:,:)      ! beta_k coefficients [SI]
     19  real,save :: mu_PEM                           ! mu coefficient [SI]
     20  real,parameter :: fluxgeo = 30e-3             ! Geothermal flux [W/m^2]
     21  real, save, allocatable :: co2_adsorbded_phys(:,:,:)  ! co2 that is in the regolith [kg/m^2]
     22  LOGICAL soil_pem  ! True by default, to run with the subsurface physic. Read in pem.def
    2423
    2524contains
  • trunk/LMDZ.COMMON/libf/evolution/interpolate_TIPEM_TIGCM.F90

    r2794 r2855  
    11   subroutine interpolate_TIPEM_TIGCM(ngrid,nslope,nsoil_PEM,nsoil_GCM,TI_PEM,TI_GCM)
    22
     3!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     4!!!
     5!!! Purpose: Transfer the thermal inertia from the PEM vertical  grid to the GCM vertical grid
     6!!!
     7!!!
     8!!! Author: LL
     9!!!
     10!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    311
    412
     
    1018!  ---------
    1119!  inputs:
    12       integer,intent(in) :: ngrid       ! # of horizontal grid points
    13       integer,intent(in) :: nslope      ! # of subslope wihtin the mesh
    14       integer,intent(in) :: nsoil_PEM   ! # of soil layers in the PEM
    15       integer,intent(in) :: nsoil_GCM   ! # of soil layers in the GCM
    16       real,intent(in) :: TI_PEM(ngrid,nsoil_PEM,nslope) ! # of soil layers in the PEM
    17       real,intent(inout) :: TI_GCM(ngrid,nsoil_GCM,nslope)      ! # of soil layers in the PEM
     20      integer,intent(in) :: ngrid                                ! # of horizontal grid points
     21      integer,intent(in) :: nslope                               ! # of subslope wihtin the mesh
     22      integer,intent(in) :: nsoil_PEM                            ! # of soil layers in the PEM
     23      integer,intent(in) :: nsoil_GCM                            ! # of soil layers in the GCM
     24      real,intent(in) :: TI_PEM(ngrid,nsoil_PEM,nslope)          !  Thermal inertia in the PEM vertical grid [J/m^2/K/s^{1/2}]
     25      real,intent(inout) :: TI_GCM(ngrid,nsoil_GCM,nslope)       ! Thermal inertia in the GCM vertical grid [J/m^2/K/s^{1/2}]
    1826
    1927!local variable
    20       integer :: ig,islope,iloop
    21 
    22 
    23 
     28      integer :: ig,islope,iloop                                 ! loop variables
    2429     do ig = 1,ngrid
    2530       do islope = 1,nslope
  • trunk/LMDZ.COMMON/libf/evolution/iostart_PEM.F90

    r2842 r2855  
    11MODULE iostart_PEM
     2!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3!!!
     4!!! Purpose: Read and write specific netcdf for the PEM
     5!!!
     6!!!
     7!!! Author: LL, inspired by iostart from the GCM
     8!!!
     9!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    210
    311    IMPLICIT NONE
  • trunk/LMDZ.COMMON/libf/evolution/lask_param_mod.F90

    r2835 r2855  
    11module lask_param_mod
     2!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3!!!
     4!!! Purpose: Define parameters from Laskar et al., 2004 evolution
     5!!!
     6!!! Author: RV
     7!!!
     8!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    29
    310implicit none
    411
    5   real,save,allocatable :: yearlask(:)
    6   real,save,allocatable :: oblask(:)
    7   real,save,allocatable :: exlask(:)
    8   real,save,allocatable :: lsplask(:)
    9   integer, save :: last_ilask
     12  real,save,allocatable :: yearlask(:)          ! year before present from Laskar et al. Tab
     13  real,save,allocatable :: oblask(:)            ! obliquity [deg]
     14  real,save,allocatable :: exlask(:)            ! excentricity[deg]
     15  real,save,allocatable :: lsplask(:)           ! ls perihelie [deg]
     16  integer, save :: last_ilask                   ! Index of the line in the file year_obl_lask.asc corresponding to the closest lower year to the current year
    1017
    1118contains
  • trunk/LMDZ.COMMON/libf/evolution/nb_time_step_GCM.F90

    r2841 r2855  
    1212!=======================================================================
    1313!
    14 ! Read initial confitions file
     14! Purpose: Read in the data_GCM_Yr*.nc the number of time step
    1515!
     16! Author: RV
    1617!=======================================================================
    1718
  • trunk/LMDZ.COMMON/libf/evolution/pem.F90

    r2849 r2855  
    8787
    8888!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! SOIL
    89       USE soil_evolution_mod, ONLY: soil_pem, soil_pem_CN
    90       use comsoil_h_PEM, only: ini_comsoil_h_PEM,end_comsoil_h_PEM,nsoilmx_PEM, &
     89      use comsoil_h_PEM, only: soil_pem,ini_comsoil_h_PEM,end_comsoil_h_PEM,nsoilmx_PEM, &
    9190                              TI_PEM,inertiedat_PEM,alph_PEM, beta_PEM,         & ! soil thermal inertia         
    9291                              tsoil_PEM, mlayer_PEM,layer_PEM,                  & !number of subsurface layers, soil mid layer depths
     
    175174      REAL , dimension(:,:), allocatable :: min_co2_ice_s_2        ! LON x LAT field : minimum of water ice at each point for the second year
    176175
    177       REAL :: global_ave_press_GCM
    178       REAL :: global_ave_press_old           ! physical point field : Global average pressure of initial/previous time step
    179       REAL :: global_ave_press_new           ! physical point field : Global average pressure of current time step
    180 
    181       REAL , dimension(:,:), allocatable ::  zplev_new
    182       REAL , dimension(:,:), allocatable ::  zplev_gcm
    183       REAL , dimension(:,:,:), allocatable ::  zplev_new_timeseries  ! same but with the time series
    184       REAL , dimension(:,:,:), allocatable :: zplev_old_timeseries   ! same but with the time series
     176      REAL :: global_ave_press_GCM           ! constant: global average pressure retrieved in the GCM [Pa]
     177      REAL :: global_ave_press_old           ! constant: Global average pressure of initial/previous time step
     178      REAL :: global_ave_press_new           ! constant: Global average pressure of current time step
     179
     180      REAL , dimension(:,:), allocatable ::  zplev_new ! Physical x Atmospheric field : mass of the atmospheric  layers in the pem at current time step [kg/m^2]
     181      REAL , dimension(:,:), allocatable ::  zplev_gcm ! same but retrieved from the gcm [kg/m^2]
     182      REAL , dimension(:,:,:), allocatable ::  zplev_new_timeseries  ! Physical x Atmospheric x Time: same as zplev_new, but in times series [kg/m ^2]
     183      REAL , dimension(:,:,:), allocatable :: zplev_old_timeseries   ! same but with the time series, for oldest time step
    185184
    186185      LOGICAL :: STOPPING_water   ! Logical : is the criterion (% of change in the surface of sublimating water ice) reached?
     
    189188      LOGICAL :: STOPPING_1_co2   ! Logical : is there still water ice to sublimate?
    190189
    191       REAL, dimension(:,:,:),allocatable  :: q_co2_GCM ! Initial amount of co2 in the first layer
    192 
    193       real,save :: m_co2, m_noco2, A , B, mmean
    194       real ,allocatable :: vmr_co2_gcm_phys(:,:) !(ngrid) ! co2 volume mixing ratio
    195       real ,allocatable :: vmr_co2_pem_phys(:,:) !(ngrid) ! co2 volume mixing ratio
    196       real ,allocatable :: q_h2o_GCM_phys(:,:)
    197       real ,allocatable :: q_co2_GCM_phys(:,:)
    198       real ,allocatable :: q_co2_PEM_phys(:,:)
    199       REAL, ALLOCATABLE ::  ps_GCM(:,:,:)
    200       REAL, ALLOCATABLE :: ps_GCM_yr1(:,:,:)
    201       REAL, ALLOCATABLE ::  vmr_co2_gcm(:,:,:)
    202       REAL, ALLOCATABLE ::  q_h2o_GCM(:,:,:)
    203       REAL ,allocatable ::  q_h2o_PEM_phys(:,:)
    204       integer :: timelen
    205       REAL :: ave
    206 
    207       REAL, ALLOCATABLE :: p(:,:)  !(ngrid,llmp1)
    208       REAL :: extra_mass           ! Extra mass of a tracer if it is greater than 1
     190      REAL, dimension(:,:,:),allocatable  :: q_co2_GCM ! Lon x Lat x Time : mass mixing ratio of co2 in the first layer [kg/kg]
     191
     192      real,save :: m_co2, m_noco2, A , B, mmean        ! Molar mass of co2, no co2 (Ar, ...), intermediate A, B for computations, mean molar mass of the layer [mol/kg]
     193      real ,allocatable :: vmr_co2_gcm_phys(:,:) ! Physics x Times  co2 volume mixing ratio retrieve from the gcm [m^3/m^3]
     194      real ,allocatable :: vmr_co2_pem_phys(:,:) ! Physics x Times  co2 volume mixing ratio used in the PEM
     195      real ,allocatable :: q_h2o_GCM_phys(:,:)   ! Physics x Times h2o mass mixing ratio in the first layer from the GCM  [kg/kg]
     196      real ,allocatable :: q_co2_GCM_phys(:,:)   ! Physics x Times co2 mass mixing ratio in the first layer from the GCM  [kg/kg]
     197      real ,allocatable :: q_co2_PEM_phys(:,:)   ! Physics x Times co2 mass mixing ratio in the first layer computed in the PEM  [kg/kg]
     198      REAL, ALLOCATABLE ::  ps_GCM(:,:,:)        ! Lon x Lat x Times: surface pressure from the GCM [Pa]
     199      REAL, ALLOCATABLE :: ps_GCM_yr1(:,:,:)     ! Lon x Lat x Times: surface pressure from the 1st year of the GCM [Pa]
     200      REAL, ALLOCATABLE ::  vmr_co2_gcm(:,:,:)   ! Lon x Lat x Times: co2 volumemixing ratio retrieve from the gcm [m^3/m^3]
     201      REAL, ALLOCATABLE ::  q_h2o_GCM(:,:,:)     ! Lon x Lat x Times: h2o volume mixing ratio retrieved from the GCM
     202      REAL ,allocatable ::  q_h2o_PEM_phys(:,:)  ! Physics x Times: h2o mass mixing ratio computed in the PEM
     203      integer :: timelen                         ! # time samples
     204      REAL :: ave                                ! intermediate varibale to compute average
     205
     206      REAL, ALLOCATABLE :: p(:,:)  ! Physics x Atmosphere: pressure to recompute and write in restart (ngrid,llmp1)
     207      REAL :: extra_mass           ! Intermediate variables Extra mass of a tracer if it is greater than 1
     208     REAL :: beta_clap_co2 = 3182.48 ! clapeyron's law for CO2
     209     REAL :: alpha_clap_co2 = 23.3494 ! Clapeyron's law for CO2
     210
    209211
    210212!!!!!!!!!!!!!!!!!!!!!!!! SLOPE
    211       character*2 :: str2
    212       REAL ,allocatable :: watercap_slope(:,:)                           !(ngrid,nslope)
    213       REAL , dimension(:,:,:), allocatable :: min_co2_ice_slope_1        ! LON x LAT field : minimum of water ice at each point for the first year
    214       REAL , dimension(:,:,:), allocatable :: min_co2_ice_slope_2        ! LON x LAT field : minimum of water ice at each point for the second year
    215       REAL , dimension(:,:,:), allocatable :: min_h2o_ice_slope_1        ! LON x LAT field : minimum of water ice at each point for the first year
    216       REAL , dimension(:,:,:), allocatable :: min_h2o_ice_slope_2        ! LON x LAT field : minimum of water ice at each point for the second year
    217 
    218 
    219       REAL , dimension(:,:,:,:), allocatable :: co2_ice_GCM_slope        ! LON x LATX NSLOPE x Times field : co2 ice given by the GCM
    220       REAL , dimension(:,:,:), allocatable :: co2_ice_GCM_phys_slope     ! LON x LATX NSLOPE x Times field : co2 ice given by the GCM
    221 
    222 
     213      REAL ,allocatable :: watercap_slope(:,:)                           ! Physics x Nslope: watercap per slope
     214      REAL , dimension(:,:,:), allocatable :: min_co2_ice_slope_1        ! LON x LAT field : minimum of co2 ice at each point for the first year [kg/m^2]
     215      REAL , dimension(:,:,:), allocatable :: min_co2_ice_slope_2        ! LON x LAT field : minimum of co2 ice at each point for the second year [kg/m^2]
     216      REAL , dimension(:,:,:), allocatable :: min_h2o_ice_slope_1        ! LON x LAT field : minimum of water ice at each point for the first year [kg/m^2]
     217      REAL , dimension(:,:,:), allocatable :: min_h2o_ice_slope_2        ! LON x LAT field : minimum of water ice at each point for the second year [kg/m^2]
     218      REAL , dimension(:,:,:,:), allocatable :: co2_ice_GCM_slope        ! LON x LATX NSLOPE x Times field : co2 ice given by the GCM [kg/m^2]
     219      REAL , dimension(:,:,:), allocatable :: co2_ice_GCM_phys_slope     ! Physics x NSLOPE x Times field : co2 ice given by the GCM  [kg/m^2]
    223220      REAL, dimension(:,:),allocatable  :: initial_co2_ice_sublim_slope    ! physical point field : Logical array indicating sublimating point of co2 ice
    224       REAL, dimension(:,:),allocatable  :: initial_h2o_ice_slope           ! physical point field : Logical array indicating sublimating point of h2o ice
    225       REAL, dimension(:,:),allocatable  :: initial_co2_ice_slope           ! physical point field : Logical array indicating sublimating point of co2 ice
    226 
    227       REAL , dimension(:,:,:), allocatable :: tendencies_co2_ice_slope     ! LON x LAT field : Tendency of evolution of perenial co2 ice over a year
    228       REAL , dimension(:,:,:), allocatable :: tendencies_h2o_ice_slope     ! LON x LAT field : Tendency of evolution of perenial co2 ice over a year
    229       REAL, dimension(:,:),allocatable  :: tendencies_co2_ice_phys_slope   ! physical point field : Tendency of evolution of perenial co2 ice over a year
    230       REAL, dimension(:,:),allocatable  :: tendencies_co2_ice_phys_slope_ini ! physical point field x nslope: Tendency of evolution of perenial co2 ice over a year in the GCM
    231       REAL, dimension(:,:),allocatable  :: tendencies_h2o_ice_phys_slope   ! physical point field : Tendency of evolution of perenial co2 ice over a year
    232 
     221      REAL, dimension(:,:),allocatable  :: initial_h2o_ice_slope           ! physical point field : Logical array indicating if there is water ice at initial state
     222      REAL, dimension(:,:),allocatable  :: initial_co2_ice_slope           ! physical point field : Logical array indicating if there is co2 ice at initial state
     223      REAL , dimension(:,:,:), allocatable :: tendencies_co2_ice_slope     ! LON x LAT x nslope field : Tendency of evolution of perenial co2 ice over a year
     224      REAL , dimension(:,:,:), allocatable :: tendencies_h2o_ice_slope     ! LON x LAT x slope field : Tendency of evolution of perenial water ice over a year
     225      REAL, dimension(:,:),allocatable  :: tendencies_co2_ice_phys_slope   ! physical point xslope field : Tendency of evolution of perenial co2 ice over a year
     226      REAL, dimension(:,:),allocatable  :: tendencies_co2_ice_phys_slope_ini ! physical point x slope field x nslope: Tendency of evolution of perenial co2 ice over a year in the GCM
     227      REAL, dimension(:,:),allocatable  :: tendencies_h2o_ice_phys_slope   ! physical pointx slope  field : Tendency of evolution of perenial co2 ice over a year
    233228      REAL,SAVE,ALLOCATABLE,DIMENSION(:) ::  co2_hmax                      ! Maximum height  for CO2 deposit on slopes (m)
    234 
    235229      REAL, PARAMETER :: rho_co2 = 1600           ! CO2 ice density (kg/m^3)
    236230      INTEGER :: iaval                            ! Index of the neighboord slope ()
     
    241235
    242236     REAL, ALLOCATABLE :: tsurf_ave(:,:,:)          ! LON x LAT x SLOPE field : Averaged Surface Temperature [K]
    243      REAL, ALLOCATABLE  :: tsurf_ave_phys(:,:)      ! IG x LAT x SLOPE field : Averaged Surface Temperature [K]
     237     REAL, ALLOCATABLE  :: tsurf_ave_phys(:,:)      ! Physic x LAT x SLOPE field : Averaged Surface Temperature [K]
    244238     REAL, ALLOCATABLE :: tsoil_ave(:,:,:,:)        ! LON x LAT x SLOPE field : Averaged Soil Temperature [K]
    245      REAL, ALLOCATABLE :: tsoil_ave_yr1(:,:,:,:)    ! LON x LAT x SLOPE field : Averaged Soil Temperature [K]
    246      REAL, ALLOCATABLE :: tsoil_ave_phys_yr1(:,:,:) !IG x SLOPE field : Averaged Soil Temperature [K]
    247 
     239     REAL, ALLOCATABLE :: tsoil_ave_yr1(:,:,:,:)    ! LON x LAT x SLOPE field : Averaged Soil Temperature  during 1st year of the GCM [K]
     240     REAL, ALLOCATABLE :: tsoil_ave_phys_yr1(:,:,:) ! Physics x SLOPE field : Averaged Soil Temperature  during 1st year [K]
    248241     REAL, ALLOCATABLE :: TI_GCM(:,:,:,:)                  ! LON x LAT x SLOPE field : Averaged Thermal Inertia  [SI]
    249      REAL, ALLOCATABLE :: tsurf_GCM_timeseries(:,:,:,:)    !LONX LAT x SLOPE XTULES field : NOn averaged Surf Temperature [K]
    250      REAL, ALLOCATABLE :: tsurf_phys_GCM_timeseries(:,:,:) !IG x SLOPE XTULES field : NOn averaged Surf Temperature [K]
     242     REAL, ALLOCATABLE :: tsurf_GCM_timeseries(:,:,:,:)    ! LON X LAT x SLOPE XTULES field : Surface Temperature in timeseries [K]
     243     REAL, ALLOCATABLE :: tsurf_phys_GCM_timeseries(:,:,:) ! Physic x SLOPE XTULES field : NOn averaged Surf Temperature [K]
    251244
    252245     REAL, ALLOCATABLE :: tsoil_phys_PEM_timeseries(:,:,:,:) !IG x SLOPE XTULES field : NOn averaged Soil Temperature [K]
    253246     REAL, ALLOCATABLE :: tsoil_GCM_timeseries(:,:,:,:,:)    !IG x SLOPE XTULES field : NOn averaged Soil Temperature [K]
    254 
    255247     REAL, ALLOCATABLE :: tsurf_ave_yr1(:,:,:)     ! LON x LAT x SLOPE field : Averaged Surface Temperature of the first year of the gcm [K]
    256      REAL, ALLOCATABLE  :: tsurf_ave_phys_yr1(:,:) ! IG x LAT x SLOPE field : Averaged Surface Temperature of first call of the gcm [K]
    257 
    258      REAL, PARAMETER :: year_step = 1
    259      INTEGER :: year_iter        !Counter for the number of PEM iteration
    260      INTEGER :: year_iter_max
    261      REAL :: timestep
    262 
    263      REAL, ALLOCATABLE :: inertiesoil(:,:)    ! Thermal inertia of the mesh for restart
    264 
    265      REAL, ALLOCATABLE :: TI_GCM_phys(:,:,:)  ! Averaged GCM Thermal Inertia  [SI]
    266      REAL, ALLOCATABLE :: TI_GCM_start(:,:,:) ! Averaged GCM Thermal Inertia  [SI]
    267 
    268      REAL,ALLOCATABLE  :: ice_depth(:,:)
    269      REAL,ALLOCATABLE  :: TI_locslope(:,:)
    270      REAL,ALLOCATABLE  :: Tsoil_locslope(:,:)
    271      REAL,ALLOCATABLE  :: Tsurf_locslope(:)
    272      REAL,ALLOCATABLE  :: alph_locslope(:,:)
    273      REAL,ALLOCATABLE  :: beta_locslope(:,:)   
    274      REAL,ALLOCATABLE  :: watersurf_density_timeseries(:,:,:,:)
    275      REAL,ALLOCATABLE  :: watersoil_density_timeseries(:,:,:,:,:)
    276      REAL,ALLOCATABLE  :: watersurf_density_phys_timeseries(:,:,:)
    277      REAL,ALLOCATABLE  :: watersurf_density_phys_ave(:,:)
    278      REAL,ALLOCATABLE  :: watersoil_density_phys_PEM_timeseries(:,:,:,:)
    279      REAL,ALLOCATABLE  :: watersoil_density_phys_PEM_ave(:,:,:)
    280      REAL,ALLOCATABLE  :: Tsurfave_before_saved(:,:)
    281      REAL, ALLOCATABLE :: delta_co2_adsorbded(:)
    282      REAL :: totmass_adsorbded
    283    real :: alpha_clap_h2o = -6143.7
    284    real :: beta_clap_h2o = 28.9074
    285 
     248     REAL, ALLOCATABLE  :: tsurf_ave_phys_yr1(:,:) ! Physic SLOPE field : Averaged Surface Temperature of first call of the gcm [K]
     249     REAL, ALLOCATABLE :: inertiesoil(:,:)    !Physic x Depth  Thermal inertia of the mesh for restart [SI]
     250
     251     REAL, ALLOCATABLE :: TI_GCM_phys(:,:,:)  ! Physic x Depth x Slope Averaged GCM Thermal Inertia per slope  [SI]
     252     REAL, ALLOCATABLE :: TI_GCM_start(:,:,:) ! Same but for the start
     253
     254     REAL,ALLOCATABLE  :: ice_depth(:,:)      ! Physic x SLope: Ice table depth [m]
     255     REAL,ALLOCATABLE  :: TI_locslope(:,:)    ! Physic x Soil: Intermediate thermal inertia  to compute Tsoil [SI]
     256     REAL,ALLOCATABLE  :: Tsoil_locslope(:,:) ! Physic x Soil: intermediate when computing Tsoil [K]
     257     REAL,ALLOCATABLE  :: Tsurf_locslope(:)   ! Physic x Soil: Intermediate surface temperature  to compute Tsoil [K]
     258     REAL,ALLOCATABLE  :: alph_locslope(:,:)  ! Physic x Soil: Intermediate to compute Tsoil [1]
     259     REAL,ALLOCATABLE  :: beta_locslope(:,:)  ! Physic x Soil : Intermediate tocompute Tsoil [K]
     260     REAL,ALLOCATABLE  :: watersurf_density_timeseries(:,:,:,:) ! Lon x Lat x Slope x Times: water surface density, time series [kg/m^3]
     261     REAL,ALLOCATABLE  :: watersoil_density_timeseries(:,:,:,:,:) ! Lon x Lat x Soil x Slope x Times water soil density, time series [kg /m^3]
     262     REAL,ALLOCATABLE  :: watersurf_density_phys_timeseries(:,:,:) ! Physic  x Slope x Times, water surface density, time series [kg/m^3]
     263     REAL,ALLOCATABLE  :: watersurf_density_phys_ave(:,:)          ! Physic  x Slope, water surface density, yearly averaged [kg/m^3]
     264     REAL,ALLOCATABLE  :: watersoil_density_phys_PEM_timeseries(:,:,:,:) ! Physic x Soil x Slope x Times, water soil density, time series [kg/m^3]
     265     REAL,ALLOCATABLE  :: watersoil_density_phys_PEM_ave(:,:,:)          ! Physic x Soil x SLopes, water soil density, yearly averaged [kg/m^3]
     266     REAL,ALLOCATABLE  :: Tsurfave_before_saved(:,:)                     ! Surface temperature saved from previous time step [K]
     267     REAL, ALLOCATABLE :: delta_co2_adsorbded(:)                         ! Physics: quantity of CO2 that is exchanged because of adsorption / desorption [kg/m^2]
     268     REAL :: totmass_adsorbded                                           ! Total mass of CO2 that is exchanged because of adsorption / desoprtion over the planets [kg]
     269     REAL :: alpha_clap_h2o = -6143.7                                    ! coeffcient to compute psat, from Murphie et Kood 2005 [K]
     270     REAL :: beta_clap_h2o = 28.9074                                     ! coefficient to compute psat, from Murphie et Kood 2005 [1]
     271     LOGICAL :: bool_sublim                                              ! logical to check if there is sublimation or not
     272
     273!! Some parameters for the PEM run
     274    REAL, PARAMETER :: year_step = 1     ! timestep for the pem
     275    INTEGER ::  year_iter                ! number of iteration
     276    INTEGER :: year_iter_max             ! maximum number of iterations before stopping
     277    REAL :: timestep                     ! timestep [s]
    286278#ifdef CPP_STD
    287279!     INTEGER :: nsplit_phys=1
    288280!     LOGICAL :: iflag_phys=.true.
    289      REAL :: frost_albedo_threshold=0.05
    290      REAL :: albedo_h2o_frost
    291      REAL,ALLOCATABLE :: co2ice(:)
     281     REAL :: frost_albedo_threshold=0.05             ! frost albedo threeshold to convert fresh frost to old ice
     282     REAL :: albedo_h2o_frost                        ! albedo of h2o frost
     283     REAL,ALLOCATABLE :: co2ice(:)                   ! Physics: co2 ice mesh averaged [kg/m^2]
    292284#endif
    293285
     
    295287
    296288! Loop variable
    297      LOGICAL :: bool_sublim
    298      INTEGER :: i,j,ig0,l,ig,nnq,t,islope,ig_loop,islope_loop,iloop,isoil
    299      REAL :: beta = 3182.48
    300      REAL :: alpha = 23.3494
    301  
     289     INTEGER :: i,j,ig0,l,ig,nnq,t,islope,ig_loop,islope_loop,iloop,isoil
    302290#ifndef CPP_STD
    303291! Parallel variables
     
    566554     print *, "Downloading data Y1..."
    567555
    568      call read_data_GCM("data_GCM_Y1.nc", min_h2o_ice_s_1,min_co2_ice_s_1,iim,jjm,nlayer,vmr_co2_gcm,ps_GCM_yr1,timelen,min_co2_ice_slope_1,min_h2o_ice_slope_1,&   
     556     call read_data_GCM("data_GCM_Y1.nc",timelen, iim,jjm, min_h2o_ice_s_1,min_co2_ice_s_1,vmr_co2_gcm,ps_GCM_yr1,min_co2_ice_slope_1,min_h2o_ice_slope_1,&   
    569557                       nslope,tsurf_ave_yr1,tsoil_ave_yr1, tsurf_GCM_timeseries,tsoil_GCM_timeseries,TI_GCM,q_co2_GCM,q_h2o_GCM,co2_ice_GCM_slope, &     
    570558                       watersurf_density_timeseries,watersoil_density_timeseries)
     
    581569     print *, "Downloading data Y2"
    582570
    583      call read_data_GCM("data_GCM_Y2.nc", min_h2o_ice_s_2,min_co2_ice_s_2,iim,jjm,nlayer,vmr_co2_gcm,ps_GCM,timelen,min_co2_ice_slope_2,min_h2o_ice_slope_2, &
     571     call read_data_GCM("data_GCM_Y2.nc",timelen,iim,jjm ,min_h2o_ice_s_2,min_co2_ice_s_2,vmr_co2_gcm,ps_GCM,min_co2_ice_slope_2,min_h2o_ice_slope_2, &
    584572                  nslope,tsurf_ave,tsoil_ave, tsurf_GCM_timeseries,tsoil_GCM_timeseries,TI_GCM,q_co2_GCM,q_h2o_GCM,co2_ice_GCM_slope, &     
    585573                  watersurf_density_timeseries,watersoil_density_timeseries)
     
    815803!---------------------------- Read the PEMstart ---------------------
    816804
    817       call pemetat0(.true.,"startfi_PEM.nc",ngrid,nsoilmx,nsoilmx_PEM,nslope,timelen,timestep,TI_PEM,tsoil_ave_phys_yr1,tsoil_PEM,ice_depth, &
    818       tsurf_ave_phys_yr1, tsurf_ave_phys, q_co2_PEM_phys, q_h2o_PEM_phys,ps_phys_timeseries,tsurf_phys_GCM_timeseries,tsoil_phys_PEM_timeseries,&
     805      call pemetat0(.false.,"startfi_PEM.nc",ngrid,nsoilmx,nsoilmx_PEM,nslope,timelen,timestep,TI_PEM,tsoil_ave_phys_yr1,tsoil_PEM,ice_depth, &
     806      tsurf_ave_phys_yr1, tsurf_ave_phys, q_co2_PEM_phys, q_h2o_PEM_phys,ps_phys_timeseries,tsoil_phys_PEM_timeseries,&
    819807      tendencies_h2o_ice_phys_slope,tendencies_co2_ice_phys_slope,co2ice_slope,qsurf_slope(:,igcm_h2o_ice,:),global_ave_press_GCM, co2_adsorbded_phys,delta_co2_adsorbded,&
    820808      watersurf_density_phys_ave,watersoil_density_phys_PEM_ave)
     
    10711059            do t=1,timelen
    10721060              if(co2_ice_GCM_phys_slope(ig,islope,t).gt.1e-3) then
    1073                 ave = ave + beta/(alpha-log(vmr_co2_pem_phys(ig,t)*ps_phys_timeseries(ig,t)/100.))
     1061                ave = ave + beta_clap_co2/(alpha_clap_co2-log(vmr_co2_pem_phys(ig,t)*ps_phys_timeseries(ig,t)/100.))
    10741062              else
    10751063                ave = ave + tsurf_phys_GCM_timeseries(ig,islope,t)
  • trunk/LMDZ.COMMON/libf/evolution/pemetat0.F90

    r2849 r2855  
    11   SUBROUTINE pemetat0(startpem_file,filename,ngrid,nsoil_GCM,nsoil_PEM,nslope,timelen,timestep,TI_PEM,tsoil_PEM_yr1,tsoil_PEM,ice_table, &
    2                       tsurf_ave_yr1,tsurf_ave,q_co2,q_h2o,ps_inst,tsurf_inst,tsoil_inst,tend_h2oglaciers,tend_co2glaciers,co2ice,waterice, &
     2                      tsurf_ave_yr1,tsurf_ave_yr2,q_co2,q_h2o,ps_inst,tsoil_inst,tend_h2oglaciers,tend_co2glaciers,co2ice,waterice, &
    33                      global_ave_pressure, m_co2_regolith_phys,deltam_co2_regolith_phys, watersurf_ave,watersoil_ave)
    44   
    55   use iostart_PEM, only:  open_startphy, close_startphy, get_field, get_var
    6    use comsoil_h_PEM, only: layer_PEM, mlayer_PEM,n_1km,fluxgeo,inertiedat_PEM
     6   use comsoil_h_PEM, only: soil_pem,layer_PEM, mlayer_PEM,n_1km,fluxgeo,inertiedat_PEM
    77   use comsoil_h, only:  volcapa,inertiedat
    8    use soil_evolution_mod, only: soil_pem,soil_pem_CN
    98   use adsorption_mod, only : regolith_co2adsorption
    109   USE temps_mod_evol, ONLY: year_PEM
     
    1211  implicit none
    1312 
    14   character(len=*), intent(in) :: filename
    15   LOGICAL,intent(in) :: startpem_file
    16   character*8 :: fichnom
    17   integer,intent(in) :: ngrid
    18   integer,intent(in) :: nsoil_GCM
    19   integer,intent(in) :: nsoil_PEM
    20   integer,intent(in) :: nslope
    21   integer,intent(in) :: timelen
    22   real, intent(in) :: tsurf_ave_yr1(ngrid,nslope)                          ! surface temperature at the first year of GCM call
    23   real,intent(in) :: tsurf_ave(ngrid,nslope)                  ! surface temperature at the current year
     13  character(len=*), intent(in) :: filename                    ! name of the startfi_PEM.nc
     14  LOGICAL,intent(in) :: startpem_file                         ! boolean to check if we read the startfile or not
     15  integer,intent(in) :: ngrid                                 ! # of physical grid points
     16  integer,intent(in) :: nsoil_GCM                             ! # of vertical grid points in the GCM
     17  integer,intent(in) :: nsoil_PEM                             ! # of vertical grid points in the PEM
     18  integer,intent(in) :: nslope                                ! # of sub-grid slopes
     19  integer,intent(in) :: timelen                               ! # time samples
     20  real, intent(in) :: tsurf_ave_yr1(ngrid,nslope)             ! surface temperature at the first year of GCM call [K]
     21  real,intent(in) :: tsurf_ave_yr2(ngrid,nslope)              ! surface temperature at the second  year of GCM call [K]
    2422  real,intent(in) :: q_co2(ngrid,timelen)                     ! MMR tracer co2 [kg/kg]
    2523  real,intent(in) :: q_h2o(ngrid,timelen)                     ! MMR tracer h2o [kg/kg]
    26   real,intent(in) :: ps_inst(ngrid,timelen)                        ! surface pressure [Pa]
    27   real,intent(in) :: tsurf_inst(ngrid,nslope,timelen) ! soil (mid-layer) temperature
    28   real,intent(in) :: timestep       ! time step
    29   real,intent(in) :: tend_h2oglaciers(ngrid,nslope)
    30   real,intent(in) :: tend_co2glaciers(ngrid,nslope)
    31   real,intent(in) :: co2ice(ngrid,nslope)
    32   real,intent(in) :: waterice(ngrid,nslope)
    33   real, intent(in) :: tsoil_PEM_yr1(ngrid,nsoil_PEM,nslope)
    34   real, intent(in) :: watersurf_ave(ngrid,nslope)     ! surface water ice density, yearly averaged  (kg/m^3)
    35   real, intent(inout) :: watersoil_ave(ngrid,nsoil_PEM,nslope)    ! surface water ice density, yearly averaged (kg/m^3)
    36   real, intent(in) :: global_ave_pressure                       ! mean average pressure on the planet
     24  real,intent(in) :: ps_inst(ngrid,timelen)                   ! surface pressure [Pa]
     25  real,intent(in) :: timestep                                 ! time step [s]
     26  real,intent(in) :: tend_h2oglaciers(ngrid,nslope)           ! tendencies on h2o glaciers
     27  real,intent(in) :: tend_co2glaciers(ngrid,nslope)           ! tendencies on co2 glaciers
     28  real,intent(in) :: co2ice(ngrid,nslope)                     ! co2 ice amount [kg/m^2]
     29  real,intent(in) :: waterice(ngrid,nslope)                   ! water ice amount [kg/m^2]
     30  real, intent(in) :: tsoil_PEM_yr1(ngrid,nsoil_PEM,nslope)   ! soil temperature during 1st year [K]
     31  real, intent(in) :: watersurf_ave(ngrid,nslope)             ! surface water ice density, yearly averaged  (kg/m^3)
     32  real, intent(inout) :: watersoil_ave(ngrid,nsoil_PEM,nslope)! surface water ice density, yearly averaged (kg/m^3)
     33  real, intent(in) :: global_ave_pressure                     ! mean average pressure on the planet [Pa]
    3734! outputs
    3835
    39   real,intent(inout) :: TI_PEM(ngrid,nsoil_PEM,nslope) !soil (mid-layer) thermal inertia (SI)
    40   real,intent(inout) :: tsoil_PEM(ngrid,nsoil_PEM,nslope) !soil (mid-layer) temperature (K)
    41   real,intent(inout) :: ice_table(ngrid,nslope) ! (m)
    42   real,intent(inout) :: tsoil_inst(ngrid,nsoil_PEM,nslope,timelen) ! instantaneous soil (mid-layer) temperature (k)
    43   real,intent(inout) :: m_co2_regolith_phys(ngrid,nsoil_PEM,nslope)   ! mass of co2 adsorbed (kg/m^2)
    44   real,intent(out) :: deltam_co2_regolith_phys(ngrid)                 ! mass of co2 that is exchanged due to adsorption desorption
     36  real,intent(inout) :: TI_PEM(ngrid,nsoil_PEM,nslope)                ! soil (mid-layer) thermal inertia in the PEM grid [SI]
     37  real,intent(inout) :: tsoil_PEM(ngrid,nsoil_PEM,nslope)             ! soil (mid-layer) temperature [K]
     38  real,intent(inout) :: ice_table(ngrid,nslope)                       ! Ice table depth [m]
     39  real,intent(inout) :: tsoil_inst(ngrid,nsoil_PEM,nslope,timelen)    ! instantaneous soil (mid-layer) temperature [K]
     40  real,intent(inout) :: m_co2_regolith_phys(ngrid,nsoil_PEM,nslope)   ! mass of co2 adsorbed [kg/m^2]
     41  real,intent(out) :: deltam_co2_regolith_phys(ngrid)                 ! mass of co2 that is exchanged due to adsorption desorption [kg/m^2]
    4542! local
    46    real :: tsoil_startPEM(ngrid,nsoil_PEM,nslope)
    47    real :: TI_startPEM(ngrid,nsoil_PEM,nslope)
    48    LOGICAL :: found
    49    integer :: iloop,ig,islope,it,isoil
    50    REAL :: TI_breccia = 750.
    51    REAL :: TI_bedrock = 2300.
    52    real :: kcond
    53    real :: delta
    54    CHARACTER*2 :: num
    55    real :: tsoil_saved(ngrid,nsoil_PEM)
    56    real :: tsoil_tmp_yr1(ngrid,nsoil_PEM,nslope)
    57    real :: tsoil_tmp_yr2(ngrid,nsoil_PEM,nslope)
    58    real :: tsoil_tmp(ngrid,nsoil_PEM,nslope)
    59    real :: alph_tmp(ngrid,nsoil_PEM-1)
    60    real :: beta_tmp(ngrid,nsoil_PEM-1)
    61    real :: co2_ads_prev(ngrid)
    62    real :: year_PEM_read
    63    real :: alpha_clap_h2o = -6143.7
    64    real :: beta_clap_h2o = 28.9074
     43   real :: tsoil_startPEM(ngrid,nsoil_PEM,nslope)                     ! soil temperature saved in the start [K]
     44   real :: TI_startPEM(ngrid,nsoil_PEM,nslope)                        ! soil thermal inertia saved in the start [SI]
     45   LOGICAL :: found                                                   ! check if variables are found in the start
     46   integer :: iloop,ig,islope,it,isoil                                ! index for loops
     47   REAL :: TI_breccia = 750.                                          ! Thermal inertia of Breccia following Wood 2009 [SI]
     48   REAL :: TI_bedrock = 2300.                                         ! Thermal inertia of Bedrock following Wood 2009 [SI]
     49   real :: kcond                                                      ! Thermal conductivity, intermediate variable [SI]
     50   real :: delta                                                      ! Depth of the interface regolith-breccia, breccia -bedrock [m]
     51   CHARACTER*2 :: num                                                 ! intermediate string to read PEM start sloped variables
     52   real :: tsoil_saved(ngrid,nsoil_PEM)                               ! saved soil temperature [K]
     53   real :: tsoil_tmp_yr1(ngrid,nsoil_PEM,nslope)                      ! intermediate soil temperature during yr1[K]
     54   real :: tsoil_tmp_yr2(ngrid,nsoil_PEM,nslope)                      ! intermediate soil temperature during yr 2 [K]
     55   real :: alph_tmp(ngrid,nsoil_PEM-1)                                ! Intermediate for tsoil computation []
     56   real :: beta_tmp(ngrid,nsoil_PEM-1)                                ! Intermediate for tsoil computatio []                               
     57   real :: year_PEM_read                                              ! Year of the PEM previous run
     58   real :: alpha_clap_h2o = -6143.7                                   ! Intermediate coefficient to compute psat using clapeyron law, Murphie et al. 2005 [K^-1]
     59   real :: beta_clap_h2o = 28.9074                                    ! Intermediate coefficient to compute psat using clapeyron law, Murphie et al. 2005 [1]
     60
    6561
    6662!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     
    6864!!! Purpose: read start_pem. Need a specific iostart_PEM
    6965!!!
    70 !!! Order: 1. Thermal Inertia
    71 !!!        2. Soil Temperature
    72 !!!        3. Ice table
    73 !!!        4. Mass of CO2 adsorbed
     66!!! Order: 1. Previous year of the PEM run
     67!!!        2. Thermal Inertia
     68!!!        3. Soil Temperature
     69!!!        4. Ice table
     70!!!        5. Mass of CO2 adsorbed
    7471!!!
    7572!!! /!\ This order must be respected !
     
    212209    call soil_pem_routine(ngrid,nsoil_PEM,.false.,TI_PEM(:,:,islope),timestep,tsurf_ave_yr1(:,islope),tsoil_tmp_yr1(:,:,islope),alph_tmp,beta_tmp)
    213210    tsoil_tmp_yr2(:,:,islope) = tsoil_tmp_yr1(:,:,islope)     
    214     call soil_pem_routine(ngrid,nsoil_PEM,.true.,TI_PEM(:,:,islope),timestep,tsurf_ave(:,islope),tsoil_tmp_yr2(:,:,islope),alph_tmp,beta_tmp)
     211    call soil_pem_routine(ngrid,nsoil_PEM,.true.,TI_PEM(:,:,islope),timestep,tsurf_ave_yr2(:,islope),tsoil_tmp_yr2(:,:,islope),alph_tmp,beta_tmp)
    215212
    216213
     
    419416     if(soil_pem) then
    420417
    421     DO ig = 1,ngrid
    422       DO islope = 1,nslope
    423 
    424      write(*,*) 'ig,islope ,ice table=',ig,islope,ice_table(ig,islope)
    425 
    426       ENDDO
    427     ENDDO
    428 
    429 !! small test
     418   
     419
     420!! Sanity check
    430421    DO ig = 1,ngrid
    431422      DO islope = 1,nslope
     
    441432     endif!soil_pem
    442433
    443       write(*,*) "construction ok, no nan"
    444      
     434
    445435
    446436   END SUBROUTINE
  • trunk/LMDZ.COMMON/libf/evolution/pemredem.F90

    r2842 r2855  
    11module pemredem
     2!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3!!!
     4!!! Purpose: Write specific netcdf restart for the PEM
     5!!!
     6!!!
     7!!! Author: LL, inspired by phyredem  from the GCM
     8!!!
     9!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    210
    311implicit none
     
    816                         day_ini,time,nslope,def_slope,subslope_dist)
    917! create physics restart file and write time-independent variables
    10   use comsoil_h_PEM, only: mlayer_PEM,fluxgeo,inertiedat_PEM
     18  use comsoil_h_PEM, only: soil_pem,mlayer_PEM,fluxgeo,inertiedat_PEM
    1119  use iostart_PEM, only : open_restartphy, close_restartphy, &
    1220                      put_var, put_field, length
     
    6977                      put_var, put_field
    7078
    71   use comsoil_h_PEM, only: mlayer_PEM,fluxgeo,inertiedat_PEM
     79  use comsoil_h_PEM, only: mlayer_PEM,fluxgeo,inertiedat_PEM,soil_pem
    7280  USE temps_mod_evol, ONLY: year_PEM
    73   use soil_evolution_mod, only: soil_pem
    7481               
    7582  implicit none
  • trunk/LMDZ.COMMON/libf/evolution/read_data_GCM.F90

    r2849 r2855  
    22! $Id $
    33!
    4 SUBROUTINE read_data_GCM(fichnom,min_h2o_ice_s,min_co2_ice_s,iim_input,jjm_input,nlayer,vmr_co2_gcm,ps_GCM,timelen, &
     4SUBROUTINE read_data_GCM(fichnom,timelen, iim_input,jjm_input,min_h2o_ice_s,min_co2_ice_s,vmr_co2_gcm,ps_GCM, &
    55             min_co2_ice_slope,min_h2o_ice_slope,nslope,tsurf_ave,tsoil_ave,tsurf_gcm,tsoil_gcm,TI_ave,q_co2_GCM,q_h2o_GCM,co2_ice_slope, &
    66             watersurf_density,watersoil_density)
     
    1010                        nf90_inquire_dimension,nf90_close
    1111      use comsoil_h, only: nsoilmx
    12       USE soil_evolution_mod, ONLY: soil_pem
     12      USE comsoil_h_PEM, ONLY: soil_pem
    1313
    1414      IMPLICIT NONE
     
    1616!=======================================================================
    1717!
    18 ! Read initial confitions file
     18! Purpose: Read initial confitions file from the GCM
    1919!
     20! Authors: RV & LL
    2021!=======================================================================
    2122
     
    2425!===============================================================================
    2526! Arguments:
     27
    2628  CHARACTER(LEN=*), INTENT(IN) :: fichnom          !--- FILE NAME
    2729  INTEGER, INTENT(IN) :: timelen                   ! number of times stored in the file
    28 
    29   INTEGER :: iim_input,jjm_input,nlayer,nslope
    30   REAL, ALLOCATABLE ::  h2o_ice_s(:,:,:)                       ! h2o_ice_s of the concatenated file
    31   REAL, ALLOCATABLE ::  co2_ice_s(:,:,:)                       ! co2_ice_s of the concatenated file
    32 
    33   REAL, ALLOCATABLE ::  h2o_ice_s_slope(:,:,:,:)                       ! co2_ice_s of the concatenated file
    34 
    35   REAL, INTENT(OUT) ::  min_h2o_ice_s(iim_input+1,jjm_input+1) ! Minimum of h2o_ice_s of the year
    36   REAL, INTENT(OUT) ::  min_co2_ice_s(iim_input+1,jjm_input+1) ! Minimum of co2_ice_s of the year
    37   REAL, INTENT(OUT) ::  min_co2_ice_slope(iim_input+1,jjm_input+1,nslope) ! Minimum of co2_ice slope of the year
    38   REAL, INTENT(OUT) ::  min_h2o_ice_slope(iim_input+1,jjm_input+1,nslope) ! Minimum of co2_ice slope of the year
    39   REAL, INTENT(OUT) ::  vmr_co2_gcm(iim_input+1,jjm_input+1,timelen)      !!!!vmr_co2_phys_gcm(iim_input+1,jjm_input+1,timelen)
    40   REAL, INTENT(OUT) ::  q_h2o_GCM(iim_input+1,jjm_input+1,timelen)
    41   REAL, INTENT(OUT) ::  q_co2_GCM(iim_input+1,jjm_input+1,timelen)
    42   REAL,  INTENT(OUT) ::  ps_GCM(iim_input+1,jjm_input+1,timelen)
     30  INTEGER :: iim_input,jjm_input,nslope            ! number of points in the lat x lon dynamical grid, number of subgrid slopes
     31
     32! Ouputs
     33  REAL, INTENT(OUT) ::  min_h2o_ice_s(iim_input+1,jjm_input+1) ! Minimum of h2o ice, mesh averaged of the year  [kg/m^2]
     34  REAL, INTENT(OUT) ::  min_co2_ice_s(iim_input+1,jjm_input+1) ! Minimum of co2 ice, mesh averaged  of the year [kg/m^2]
     35  REAL, INTENT(OUT) ::  min_co2_ice_slope(iim_input+1,jjm_input+1,nslope) ! Minimum of co2 ice  per slope of the year [kg/m^2]
     36  REAL, INTENT(OUT) ::  min_h2o_ice_slope(iim_input+1,jjm_input+1,nslope) ! Minimum of h2o ice per slope of the year [kg/m^2]
     37  REAL, INTENT(OUT) ::  vmr_co2_gcm(iim_input+1,jjm_input+1,timelen)      ! CO2 volume mixing ratio in the first layer  [mol/m^3]
     38  REAL, INTENT(OUT) ::  q_h2o_GCM(iim_input+1,jjm_input+1,timelen)        ! H2O mass mixing ratio in the first layer [kg/m^3]
     39  REAL, INTENT(OUT) ::  q_co2_GCM(iim_input+1,jjm_input+1,timelen)        ! CO2 mass mixing ratio in the first layer [kg/m^3]
     40  REAL,  INTENT(OUT) ::  ps_GCM(iim_input+1,jjm_input+1,timelen)          ! Surface Pressure [Pa]
     41  REAL, INTENT(OUT) ::  co2_ice_slope(iim_input+1,jjm_input+1,nslope,timelen) ! co2 ice amount per  slope of the year [kg/m^2]
    4342
    4443!SOIL
    45   REAL, INTENT(OUT) ::  tsurf_ave(iim_input+1,jjm_input+1,nslope) ! Average surface temperature of the concatenated file
    46   REAL, INTENT(OUT) ::  tsoil_ave(iim_input+1,jjm_input+1,nsoilmx,nslope) ! Average soil temperature of the concatenated file
    47 
    48   REAL ,INTENT(OUT) ::  tsurf_gcm(iim_input+1,jjm_input+1,nslope,timelen) ! Surface temperature of the concatenated file
    49   REAL , INTENT(OUT) ::  tsoil_gcm(iim_input+1,jjm_input+1,nsoilmx,nslope,timelen) ! Soil temperature of the concatenated file
    50   REAL , INTENT(OUT) ::  watersurf_density(iim_input+1,jjm_input+1,nslope,timelen) ! Soil temperature of the concatenated file
    51   REAL , INTENT(OUT) ::  watersoil_density(iim_input+1,jjm_input+1,nsoilmx,nslope,timelen) ! Soil temperature of the concatenated file
    52 
    53   REAL ::  TI_gcm(iim_input+1,jjm_input+1,nsoilmx,nslope,timelen) ! Thermal Inertia  of the concatenated file
    54   REAL, INTENT(OUT) ::  TI_ave(iim_input+1,jjm_input+1,nsoilmx,nslope) ! Average Thermal Inertia  of the concatenated file
    55   REAL, INTENT(OUT) ::  co2_ice_slope(iim_input+1,jjm_input+1,nslope,timelen) ! Minimum of co2_ice slope of the year
     44  REAL, INTENT(OUT) ::  tsurf_ave(iim_input+1,jjm_input+1,nslope)         ! Average surface temperature of the concatenated file [K]
     45  REAL, INTENT(OUT) ::  tsoil_ave(iim_input+1,jjm_input+1,nsoilmx,nslope) ! Average soil temperature of the concatenated file [K]
     46  REAL ,INTENT(OUT) ::  tsurf_gcm(iim_input+1,jjm_input+1,nslope,timelen)                  ! Surface temperature of the concatenated file, time series [K]
     47  REAL , INTENT(OUT) ::  tsoil_gcm(iim_input+1,jjm_input+1,nsoilmx,nslope,timelen)         ! Soil temperature of the concatenated file, time series [K]
     48  REAL , INTENT(OUT) ::  watersurf_density(iim_input+1,jjm_input+1,nslope,timelen)         ! Water density at the surface, time series [kg/m^3]
     49  REAL , INTENT(OUT) ::  watersoil_density(iim_input+1,jjm_input+1,nsoilmx,nslope,timelen) ! Water density in the soil layer, time series [kg/m^3]
     50  REAL ::  TI_gcm(iim_input+1,jjm_input+1,nsoilmx,nslope,timelen)                          ! Thermal Inertia  of the concatenated file, times series [SI]
     51  REAL, INTENT(OUT) ::  TI_ave(iim_input+1,jjm_input+1,nsoilmx,nslope)                     ! Average Thermal Inertia  of the concatenated file [SI]
    5652!===============================================================================
    5753!   Local Variables
    58   CHARACTER(LEN=256) :: msg, var, modname
    59   INTEGER,PARAMETER :: length=100
    60   INTEGER :: iq, fID, vID, idecal
    61   INTEGER :: ierr
     54  CHARACTER(LEN=256) :: msg, var, modname               ! for reading
     55  INTEGER :: iq, fID, vID, idecal                       ! for reading
     56  INTEGER :: ierr                                       ! for reading
    6257  CHARACTER(len=12) :: start_file_type="earth" ! default start file type
    6358
     
    6661
    6762  INTEGER :: edges(4),corner(4)
    68   INTEGER :: i,j,t
    69   real,save :: m_co2, m_noco2, A , B, mmean
    70 
    71   INTEGER :: islope
    72   CHARACTER*2 :: num
     63  INTEGER :: i,j,t                                                     ! loop variables
     64  real,save :: m_co2, m_noco2, A , B, mmean                            ! Molar Mass of co2 and no co2, A;B intermediate variables to compute the mean molar mass of the layer
     65
     66  INTEGER :: islope                                                    ! loop for variables
     67  CHARACTER*2 :: num                                                   ! for reading sloped variables
     68  REAL, ALLOCATABLE ::  h2o_ice_s(:,:,:)                               ! h2o ice, mesh averaged, of the concatenated file [kg/m^2]
     69  REAL, ALLOCATABLE ::  co2_ice_s(:,:,:)                               ! co2 ice, mesh averaged, of the concatenated file [kg/m^2]
     70  REAL, ALLOCATABLE ::  h2o_ice_s_slope(:,:,:,:)                       ! h2o ice per slope of the concatenated file [kg/m^2]
    7371
    7472!-----------------------------------------------------------------------
    75   modname="pemetat0"
     73  modname="read_data_gcm"
    7674
    7775      m_co2 = 44.01E-3  ! CO2 molecular mass (kg/mol)   
  • trunk/LMDZ.COMMON/libf/evolution/recomp_orb_param_mod.F90

    r2849 r2855  
    11      MODULE recomp_orb_param_mod
    2 
     2!=======================================================================
     3!
     4! Purpose: Recompute orbit parameters based on values by Laskar et al.,
     5! 2004
     6!
     7! Authors: RV
     8!=======================================================================
    39      IMPLICIT NONE
    410
  • trunk/LMDZ.COMMON/libf/evolution/reshape_XIOS_output.F90

    r2850 r2855  
    11program reshape_XIOS_output
     2
     3!=======================================================================
     4!
     5! Purpose: Read XIOS files, and convert them into the correct GCM grid
     6!          XIOS  longitudes start at -180 but stop before -180 (not duplicated)
     7!          We basically add the last point, and complete the XIOS file. Looped
     8!          over the two GCM runs
     9!
     10! Authors: RV & LL
     11!=======================================================================
    212    use netcdf
    313    implicit none
  • trunk/LMDZ.COMMON/libf/evolution/soil_TIfeedback_PEM.F90

    r2835 r2855  
    3333      INTEGER, INTENT(IN) :: nsoil                  ! Number of soil layers
    3434      REAL :: icedepth                  ! Ice cover thickness (m)
    35       REAL :: inert_h2o_ice = 800.
    36       REAL :: rho_ice = 920.
    37       REAL :: prev_thermi(ngrid,nsoil)
     35      REAL :: inert_h2o_ice = 800.      ! surface water ice thermal inertia [SI]
     36      REAL :: rho_ice = 920.            ! density of water ice [kg/m^3]
     37      REAL :: prev_thermi(ngrid,nsoil)  ! previous thermal inertia [SI]
    3838!Inputs
    3939!------
  • trunk/LMDZ.COMMON/libf/evolution/soil_pem.F90

    r2835 r2855  
    2525!  ---------
    2626!  inputs:
    27       integer,intent(in) :: ngrid       ! number of (horizontal) grid-points
    28       integer,intent(in) :: nsoil       ! number of soil layers 
    29       logical,intent(in) :: firstcall ! identifier for initialization call
    30       real,intent(in) :: therm_i(ngrid,nsoil) ! thermal inertia
    31       real,intent(in) :: timestep           ! time step
    32       real,intent(in) :: tsurf(ngrid)   ! surface temperature
     27      integer,intent(in) :: ngrid       ! number of (horizontal) grid-points
     28      integer,intent(in) :: nsoil       ! number of soil layers 
     29      logical,intent(in) :: firstcall    ! identifier for initialization call
     30      real,intent(in) :: therm_i(ngrid,nsoil) ! thermal inertia [SI]
     31      real,intent(in) :: timestep       ! time step [s]
     32      real,intent(in) :: tsurf(ngrid)   ! surface temperature [K]
    3333 
    3434! outputs:
    35       real,intent(inout) :: tsoil(ngrid,nsoil) ! soil (mid-layer) temperature
    36       real,intent(inout) :: alph_PEM(ngrid,nsoil-1)
    37       real,intent(inout) :: beta_PEM(ngrid,nsoil-1)
     35      real,intent(inout) :: tsoil(ngrid,nsoil) ! soil (mid-layer) temperature [K]
     36      real,intent(inout) :: alph_PEM(ngrid,nsoil-1) ! intermediate for computations [1]
     37      real,intent(inout) :: beta_PEM(ngrid,nsoil-1) ! intermediate for computations [K]
    3838     
    3939! local variables:
  • trunk/LMDZ.COMMON/libf/evolution/soil_settings_PEM.F

    r2835 r2855  
    1212
    1313!======================================================================
    14 !  Author: Ehouarn Millour (07/2006)
     14!  Author: LL, based on work by  Ehouarn Millour (07/2006)
    1515!
    1616!  Purpose: Read and/or initialise soil depths and properties
     
    2222!
    2323!
    24 !  This subroutine reads from a NetCDF file (opened by the caller)
    25 !  of "startfi.nc" format.
    2624!  The various actions and variable read/initialized are:
    27 !  1. Check out the number of soil layers (in datafile); if it isn't equal
    28 !     to nsoil, then some interpolation will be required
    29 !     Also check if data in file "startfi.nc" is in older format (ie:
    30 !     thermal inertia was depth-independent; and there was no "depth"
    31 !     coordinate.
    32 !     Read/build layer (and midlayer) depths
    33 !  2. Read volumetric specific heat (or initialise it to default value)
    34 !  3. Read Thermal inertia
    35 !  4. Read soil temperatures
    36 !  5. Interpolate thermal inertia and temperature on the new grid, if
     25!  1. Read/build layer (and midlayer) depth
     26!  2. Interpolate thermal inertia and temperature on the new grid, if
    3727!     necessary
    3828!======================================================================
     
    4636      integer,intent(in) :: nsoil_PEM   ! # of soil layers in the PEM
    4737      integer,intent(in) :: nsoil_GCM   ! # of soil layers in the GCM
    48       real,intent(in) :: TI_GCM(ngrid,nsoil_GCM,nslope) ! # of soil layers in the GCM
    49       real,intent(inout) :: TI_PEM(ngrid,nsoil_PEM,nslope)      ! # of soil layers in the PEM
     38      real,intent(in) :: TI_GCM(ngrid,nsoil_GCM,nslope) !  Thermal inertia  in the GCM [SI]
     39      real,intent(inout) :: TI_PEM(ngrid,nsoil_PEM,nslope)      ! Thermal inertia   in the PEM [SI]
    5040
    5141!======================================================================
Note: See TracChangeset for help on using the changeset viewer.