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

Timestamp:

Jan 4, 2024, 11:37:33 AM (11 months ago)

Author:

llange

Message:

Mars PEM
Fixing a small bug: the subroutine compute_icetable was always called, even if tthe option 'icetable_equilibrium' was set to false in the run_PEM.def. It is now fixed by adding a flag before the call.
LL

File:

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

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

@@ -16 +16 @@
 use comsoil_h,                  only: volcapa, inertiedat
 use adsorption_mod,             only: regolith_adsorption, adsorption_pem
-use ice_table_mod,              only: computeice_table_equilibrium
+use ice_table_mod,              only: computeice_table_equilibrium, icetable_equilibrium
 use constants_marspem_mod,      only: alpha_clap_h2o, beta_clap_h2o, TI_breccia, TI_bedrock
 use soil_thermalproperties_mod, only: update_soil_thermalproperties
@@ -102 +102 @@
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 
-!0. Check if the start_PEM exist.
+!0.1 Check if the start_PEM exist.
 
 inquire(file = filename,exist = startpem_file)
@@ -108 +108 @@
 write(*,*)'Is start PEM?',startpem_file
 
+!0.2 Set to default values
+ice_table = -1.  ! by default, no ice table
+ice_table_thickness = -1.
 !1. Run
 if (startpem_file) then
@@ -272 +275 @@
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 !3. Ice Table
-        call get_field("ice_table",ice_table,found)
-        if (.not. found) then
-            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_avg,watersoil_avg, TI_PEM(:,1,:),ice_table,ice_table_thickness)
-            call update_soil_thermalproperties(ngrid,nslope,nsoil_PEM,tend_h2o_ice,h2o_ice,global_avg_pressure,ice_table,ice_table_thickness,TI_PEM)
-            do islope = 1,nslope
-                call soil_pem_ini(ngrid,nsoil_PEM,TI_PEM(:,:,islope),tsurf_avg_yr2(:,islope),tsoil_PEM(:,:,islope))
-            enddo
+        if(icetable_equilibrium) then
+            call get_field("ice_table",ice_table,found)
+            if (.not. found) then
+                write(*,*)'PEM settings: failed loading <ice_table>'
+                write(*,*)'will reconstruct the values of the ice table given the current state'
+                call computeice_table_equilibrium(ngrid,nslope,nsoil_PEM,watercaptag,watersurf_avg,watersoil_avg, TI_PEM(:,1,:),ice_table,ice_table_thickness)
+                call update_soil_thermalproperties(ngrid,nslope,nsoil_PEM,tend_h2o_ice,h2o_ice,global_avg_pressure,ice_table,ice_table_thickness,TI_PEM)
+                do islope = 1,nslope
+                    call soil_pem_ini(ngrid,nsoil_PEM,TI_PEM(:,:,islope),tsurf_avg_yr2(:,islope),tsoil_PEM(:,:,islope))
+                enddo
+            endif
+            write(*,*) 'PEMETAT0: ICE TABLE done'
         endif
-
-        write(*,*) 'PEMETAT0: ICE TABLE done'
 
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
@@ -442 +444 @@
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 !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_avg,watersoil_avg,TI_PEM(:,1,:),ice_table,ice_table_thickness)
-        call update_soil_thermalproperties(ngrid,nslope,nsoil_PEM,tend_h2o_ice,h2o_ice,global_avg_pressure,ice_table,ice_table_thickness,TI_PEM)
-        do islope = 1,nslope
-            call soil_pem_ini(ngrid,nsoil_PEM,TI_PEM(:,:,islope),tsurf_avg_yr2(:,islope),tsoil_PEM(:,:,islope))
-        enddo
-        write(*,*) 'PEMETAT0: Ice table done'
-
+        if(icetable_equilibrium) then
+            call computeice_table_equilibrium(ngrid,nslope,nsoil_PEM,watercaptag,watersurf_avg,watersoil_avg,TI_PEM(:,1,:),ice_table,ice_table_thickness)
+            call update_soil_thermalproperties(ngrid,nslope,nsoil_PEM,tend_h2o_ice,h2o_ice,global_avg_pressure,ice_table,ice_table_thickness,TI_PEM)
+            do islope = 1,nslope
+                call soil_pem_ini(ngrid,nsoil_PEM,TI_PEM(:,:,islope),tsurf_avg_yr2(:,islope),tsoil_PEM(:,:,islope))
+            enddo
+            write(*,*) 'PEMETAT0: Ice table done'
+        endif
+        
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 !d) Regolith adsorbed