Changeset 3123 for trunk/LMDZ.COMMON/libf/evolution/pemetat0.F90

Timestamp:

Nov 11, 2023, 5:34:53 PM (13 months ago)

Author:

llange

Message:

MARS PEM
1) Adapting the Mars PEM soil grid to the one of the PCM. The first layers in the PEM follow those from the PCM (PYM grid), and then, for layers at depth, we use the classical power low grid.
2) Correction when reading the soil temperature, water density at the surface, and initialization of the ice table.
LL

File:

• trunk/LMDZ.COMMON/libf/evolution/pemetat0.F90 (modified) (3 diffs)

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

@@ -262 +262 @@
       write(*,*)'PEM settings: failed loading <ice_table>'
       write(*,*)'will reconstruct the values of the ice table given the current state'
+      ice_table(:,:) = -1  ! by default, no ice table  
+      ice_table_thickness(:,:) = -1  ! by default, no ice table
      call computeice_table_equilibrium(ngrid,nslope,nsoil_PEM,watercaptag,watersurf_ave,watersoil_ave, TI_PEM(:,1,:),ice_table,ice_table_thickness)
      call update_soil_thermalproperties(ngrid,nslope,nsoil_PEM,tend_h2oglaciers,waterice,global_ave_pressure,ice_table,ice_table_thickness,TI_PEM)
@@ -427 +429 @@
       call soil_pem_ini(ngrid,nsoil_PEM,TI_PEM(:,:,islope),tsurf_ave_yr2(:,islope),tsoil_PEM(:,:,islope))
       call soil_pem_compute(ngrid,nsoil_PEM,.true.,TI_PEM(:,:,islope),timestep,tsurf_ave_yr2(:,islope),tsoil_PEM(:,:,islope))
-
-
-    do it = 1,timelen
+  
+! First raw initialization 
+      do it = 1,timelen
         do isoil = nsoil_GCM+1,nsoil_PEM
-      call soil_pem_ini(ngrid,nsoil_PEM,TI_PEM(:,:,islope),tsurf_ave_yr2(:,islope),tsoil_inst(:,:,islope,it))
+          tsoil_inst(:,isoil,islope,it) = tsoil_PEM(:,isoil,islope)
         enddo
-     enddo
-
+      enddo
+
+      do it = 1,timelen
         do isoil = nsoil_GCM+1,nsoil_PEM
-          do ig = 1,ngrid
-            watersoil_ave(ig,isoil,islope) = exp(beta_clap_h2o/tsoil_PEM(ig,isoil,islope) + alpha_clap_h2o)/tsoil_PEM(ig,isoil,islope)*mmol(igcm_h2o_vap)/(mugaz*r)
-          enddo
+          call soil_pem_ini(ngrid,nsoil_PEM,TI_PEM(:,:,islope),tsurf_ave_yr2(:,islope),tsoil_inst(:,:,islope,it))
         enddo
+       enddo
+
+       do isoil = nsoil_GCM+1,nsoil_PEM
+         do ig = 1,ngrid
+           watersoil_ave(ig,isoil,islope) = exp(beta_clap_h2o/tsoil_PEM(ig,isoil,islope) + alpha_clap_h2o)/tsoil_PEM(ig,isoil,islope)*mmol(igcm_h2o_vap)/(mugaz*r)
+         enddo
+       enddo
 enddo !islope
 write(*,*) 'PEMETAT0: TSOIL done'
@@ -445 +453 @@
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 !c) Ice table
+       ice_table(:,:) = -1.  ! by default, no ice table
+       ice_table_thickness(:,:) = -1.
        call computeice_table_equilibrium(ngrid,nslope,nsoil_PEM,watercaptag,watersurf_ave,watersoil_ave,TI_PEM(:,1,:),ice_table,ice_table_thickness)
        call update_soil_thermalproperties(ngrid,nslope,nsoil_PEM,tend_h2oglaciers,waterice,global_ave_pressure,ice_table,ice_table_thickness,TI_PEM)