Changeset 2885


Ignore:
Timestamp:
Jan 31, 2023, 9:52:30 PM (23 months ago)
Author:
romain.vande
Message:

Mars PCM:
Move a endif misplaced
RV

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

Legend:

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

    r2863 r2885  
    2222  real, save, allocatable :: h2o_adsorbded_phys(:,:,:)  ! h2o that is in the regolith [kg/m^2]
    2323  LOGICAL soil_pem  ! True by default, to run with the subsurface physic. Read in pem.def
     24  real,save,allocatable,dimension(:) :: water_reservoir      ! Large reserve of water   [kg/m^2]
    2425
    2526contains
     
    4445    allocate(co2_adsorbded_phys(ngrid,nsoilmx_PEM,nslope))
    4546    allocate(h2o_adsorbded_phys(ngrid,nsoilmx_PEM,nslope))
     47    allocate(water_reservoir(ngrid))
    4648  end subroutine ini_comsoil_h_PEM
    4749
     
    6466    if (allocated(co2_adsorbded_phys)) deallocate(co2_adsorbded_phys)
    6567    if (allocated(h2o_adsorbded_phys)) deallocate(h2o_adsorbded_phys)
     68    if (allocated(water_reservoir)) deallocate(water_reservoir)
    6669  end subroutine end_comsoil_h_PEM
    6770
  • trunk/LMDZ.COMMON/libf/evolution/criterion_ice_stop_mod_slope.F90

    r2835 r2885  
    22! $Id $
    33!
    4 SUBROUTINE criterion_ice_stop_slope(cell_area,ini_surf,qsurf,STOPPING,ngrid,initial_h2o_ice,global_ave_press_GCM,global_ave_press_new,nslope)
     4SUBROUTINE criterion_ice_stop_slope(cell_area,ini_surf,qsurf,STOPPING,STOPPING_ps,ngrid,initial_h2o_ice,global_ave_press_GCM,global_ave_press_new,nslope)
    55
    66  USE temps_mod_evol, ONLY: alpha_criterion
     
    2828
    2929!   OUTPUT
    30   LOGICAL, intent(out) :: STOPPING              ! Logical : is the criterion reached?
     30  LOGICAL, intent(out) :: STOPPING,STOPPING_ps              ! Logical : is the criterion reached?
    3131
    3232!   local:
     
    3939!   initialisation to false
    4040    STOPPING=.FALSE.
     41    STOPPING_ps=.FALSE.
    4142
    4243!   computation of the actual surface
     
    6768  if(global_ave_press_new.LT.global_ave_press_GCM*(0.9) .OR. &
    6869     global_ave_press_new.GT.global_ave_press_GCM*(1.1)) then
    69     STOPPING=.TRUE.
     70    STOPPING_ps=.TRUE.
    7071    print *, "Reason of stopping : The global pressure reach the threshold:"
    7172    print *, "Current global pressure=", global_ave_press_new
  • trunk/LMDZ.COMMON/libf/evolution/pem.F90

    r2863 r2885  
    3737                           albedodat, zmea, zstd, zsig, zgam, zthe, &
    3838                           hmons, summit, base,albedo_h2o_frost, &
    39                            frost_albedo_threshold,emissiv
     39                           frost_albedo_threshold,emissiv,watercaptag
    4040      use dimradmars_mod, only: totcloudfrac, albedo
    4141      use dust_param_mod, only: tauscaling
     
    9090                              TI_PEM,inertiedat_PEM,alph_PEM, beta_PEM,         & ! soil thermal inertia         
    9191                              tsoil_PEM, mlayer_PEM,layer_PEM,                  & !number of subsurface layers, soil mid layer depths
    92                               fluxgeo, co2_adsorbded_phys, h2o_adsorbded_phys     ! geothermal flux, mass of co2 & h2o in the regolith
     92                              fluxgeo, co2_adsorbded_phys, h2o_adsorbded_phys,  & ! geothermal flux, mass of co2 & h2o in the regolith
     93                              water_reservoir
    9394      use adsorption_mod, only : regolith_adsorption
    9495
     
    187188      LOGICAL :: STOPPING_co2     ! Logical : is the criterion (% of change in the surface of sublimating water ice) reached?
    188189      LOGICAL :: STOPPING_1_co2   ! Logical : is there still water ice to sublimate?
     190      LOGICAL :: STOPPING_pressure
     191      INTEGER :: criterion_stop            !Which criterion is reached : 1=surf h20, 2=surf_co2, 3=ps, 4=orb.param
     192
    189193
    190194      REAL, dimension(:,:,:),allocatable  :: q_co2_GCM ! Lon x Lat x Time : mass mixing ratio of co2 in the first layer [kg/kg]
     
    275279    INTEGER :: year_iter_max             ! maximum number of iterations before stopping
    276280    REAL :: timestep                     ! timestep [s]
     281    REAL WC_sum
     282
    277283#ifdef CPP_STD
    278284!     INTEGER :: nsplit_phys=1
     
    418424       enddo
    419425     enddo
     426
     427     call surfini(ngrid,qsurf)
    420428
    421429#else
     
    793801      watersurf_density_phys_ave,watersoil_density_phys_PEM_ave, &
    794802      co2_adsorbded_phys,delta_co2_adsorbded,h2o_adsorbded_phys,delta_h2o_adsorbded)
     803
     804      do ig=1,ngrid
     805         do islope=1,nslope
     806           qsurf_slope(ig,igcm_h2o_ice,islope)=qsurf_slope(ig,igcm_h2o_ice,islope)+watercap_slope(ig,islope)+water_reservoir(ig)
     807         enddo
     808      enddo
    795809
    796810    if(soil_pem) then
     
    11391153      call criterion_ice_stop_water_slope(cell_area,ini_surf_h2o,qsurf_slope(:,igcm_h2o_ice,:),STOPPING_water,ngrid,initial_h2o_ice_slope)
    11401154
    1141       call criterion_ice_stop_slope(cell_area,ini_surf_co2,co2ice_slope,STOPPING_co2,ngrid,initial_co2_ice_sublim_slope,global_ave_press_GCM,global_ave_press_new,nslope)
     1155      call criterion_ice_stop_slope(cell_area,ini_surf_co2,co2ice_slope,STOPPING_co2,STOPPING_pressure,ngrid,initial_co2_ice_sublim_slope,global_ave_press_GCM,global_ave_press_new,nslope)
    11421156
    11431157      year_iter=year_iter+dt_pem
     
    11461160      if (STOPPING_water) then
    11471161        print *, "STOPPING because surface of water ice sublimating is too low, see message above", STOPPING_water
     1162        criterion_stop=1
     1163      endif
     1164      if (STOPPING_1_water) then
     1165        print *, "STOPPING because tendencies on water ice=0, see message above", STOPPING_1_water
     1166        criterion_stop=1
     1167      endif
     1168      if (STOPPING_co2) then
     1169        print *, "STOPPING because surface of co2 ice sublimating is too low, see message above", STOPPING_co2
     1170        criterion_stop=2
     1171      endif
     1172      if (STOPPING_1_co2) then
     1173        print *, "STOPPING because tendencies on co2 ice=0, see message above", STOPPING_1_co2
     1174        criterion_stop=2
     1175      endif
     1176      if (STOPPING_pressure) then
     1177        print *, "STOPPING because surface global pressure changed too much, see message above", STOPPING_pressure
     1178        criterion_stop=3
    11481179      endif
    11491180      if (year_iter.ge.year_iter_max) then
    11501181        print *, "STOPPING because maximum number of iterations reached"
     1182        criterion_stop=4
    11511183      endif
    1152       if (STOPPING_1_water) then
    1153         print *, "STOPPING because tendencies on water ice=0, see message above", STOPPING_1_water
    1154       endif
    1155       if (STOPPING_co2) then
    1156         print *, "STOPPING because surface of co2 ice sublimating is too low or global pressure changed too much, see message above", STOPPING_co2
    1157       endif
    1158       if (STOPPING_1_co2) then
    1159         print *, "STOPPING because tendencies on co2 ice=0, see message above", STOPPING_1_co2
    1160       endif
    1161 
    1162 
    1163 
    1164       if (STOPPING_water .or. STOPPING_1_water .or. STOPPING_co2 .or. STOPPING_1_co2)  then
     1184
     1185
     1186
     1187      if (STOPPING_water .or. STOPPING_1_water .or. STOPPING_co2 .or. STOPPING_1_co2 .or. STOPPING_pressure)  then
    11651188        exit
    11661189      else
     
    11971220        qsurf(ig,igcm_co2_ice)=co2ice(ig)
    11981221#endif
    1199       ENDDO ! of DO ig=1,ngrid
    1200 ! H2o ice
    1201       DO ig = 1,ngrid
    1202         qsurf(ig,igcm_h2o_ice) = 0.
    1203         DO islope = 1,nslope
    1204           qsurf(ig,igcm_h2o_ice) = qsurf(ig,igcm_h2o_ice) + qsurf_slope(ig,igcm_h2o_ice,islope) &
    1205                                    * subslope_dist(ig,islope) /          &
    1206                                   cos(pi*def_slope_mean(islope)/180.)
    1207         ENDDO
    12081222      ENDDO ! of DO ig=1,ngrid
    12091223
     
    13191333       enddo
    13201334     enddo
     1335
     1336     DO ig=1,ngrid
     1337       if(watercaptag(ig)) then
     1338         print *, "INNN"
     1339         WC_sum=0.
     1340         DO islope=1,nslope
     1341           if(qsurf_slope(ig,igcm_h2o_ice,islope).GT. (watercap_slope(ig,islope)+water_reservoir(ig))) then
     1342             qsurf_slope(ig,igcm_h2o_ice,islope)=qsurf_slope(ig,igcm_h2o_ice,islope)-(watercap_slope(ig,islope)+water_reservoir(ig))
     1343           else
     1344             watercap_slope(ig,islope)=watercap_slope(ig,islope)+qsurf_slope(ig,igcm_h2o_ice,islope)-(watercap_slope(ig,islope)+water_reservoir(ig)/nslope)
     1345             qsurf_slope(ig,igcm_h2o_ice,islope)=0.
     1346           endif
     1347           WC_sum=WC_sum+watercap_slope(ig,islope)*subslope_dist(ig,islope)/cos(pi*def_slope_mean(islope)/180.)
     1348           watercap_slope(ig,islope)=0.
     1349         enddo 
     1350           water_reservoir(ig)=water_reservoir(ig)+WC_sum
     1351       endif
     1352     enddo
     1353
     1354! H2o ice
     1355      DO ig = 1,ngrid
     1356        qsurf(ig,igcm_h2o_ice) = 0.
     1357        DO islope = 1,nslope
     1358          qsurf(ig,igcm_h2o_ice) = qsurf(ig,igcm_h2o_ice) + qsurf_slope(ig,igcm_h2o_ice,islope) &
     1359                                   * subslope_dist(ig,islope) /          &
     1360                                  cos(pi*def_slope_mean(islope)/180.)
     1361        ENDDO
     1362      ENDDO ! of DO ig=1,ngrid
     1363
     1364     DO ig=1,ngrid
     1365!       DO islope=1,nslope
     1366         if(qsurf(ig,igcm_h2o_ice).GT.500) then
     1367            watercaptag(ig)=.true.
     1368       DO islope=1,nslope
     1369            qsurf_slope(ig,igcm_h2o_ice,islope)=qsurf_slope(ig,igcm_h2o_ice,islope)-250
     1370            water_reservoir(ig)=water_reservoir(ig)+250
     1371       ENDDO           
     1372         endif
     1373!       enddo
     1374     enddo
     1375
     1376     DO ig=1,ngrid
     1377!       DO islope=1,nslope
     1378         if(water_reservoir(ig).LT. 10) then
     1379            watercaptag(ig)=.false.
     1380            qsurf(ig,igcm_h2o_ice)=qsurf(ig,igcm_h2o_ice)+water_reservoir(ig)
     1381       DO islope=1,nslope
     1382            qsurf_slope(ig,igcm_h2o_ice,islope)=qsurf_slope(ig,igcm_h2o_ice,islope)+water_reservoir(ig)
     1383       ENDDO
     1384         endif
     1385!       enddo
     1386     enddo
     1387
     1388      DO ig = 1,ngrid
     1389        watercap(ig) = 0.
     1390        DO islope = 1,nslope
     1391          watercap(ig) = watercap(ig) + watercap_slope(ig,islope) &
     1392                                   * subslope_dist(ig,islope) /          &
     1393                                  cos(pi*def_slope_mean(islope)/180.)
     1394        ENDDO
     1395      ENDDO ! of DO ig=1,ngrid
     1396
    13211397     
    13221398!------------------------
     
    13861462
    13871463         call pemdem1("restartfi_PEM.nc",year_iter,nsoilmx_PEM,ngrid,nslope , &
    1388                       tsoil_PEM, TI_PEM, ice_depth,co2_adsorbded_phys,h2o_adsorbded_phys)
     1464                      tsoil_PEM, TI_PEM, ice_depth,co2_adsorbded_phys,h2o_adsorbded_phys,water_reservoir)
     1465
     1466     call info_run_PEM(year_iter, criterion_stop)
    13891467
    13901468      print *, "restartfi_PEM.nc has been written"
  • trunk/LMDZ.COMMON/libf/evolution/pemetat0.F90

    r2863 r2885  
    44   
    55   use iostart_PEM, only:  open_startphy, close_startphy, get_field, get_var
    6    use comsoil_h_PEM, only: soil_pem,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, water_reservoir
    77   use comsoil_h, only:  volcapa,inertiedat
    88   use adsorption_mod, only : regolith_adsorption
    99   USE temps_mod_evol, ONLY: year_PEM
     10#ifndef CPP_STD   
     11   use surfdat_h, only: watercaptag
     12#endif
    1013
    1114
  • trunk/LMDZ.COMMON/libf/evolution/pemredem.F90

    r2863 r2885  
    7373subroutine pemdem1(filename,year_iter,nsoil_PEM,ngrid,nslope, &
    7474                    tsoil_slope_PEM,inertiesoil_slope_PEM,ice_table, &
    75                     m_co2_regolith,m_h2o_regolith)
     75                    m_co2_regolith,m_h2o_regolith,water_reservoir)
    7676  ! write time-dependent variable to restart file
    7777  use iostart_PEM, only : open_restartphy, close_restartphy, &
     
    9696  real,intent(in) :: m_co2_regolith(ngrid,nsoil_PEM,nslope)
    9797  real,intent(in) :: m_h2o_regolith(ngrid,nsoil_PEM,nslope)
     98  real,intent(IN) :: water_reservoir(ngrid)
    9899  integer :: islope
    99100  CHARACTER*2 :: num 
     
    115116  ! set time counter in file
    116117  call put_var("Time","Year of simulation",year_tot)
     118
     119  call put_field("water_reservoir","water_reservoir", &
     120                 water_reservoir,year_tot)
    117121
    118122    if(soil_pem) then
  • trunk/LMDZ.COMMON/libf/evolution/read_data_GCM.F90

    r2855 r2885  
    6969  REAL, ALLOCATABLE ::  co2_ice_s(:,:,:)                               ! co2 ice, mesh averaged, of the concatenated file [kg/m^2]
    7070  REAL, ALLOCATABLE ::  h2o_ice_s_slope(:,:,:,:)                       ! h2o ice per slope of the concatenated file [kg/m^2]
     71  REAL, ALLOCATABLE :: watercap_slope(:,:,:,:)
    7172
    7273!-----------------------------------------------------------------------
     
    8081      allocate(co2_ice_s(iim+1,jjm+1,timelen))
    8182      allocate(h2o_ice_s_slope(iim+1,jjm+1,nslope,timelen))
     83      allocate(watercap_slope(iim_input+1,jjm_input+1,nslope,timelen))
    8284      allocate(h2o_ice_s(iim+1,jjm+1,timelen))
    8385
     
    133135
    134136     print *, "Downloading data for h2o_ice_s_slope done"
     137     print *, "Downloading data for watercap_slope ..."
     138
     139DO islope=1,nslope
     140  write(num,fmt='(i2.2)') islope
     141  call get_var3("watercap_slope"//num,watercap_slope(:,:,islope,:))
     142ENDDO
     143
     144     print *, "Downloading data for watercap_slope done"
    135145     print *, "Downloading data for tsurf_slope ..."
    136146
     
    198208
    199209  print *, "Computing the min of h2o_ice_slope"
    200   min_h2o_ice_slope(:,:,:)=minval(h2o_ice_s_slope,4)
     210  min_h2o_ice_slope(:,:,:)=minval(h2o_ice_s_slope+watercap_slope,4)
    201211  print *, "Computing the min of co2_ice_slope"
    202212  min_co2_ice_slope(:,:,:)=minval(co2_ice_slope,4)
     
    225235            min_co2_ice_slope(i,j,islope)  = 0.
    226236          endif
    227           if (min_h2o_ice_slope(i,j,islope).LT.0) then
    228             min_h2o_ice_slope(i,j,islope)  = 0.
    229           endif
     237!          if (min_h2o_ice_slope(i,j,islope).LT.0) then
     238!            min_h2o_ice_slope(i,j,islope)  = 0.
     239!          endif
    230240       ENDDO
    231241    ENDDO
Note: See TracChangeset for help on using the changeset viewer.