Changeset 3456


Ignore:
Timestamp:
Oct 11, 2024, 10:37:48 AM (6 weeks ago)
Author:
jbclement
Message:

PEM:
Cleaning of the adsorption module to make the debugging easier.
JBC

Location:
trunk/LMDZ.COMMON/libf/evolution
Files:
4 edited

Legend:

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

    r3264 r3456  
    1 module adsorption_mod
    2 
    3   implicit none
    4 
    5   LOGICAL adsorption_pem ! True by default, to compute adsorption/desorption. Read in  pem.def
    6   real, save, allocatable :: co2_adsorbded_phys(:,:,:)  ! co2 that is in the regolith [kg/m^2]
    7   real, save, allocatable :: h2o_adsorbded_phys(:,:,:)  ! h2o that is in the regolith [kg/m^2]
    8 
    9   contains
     1MODULE adsorption_mod
     2
     3implicit none
     4
     5logical                                   :: adsorption_pem     ! True by default, to compute adsorption/desorption. Read in pem.def
     6real, save, allocatable, dimension(:,:,:) :: co2_adsorbded_phys ! co2 that is in the regolith [kg/m^2]
     7real, save, allocatable, dimension(:,:,:) :: h2o_adsorbded_phys ! h2o that is in the regolith [kg/m^2]
     8
     9!=======================================================================
     10contains
     11!=======================================================================
    1012
    1113!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     
    1618!!!
    1719!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    18   subroutine ini_adsorption_h_PEM(ngrid,nslope,nsoilmx_PEM)
    19 
    20   implicit none
    21   integer,intent(in) :: ngrid ! number of atmospheric columns
    22   integer,intent(in) :: nslope ! number of slope within a mesh
    23   integer,intent(in) :: nsoilmx_PEM ! number of soil layer in the PEM
    24     allocate(co2_adsorbded_phys(ngrid,nsoilmx_PEM,nslope))
    25     allocate(h2o_adsorbded_phys(ngrid,nsoilmx_PEM,nslope))
    26   end subroutine ini_adsorption_h_PEM
    27 
    28 !!! -----------------------------------------------
    29 
    30   subroutine end_adsorption_h_PEM
    31 
    32   implicit none
    33     if (allocated(co2_adsorbded_phys)) deallocate(co2_adsorbded_phys)
    34     if (allocated(h2o_adsorbded_phys)) deallocate(h2o_adsorbded_phys)
    35   end subroutine end_adsorption_h_PEM
    36 
    37 !!! -----------------------------------------------
    38 
    39   subroutine regolith_adsorption(ngrid,nslope,nsoil_PEM,timelen,tend_h2oglaciers,tend_co2glaciers,waterice,co2ice,tsoil_PEM,TI_PEM,ps,q_co2,q_h2o, &
    40                                 m_h2o_completesoil,delta_mh2oreg, m_co2_completesoil,delta_mco2reg)
    41 
    42 ! inputs
    43  INTEGER,INTENT(IN) :: ngrid, nslope, nsoil_PEM,timelen ! size dimension: physics x subslope x soil x timeseries
    44  REAL,INTENT(IN) :: tend_h2oglaciers(ngrid,nslope),tend_co2glaciers(ngrid,nslope) !tendancies on the glaciers [1]
    45  REAL,INTENT(IN) :: waterice(ngrid,nslope)              ! water ice at the surface [kg/m^2]
    46  REAL,INTENT(IN) :: co2ice(ngrid,nslope)                ! co2 ice at the surface [kg/m^2]
    47  REAL,INTENT(IN) :: TI_PEM(ngrid,nsoil_PEM,nslope)      ! Soil Thermal inertia (J/K/^2/s^1/2)
    48  REAL,INTENT(IN) :: tsoil_PEM(ngrid,nsoil_PEM,nslope)   ! Soil temperature (K)
    49  REAL,INTENT(IN) :: ps(ngrid,timelen)                   ! Average surface pressure [Pa]
    50  REAL,INTENT(IN) :: q_co2(ngrid,timelen)                ! Mass mixing ratio of co2 in the first layer (kg/kg)
    51  REAL,INTENT(IN) :: q_h2o(ngrid,timelen)                ! Mass mixing ratio of H2o in the first layer (kg/kg)
    52 
    53 ! outputs
    54  REAL,INTENT(INOUT) :: m_h2o_completesoil(ngrid,nsoil_PEM,nslope) ! Density of h2o adsorbed (kg/m^3)(ngrid,nsoil_PEM,nslope)     
    55  REAL,INTENT(OUT) :: delta_mh2oreg(ngrid)                         ! Difference density of h2o adsorbed (kg/m^3)
    56 
    57  REAL,INTENT(INOUT) :: m_co2_completesoil(ngrid,nsoil_PEM,nslope)   ! Density of co2 adsorbed (kg/m^3)
    58  REAL,INTENT(OUT) :: delta_mco2reg(ngrid)                            ! Difference density of co2 adsorbed (kg/m^3)
    59  
    60 ! local variables
    61  REAL :: theta_h2o_adsorbded(ngrid,nsoil_PEM,nslope) ! Fraction of the pores occupied by H2O molecules
     20SUBROUTINE ini_adsorption_h_PEM(ngrid,nslope,nsoilmx_PEM)
     21
     22implicit none
     23
     24integer, intent(in) :: ngrid       ! number of atmospheric columns
     25integer, intent(in) :: nslope      ! number of slope within a mesh
     26integer, intent(in) :: nsoilmx_PEM ! number of soil layer in the PEM
     27
     28allocate(co2_adsorbded_phys(ngrid,nsoilmx_PEM,nslope))
     29allocate(h2o_adsorbded_phys(ngrid,nsoilmx_PEM,nslope))
     30
     31END SUBROUTINE ini_adsorption_h_PEM
     32
     33!=======================================================================
     34SUBROUTINE end_adsorption_h_PEM
     35
     36if (allocated(co2_adsorbded_phys)) deallocate(co2_adsorbded_phys)
     37if (allocated(h2o_adsorbded_phys)) deallocate(h2o_adsorbded_phys)
     38
     39END SUBROUTINE end_adsorption_h_PEM
     40
     41!=======================================================================
     42SUBROUTINE regolith_adsorption(ngrid,nslope,nsoil_PEM,timelen,tend_h2oglaciers,tend_co2glaciers,waterice,co2ice,tsoil_PEM,TI_PEM,ps,q_co2,q_h2o, &
     43                               m_h2o_completesoil,delta_mh2oreg, m_co2_completesoil,delta_mco2reg)
     44
     45implicit none
     46
     47! Inputs
     48integer,                                    intent(in) :: ngrid, nslope, nsoil_PEM, timelen  ! size dimension: physics x subslope x soil x timeseries
     49real,    dimension(ngrid,nslope),           intent(in) :: tend_h2oglaciers, tend_co2glaciers ! tendancies on the glaciers [1]
     50real,    dimension(ngrid,nslope),           intent(in) :: waterice                           ! water ice at the surface [kg/m^2]
     51real,    dimension(ngrid,nslope),           intent(in) :: co2ice                             ! co2 ice at the surface [kg/m^2]
     52real,    dimension(ngrid,nsoil_PEM,nslope), intent(in) :: TI_PEM                             ! Soil Thermal inertia (J/K/^2/s^1/2)
     53real,    dimension(ngrid,nsoil_PEM,nslope), intent(in) :: tsoil_PEM                          ! Soil temperature (K)
     54real,    dimension(ngrid,timelen),          intent(in) :: ps                                 ! Average surface pressure [Pa]
     55real,    dimension(ngrid,timelen),          intent(in) :: q_co2                              ! Mass mixing ratio of co2 in the first layer (kg/kg)
     56real,    dimension(ngrid,timelen),          intent(in) :: q_h2o                              ! Mass mixing ratio of H2o in the first layer (kg/kg)
     57
     58! Outputs
     59real, dimension(ngrid), intent(out) :: delta_mh2oreg ! Difference density of h2o adsorbed (kg/m^3)
     60real, dimension(ngrid), intent(out) :: delta_mco2reg ! Difference density of co2 adsorbed (kg/m^3)
     61real, dimension(ngrid,nsoil_PEM,nslope), intent(inout) :: m_co2_completesoil ! Density of co2 adsorbed (kg/m^3)
     62real, dimension(ngrid,nsoil_PEM,nslope), intent(inout) :: m_h2o_completesoil ! Density of h2o adsorbed (kg/m^3)
     63
     64! Local variables
     65real, dimension(ngrid,nsoil_PEM,nslope) :: theta_h2o_adsorbded ! Fraction of the pores occupied by H2O molecules
    6266! -------------
    63 
    64 
    6567! Compute H2O adsorption, then CO2 adsorption
    66 
    67   call regolith_h2oadsorption(ngrid,nslope,nsoil_PEM,timelen,tend_h2oglaciers,tend_co2glaciers,waterice,co2ice,ps,q_co2,q_h2o,tsoil_PEM,TI_PEM, &
    68                                    theta_h2o_adsorbded,m_h2o_completesoil,delta_mh2oreg)
    69 
    70 
    71   call regolith_co2adsorption(ngrid,nslope,nsoil_PEM,timelen,tend_h2oglaciers,tend_co2glaciers,waterice,co2ice,ps,q_co2,q_h2o, &
    72                                                           tsoil_PEM,TI_PEM,m_co2_completesoil,delta_mco2reg)
    73 
    74    RETURN
    75   end subroutine
    76 
    77 !------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
    78 
    79   subroutine regolith_h2oadsorption(ngrid,nslope,nsoil_PEM,timelen,tend_h2oglaciers,tend_co2glaciers,waterice,co2ice,ps,q_co2,q_h2o,tsoil_PEM,TI_PEM, &
    80                                    theta_h2o_adsorbded,m_h2o_completesoil,delta_mreg)
     68call regolith_h2oadsorption(ngrid,nslope,nsoil_PEM,timelen,tend_h2oglaciers,tend_co2glaciers,waterice,co2ice,ps,q_co2,q_h2o,tsoil_PEM,TI_PEM, &
     69                            theta_h2o_adsorbded,m_h2o_completesoil,delta_mh2oreg)
     70call regolith_co2adsorption(ngrid,nslope,nsoil_PEM,timelen,tend_h2oglaciers,tend_co2glaciers,waterice,co2ice,ps,q_co2,q_h2o, &
     71                            tsoil_PEM,TI_PEM,m_co2_completesoil,delta_mco2reg)
     72
     73END SUBROUTINE
     74
     75!=======================================================================
     76SUBROUTINE regolith_h2oadsorption(ngrid,nslope,nsoil_PEM,timelen,tend_h2oglaciers,tend_co2glaciers,waterice,co2ice,ps,q_co2,q_h2o,tsoil_PEM,TI_PEM, &
     77                                  theta_h2o_adsorbded,m_h2o_completesoil,delta_mreg)
    8178
    8279!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     
    9087use comsoil_h_PEM,         only: layer_PEM, index_breccia
    9188use comslope_mod,          only: subslope_dist, def_slope_mean
    92 use vertical_layers_mod,   only: ap,bp
     89use vertical_layers_mod,   only: ap, bp
    9390use constants_marspem_mod, only: alpha_clap_h2o, beta_clap_h2o, m_h2o, m_co2,m_noco2, rho_regolith
    9491
     
    9996#endif
    10097
    101  implicit none
    102 ! inputs
    103  INTEGER,INTENT(IN) :: ngrid, nslope, nsoil_PEM,timelen ! Size dimension
    104  REAL,INTENT(IN) :: tend_h2oglaciers(ngrid,nslope),tend_co2glaciers(ngrid,nslope) ! Tendencies on the glaciers ()
    105  REAL,INTENT(IN) :: waterice(ngrid,nslope)              ! Water ice at the surface [kg/m^2]
    106  REAL,INTENT(IN) :: co2ice(ngrid,nslope)                ! CO2 ice at the surface [kg/m^2]
    107  REAL,INTENT(IN) :: ps(ngrid,timelen)                   ! Surface pressure (Pa)         
    108  REAL,INTENT(IN) :: q_co2(ngrid,timelen)                ! Mass mixing ratio of co2 in the first layer (kg/kg)
    109  REAL,INTENT(IN) :: q_h2o(ngrid,timelen)                ! Mass mixing ratio of H2o in the first layer (kg/kg)
    110  REAL,INTENT(IN) :: TI_PEM(ngrid,nsoil_PEM,nslope)      ! Soil Thermal inertia (J/K/^2/s^1/2)
    111  REAL,INTENT(IN) :: tsoil_PEM(ngrid,nsoil_PEM,nslope)   ! Soil temperature (K)
    112 
    113 ! outputs
    114  REAL,INTENT(INOUT) :: m_h2o_completesoil(ngrid,nsoil_PEM,nslope) ! Density of h2o adsorbed (kg/m^3)(ngrid,nsoil_PEM,nslope)     
    115  REAL,INTENT(OUT) :: theta_h2o_adsorbded(ngrid,nsoil_PEM,nslope)  ! Fraction of the pores occupied by H2O molecules
    116  REAL,INTENT(OUT) :: delta_mreg(ngrid)                            ! Difference density of h2o adsorbed (kg/m^3)
    117 
    118 ! constants
    119  REAL :: Ko =  1.57e-8            ! Jackosky et al. 1997
    120  REAL :: e = 2573.9               ! Jackosky et al. 1997
    121  REAL :: mu = 0.48                ! Jackosky et al. 1997
    122  real :: m_theta = 2.84e-7        ! Mass of h2o per m^2 absorbed Jackosky et al. 1997
    123 ! real :: as = 18.9e3             ! Specific area, Buhler & Piqueux 2021
    124  real :: as = 9.48e4              ! Specific area, Zent
    125  real ::  inertie_thresold = 800. ! TI > 800 means cementation
    126 
    127  ! local variables
    128  REAL :: deltam_reg_complete(ngrid,index_breccia,nslope) ! Difference in the mass per slope and soil layer (kg/m^3)
    129  real :: K                        ! Used to compute theta
    130  integer ig, iloop, islope, it    ! For loops
    131  INTEGER :: ispermanent_co2glaciers(ngrid,nslope)      ! Check if the co2 glacier is permanent
    132  INTEGER :: ispermanent_h2oglaciers(ngrid,nslope)      ! Check if the h2o glacier is permanent
    133  REAL :: deltam_reg_slope(ngrid,nslope)                ! Difference density of h2o adsorbed per slope (kg/m^3)
    134  REAL :: dm_h2o_regolith_slope(ngrid,nsoil_PEM,nslope) ! Elementary h2o mass adsorded per mesh per slope
    135  real :: A,B                                           ! Used to compute the mean mass above the surface
    136  real :: p_sat                                         ! Saturated vapor pressure of ice
    137  real,allocatable :: mass_mean(:,:)                    ! Mean mass above the surface
    138  real,allocatable :: zplev_mean(:,:)                   ! Pressure above the surface
    139  real,allocatable :: pvapor(:,:)                       ! Partial pressure above the surface
    140  real, allocatable :: pvapor_avg(:)                    ! Yearly averaged
     98implicit none
     99
     100! Inputs
     101integer,                                    intent(in) :: ngrid, nslope, nsoil_PEM,timelen   ! Size dimension
     102real,    dimension(ngrid,nslope),           intent(in) :: tend_h2oglaciers, tend_co2glaciers ! Tendencies on the glaciers ()
     103real,    dimension(ngrid,nslope),           intent(in) :: waterice                           ! Water ice at the surface [kg/m^2]
     104real,    dimension(ngrid,nslope),           intent(in) :: co2ice                             ! CO2 ice at the surface [kg/m^2]
     105real,    dimension(ngrid,timelen),          intent(in) :: ps                                 ! Surface pressure (Pa)
     106real,    dimension(ngrid,timelen),          intent(in) :: q_co2                              ! Mass mixing ratio of co2 in the first layer (kg/kg)
     107real,    dimension(ngrid,timelen),          intent(in) :: q_h2o                              ! Mass mixing ratio of H2o in the first layer (kg/kg)
     108real,    dimension(ngrid,nsoil_PEM,nslope), intent(in) :: TI_PEM                             ! Soil Thermal inertia (J/K/^2/s^1/2)
     109real,    dimension(ngrid,nsoil_PEM,nslope), intent(in) :: tsoil_PEM                          ! Soil temperature (K)
     110
     111! Outputs
     112real, dimension(ngrid,nsoil_PEM,nslope), intent(inout) :: m_h2o_completesoil ! Density of h2o adsorbed (kg/m^3)(ngrid,nsoil_PEM,nslope)
     113real, dimension(ngrid,nsoil_PEM,nslope), intent(out) :: theta_h2o_adsorbded ! Fraction of the pores occupied by H2O molecules
     114real, dimension(ngrid),                  intent(out) :: delta_mreg          ! Difference density of h2o adsorbed (kg/m^3)
     115
     116! Constants
     117real :: Ko =  1.57e-8            ! Jackosky et al. 1997
     118real :: e = 2573.9               ! Jackosky et al. 1997
     119real :: mu = 0.48                ! Jackosky et al. 1997
     120real :: m_theta = 2.84e-7        ! Mass of h2o per m^2 absorbed Jackosky et al. 1997
     121! real :: as = 18.9e3              ! Specific area, Buhler & Piqueux 2021
     122real :: as = 9.48e4              ! Specific area, Zent
     123real ::  inertie_thresold = 800. ! TI > 800 means cementation
     124
     125! Local variables
     126real,    dimension(ngrid,index_breccia,nslope) :: deltam_reg_complete     ! Difference in the mass per slope and soil layer (kg/m^3)
     127real                                           :: K                       ! Used to compute theta
     128integer                                        :: ig, iloop, islope, it   ! For loops
     129logical, dimension(ngrid,nslope)               :: ispermanent_co2glaciers ! Check if the co2 glacier is permanent
     130logical, dimension(ngrid,nslope)               :: ispermanent_h2oglaciers ! Check if the h2o glacier is permanent
     131real,    dimension(ngrid,nslope)               :: deltam_reg_slope        ! Difference density of h2o adsorbed per slope (kg/m^3)
     132real,    dimension(ngrid,nsoil_PEM,nslope)     :: dm_h2o_regolith_slope   ! Elementary h2o mass adsorded per mesh per slope
     133real                                           :: A, B                    ! Used to compute the mean mass above the surface
     134real                                           :: p_sat                   ! Saturated vapor pressure of ice
     135real,    dimension(:,:), allocatable           :: mass_mean               ! Mean mass above the surface
     136real,    dimension(:,:), allocatable           :: zplev_mean              ! Pressure above the surface
     137real,    dimension(:,:), allocatable           :: pvapor                  ! Partial pressure above the surface
     138real,    dimension(:)  , allocatable           :: pvapor_avg              ! Yearly averaged
    141139
    142140! 0. Some initializations
    143 
    144 
    145      allocate(mass_mean(ngrid,timelen))
    146      allocate(zplev_mean(ngrid,timelen))
    147      allocate(pvapor(ngrid,timelen))
    148      allocate(pvapor_avg(ngrid))
    149      A =(1/m_co2 - 1/m_noco2)
    150      B=1/m_noco2
    151      theta_h2o_adsorbded(:,:,:) = 0.
    152      dm_h2o_regolith_slope(:,:,:) = 0.
    153 
    154 #ifndef CPP_STD
    155 
    156 !0.1 Look at perennial ice
    157   do ig = 1,ngrid
    158     do islope = 1,nslope
    159      if((abs(tend_h2oglaciers(ig,islope)).gt.1e-5).and.(abs(waterice(ig,islope)).gt.0)) then
    160         ispermanent_h2oglaciers(ig,islope) = 1
    161      else
    162         ispermanent_h2oglaciers(ig,islope) = 0
    163      endif
    164 
    165      if((abs(tend_co2glaciers(ig,islope)).gt.1e-5).and.(abs(co2ice(ig,islope)).gt.0)) then
    166         ispermanent_co2glaciers(ig,islope) = 1
    167      else
    168         ispermanent_co2glaciers(ig,islope) = 0
    169      endif
    170     enddo
    171    enddo
    172 
    173 !   0.2 Compute the partial pressure of vapor
    174 !a. the molecular mass into the column
    175      do ig = 1,ngrid
    176        mass_mean(ig,:) = 1/(A*q_co2(ig,:) +B)
    177      enddo
    178 
    179 
    180 ! b. pressure level
    181      do it = 1,timelen
    182        do ig = 1,ngrid
    183          zplev_mean(ig,it) = ap(1) + bp(1)*ps(ig,it)
    184        enddo
    185      enddo
    186 ! c. Vapor pressure
    187      pvapor(:,:) = mass_mean(:,:)/m_h2o*q_h2o(:,:)*zplev_mean(:,:)
    188      pvapor_avg(:) = sum(pvapor(:,:),2)/timelen
    189 #endif
    190      deallocate(pvapor)
    191      deallocate(zplev_mean)
    192      deallocate(mass_mean)
    193 #ifndef CPP_STD
    194 
    195 ! 1. we compute the mass of H2O adsorded in each layer of the meshes 
    196 
    197  do ig = 1,ngrid
    198   do islope = 1,nslope
    199     do iloop = 1,index_breccia
    200       K = Ko*exp(e/tsoil_PEM(ig,iloop,islope))
    201       if(TI_PEM(ig,iloop,islope).lt.inertie_thresold)  then
    202         theta_h2o_adsorbded(ig,iloop,islope) = (K*pvapor_avg(ig)/(1+K*pvapor_avg(ig)))**mu
    203       else
    204         p_sat =exp(beta_clap_h2o/tsoil_PEM(ig,iloop,islope) +alpha_clap_h2o) ! we assume fixed temperature in the ice ... not really good but ...
    205         theta_h2o_adsorbded(ig,iloop,islope) = (K*p_sat/(1+K*p_sat))**mu
    206       endif
    207         dm_h2o_regolith_slope(ig,iloop,islope) = as*theta_h2o_adsorbded(ig,iloop,islope)*m_theta*rho_regolith
    208    enddo
    209   enddo
    210  enddo
    211 
    212 ! 2. Check the exchange between the atmosphere and the regolith
    213 
    214   do ig = 1,ngrid
    215    delta_mreg(ig) = 0.
    216    do islope = 1,nslope
    217     deltam_reg_slope(ig,islope) = 0.
    218     do iloop = 1,index_breccia
    219        if((TI_PEM(ig,iloop,islope).lt.inertie_thresold).and.(ispermanent_h2oglaciers(ig,islope).eq.0).and.(ispermanent_co2glaciers(ig,islope).eq.0)) then
    220              if(iloop==1) then
    221                  deltam_reg_complete(ig,iloop,islope) = (dm_h2o_regolith_slope(ig,iloop,islope) - m_h2o_completesoil(ig,iloop,islope)) &
    222                                                    *(layer_PEM(iloop))
    223              else
    224                  deltam_reg_complete(ig,iloop,islope) = (dm_h2o_regolith_slope(ig,iloop,islope) - m_h2o_completesoil(ig,iloop,islope)) &
    225                                                    *(layer_PEM(iloop) - layer_PEM(iloop-1))
    226              endif
    227        else ! NO EXCHANGE AS ICE BLOCK THE DYNAMIC!
    228              deltam_reg_complete(ig,iloop,islope) = 0.
    229        endif
    230        deltam_reg_slope(ig,islope) = deltam_reg_slope(ig,islope) + deltam_reg_complete(ig,iloop,islope)
    231     enddo
    232    delta_mreg(ig) = delta_mreg(ig) + deltam_reg_slope(ig,islope)*subslope_dist(ig,islope)/cos(pi*def_slope_mean(islope)/180.)
    233    enddo
    234   enddo
    235    m_h2o_completesoil(:,:,:) = dm_h2o_regolith_slope(:,:,:)
    236 
    237  RETURN
    238 #endif
    239   end subroutine
    240 
    241 !------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
     141allocate(mass_mean(ngrid,timelen),zplev_mean(ngrid,timelen),pvapor(ngrid,timelen),pvapor_avg(ngrid))
     142A = 1./m_co2 - 1./m_noco2
     143B = 1./m_noco2
     144theta_h2o_adsorbded = 0.
     145dm_h2o_regolith_slope = 0.
     146ispermanent_h2oglaciers = .false.
     147ispermanent_co2glaciers = .false.
     148
     149#ifndef CPP_STD
     150    ! 0.1 Look at perennial ice
     151    do ig = 1,ngrid
     152        do islope = 1,nslope
     153            if (abs(tend_h2oglaciers(ig,islope)) > 1.e-5 .and. abs(waterice(ig,islope)) > 0.) ispermanent_h2oglaciers(ig,islope) = .true.
     154            if (abs(tend_co2glaciers(ig,islope)) > 1.e-5 .and. abs(co2ice(ig,islope)) > 0.) ispermanent_co2glaciers(ig,islope) = .true.
     155        enddo
     156    enddo
     157
     158    ! 0.2 Compute the partial pressure of vapor
     159    ! a. the molecular mass into the column
     160    do ig = 1,ngrid
     161        mass_mean(ig,:) = 1/(A*q_co2(ig,:) + B)
     162    enddo
     163
     164    ! b. pressure level
     165    do it = 1,timelen
     166        do ig = 1,ngrid
     167            zplev_mean(ig,it) = ap(1) + bp(1)*ps(ig,it)
     168        enddo
     169    enddo
     170
     171    ! c. Vapor pressure
     172    pvapor = mass_mean/m_h2o*q_h2o*zplev_mean
     173    pvapor_avg = sum(pvapor,2)/timelen
     174#endif
     175deallocate(pvapor,zplev_mean,mass_mean)
     176
     177#ifndef CPP_STD
     178    ! 1. we compute the mass of H2O adsorded in each layer of the meshes
     179    do ig = 1,ngrid
     180        do islope = 1,nslope
     181            do iloop = 1,index_breccia
     182                K = Ko*exp(e/tsoil_PEM(ig,iloop,islope))
     183                if (TI_PEM(ig,iloop,islope) < inertie_thresold) then
     184                    theta_h2o_adsorbded(ig,iloop,islope) = (K*pvapor_avg(ig)/(1. + K*pvapor_avg(ig)))**mu
     185                else
     186                    p_sat = exp(beta_clap_h2o/tsoil_PEM(ig,iloop,islope) + alpha_clap_h2o) ! we assume fixed temperature in the ice ... not really good but ...
     187                    theta_h2o_adsorbded(ig,iloop,islope) = (K*p_sat/(1. + K*p_sat))**mu
     188                endif
     189                dm_h2o_regolith_slope(ig,iloop,islope) = as*theta_h2o_adsorbded(ig,iloop,islope)*m_theta*rho_regolith
     190            enddo
     191        enddo
     192    enddo
     193
     194    ! 2. Check the exchange between the atmosphere and the regolith
     195    do ig = 1,ngrid
     196        delta_mreg(ig) = 0.
     197        do islope = 1,nslope
     198            deltam_reg_slope(ig,islope) = 0.
     199            do iloop = 1,index_breccia
     200                if (TI_PEM(ig,iloop,islope) < inertie_thresold .and. .not. ispermanent_h2oglaciers(ig,islope) .and. .not. ispermanent_co2glaciers(ig,islope)) then
     201                    if (iloop == 1) then
     202                        deltam_reg_complete(ig,iloop,islope) = (dm_h2o_regolith_slope(ig,iloop,islope) - m_h2o_completesoil(ig,iloop,islope))*(layer_PEM(iloop))
     203                    else
     204                        deltam_reg_complete(ig,iloop,islope) = (dm_h2o_regolith_slope(ig,iloop,islope) - m_h2o_completesoil(ig,iloop,islope))*(layer_PEM(iloop) - layer_PEM(iloop - 1))
     205                    endif
     206                else ! NO EXCHANGE AS ICE BLOCK THE DYNAMIC!
     207                    deltam_reg_complete(ig,iloop,islope) = 0.
     208                endif
     209                deltam_reg_slope(ig,islope) = deltam_reg_slope(ig,islope) + deltam_reg_complete(ig,iloop,islope)
     210            enddo
     211            delta_mreg(ig) = delta_mreg(ig) + deltam_reg_slope(ig,islope)*subslope_dist(ig,islope)/cos(pi*def_slope_mean(islope)/180.)
     212        enddo
     213    enddo
     214    m_h2o_completesoil = dm_h2o_regolith_slope
     215#endif
     216END SUBROUTINE
     217
     218!=======================================================================
    242219!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    243220!!!
     
    247224!!!
    248225!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    249 
    250   SUBROUTINE regolith_co2adsorption(ngrid,nslope,nsoil_PEM,timelen,tend_h2oglaciers,tend_co2glaciers,waterice,co2ice,ps,q_co2,q_h2o,&
    251                                    tsoil_PEM,TI_PEM,m_co2_completesoil,delta_mreg)
     226SUBROUTINE regolith_co2adsorption(ngrid,nslope,nsoil_PEM,timelen,tend_h2oglaciers,tend_co2glaciers,waterice,co2ice,ps,q_co2,q_h2o,tsoil_PEM,TI_PEM,m_co2_completesoil,delta_mreg)
    252227
    253228use comsoil_h_PEM,         only: layer_PEM, index_breccia, index_breccia
     
    262237#endif
    263238
    264       IMPLICIT NONE
    265 ! Inputs: 
    266  INTEGER,INTENT(IN) :: ngrid, nslope, nsoil_PEM,timelen             ! Size dimension
    267  REAL,INTENT(IN) :: ps(ngrid,timelen)                               ! Average surface pressure [Pa]
    268  REAL,INTENT(IN) :: tsoil_PEM(ngrid,nsoil_PEM,nslope)               ! Mean Soil Temperature [K]
    269  REAL,INTENT(IN) :: TI_PEM(ngrid,nsoil_PEM,nslope)                  ! Mean Thermal Inertia [USI]
    270  REAL,INTENT(IN) :: tend_h2oglaciers(ngrid,nslope),tend_co2glaciers(ngrid,nslope) ! Tendencies on the glaciers ()
    271  REAL,INTENT(IN) :: q_co2(ngrid,timelen),q_h2o(ngrid,timelen)       ! Mass mixing ratio of co2 and h2o in the first layer (kg/kg)
    272  REAL,INTENT(IN) :: waterice(ngrid,nslope)                          ! Water ice at the surface [kg/m^2]
    273  REAL,INTENT(IN) :: co2ice(ngrid,nslope)                            ! CO2 ice at the surface [kg/m^2]
     239implicit none
     240
     241! Inputs:
     242integer,                                    intent(in) :: ngrid, nslope, nsoil_PEM,timelen   ! Size dimension
     243real,    dimension(ngrid,timelen),          intent(in) :: ps                                 ! Average surface pressure [Pa]
     244real,    dimension(ngrid,nsoil_PEM,nslope), intent(in) :: tsoil_PEM                          ! Mean Soil Temperature [K]
     245real,    dimension(ngrid,nsoil_PEM,nslope), intent(in) :: TI_PEM                             ! Mean Thermal Inertia [USI]
     246real,    dimension(ngrid,nslope),           intent(in) :: tend_h2oglaciers, tend_co2glaciers ! Tendencies on the glaciers ()
     247real,    dimension(ngrid,timelen),          intent(in) :: q_co2, q_h2o                       ! Mass mixing ratio of co2 and h2o in the first layer (kg/kg)
     248real,    dimension(ngrid,nslope),           intent(in) :: waterice                           ! Water ice at the surface [kg/m^2]
     249real,    dimension(ngrid,nslope),           intent(in) :: co2ice                             ! CO2 ice at the surface [kg/m^2]
    274250
    275251! Outputs:
    276  REAL,INTENT(INOUT) :: m_co2_completesoil(ngrid,nsoil_PEM,nslope)  ! Density of co2 adsorbed (kg/m^3)
    277  REAL,INTENT(OUT) :: delta_mreg(ngrid)                              ! Difference density of co2 adsorbed (kg/m^3)
    278  
     252real, dimension(ngrid,nsoil_PEM,nslope), intent(inout) :: m_co2_completesoil ! Density of co2 adsorbed (kg/m^3)
     253real, dimension(ngrid), intent(out) :: delta_mreg ! Difference density of co2 adsorbed (kg/m^3)
     254
    279255! Constants:
    280 
    281  REAL :: alpha = 7.512e-6         ! Zent & Quinn 1995
    282  REAL :: beta =  -1541.5          ! Zent & Quinn 1995
    283  REAL ::  inertie_thresold = 800. ! TI > 800 means cementation
    284  real :: m_theta = 4.27e-7        ! Mass of co2 per m^2 absorbed
     256real :: alpha = 7.512e-6        ! Zent & Quinn 1995
     257real :: beta = -1541.5          ! Zent & Quinn 1995
     258real :: inertie_thresold = 800. ! TI > 800 means cementation
     259real :: m_theta = 4.27e-7       ! Mass of co2 per m^2 absorbed
    285260! real :: as = 18.9e3             ! Specific area, Buhler & Piqueux 2021
    286  real :: as = 9.48e4              ! Same as previous but from zent
    287 ! Local         
    288  real :: A,B                                             ! Used to compute the mean mass above the surface
    289  INTEGER :: ig,islope,iloop,it                           ! For loops
    290  REAL :: dm_co2_regolith_slope(ngrid,nsoil_PEM,nslope)   ! Elementary mass adsorded per mesh per slope
    291  INTEGER :: ispermanent_co2glaciers(ngrid,nslope)        ! Check if the co2 glacier is permanent
    292  INTEGER :: ispermanent_h2oglaciers(ngrid,nslope)        ! Check if the h2o glacier is permanent
    293  REAL :: deltam_reg_complete(ngrid,index_breccia,nslope) ! Difference in the mass per slope and soil layer (kg/m^3)
    294  REAL :: deltam_reg_slope(ngrid,nslope)                  ! Difference in the mass per slope  (kg/m^3)
    295  REAL :: m_h2o_adsorbed(ngrid,nsoil_PEM,nslope)          ! Density of CO2 adsorbed (kg/m^3)
    296  REAL :: theta_h2o_adsorbed(ngrid,nsoil_PEM,nslope)      ! Fraction of the pores occupied by H2O molecules
    297  REAL :: delta_mh2o(ngrid)                               ! Difference density of h2o adsorbed (kg/m^3)
    298 !timelen array are allocated because heavy ...
    299  real,allocatable :: mass_mean(:,:)                      ! Mean mass above the surface
    300  real,allocatable :: zplev_mean(:,:)                     ! Pressure above the surface
    301  real,allocatable :: pco2(:,:)                           ! Partial pressure above the surface
    302  real, allocatable :: pco2_avg(:)                        ! Yearly averaged
     261real :: as = 9.48e4             ! Same as previous but from zent
     262
     263! Local
     264real                                           :: A, B                    ! Used to compute the mean mass above the surface
     265integer                                        :: ig, islope, iloop, it   ! For loops
     266real,    dimension(ngrid,nsoil_PEM,nslope)     :: dm_co2_regolith_slope   ! Elementary mass adsorded per mesh per slope
     267logical, dimension(ngrid,nslope)               :: ispermanent_co2glaciers ! Check if the co2 glacier is permanent
     268logical, dimension(ngrid,nslope)               :: ispermanent_h2oglaciers ! Check if the h2o glacier is permanent
     269real,    dimension(ngrid,index_breccia,nslope) :: deltam_reg_complete     ! Difference in the mass per slope and soil layer (kg/m^3)
     270real,    dimension(ngrid,nslope)               :: deltam_reg_slope        ! Difference in the mass per slope  (kg/m^3)
     271real,    dimension(ngrid,nsoil_PEM,nslope)     :: m_h2o_adsorbed          ! Density of CO2 adsorbed (kg/m^3)
     272real,    dimension(ngrid,nsoil_PEM,nslope)     :: theta_h2o_adsorbed      ! Fraction of the pores occupied by H2O molecules
     273real,    dimension(ngrid)                      :: delta_mh2o              ! Difference density of h2o adsorbed (kg/m^3)
     274! timelen array are allocated because heavy...
     275real,    dimension(:,:), allocatable           :: mass_mean               ! Mean mass above the surface
     276real,    dimension(:,:), allocatable           :: zplev_mean              ! Pressure above the surface
     277real,    dimension(:,:), allocatable           :: pco2                    ! Partial pressure above the surface
     278real,    dimension(:),   allocatable           :: pco2_avg                ! Yearly averaged
    303279
    304280! 0. Some initializations
    305 
    306      allocate(mass_mean(ngrid,timelen))
    307      allocate(zplev_mean(ngrid,timelen))
    308      allocate(pco2(ngrid,timelen))
    309      allocate(pco2_avg(ngrid))
    310 
    311 
    312  
    313       m_h2o_adsorbed(:,:,:) = 0.
    314       A =(1/m_co2 - 1/m_noco2)
    315       B=1/m_noco2
    316 
    317      dm_co2_regolith_slope(:,:,:) = 0
    318      delta_mreg(:) = 0.
    319 
    320 #ifndef CPP_STD
    321 
    322 !0.1 Look at perennial ice
    323   do ig = 1,ngrid
    324     do islope = 1,nslope
    325      if((abs(tend_h2oglaciers(ig,islope)).gt.1e-5).and.(abs(waterice(ig,islope)).gt.0)) then
    326         ispermanent_h2oglaciers(ig,islope) = 1
    327      else
    328         ispermanent_h2oglaciers(ig,islope) = 0
    329      endif
    330 
    331      if((abs(tend_co2glaciers(ig,islope)).gt.1e-5).and.(abs(co2ice(ig,islope)).gt.0)) then
    332         ispermanent_co2glaciers(ig,islope) = 1
    333      else
    334         ispermanent_co2glaciers(ig,islope) = 0
    335      endif
    336     enddo
    337    enddo
    338 
    339 !   0.2  Compute the partial pressure of CO2
    340 !a. the molecular mass into the column
    341      do ig = 1,ngrid
    342        mass_mean(ig,:) = 1./(A*q_co2(ig,:) +B)
    343      enddo
    344 
    345 ! b. pressure level
    346      do it = 1,timelen
    347        do ig = 1,ngrid
    348          zplev_mean(ig,it) = ap(1) + bp(1)*ps(ig,it)
    349        enddo
    350      enddo
    351 
    352 ! c. Vapor pressure
    353      pco2(:,:) = mass_mean(:,:)/m_co2*q_co2(:,:)*zplev_mean(:,:)
    354      pco2_avg(:) = sum(pco2(:,:),2)/timelen
    355 
    356      deallocate(zplev_mean)
    357      deallocate(mass_mean)
    358      deallocate(pco2)
    359 
    360 
    361 ! 1. Compute the fraction of the pores occupied by H2O
    362 
    363  call regolith_h2oadsorption(ngrid,nslope,nsoil_PEM,timelen,tend_h2oglaciers,tend_co2glaciers,waterice,co2ice,ps,q_co2,q_h2o,tsoil_PEM,TI_PEM, &
    364                                    theta_h2o_adsorbed, m_h2o_adsorbed,delta_mh2o)
    365 
    366 ! 2.  we compute the mass of co2 adsorded in each layer of the meshes 
    367 
    368  do ig = 1,ngrid
    369   do islope = 1,nslope
    370     do iloop = 1,index_breccia
    371      if((TI_PEM(ig,iloop,islope).lt.inertie_thresold).and.(ispermanent_h2oglaciers(ig,islope).eq.0).and.(ispermanent_co2glaciers(ig,islope).eq.0)) then
    372      dm_co2_regolith_slope(ig,iloop,islope) = as*rho_regolith*m_theta*(1-theta_h2o_adsorbed(ig,iloop,islope))*alpha*pco2_avg(ig)/ &
    373                                              (alpha*pco2_avg(ig)+sqrt(tsoil_PEM(ig,iloop,islope))*exp(beta/tsoil_PEM(ig,iloop,islope)))
    374      else
    375         if(abs(m_co2_completesoil(ig,iloop,islope)).lt.(1e-10)) then !!! we are at first call
    376           dm_co2_regolith_slope(ig,iloop,islope) = as*rho_regolith*m_theta*(1-theta_h2o_adsorbed(ig,iloop,islope))*alpha*pco2_avg(ig) &
    377                                                  /(alpha*pco2_avg(ig)+sqrt(tsoil_PEM(ig,iloop,islope))*exp(beta/tsoil_PEM(ig,iloop,islope)))
    378         else ! no change: permanent ice stick the atoms of CO2
    379           dm_co2_regolith_slope(ig,iloop,islope) = m_co2_completesoil(ig,iloop,islope)
    380         endif
    381      endif
    382     enddo
    383  enddo
    384 enddo
    385 
    386 ! 3. Check the exchange between the atmosphere and the regolith
    387 
    388   do ig = 1,ngrid
    389    delta_mreg(ig) = 0.
    390    do islope = 1,nslope
    391     deltam_reg_slope(ig,islope) = 0.
    392     do iloop = 1,index_breccia
    393        if((TI_PEM(ig,iloop,islope).lt.inertie_thresold).and.(ispermanent_h2oglaciers(ig,islope).eq.0).and.(ispermanent_co2glaciers(ig,islope).eq.0)) then
    394              if(iloop == 1) then
    395                  deltam_reg_complete(ig,iloop,islope) = (dm_co2_regolith_slope(ig,iloop,islope) - m_co2_completesoil(ig,iloop,islope)) &
    396                                                    *(layer_PEM(iloop))               
    397              else
    398                  deltam_reg_complete(ig,iloop,islope) = (dm_co2_regolith_slope(ig,iloop,islope) - m_co2_completesoil(ig,iloop,islope)) &
    399                                                    *(layer_PEM(iloop) - layer_PEM(iloop-1))
    400              endif
    401        else ! NO EXCHANGE AS ICE BLOCK THE DYNAMIC!
    402              deltam_reg_complete(ig,iloop,islope) = 0.
    403        endif
    404        deltam_reg_slope(ig,islope) = deltam_reg_slope(ig,islope) + deltam_reg_complete(ig,iloop,islope)
    405     enddo
    406    delta_mreg(ig) = delta_mreg(ig) + deltam_reg_slope(ig,islope)*subslope_dist(ig,islope)/cos(pi*def_slope_mean(islope)/180.)
    407    enddo
    408   enddo
    409   m_co2_completesoil(:,:,:) = dm_co2_regolith_slope(:,:,:)
    410 
    411 !=======================================================================
    412       RETURN
    413 #endif
    414       END
    415 
    416 
    417    end module
     281allocate(mass_mean(ngrid,timelen),zplev_mean(ngrid,timelen),pco2(ngrid,timelen),pco2_avg(ngrid))
     282m_h2o_adsorbed = 0.
     283A = 1./m_co2 - 1./m_noco2
     284B = 1./m_noco2
     285dm_co2_regolith_slope = 0.
     286delta_mreg = 0.
     287ispermanent_h2oglaciers = .false.
     288ispermanent_co2glaciers = .false.
     289
     290#ifndef CPP_STD
     291    ! 0.1 Look at perennial ice
     292    do ig = 1,ngrid
     293        do islope = 1,nslope
     294            if (abs(tend_h2oglaciers(ig,islope)) > 1.e-5 .and. abs(waterice(ig,islope)) > 0.) ispermanent_h2oglaciers(ig,islope) = .true.
     295            if (abs(tend_co2glaciers(ig,islope)) > 1.e-5 .and. abs(co2ice(ig,islope)) > 0.) ispermanent_co2glaciers(ig,islope) = .true.
     296        enddo
     297    enddo
     298
     299    ! 0.2  Compute the partial pressure of CO2
     300    ! a. the molecular mass into the column
     301    do ig = 1,ngrid
     302        mass_mean(ig,:) = 1./(A*q_co2(ig,:) + B)
     303    enddo
     304
     305    ! b. pressure level
     306    do it = 1,timelen
     307        do ig = 1,ngrid
     308            zplev_mean(ig,it) = ap(1) + bp(1)*ps(ig,it)
     309        enddo
     310    enddo
     311
     312    ! c. Vapor pressure
     313    pco2 = mass_mean/m_co2*q_co2*zplev_mean
     314    pco2_avg(:) = sum(pco2(:,:),2)/timelen
     315
     316    deallocate(zplev_mean,mass_mean,pco2)
     317
     318    ! 1. Compute the fraction of the pores occupied by H2O
     319    call regolith_h2oadsorption(ngrid,nslope,nsoil_PEM,timelen,tend_h2oglaciers,tend_co2glaciers,waterice,co2ice,ps,q_co2,q_h2o,tsoil_PEM,TI_PEM, &
     320                                theta_h2o_adsorbed, m_h2o_adsorbed,delta_mh2o)
     321
     322    ! 2. we compute the mass of co2 adsorded in each layer of the meshes
     323    do ig = 1,ngrid
     324        do islope = 1,nslope
     325            do iloop = 1,index_breccia
     326                if (TI_PEM(ig,iloop,islope) < inertie_thresold .and. .not. ispermanent_h2oglaciers(ig,islope) .and. .not. ispermanent_co2glaciers(ig,islope)) then
     327                    dm_co2_regolith_slope(ig,iloop,islope) = as*rho_regolith*m_theta*(1. - theta_h2o_adsorbed(ig,iloop,islope))*alpha*pco2_avg(ig)/ &
     328                                                             (alpha*pco2_avg(ig) + sqrt(tsoil_PEM(ig,iloop,islope))*exp(beta/tsoil_PEM(ig,iloop,islope)))
     329                else
     330                    if (abs(m_co2_completesoil(ig,iloop,islope)) < 1.e-10) then !!! we are at first call
     331                        dm_co2_regolith_slope(ig,iloop,islope) = as*rho_regolith*m_theta*(1. - theta_h2o_adsorbed(ig,iloop,islope))*alpha*pco2_avg(ig) &
     332                                                                 /(alpha*pco2_avg(ig)+sqrt(tsoil_PEM(ig,iloop,islope))*exp(beta/tsoil_PEM(ig,iloop,islope)))
     333                    else ! no change: permanent ice stick the atoms of CO2
     334                        dm_co2_regolith_slope(ig,iloop,islope) = m_co2_completesoil(ig,iloop,islope)
     335                    endif
     336                endif
     337            enddo
     338        enddo
     339    enddo
     340
     341    ! 3. Check the exchange between the atmosphere and the regolith
     342    do ig = 1,ngrid
     343        delta_mreg(ig) = 0.
     344        do islope = 1,nslope
     345            deltam_reg_slope(ig,islope) = 0.
     346            do iloop = 1,index_breccia
     347                if (TI_PEM(ig,iloop,islope) < inertie_thresold .and. .not. ispermanent_h2oglaciers(ig,islope) .and. .not. ispermanent_co2glaciers(ig,islope)) then
     348                    if (iloop == 1) then
     349                        deltam_reg_complete(ig,iloop,islope) = (dm_co2_regolith_slope(ig,iloop,islope) - m_co2_completesoil(ig,iloop,islope))*(layer_PEM(iloop))
     350                    else
     351                        deltam_reg_complete(ig,iloop,islope) = (dm_co2_regolith_slope(ig,iloop,islope) - m_co2_completesoil(ig,iloop,islope))*(layer_PEM(iloop) - layer_PEM(iloop - 1))
     352                    endif
     353                else ! NO EXCHANGE AS ICE BLOCK THE DYNAMIC!
     354                    deltam_reg_complete(ig,iloop,islope) = 0.
     355                endif
     356                deltam_reg_slope(ig,islope) = deltam_reg_slope(ig,islope) + deltam_reg_complete(ig,iloop,islope)
     357            enddo
     358            delta_mreg(ig) = delta_mreg(ig) + deltam_reg_slope(ig,islope)*subslope_dist(ig,islope)/cos(pi*def_slope_mean(islope)/180.)
     359        enddo
     360    enddo
     361    m_co2_completesoil = dm_co2_regolith_slope
     362#endif
     363
     364END SUBROUTINE regolith_co2adsorption
     365
     366END MODULE adsorption_mod
     367
  • trunk/LMDZ.COMMON/libf/evolution/changelog.txt

    r3446 r3456  
    436436== 30/09/2024 == JBC
    437437Correction of the launching script for a relaunch situation.
     438
     439== 11/10/2024 == JBC
     440Cleaning of the adsorption module to make the debugging easier.
  • trunk/LMDZ.COMMON/libf/evolution/deftank/multiple_exec.sh

    r3354 r3456  
    11#!/bin/bash
    2 #############################################################
    3 ### Script to exectute mutliple scripts in subdirectories ###
    4 #############################################################
     2############################################################
     3### Script to execute multiple scripts in subdirectories ###
     4############################################################
    55
    66# Name of the file to execute in all subdirectories
  • trunk/LMDZ.COMMON/libf/evolution/pemetat0.F90

    r3430 r3456  
    217217! b. Special case for inertiedat, inertiedat_PEM
    218218        call get_field("inertiedat_PEM",inertiedat_PEM,found)
    219         if (.not.found) then
     219        if (.not. found) then
    220220            do iloop = 1,nsoil_PCM
    221221                inertiedat_PEM(:,iloop) = inertiedat(:,iloop)
     
    259259            write(num,fmt='(i2.2)') islope
    260260            call get_field("tsoil_PEM_slope"//num,tsoil_startPEM(:,:,islope),found)
    261             if (.not.found) then
     261            if (.not. found) then
    262262                write(*,*)'PEM settings: failed loading <tsoil_PEM_slope'//num//'>'
    263263                write(*,*)'will reconstruct the values of Tsoil'
     
    346346
    347347            if (.not. found) deltam_co2_regolith_phys = 0.
    348             if (.not.found2) deltam_h2o_regolith_phys = 0.
     348            if (.not. found2) deltam_h2o_regolith_phys = 0.
    349349
    350350            write(*,*) 'PEMETAT0: CO2 & H2O adsorption done'
Note: See TracChangeset for help on using the changeset viewer.