Index: palm/trunk/SOURCE/check_parameters.f90
===================================================================
--- palm/trunk/SOURCE/check_parameters.f90	(revision 1750)
+++ palm/trunk/SOURCE/check_parameters.f90	(revision 1757)
@@ -4284,10 +4284,10 @@
 !
 !-- Check nudging and large scale forcing from external file
-    IF ( nudging .AND. ( .NOT. large_scale_forcing ) )  THEN
-       message_string = 'Nudging requires large_scale_forcing = .T.. &'//      &
-                        'Surface fluxes and geostrophic wind should be &'//    &
-                        'prescribed in file LSF_DATA'
-       CALL message( 'check_parameters', 'PA0374', 1, 2, 0, 6, 0 )
-    ENDIF
+!!!    IF ( nudging .AND. ( .NOT. large_scale_forcing ) )  THEN
+!!!       message_string = 'Nudging requires large_scale_forcing = .T.. &'//      &
+!!!                        'Surface fluxes and geostrophic wind should be &'//    &
+!!!                        'prescribed in file LSF_DATA'
+!!!       CALL message( 'check_parameters', 'PA0374', 1, 2, 0, 6, 0 )
+!!    ENDIF
 
     IF ( large_scale_forcing .AND. ( bc_lr /= 'cyclic'  .OR.                   &
Index: palm/trunk/SOURCE/land_surface_model.f90
===================================================================
--- palm/trunk/SOURCE/land_surface_model.f90	(revision 1750)
+++ palm/trunk/SOURCE/land_surface_model.f90	(revision 1757)
@@ -19,5 +19,9 @@
 ! Current revisions:
 ! -----------------
-! 
+! Bugfix: set tm_soil_m to zero after allocation. Added parameter 
+! unscheduled_radiation_calls to control calls of the radiation model based on
+! the skin temperature change during one time step (preliminary version). Set 
+! qsws_soil_eb to zero at model start (previously set to qsws_eb). Removed MAX
+! function as it cannot be vectorized.
 ! 
 ! Former revisions:
@@ -132,5 +136,6 @@
     USE radiation_model_mod,                                                   &
         ONLY:  force_radiation_call, radiation_scheme, rad_net, rad_sw_in,     &
-               rad_lw_out, rad_lw_out_change_0, sigma_sb
+               rad_lw_out, rad_lw_out_change_0, sigma_sb,                      &
+               unscheduled_radiation_calls
         
 #if defined ( __rrtmg )
@@ -320,6 +325,6 @@
 
     REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::                                 &
-              tt_soil_m, & !< t_soil storage array 
-              tm_soil_m, & !< m_soil storage array 
+              tt_soil_m, & !< t_soil storage array
+              tm_soil_m, & !< m_soil storage array
               root_fr      !< root fraction (sum=1) 
 
@@ -668,4 +673,5 @@
        tt_surface_m = 0.0_wp
        tt_soil_m    = 0.0_wp
+       tm_soil_m    = 0.0_wp
        tm_liq_eb_m  = 0.0_wp
        c_liq        = 0.0_wp
@@ -681,5 +687,5 @@
 
        qsws_liq_eb  = 0.0_wp
-       qsws_soil_eb = qsws_eb
+       qsws_soil_eb = 0.0_wp
        qsws_veg_eb  = 0.0_wp
 
@@ -1235,7 +1241,14 @@
 
 !
-!--          In case of fast changes in the skin temperature, it is required to
-!--          update the radiative fluxes in order to keep the solution stable
-             IF ( ABS( t_surface_p(j,i) - t_surface(j,i) ) > 0.2_wp )  THEN
+!--          In case of fast changes in the skin temperature, it is possible to
+!--          update the radiative fluxes independently from the prescribed
+!--          radiation call frequency. This effectively prevents oscillations,
+!--          especially when setting skip_time_do_radiation /= 0. The threshold
+!--          value of 0.2 used here is just a first guess. This method should be
+!--          revised in the future as tests have shown that the threshold is 
+!--          often reached, when no oscillations would occur (causes immense
+!--          computing time for the radiation code).
+             IF ( ABS( t_surface_p(j,i) - t_surface(j,i) ) > 0.2_wp .AND.      &
+                  unscheduled_radiation_calls )  THEN
                 force_radiation_call_l = .TRUE.
              ENDIF
@@ -1388,7 +1401,7 @@
 !
 !--    Make a logical OR for all processes. Force radiation call if at
-!--    least one processor
-       IF ( intermediate_timestep_count == intermediate_timestep_count_max-1 ) &
-       THEN
+!--    least one processor reached the threshold change in skin temperature
+       IF ( unscheduled_radiation_calls .AND. intermediate_timestep_count      &
+            == intermediate_timestep_count_max-1 )  THEN
 #if defined( __parallel )
           IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
@@ -1628,5 +1641,10 @@
 !
 !--             Account for dry soils (find a better solution here!)
-                m_soil_p(:,j,i) = MAX(m_soil_p(:,j,i),0.0_wp)
+                DO  k = nzb_soil, nzt_soil
+                   IF ( m_soil_p(k,j,i) < 0.0_wp )  m_soil_p(k,j,i) = 0.0_wp
+                ENDDO
+
+
+
 
 !
Index: palm/trunk/SOURCE/nudging.f90
===================================================================
--- palm/trunk/SOURCE/nudging.f90	(revision 1750)
+++ palm/trunk/SOURCE/nudging.f90	(revision 1757)
@@ -19,5 +19,6 @@
 ! Current revisions:
 ! ------------------
-! 
+! Bugfix: allow for using higher vertical resolution in nudging file than grid
+! spacing in the LES model
 ! 
 ! Former revisions:
@@ -208,5 +209,5 @@
 
           DO  k = nzb, nzt+1
-             IF ( highheight < zu(k) )  THEN
+             DO WHILE ( highheight < zu(k) )
                 lowheight  = highheight
                 lowtnudge  = hightnudge
@@ -229,5 +230,5 @@
                    CALL message( 'nudging', 'PA0364', 1, 2, 0, 6, 0 )
                 ENDIF
-             ENDIF
+             ENDDO
 
 !
Index: palm/trunk/SOURCE/package_parin.f90
===================================================================
--- palm/trunk/SOURCE/package_parin.f90	(revision 1750)
+++ palm/trunk/SOURCE/package_parin.f90	(revision 1757)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! -----------------
-! 
+! Added parameter unscheduled_radiation_calls
 ! 
 ! Former revisions:
@@ -161,5 +161,5 @@
               emissivity, lambda, lw_radiation, net_radiation, radiation,      &
               radiation_scheme, skip_time_do_radiation, sw_radiation,          &
-              time_utc_init
+              time_utc_init, unscheduled_radiation_calls
                 
 
@@ -242,5 +242,6 @@
                                   lambda, lw_radiation, net_radiation,         &
                                   radiation_scheme, skip_time_do_radiation,    &
-                                  sw_radiation, time_utc_init 
+                                  sw_radiation, time_utc_init,                 &
+                                  unscheduled_radiation_calls
 
     NAMELIST /spectra_par/        averaging_interval_sp, comp_spectra_level,   &
Index: palm/trunk/SOURCE/prognostic_equations.f90
===================================================================
--- palm/trunk/SOURCE/prognostic_equations.f90	(revision 1750)
+++ palm/trunk/SOURCE/prognostic_equations.f90	(revision 1757)
Index: palm/trunk/SOURCE/radiation_model.f90
===================================================================
--- palm/trunk/SOURCE/radiation_model.f90	(revision 1750)
+++ palm/trunk/SOURCE/radiation_model.f90	(revision 1757)
@@ -19,5 +19,6 @@
 ! Current revisions:
 ! -----------------
-! 
+! Added parameter unscheduled_radiation_calls. Bugfix: interpolation of sounding
+! profiles for pressure and temperature above the LES domain.
 ! 
 ! Former revisions:
@@ -164,10 +165,11 @@
 
 
-    LOGICAL ::  constant_albedo = .FALSE.,       & !< flag parameter indicating whether the albedo may change depending on zenith
-                force_radiation_call = .FALSE.,  & !< flag parameter for unscheduled radiation calls
-                lw_radiation = .TRUE.,           & !< flag parameter indicing whether longwave radiation shall be calculated
-                radiation = .FALSE.,             & !< flag parameter indicating whether the radiation model is used
-                sun_up    = .TRUE.,              & !< flag parameter indicating whether the sun is up or down
-                sw_radiation = .TRUE.              !< flag parameter indicing whether shortwave radiation shall be calculated
+    LOGICAL ::  unscheduled_radiation_calls = .TRUE., & !< flag parameter indicating whether additional calls of the radiation code are allowed
+                constant_albedo = .FALSE.,            & !< flag parameter indicating whether the albedo may change depending on zenith
+                force_radiation_call = .FALSE.,       & !< flag parameter for unscheduled radiation calls
+                lw_radiation = .TRUE.,                & !< flag parameter indicating whether longwave radiation shall be calculated
+                radiation = .FALSE.,                  & !< flag parameter indicating whether the radiation model is used
+                sun_up    = .TRUE.,                   & !< flag parameter indicating whether the sun is up or down
+                sw_radiation = .TRUE.                   !< flag parameter indicing whether shortwave radiation shall be calculated
 
 
@@ -372,5 +374,5 @@
            rad_sw_out_av, rad_sw_cs_hr, rad_sw_cs_hr_av, rad_sw_hr,            &
            rad_sw_hr_av, sigma_sb, skip_time_do_radiation, sw_radiation,       &
-           time_radiation, time_utc_init
+           time_radiation, time_utc_init, unscheduled_radiation_calls
 
 
@@ -1229,9 +1231,6 @@
        q_snd = 0.0_wp
 
-       hyp_snd(nzt+2:nzt_rad) = hyp_snd_tmp(nz_snd_start:nz_snd_end)
-       t_snd(nzt+2:nzt_rad)   = t_snd_tmp(nz_snd_start:nz_snd_end)
-
-       DEALLOCATE ( hyp_snd_tmp )
-       DEALLOCATE ( t_snd_tmp )
+       hyp_snd(nzt+2:nzt_rad) = hyp_snd_tmp(nz_snd_start+1:nz_snd_end)
+       t_snd(nzt+2:nzt_rad)   = t_snd_tmp(nz_snd_start+1:nz_snd_end)
 
        nc_stat = NF90_CLOSE( id )
Index: palm/trunk/SOURCE/read_3d_binary.f90
===================================================================
--- palm/trunk/SOURCE/read_3d_binary.f90	(revision 1750)
+++ palm/trunk/SOURCE/read_3d_binary.f90	(revision 1757)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! -----------------
-! 
+! Bugfix in allocation of radiative heating rates
 ! 
 ! Former revisions:
@@ -898,5 +898,5 @@
 
                 CASE ( 'rad_lw_cs_hr_av' )
-                   IF ( .NOT. ALLOCATED( rad_lw_out_av ) )  THEN
+                   IF ( .NOT. ALLOCATED( rad_lw_cs_hr_av ) )  THEN
                       ALLOCATE( rad_lw_cs_hr_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
                    ENDIF
@@ -962,5 +962,5 @@
 
                 CASE ( 'rad_sw_cs_hr_av' )
-                   IF ( .NOT. ALLOCATED( rad_sw_out_av ) )  THEN
+                   IF ( .NOT. ALLOCATED( rad_sw_cs_hr_av ) )  THEN
                       ALLOCATE( rad_sw_cs_hr_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
                    ENDIF
Index: palm/trunk/SOURCE/surface_layer_fluxes.f90
===================================================================
--- palm/trunk/SOURCE/surface_layer_fluxes.f90	(revision 1750)
+++ palm/trunk/SOURCE/surface_layer_fluxes.f90	(revision 1757)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! ------------------
-! 
+! Minor fixes.
 !
 ! Former revisions:
@@ -358,8 +358,7 @@
 
 !
-!--       Invert array and switch sign, then calculate Obukhov length and bulk
-!--       Richardson number
-          z0h_min = MINVAL(z0h)
-          z0_min  = MINVAL(z0)
+!--       Save roughness lengths to temporary variables
+          z0h_min = z0h(nys,nxl)
+          z0_min  = z0(nys,nxl)
           
 !
@@ -641,6 +640,4 @@
              DO  j = nys, nyn
 
-                k   = nzb_s_inner(j,i)
-
 !
 !--             If the bulk Richardson number is outside the range of the lookup
@@ -1055,9 +1052,10 @@
              !$acc loop independent
              DO  j = nysg, nyng
-                e(nzb_s_inner(j,i)+1,j,i) = ( us(j,i) / 0.1_wp )**2
+                k = nzb_s_inner(j,i)
+                e(k+1,j,i) = ( us(j,i) / 0.1_wp )**2
 !
 !--             As a test: cm = 0.4
-!               e(nzb_s_inner(j,i)+1,j,i) = ( us(j,i) / 0.4_wp )**2
-                e(nzb_s_inner(j,i),j,i)   = e(nzb_s_inner(j,i)+1,j,i)
+!               e(k+1,j,i) = ( us(j,i) / 0.4_wp )**2
+                e(k,j,i)   = e(k+1,j,i)
              ENDDO
           ENDDO
