Changeset 3250 for trunk/LMDZ.MARS/libf/phymars

Timestamp:

Mar 4, 2024, 11:03:35 AM (10 months ago)

Author:

llange

Message:

Mars PCM

• Fixing bug in initestphys and waterice_tifeedback: when changing the soil

thermal properties to account for the high thermal inertia of subsurface ice,
the wrong layer was considered as 'mixed' (it was the layer after the correct
one which was modified).

LL

Location:

trunk/LMDZ.MARS/libf/phymars

Files:

• dyn1d/init_testphys1d_mod.F90 (modified) (2 diffs)
• waterice_tifeedback_mod.F90 (modified) (2 diffs)

trunk/LMDZ.MARS/libf/phymars/dyn1d/init_testphys1d_mod.F90

@@ -678 +678 @@
             do isoil = 1,nsoil
                 if ((ice_depth >= layer(isoil)) .and. (ice_depth < layer(isoil + 1))) then
-                    iref = isoil + 1
+                    iref = isoil 
                     exit
                 endif
@@ -685 +685 @@
             inertiedat(1,:iref - 1) = inertiedat(1,1)
             ! We compute the transition in layer(iref)
-            inertiedat(1,iref) = sqrt((layer(iref) - layer(iref - 1))/(((ice_depth - layer(iref - 1))/inertiedat(1,1)**2) + ((layer(iref) - ice_depth)/inertieice**2)))
+            inertiedat(1,iref) = sqrt((layer(iref+1) - layer(iref))/(((ice_depth - layer(iref))/inertiedat(1,1)**2) + ((layer(iref+1) - ice_depth)/inertieice**2)))
             ! Finally, we compute the underlying ice:
             inertiedat(1,iref + 1:) = inertieice

trunk/LMDZ.MARS/libf/phymars/waterice_tifeedback_mod.F90

@@ -107 +107 @@
               DO ik=1,nsoil-1
                  IF ((icedepth.ge.layer(ik)).and.icedepth.lt.layer(ik+1)) THEN
-                    iref=ik+1
+                    iref=ik
                  EXIT
                  ENDIF
@@ -116 +116 @@
               ENDDO
 !             Transition (based on the equations of thermal conduction):
-              newtherm_i(ig,iref,islope)=sqrt( (layer(iref)-layer(iref-1))/(((icedepth-layer(iref-1))/inert_h2o_ice**2) + &
-              ((layer(iref)-icedepth)/inertiedat(ig,ik)**2) ) )
+              newtherm_i(ig,iref,islope)=sqrt( (layer(iref+1)-layer(iref))/(((icedepth-layer(iref))/inert_h2o_ice**2) + &
+              ((layer(iref+1)-icedepth)/inertiedat(ig,ik)**2) ) )
 !              Underlying regolith:
               DO ik=iref+1,nsoil