Index: palm/trunk/SOURCE/land_surface_model_mod.f90
===================================================================
--- palm/trunk/SOURCE/land_surface_model_mod.f90	(revision 2237)
+++ palm/trunk/SOURCE/land_surface_model_mod.f90	(revision 2242)
@@ -25,4 +25,8 @@
 ! -----------------
 ! $Id$
+! Changed soil configuration to 8 layers. The number of soil layers is now
+! freely adjustable via the NAMELIST.
+! 
+! 2237 2017-05-31 10:34:53Z suehring
 ! Bugfix in write restart data
 ! 
@@ -225,11 +229,9 @@
        REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  var_2D !< 2D prognostic variable
     END TYPE surf_type_lsm
+
 !
 !-- LSM model constants
-    INTEGER(iwp), PARAMETER :: nzb_soil = 0, & !< bottom of the soil model (to be switched)
-                               nzt_soil = 3, & !< top of the soil model (to be switched)
-                               nzs = 4         !< number of soil layers (fixed for now)
-
-    REAL(wp), PARAMETER ::                     &
+
+    REAL(wp), PARAMETER  ::                    &
               b_ch               = 6.04_wp,    & ! Clapp & Hornberger exponent
               lambda_h_dry       = 0.19_wp,    & ! heat conductivity for dry soil   
@@ -242,9 +244,19 @@
 
 
+    REAL(wp), DIMENSION(0:7), PARAMETER  :: zs_default =                       & ! default soil layer configuration
+                                            (/ 0.005_wp, 0.02_wp, 0.04_wp,     &
+                                               0.07_wp, 0.15_wp, 0.28_wp,      &
+                                               1.00_wp, 2.89_wp /)
+
+
 !
 !-- LSM variables
-    INTEGER(iwp) :: veg_type  = 2, & !< NAMELIST veg_type_2d
-                    soil_type = 3    !< NAMELIST soil_type_2d
-
+    INTEGER(iwp) :: nzb_soil = 0,  & !< bottom of the soil model (Earth's surface)
+                    nzt_soil = 0,  & !< top of the soil model
+                    nzs = 0,       & !< number of soil layers
+                    veg_type  = 2, & !< default NAMELIST veg_type_2d
+                    soil_type = 3    !< default NAMELIST soil_type_2d
+   
+        
     LOGICAL :: conserve_water_content = .TRUE.,  & !< open or closed bottom surface for the soil model
                force_radiation_call_l = .FALSE., & !< flag to force calling of radiation routine
@@ -288,20 +300,17 @@
                 z0q_eb = 9999999.9_wp                     !< NAMELIST z0q (lsm_par)
 
-    REAL(wp), DIMENSION(nzb_soil:nzt_soil) :: &
-              ddz_soil_stag,                  & !< 1/dz_soil_stag
-              dz_soil_stag,                   & !< soil grid spacing (center-center)
-              root_extr = 0.0_wp,             & !< root extraction
-              root_fraction = (/9999999.9_wp, 9999999.9_wp,    &
-                                9999999.9_wp, 9999999.9_wp /), & !< distribution of root surface area to the individual soil layers
-              zs = (/0.07_wp, 0.28_wp, 1.00_wp,  2.89_wp/),    & !< soil layer depths (m)
-              soil_moisture = 0.0_wp          !< soil moisture content (m3/m3)
-
-    REAL(wp), DIMENSION(nzb_soil:nzt_soil+1) ::   &
-              soil_temperature = (/290.0_wp, 287.0_wp, 285.0_wp,  283.0_wp,    & !< soil temperature (K)
-                                   283.0_wp /),                                &                                   
-              ddz_soil,                                                        & !< 1/dz_soil
-              dz_soil                                                            !< soil grid spacing (edge-edge)
-
-
+    REAL(wp), DIMENSION(:), ALLOCATABLE  :: ddz_soil,         & !< 1/dz_soil
+                                            ddz_soil_stag,    & !< 1/dz_soil_stag
+                                            dz_soil,          & !< soil grid spacing (edge-edge)
+                                            dz_soil_stag,     & !< soil grid spacing (center-center)
+                                            root_extr           !< root extraction
+
+
+                                            
+    REAL(wp), DIMENSION(0:20)  ::  root_fraction = 9999999.9_wp,               & !< distribution of root surface area to the individual soil layers
+                                   soil_moisture = 0.0_wp,                     & !< NAMELIST soil moisture content (m3/m3)
+                                   soil_temperature = 300.0_wp,                & !< NAMELIST soil temperature (K) +1
+                                   zs = 9999999.9_wp                             !< soil layer depths
+                                   
 #if defined( __nopointer )
     TYPE(surf_type_lsm), TARGET  ::  t_soil_h,    & !< Soil temperature (K), horizontal surface elements
@@ -533,27 +542,49 @@
 
 !
-!-- Root distribution (sum = 1)  level 1, level 2, level 3, level 4,
-    REAL(wp), DIMENSION(0:3,1:20), PARAMETER :: root_distribution = RESHAPE( (/ &
-                                 0.24_wp, 0.41_wp, 0.31_wp, 0.04_wp,            & !  1
-                                 0.35_wp, 0.38_wp, 0.23_wp, 0.04_wp,            & !  2
-                                 0.26_wp, 0.39_wp, 0.29_wp, 0.06_wp,            & !  3
-                                 0.26_wp, 0.38_wp, 0.29_wp, 0.07_wp,            & !  4
-                                 0.24_wp, 0.38_wp, 0.31_wp, 0.07_wp,            & !  5
-                                 0.25_wp, 0.34_wp, 0.27_wp, 0.14_wp,            & !  6
-                                 0.27_wp, 0.27_wp, 0.27_wp, 0.09_wp,            & !  7
-                                 1.00_wp, 0.00_wp, 0.00_wp, 0.00_wp,            & !  8
-                                 0.47_wp, 0.45_wp, 0.08_wp, 0.00_wp,            & !  9
-                                 0.24_wp, 0.41_wp, 0.31_wp, 0.04_wp,            & ! 10
-                                 0.17_wp, 0.31_wp, 0.33_wp, 0.19_wp,            & ! 11
-                                 0.00_wp, 0.00_wp, 0.00_wp, 0.00_wp,            & ! 12
-                                 0.25_wp, 0.34_wp, 0.27_wp, 0.11_wp,            & ! 13
-                                 0.00_wp, 0.00_wp, 0.00_wp, 0.00_wp,            & ! 14
-                                 0.00_wp, 0.00_wp, 0.00_wp, 0.00_wp,            & ! 15
-                                 0.23_wp, 0.36_wp, 0.30_wp, 0.11_wp,            & ! 16 
-                                 0.23_wp, 0.36_wp, 0.30_wp, 0.11_wp,            & ! 17 
-                                 0.19_wp, 0.35_wp, 0.36_wp, 0.10_wp,            & ! 18
-                                 0.19_wp, 0.35_wp, 0.36_wp, 0.10_wp,            & ! 19
-                                 0.00_wp, 0.00_wp, 0.00_wp, 0.00_wp             & ! 20
-                                 /), (/ 4, 20 /) )
+!-- Root distribution (sum = 1)  level 1-8
+!--                               
+    REAL(wp), DIMENSION(0:7,1:20), PARAMETER :: root_distribution = RESHAPE( (/&
+                                   0.035_wp, 0.069_wp, 0.069_wp, 0.108_wp,     &
+                                   0.195_wp, 0.214_wp, 0.284_wp, 0.026_wp,     & ! 1
+                                   0.050_wp, 0.100_wp, 0.100_wp, 0.136_wp,     &
+                                   0.181_wp, 0.192_wp, 0.215_wp, 0.026_wp,     & ! 2
+                                   0.038_wp, 0.075_wp, 0.075_wp, 0.111_wp,     &
+                                   0.185_wp, 0.203_wp, 0.273_wp, 0.040_wp,     & ! 3
+                                   0.038_wp, 0.075_wp, 0.075_wp, 0.110_wp,     &
+                                   0.180_wp, 0.199_wp, 0.277_wp, 0.046_wp,     & ! 4
+                                   0.035_wp, 0.069_wp, 0.069_wp, 0.105_wp,     &
+                                   0.180_wp, 0.201_wp, 0.295_wp, 0.046_wp,     & ! 5
+                                   0.035_wp, 0.072_wp, 0.072_wp, 0.105_wp,     &
+                                   0.161_wp, 0.180_wp, 0.282_wp, 0.093_wp,     & ! 6
+                                   0.040_wp, 0.077_wp, 0.077_wp, 0.112_wp,     &
+                                   0.176_wp, 0.192_wp, 0.266_wp, 0.060_wp,     & ! 7
+                                   0.142_wp, 0.286_wp, 0.286_wp, 0.286_wp,     &
+                                   0.000_wp, 0.000_wp, 0.000_wp, 0.000_wp,     & ! 8
+                                   0.068_wp, 0.134_wp, 0.134_wp, 0.177_wp,     &
+                                   0.214_wp, 0.203_wp, 0.070_wp, 0.000_wp,     & ! 9
+                                   0.035_wp, 0.068_wp, 0.068_wp, 0.108_wp,     &
+                                   0.195_wp, 0.215_wp, 0.285_wp, 0.026_wp,     & ! 10
+                                   0.025_wp, 0.048_wp, 0.048_wp, 0.078_wp,     &
+                                   0.147_wp, 0.175_wp, 0.353_wp, 0.126_wp,     & ! 11
+                                   0.000_wp, 0.000_wp, 0.000_wp, 0.000_wp,     &
+                                   0.000_wp, 0.000_wp, 0.000_wp, 0.000_wp,     & ! 12
+                                   0.036_wp, 0.072_wp, 0.072_wp, 0.103_wp,     &
+                                   0.163_wp, 0.180_wp, 0.273_wp, 0.074_wp,     & ! 13
+                                   0.000_wp, 0.000_wp, 0.000_wp, 0.000_wp,     &
+                                   0.000_wp, 0.000_wp, 0.000_wp, 0.000_wp,     & ! 14
+                                   0.000_wp, 0.000_wp, 0.000_wp, 0.000_wp,     &
+                                   0.000_wp, 0.000_wp, 0.000_wp, 0.000_wp,     & ! 15
+                                   0.032_wp, 0.066_wp, 0.066_wp, 0.100_wp,     &
+                                   0.172_wp, 0.192_wp, 0.299_wp, 0.073_wp,     & ! 16
+                                   0.032_wp, 0.066_wp, 0.066_wp, 0.100_wp,     &
+                                   0.172_wp, 0.192_wp, 0.299_wp, 0.073_wp,     & ! 17
+                                   0.028_wp, 0.055_wp, 0.055_wp, 0.087_wp,     &
+                                   0.166_wp, 0.195_wp, 0.348_wp, 0.066_wp,     & ! 18
+                                   0.028_wp, 0.055_wp, 0.055_wp, 0.087_wp,     &
+                                   0.166_wp, 0.195_wp, 0.348_wp, 0.066_wp,     & ! 19
+                                   0.000_wp, 0.000_wp, 0.000_wp, 0.000_wp,     &
+                                   0.000_wp, 0.000_wp, 0.000_wp, 0.000_wp      & ! 20
+                                 /), (/ 8, 20 /) )
+
 
 !
@@ -930,11 +961,4 @@
 
     IF ( veg_type == 0 )  THEN
-       IF ( SUM( root_fraction ) /= 1.0_wp )  THEN
-          message_string = 'veg_type = 0 (user_defined)'//                  &
-                           'requires setting of root_fraction(0:3)'//       &
-                           '/= 9999999.9 and SUM(root_fraction) = 1'
-          CALL message( 'check_parameters', 'PA0401', 1, 2, 0, 6, 0 )
-       ENDIF
- 
        IF ( min_canopy_resistance == 9999999.9_wp )  THEN
           message_string = 'veg_type = 0 (user defined)'//                  &
@@ -1079,5 +1103,39 @@
     ENDIF
        
-       
+
+
+!
+!-- Determine number of soil layers to be used and check whether an appropriate
+!-- root fraction is prescribed
+    nzb_soil = 0
+    nzt_soil = 0
+    IF ( ALL( zs == 9999999.9_wp ) )  THEN
+       nzt_soil = 7
+       zs(nzb_soil:nzt_soil) = zs_default
+    ELSE
+       DO k = 0, 19
+          IF ( zs(k) /= 9999999.9_wp )  THEN
+             nzt_soil = nzt_soil + 1
+             IF ( root_fraction(k) == 9999999.9_wp )  THEN
+                message_string = 'manual setting of zs '//                     &
+                                 'requires adequate setting of root_fraction'//&
+                                 '/= 9999999.9 ' //                            &
+                                 'and SUM(root_fraction) = 1'
+                CALL message( 'check_parameters', 'PA0452', 1, 2, 0, 6, 0 )
+             ENDIF
+          ENDIF
+       ENDDO    
+    ENDIF
+
+    IF ( veg_type == 0 )  THEN
+       IF ( SUM( root_fraction(nzb_soil:nzt_soil) ) /= 1.0_wp )  THEN
+          message_string = 'veg_type = 0 (user_defined)'//                  &
+                           'requires setting of root_fraction'//            &
+                           '/= 9999999.9 and SUM(root_fraction) = 1'
+          CALL message( 'check_parameters', 'PA0401', 1, 2, 0, 6, 0 )
+       ENDIF
+    ENDIF
+    
+  
  END SUBROUTINE lsm_check_parameters
  
@@ -1225,5 +1283,9 @@
        ELSE
           pt1 = pt(k,j,i)
-          qv1 = q(k,j,i)
+          IF ( humidity )  THEN
+             qv1 = q(k,j,i)
+          ELSE
+             qv1 = 0.0_wp
+          ENDIF
        ENDIF
 !
@@ -1789,5 +1851,5 @@
 
 !
-!--    Set inital values for prognostic quantities
+!--    Set initial values for prognostic quantities
 !--    Horizontal surfaces
        tt_surface_h_m%var_1d = 0.0_wp
@@ -1842,4 +1904,5 @@
        ENDDO
 
+        
 !
 !--    Allocate 3D soil model arrays
@@ -2154,12 +2217,13 @@
           z0h_factor = z0h_eb / ( z0_eb + 1.0E-20_wp )
 
-          IF ( ANY( root_fraction == 9999999.9_wp ) )  THEN
+          IF ( ALL( root_fraction == 9999999.9_wp ) )  THEN
              DO  m = 1, surf_lsm_h%ns
                 i   = surf_lsm_h%i(m)            
                 j   = surf_lsm_h%j(m)
 
-                DO  k = nzb_soil, nzt_soil
+
+                DO  k = nzb_soil, nzt_soil 
                    surf_lsm_h%root_fr(k,m) = root_distribution(k,veg_type)
-                   root_fraction(k)    = root_distribution(k,veg_type)
+                   root_fraction(k)        = root_distribution(k,veg_type)
                 ENDDO
              ENDDO
@@ -2302,5 +2366,15 @@
 
        INTEGER(iwp) ::  l !< index indicating facing of surface array 
-
+   
+!
+!--    Allocate global 1D arrays
+       ALLOCATE ( ddz_soil(nzb_soil:nzt_soil+1) )
+       ALLOCATE ( ddz_soil_stag(nzb_soil:nzt_soil) )
+       ALLOCATE ( dz_soil(nzb_soil:nzt_soil+1) )
+       ALLOCATE ( dz_soil_stag(nzb_soil:nzt_soil) )
+       ALLOCATE ( root_extr(nzb_soil:nzt_soil) )
+  
+       root_extr = 0.0_wp 
+       
 !
 !--    Allocate surface and soil temperature / humidity. Please note, 
@@ -2494,5 +2568,5 @@
                                   conserve_water_content,                      &
                                   f_shortwave_incoming, field_capacity,        & 
-                                  aero_resist_kray, hydraulic_conductivity,                      &
+                                  aero_resist_kray, hydraulic_conductivity,    &
                                   lambda_surface_stable,                       &
                                   lambda_surface_unstable, leaf_area_index,    &
