Changeset 2937

Timestamp:

Apr 13, 2023, 11:04:49 AM (20 months ago)

Author:

llange

Message:

PEM
Thermal Inertia is now only read in the start and not read in the xios
file
Fixing few bug (problem of allocation and index for update_soil and some
variables in the pem)
LL

Location:

trunk/LMDZ.COMMON/libf/evolution

Files:

• evol_h2o_ice_s_mod_slope.F90 (modified) (1 diff)
• pem.F90 (modified) (10 diffs)
• pemetat0.F90 (modified) (2 diffs)
• read_data_GCM.F90 (modified) (7 diffs)
• update_soil.F90 (modified) (1 diff)

trunk/LMDZ.COMMON/libf/evolution/evol_h2o_ice_s_mod_slope.F90

@@ -80 +80 @@
     print *, "Tendencies on ice increasing=", pos_tend
     print *, "This can be due to the absence of water ice in the PCM run!!"
-      call criterion_waterice_stop(cell_area,1.,qsurf(:,:)*0.,STOPPING,ngrid,cell_area)
+    print *,cell_area,qsurf(:,:),ngrid,STOPPING
+      call criterion_waterice_stop(cell_area,1.,qsurf(:,:)*0.,STOPPING,ngrid,qsurf(:,:)*0.)
       do i=1,ngrid
          do islope=1,nslope

trunk/LMDZ.COMMON/libf/evolution/pem.F90

@@ -195 +195 @@
 !!!!!!!!!!!!!!!!!!!!!!!! SLOPE
       REAL ,allocatable :: watercap_slope(:,:)                           ! Physics x Nslope: watercap per slope 
-      REAL ,allocatable :: watercap_slope_saved                          ! Value saved at the previous time step
+      REAL ::              watercap_slope_saved                          ! Value saved at the previous time step
       REAL , dimension(:,:), allocatable :: min_co2_ice_1                ! ngrid field : minimum of co2 ice at each point for the first year [kg/m^2]
       REAL , dimension(:,:), allocatable :: min_co2_ice_2                ! ngrid field : minimum of co2 ice at each point for the second year [kg/m^2]
@@ -223 +223 @@
      REAL, ALLOCATABLE :: inertiesoil(:,:)    !Physic x Depth  Thermal inertia of the mesh for restart [SI]
 
-     REAL, ALLOCATABLE :: TI_GCM(:,:,:)       ! Physic x Depth x Slope Averaged GCM Thermal Inertia per slope  [SI]
-     REAL, ALLOCATABLE :: TI_GCM_start(:,:,:) ! Same but for the start
+     REAL, ALLOCATABLE :: TI_GCM(:,:,:) ! Same but for the start
 
      REAL,ALLOCATABLE  :: ice_depth(:,:)      ! Physic x SLope: Ice table depth [m]
@@ -369 +368 @@
 
       allocate(watercap_slope(ngrid,nslope))
-      allocate(TI_GCM_start(ngrid,nsoilmx,nslope))
+      allocate(TI_GCM(ngrid,nsoilmx,nslope))
       allocate(inertiesoil(ngrid,nsoilmx))
 
@@ -380 +379 @@
               watercap,inertiesoil,nslope,tsurf_slope,           &
               tsoil_slope,co2ice_slope,def_slope,def_slope_mean, &
-              subslope_dist,major_slope,albedo_slope,emiss_slope, TI_GCM_start,     &
+              subslope_dist,major_slope,albedo_slope,emiss_slope, TI_GCM,     &
               qsurf_slope,watercap_slope)
 
-     if(soil_pem) then
-       deallocate(TI_GCM_start) !not used then
-     endif
+
 
 ! Remove unphysical values of surface tracer
@@ -496 +493 @@
      allocate(tsurf_GCM_timeseries(ngrid,nslope,timelen))
      allocate(tsoil_GCM_timeseries(ngrid,nsoilmx,nslope,timelen))
-     allocate(TI_GCM(ngrid,nsoilmx,nslope))
      allocate(q_co2_PEM_phys(ngrid,timelen))
      allocate(q_h2o_PEM_phys(ngrid,timelen))
      allocate(co2_ice_GCM_slope(ngrid,nslope,timelen))
      allocate(watersurf_density_ave(ngrid,nslope))
-     allocate(watersoil_density_timeseries(nslope,nsoilmx,nslope,timelen))
+     allocate(watersoil_density_timeseries(ngrid,nsoilmx,nslope,timelen))
 
      allocate(tsoil_ave_PEM_yr1(ngrid,nsoilmx_PEM,nslope))
@@ -515 +511 @@
 
      call read_data_GCM("data_GCM_Y1.nc",timelen, iim,jjm,ngrid,nslope,vmr_co2_gcm,ps_timeseries,min_co2_ice_1,min_h2o_ice_1,&   
-                       tsurf_ave_yr1,tsoil_ave_yr1, tsurf_GCM_timeseries,tsoil_GCM_timeseries,TI_GCM,q_co2_PEM_phys,q_h2o_PEM_phys,co2_ice_GCM_slope, &      
+                       tsurf_ave_yr1,tsoil_ave_yr1, tsurf_GCM_timeseries,tsoil_GCM_timeseries,q_co2_PEM_phys,q_h2o_PEM_phys,co2_ice_GCM_slope, &      
                        watersurf_density_ave,watersoil_density_timeseries)
 
@@ -524 +520 @@
 
      call read_data_GCM("data_GCM_Y2.nc",timelen,iim,jjm,ngrid,nslope,vmr_co2_gcm,ps_timeseries,min_co2_ice_2,min_h2o_ice_2, &
-                  tsurf_ave,tsoil_ave, tsurf_GCM_timeseries,tsoil_GCM_timeseries,TI_GCM,q_co2_PEM_phys,q_h2o_PEM_phys,co2_ice_GCM_slope, &      
+                  tsurf_ave,tsoil_ave, tsurf_GCM_timeseries,tsoil_GCM_timeseries,q_co2_PEM_phys,q_h2o_PEM_phys,co2_ice_GCM_slope, &      
                   watersurf_density_ave,watersoil_density_timeseries)
 
@@ -835 +831 @@
 ! II.b. Evolution of the ice
       print *, "Evolution of h2o ice"
+     
      call evol_h2o_ice_s_slope(qsurf_slope(:,igcm_h2o_ice,:),tendencies_h2o_ice,iim,jjm,ngrid,cell_area,STOPPING_1_water,nslope)
-
+     
      DO islope=1, nslope
        write(str2(1:2),'(i2.2)') islope
@@ -994 +991 @@
      call computeice_table_equilibrium(ngrid,nslope,nsoilmx_PEM,watercaptag,watersurf_density_ave,watersoil_density_PEM_ave,ice_depth)
        
-      print *, "Update soil propreties"
+     print *, "Update soil propreties"
+
 ! II_d.4 Update the soil thermal properties
-      call update_soil(ngrid,nslope,nsoilmx,nsoilmx_PEM,tendencies_h2o_ice,qsurf_slope(:,igcm_h2o_ice,:),global_ave_press_new, &
+      call update_soil(ngrid,nslope,nsoilmx_PEM,tendencies_h2o_ice,qsurf_slope(:,igcm_h2o_ice,:),global_ave_press_new, &
         ice_depth,TI_PEM)
 
@@ -1165 +1163 @@
      fluxgeo_slope(:,:) = fluxgeo
      call interpolate_TIPEM_TIGCM(ngrid,nslope,nsoilmx_PEM,nsoilmx,TI_PEM,TI_GCM)
-     tsoil_slope(:,:,:) = tsoil_phys_PEM_timeseries(:,:,:,timelen)
-   else
-      TI_GCM(:,:,:)=TI_GCM_start(:,:,:)
+     tsoil_slope(:,:,:) = tsoil_phys_PEM_timeseries(:,:,:,timelen)      
    endif !soil_pem
 

trunk/LMDZ.COMMON/libf/evolution/pemetat0.F90

@@ -90 +90 @@
 
 write(*,*)'Is start PEM?',startpem_file
-
+startpem_file = .true.
 
 !1. Run
@@ -287 +287 @@
   DO islope=1,nslope
    write(num,fmt='(i2.2)') islope
-   call get_field("mco2_reg_ads_slope"//num,m_co2_regolith_phys(:,:,islope),found)
+   call get_field("mh2o_reg_ads_slope"//num,m_co2_regolith_phys(:,:,islope),found2)
     if((.not.found2)) then
        m_h2o_regolith_phys(:,:,:) = 0.

trunk/LMDZ.COMMON/libf/evolution/read_data_GCM.F90

@@ -3 +3 @@
 !
 SUBROUTINE read_data_GCM(fichnom,timelen, iim_input,jjm_input,ngrid,nslope,vmr_co2_gcm_phys,ps_timeseries, & 
-             min_co2_ice,min_h2o_ice,tsurf_ave,tsoil_ave,tsurf_gcm,tsoil_gcm,TI_ave,q_co2,q_h2o,co2_ice_slope, &
+             min_co2_ice,min_h2o_ice,tsurf_ave,tsoil_ave,tsurf_gcm,tsoil_gcm,q_co2,q_h2o,co2_ice_slope, &
              watersurf_density_ave,watersoil_density)
 
@@ -44 +44 @@
   REAL , INTENT(OUT) ::  watersurf_density_ave(ngrid,nslope)             ! Water density at the surface [kg/m^3]
   REAL , INTENT(OUT) ::  watersoil_density(ngrid,nsoilmx,nslope,timelen) ! Water density in the soil layer, time series [kg/m^3]
-  REAL, INTENT(OUT) ::  TI_ave(ngrid,nsoilmx,nslope)                     ! Average Thermal Inertia  of the concatenated file [SI]
 !===============================================================================
 !   Local Variables 
@@ -71 +70 @@
   REAL ::  tsurf_gcm_dyn(iim_input+1,jjm_input+1,nslope,timelen)       ! Surface temperature of the concatenated file, time series [K]
   REAL ::  tsoil_gcm_dyn(iim_input+1,jjm_input+1,nsoilmx,nslope,timelen)! Soil temperature of the concatenated file, time series [K]
-  REAL ::  TI_gcm(iim_input+1,jjm_input+1,nsoilmx,nslope,timelen)      ! Thermal Inertia  of the concatenated file, times series [SI]
-  REAL ::  TI_ave_dyn(iim_input+1,jjm_input+1,nsoilmx,nslope)          ! Average Thermal Inertia  of the concatenated file [SI]
   REAL ::  q_co2_dyn(iim_input+1,jjm_input+1,timelen)                  ! CO2 mass mixing ratio in the first layer [kg/m^3]
   REAL ::  q_h2o_dyn(iim_input+1,jjm_input+1,timelen)                  ! H2O mass mixing ratio in the first layer [kg/m^3]
@@ -156 +153 @@
 
      print *, "Downloading data for tsoil_slope done"
-     print *, "Downloading data for inertiesoil_slope ..."
-
-DO islope=1,nslope
-  write(num,fmt='(i2.2)') islope
-  call get_var4("inertiesoil_slope"//num,TI_gcm(:,:,:,islope,:))
-ENDDO
-
-     print *, "Downloading data for inertiesoil_slope done"
 
      print *, "Downloading data for watersoil_density ..."
@@ -196 +185 @@
     if(soil_pem) then
       call get_var4("tsoil",tsoil_gcm_dyn(:,:,:,1,:))
-      call get_var4("inertiesoil",TI_gcm(:,:,:,1,:))
     endif !soil_pem
   endif !nslope=1
@@ -221 +209 @@
     print *, "Computing average of tsoil"
     tsoil_ave_dyn(:,:,:,:)=SUM(tsoil_gcm_dyn(:,:,:,:,:),5)/timelen
-    print *, "Computing average of TI"
-    TI_ave_dyn(:,:,:,:)=SUM(TI_gcm(:,:,:,:,:),5)/timelen
     print *, "Computing average of watersurf_density"
     watersurf_density_ave(:,:) = SUM(watersurf_density(:,:,:),3)/timelen
@@ -269 +255 @@
        CALL gr_dyn_fi(1,iim_input+1,jjm_input+1,ngrid,min_h2o_ice_dyn(:,:,islope),min_h2o_ice(:,islope))
        if(soil_pem) then
-         CALL gr_dyn_fi(nsoilmx,iim_input+1,jjm_input+1,ngrid,TI_ave_dyn(:,:,:,islope),TI_ave(:,:,islope))
        DO l=1,nsoilmx
          CALL gr_dyn_fi(1,iim_input+1,jjm_input+1,ngrid,tsoil_ave_dyn(:,:,l,islope),tsoil_ave(:,l,islope))

trunk/LMDZ.COMMON/libf/evolution/update_soil.F90

@@ -44 +44 @@
    regolith_inertia(:,islope) = inertiedat_PEM(:,1)
    do ig = 1,ngrid
+      write(*,*) 'ig,islope',ig,islope,tendencies_waterice(ig,islope),waterice(ig,islope)
       if((tendencies_waterice(ig,islope).lt.-1e-5).and.(waterice(ig,islope).eq.0)) then
               regolith_inertia(ig,islope) = inertie_averaged