Changeset 3500

Timestamp:

Nov 8, 2024, 10:56:18 AM (4 weeks ago)

Author:

afalco

Message:

Pluto PCM: initialize soil layers when running without startphy_file.
AF

Location:

trunk/LMDZ.PLUTO/libf/phypluto

Files:

• physiq_mod.F90 (modified) (4 diffs)
• soil_settings.F (modified) (13 diffs)

trunk/LMDZ.PLUTO/libf/phypluto/physiq_mod.F90

@@ -527 +527 @@
       real :: annual_skin ! anuual skin depth (m)
 
+      ! when no startfi file is asked for init
+      real alpha,lay1 ! coefficients for building layers
+      integer iloop
 
       ! flags to trigger extra sanity checks
@@ -599 +602 @@
 #ifndef MESOSCALE
          if (.not.startphy_file) then
+           ! starting without startfi.nc and with callsoil
+           ! is not yet possible as soildepth default is not defined
+           if (callsoil) then
+              ! default mlayer distribution, following a power law:
+              !  mlayer(k)=lay1*alpha**(k-1/2)
+              lay1=2.e-4
+              alpha=2
+              do iloop=0,nsoilmx-1
+              mlayer(iloop)=lay1*(alpha**(iloop-0.5))
+          enddo
+              lay1=sqrt(mlayer(0)*mlayer(1))
+              alpha=mlayer(1)/mlayer(0)
+              do iloop=1,nsoilmx
+                 layer(iloop)=lay1*(alpha**(iloop-1))
+              enddo
+           endif
            ! additionnal "academic" initialization of physics
            if (is_master) write(*,*) "Physiq: initializing tsurf(:) to pt(:,1) !!"
@@ -1011 +1030 @@
       ! enddo
 
-      if (callsoil) then
+      if (firstcall.and.callsoil) then
          ! AF24 Originally in soil.F, but inertiedat is modified by surfprop
          ! Additional checks: is the vertical discretization sufficient
@@ -1019 +1038 @@
             inertia_min=minval(inertiedat(ig,:))
             inertia_max=maxval(inertiedat(ig,:))
+            print*, "inertia min max" , inertia_min, inertia_max
             ! diurnal and annual skin depth
             diurnal_skin=(inertia_min/volcapa)*sqrt(daysec/pi)

trunk/LMDZ.PLUTO/libf/phypluto/soil_settings.F

@@ -39 +39 @@
 !  ---------
 !  inputs:
-      integer,intent(in) :: nid ! Input Netcdf file ID 
+      integer,intent(in) :: nid ! Input Netcdf file ID
       integer,intent(in) :: ngrid       ! # of horizontal grid points
       integer,intent(in) :: nsoil       ! # of soil layers
@@ -55 +55 @@
       integer ndims     ! # of dimensions of read <inertiedat> data
       integer ig,iloop  ! loop counters
-      
+
       integer edges(3),corner(3) ! to read a specific time
 
@@ -89 +89 @@
         endif
 
-        ! allocate oldmlayer  
+        ! allocate oldmlayer
         allocate(oldmlayer(dimlen),stat=ierr)
         if (ierr.ne.0) then
@@ -115 +115 @@
      &              'and thermal inertia will be required!'
         endif
-        
+
 ! 1.2 Find out the # of dimensions <inertiedat> was defined as using
       ndims=inquire_field_ndims("inertiedat")
@@ -135 +135 @@
           oldmlayer(iloop)=sqrt(887.75/3.14)*((2.**(iloop-0.5))-1.)
         enddo
-      else ! Look for depth 
+      else ! Look for depth
         ! read <depth> coordinate
         if (interpol) then !put values in oldmlayer
@@ -180 +180 @@
 ! 2. Volumetric heat capacity (note: it is declared in comsoil_h)
 ! ---------------------------
-! "volcapa" is (so far) 0D and written in "controle" table of startfi file 
+! "volcapa" is (so far) 0D and written in "controle" table of startfi file
 ! volcapa is read or set when "controle" is read (see tabfi.F)
 ! Just in case, we check here that it is not zero. If it is, we
@@ -210 +210 @@
          call abort
        endif
-       
+
        write(*,*)' => Building soil thermal inertia (using reference sur
      &face thermal inertia)'
@@ -242 +242 @@
        endif ! of if (interpol)
       endif ! of if (ndims.eq.1)
-      
+
 ! 4. Read soil temperatures
 ! -------------------------
@@ -275 +275 @@
            call abort
          endif
-       endif ! of if (interpol) 
+       endif ! of if (interpol)
       endif! of if (.not.ok)
 
@@ -304 +304 @@
         deallocate(oldval)
         deallocate(newval)
-        interpol=.false. ! no need for interpolation any more      
+        interpol=.false. ! no need for interpolation any more
       endif !of if (olddepthdef)
 
@@ -325 +325 @@
           inertiedat(ig,:)=newval(:)
         enddo
-        
+
       ! soil temperature
         ! vertical coordinate
@@ -339 +339 @@
           tsoil(ig,:)=newval(:)
         enddo
-        
+
         !cleanup
         deallocate(oldgrid)
@@ -347 +347 @@
         deallocate(oldtsoil)
       endif ! of if (interpol)
-      
+
 ! 6. Report min and max values of soil temperatures and thermal inertias
 ! ----------------------------------------------------------------------