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

Timestamp:

Dec 8, 2025, 11:27:43 AM (7 days ago)

Author:

jbclement

Message:

PEM:

• Removing completely the ice metamorphism computed by a threshold at the end of the PCM (which was commented).
• Addition of a module "metamorphism" to compute the PCM frost at the PEM beginning and give it back to the PCM at the PEM end. The frost is considered as the ice given by the PCM "startfi.nc" which is above the yearly minimum. Thereby, metamorphism is performed through this operation.
• Ice reservoirs representation in the PEM is modernized.

JBC

File:

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

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

@@ -25 +25 @@
 use glaciers_mod,               only: rho_co2ice, rho_h2oice
 use comcstfi_h,                 only: r, mugaz, pi
-use surfdat_h,                  only: watercaptag, perennial_co2ice
+use surfdat_h,                  only: watercaptag, perennial_co2ice, qsurf
+use metamorphism,               only: frost4PCM, iPCM_h2ofrost, iPCM_co2frost
 
 implicit none
@@ -80 +81 @@
 !!!        4. Mass of CO2 & H2O adsorbed
 !!!
-!!! /!\ This order must be respected !
+!!! /!\ This order must be respected!
 !!! Author: LL
 !!!
@@ -102 +103 @@
 
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-    ! h2o ice
-    h2o_ice = 0.
+    ! H2O ice
+    h2o_ice(:,:) = 0.
     call get_field("h2o_ice",h2o_ice,found)
     if (.not. found) then
@@ -110 +111 @@
         write(*,*)'with default value ''ini_huge_h2oice'''
         do ig = 1,ngrid
-            if (watercaptag(ig)) h2o_ice(ig,:) = ini_huge_h2oice
+            if (watercaptag(ig)) h2o_ice(ig,:) = ini_huge_h2oice + qsurf(ig,iPCM_h2ofrost,:) - frost4PCM(ig,:)%h2o
         enddo
     endif
 
-    ! co2 ice
-    co2_ice = 0.
+    ! CO2 ice
+    co2_ice(:,:) = 0.
     call get_field("co2_ice",co2_ice,found)
     if (.not. found) then
         write(*,*)'Pemetat0: failed loading <co2_ice>'
         write(*,*)'will reconstruct the values from perennial_co2ice'
-        co2_ice = perennial_co2ice
+        co2_ice(:,:) = perennial_co2ice(:,:) + qsurf(:,iPCM_co2frost,:) - frost4PCM(:,:)%co2
     endif
 
@@ -208 +209 @@
                         delta = depth_breccia
                         TI_PEM(ig,index_breccia + 1,islope) = sqrt((layer_PEM(index_breccia + 1) - layer_PEM(index_breccia))/ &
-                                                                  (((delta - layer_PEM(index_breccia))/(TI_PEM(ig,index_breccia,islope)**2)) + &
-                                                                  ((layer_PEM(index_breccia + 1) - delta)/(TI_breccia**2))))
-                        do iloop=index_breccia + 2,index_bedrock
+                                                                   (((delta - layer_PEM(index_breccia))/(TI_PEM(ig,index_breccia,islope)**2)) + &
+                                                                   ((layer_PEM(index_breccia + 1) - delta)/(TI_breccia**2))))
+                        do iloop = index_breccia + 2,index_bedrock
                             TI_PEM(ig,iloop,islope) = TI_breccia
                         enddo
@@ -246 +247 @@
             do ig = 1,ngrid
                 if (inertiedat_PEM(ig,index_breccia) < TI_breccia) then
-                    inertiedat_PEM(ig,index_breccia + 1) = sqrt((layer_PEM(index_breccia+1)-layer_PEM(index_breccia))/ &
-                                                              (((delta-layer_PEM(index_breccia))/(inertiedat(ig,index_breccia)**2)) + &
-                                                              ((layer_PEM(index_breccia+1)-delta)/(TI_breccia**2))))
+                    inertiedat_PEM(ig,index_breccia + 1) = sqrt((layer_PEM(index_breccia + 1) - layer_PEM(index_breccia))/ &
+                                                                (((delta - layer_PEM(index_breccia))/(inertiedat(ig,index_breccia)**2)) + &
+                                                                ((layer_PEM(index_breccia + 1) - delta)/(TI_breccia**2))))
 
                     do iloop = index_breccia + 2,index_bedrock
@@ -254 +255 @@
                     enddo
                 else
-                    do iloop=index_breccia + 1,index_bedrock
+                    do iloop = index_breccia + 1,index_bedrock
                         inertiedat_PEM(ig,iloop) = inertiedat_PEM(ig,nsoil_PCM)
                     enddo
@@ -263 +264 @@
             delta = depth_bedrock
             do ig = 1,ngrid
-                inertiedat_PEM(ig,index_bedrock + 1) = sqrt((layer_PEM(index_bedrock+1)-layer_PEM(index_bedrock))/ &
-                                                           (((delta-layer_PEM(index_bedrock))/(inertiedat_PEM(ig,index_bedrock)**2))+ &
-                                                           ((layer_PEM(index_bedrock+1)-delta)/(TI_bedrock**2))))
-            enddo
-
-            do iloop = index_bedrock + 2, nsoil_PEM
+                inertiedat_PEM(ig,index_bedrock + 1) = sqrt((layer_PEM(index_bedrock + 1) - layer_PEM(index_bedrock))/ &
+                                                            (((delta - layer_PEM(index_bedrock))/(inertiedat_PEM(ig,index_bedrock)**2))+ &
+                                                            ((layer_PEM(index_bedrock + 1) - delta)/(TI_bedrock**2))))
+            enddo
+
+            do iloop = index_bedrock + 2,nsoil_PEM
                 do ig = 1,ngrid
                     inertiedat_PEM(ig,iloop) = TI_bedrock
@@ -365 +366 @@
 
     ! h2o ice
-    h2o_ice = 0.
+    h2o_ice(:,:) = 0.
     write(*,*)'So ''h2o_ice'' is initialized with default value ''ini_huge_h2oice'' where ''watercaptag'' is true.'
     do ig = 1,ngrid
-        if (watercaptag(ig)) h2o_ice(ig,:) = ini_huge_h2oice
+        if (watercaptag(ig)) h2o_ice(ig,:) = ini_huge_h2oice + qsurf(ig,iPCM_h2ofrost,:) - frost4PCM(ig,:)%h2o
     enddo
 
     ! co2 ice
     write(*,*)'So ''co2_ice'' is initialized with ''perennial_co2ice'' found in the PCM.'
-    co2_ice = perennial_co2ice
+    co2_ice(:,:) = perennial_co2ice(:,:) + qsurf(:,iPCM_co2frost,:) - frost4PCM(:,:)%co2
 
     ! Layerings