Index: /palm/trunk/SOURCE/Makefile
===================================================================
--- /palm/trunk/SOURCE/Makefile	(revision 2319)
+++ /palm/trunk/SOURCE/Makefile	(revision 2320)
@@ -25,4 +25,8 @@
 # -----------------
 # $Id$
+# -ls_forcing nudging
+# +large_scale_forcing_nudging
+# 
+# 2318 2017-07-20 17:27:44Z suehring
 # Add further dependencies on surface_mod
 # 
@@ -370,6 +374,6 @@
 	lpm_set_attributes.f90 lpm_splitting.f90 \
 	lpm_write_exchange_statistics.f90 lpm_write_restart_file.f90 \
-	ls_forcing_mod.f90 message.f90 microphysics_mod.f90 modules.f90 mod_kinds.f90 \
-	mod_particle_attributes.f90 netcdf_interface_mod.f90 nudging_mod.f90 \
+        large_scale_forcing_nudging_mod.f90 message.f90 microphysics_mod.f90 \
+        modules.f90 mod_kinds.f90 mod_particle_attributes.f90 netcdf_interface_mod.f90 \
 	outflow_turbulence.f90 package_parin.f90 palm.f90 parin.f90 plant_canopy_model_mod.f90 \
 	pmc_interface_mod.f90 pmc_child_mod.f90 pmc_general_mod.f90 pmc_handle_communicator_mod.f90 \
@@ -453,5 +457,5 @@
 check_open.o: modules.o mod_kinds.o mod_particle_attributes.o \
 	netcdf_interface_mod.o posix_calls_from_fortran.o
-check_parameters.o: modules.o mod_kinds.o land_surface_model_mod.o \
+check_parameters.o: modules.o mod_kinds.o land_surface_model_mod.o large_scale_forcing_nudging_mod.o \
 	microphysics_mod.o netcdf_interface_mod.o plant_canopy_model_mod.o pmc_interface_mod.o \
 	radiation_model_mod.o spectra_mod.o subsidence_mod.o synthetic_turbulence_generator_mod.o \
@@ -493,14 +497,14 @@
 fft_xy_mod.o: modules.o mod_kinds.o singleton_mod.o temperton_fft_mod.o
 flow_statistics.o: modules.o cpulog_mod.o mod_kinds.o land_surface_model_mod.o \
-   netcdf_interface_mod.o radiation_model_mod.o surface_mod.o
+                   large_scale_forcing_nudging_mod.o netcdf_interface_mod.o radiation_model_mod.o surface_mod.o
 global_min_max.o: modules.o mod_kinds.o
-header.o: modules.o cpulog_mod.o mod_kinds.o netcdf_interface_mod.o land_surface_model_mod.o\
-   microphysics_mod.o plant_canopy_model_mod.o pmc_handle_communicator_mod.o pmc_interface_mod.o \
-   radiation_model_mod.o spectra_mod.o subsidence_mod.o surface_mod.o virtual_flight_mod.o
+header.o: modules.o cpulog_mod.o mod_kinds.o large_scale_forcing_nudging_mod.o netcdf_interface_mod.o \
+   land_surface_model_mod.o microphysics_mod.o plant_canopy_model_mod.o pmc_handle_communicator_mod.o \
+   pmc_interface_mod.o radiation_model_mod.o spectra_mod.o subsidence_mod.o surface_mod.o virtual_flight_mod.o
 inflow_turbulence.o: modules.o cpulog_mod.o mod_kinds.o
 init_1d_model.o: modules.o mod_kinds.o
 init_3d_model.o: modules.o mod_kinds.o advec_ws.o cpulog_mod.o disturb_heatflux.o land_surface_model_mod.o \
-   lpm_init.o ls_forcing_mod.o netcdf_interface_mod.o plant_canopy_model_mod.o \
-   radiation_model_mod.o random_function_mod.o random_generator_parallel_mod.o \
+   large_scale_forcing_nudging_mod.o lpm_init.o netcdf_interface_mod.o \
+   plant_canopy_model_mod.o radiation_model_mod.o random_function_mod.o random_generator_parallel_mod.o \
    microphysics_mod.o mod_particle_attributes.o surface_layer_fluxes_mod.o \
    urban_surface_mod.o virtual_flight_mod.o surface_mod.o wind_turbine_model_mod.o 
@@ -518,4 +522,5 @@
 interaction_droplets_ptq.o: modules.o mod_kinds.o
 land_surface_model_mod.o: modules.o mod_kinds.o radiation_model_mod.o surface_mod.o
+large_scale_forcing_nudging_mod.o: modules.o cpulog_mod.o mod_kinds.o surface_mod.o
 local_stop.o: modules.o mod_kinds.o pmc_interface_mod.o
 local_tremain.o: modules.o cpulog_mod.o mod_kinds.o
@@ -551,5 +556,4 @@
 lpm_write_exchange_statistics.o: modules.o mod_kinds.o mod_particle_attributes.o
 lpm_write_restart_file.o: modules.o mod_kinds.o mod_particle_attributes.o
-ls_forcing_mod.o: modules.o cpulog_mod.o mod_kinds.o surface_mod.o
 message.o: modules.o mod_kinds.o pmc_interface_mod.o
 microphysics_mod.o: modules.o cpulog_mod.o mod_kinds.o surface_mod.o
@@ -559,8 +563,7 @@
 netcdf_interface_mod.o: netcdf_interface_mod.f90 modules.o mod_kinds.o \
    land_surface_model_mod.o radiation_model_mod.o spectra_mod.o urban_surface_mod.o
-nudging_mod.o: modules.o cpulog_mod.o mod_kinds.o
 outflow_turbulence.o: modules.o cpulog_mod.o mod_kinds.o
 package_parin.o: modules.o mod_kinds.o mod_particle_attributes.o
-palm.o: modules.o cpulog_mod.o land_surface_model_mod.o ls_forcing_mod.o mod_kinds.o nudging_mod.o\
+palm.o: modules.o cpulog_mod.o land_surface_model_mod.o mod_kinds.o \
 	pmc_interface_mod.o surface_layer_fluxes_mod.o
 parin.o: modules.o cpulog_mod.o land_surface_model_mod.o mod_kinds.o netcdf_interface_mod.o \
@@ -589,6 +592,7 @@
 	advec_ws.o buoyancy.o calc_radiation.o coriolis.o \
 	cpulog_mod.o diffusion_e.o diffusion_s.o diffusion_u.o diffusion_v.o diffusion_w.o \
-	eqn_state_seawater.o mod_kinds.o microphysics_mod.o \
-	nudging_mod.o plant_canopy_model_mod.o production_e.o radiation_model_mod.o \
+	eqn_state_seawater.o large_scale_forcing_nudging_mod.o modules.o cpulog_mod.o mod_kinds.o surface_mod.o \
+        mod_kinds.o microphysics_mod.o \
+	plant_canopy_model_mod.o production_e.o radiation_model_mod.o \
 	subsidence_mod.o surface_mod.o user_actions.o wind_turbine_model_mod.o
 progress_bar_mod.o: modules.o mod_kinds.o
@@ -620,5 +624,5 @@
 time_integration.o: modules.o advec_ws.o buoyancy.o calc_mean_profile.o \
 	cpulog_mod.o data_output_flight.o disturb_heatflux.o interaction_droplets_ptq.o land_surface_model_mod.o \
-	ls_forcing_mod.o mod_kinds.o nudging_mod.o pmc_interface_mod.o production_e.o \
+	large_scale_forcing_nudging_mod.o mod_kinds.o pmc_interface_mod.o production_e.o \
 	prognostic_equations.o progress_bar_mod.o radiation_model_mod.o \
 	spectra_mod.o user_actions.o microphysics_mod.o synthetic_turbulence_generator_mod.o \
Index: /palm/trunk/SOURCE/boundary_conds.f90
===================================================================
--- /palm/trunk/SOURCE/boundary_conds.f90	(revision 2319)
+++ /palm/trunk/SOURCE/boundary_conds.f90	(revision 2320)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Remove unused control parameter large_scale_forcing from only-list
+! 
+! 2292 2017-06-20 09:51:42Z schwenkel
 ! Implementation of new microphysic scheme: cloud_scheme = 'morrison' 
 ! includes two more prognostic equations for cloud drop concentration (nc)  
@@ -166,5 +169,5 @@
                ibc_pt_b, ibc_pt_t, ibc_q_b, ibc_q_t, ibc_s_b, ibc_s_t,         &
                ibc_sa_t, ibc_uv_b, ibc_uv_t, inflow_l, inflow_n, inflow_r,     &
-               inflow_s, intermediate_timestep_count, large_scale_forcing,     &
+               inflow_s, intermediate_timestep_count,                          &
                microphysics_morrison, microphysics_seifert, nest_domain,       &
                nest_bound_l, nest_bound_s, nudging, ocean, outflow_l,          &
Index: /palm/trunk/SOURCE/check_parameters.f90
===================================================================
--- /palm/trunk/SOURCE/check_parameters.f90	(revision 2319)
+++ /palm/trunk/SOURCE/check_parameters.f90	(revision 2320)
@@ -20,9 +20,12 @@
 ! Current revisions:
 ! -----------------
-!
-!
+! 
+! 
 ! Former revisions:
 ! -----------------
 ! $Id$
+! Modularize large-scale forcing and nudging
+!
+! 2312 2017-07-14 20:26:51Z hoffmann
 ! PA0349 and PA0420 removed.
 !
@@ -498,4 +501,8 @@
         ONLY: lsm_check_data_output, lsm_check_data_output_pr,                 &
               lsm_check_parameters
+
+    USE lsf_nudging_mod,                                                       &
+        ONLY:  lsf_nudging_check_parameters, lsf_nudging_check_data_output_pr, &
+               qsws_surf, shf_surf
 
     USE kinds
@@ -1227,4 +1234,8 @@
 !-- When radiation model is used, peform addtional checks
     IF ( radiation )  CALL radiation_check_parameters
+!
+!-- When large-scale forcing or nudging is used, peform addtional checks
+    IF ( large_scale_forcing  .OR.  nudging )  CALL lsf_nudging_check_parameters
+
 
 
@@ -2848,100 +2859,4 @@
              ENDIF
 
-          CASE ( 'td_lsa_lpt' )
-             IF (  .NOT.  large_scale_forcing )  THEN
-                message_string = 'data_output_pr = ' //                        &
-                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
-                                 'lemented for large_scale_forcing = .FALSE.'
-                CALL message( 'check_parameters', 'PA0393', 1, 2, 0, 6, 0 )
-             ELSE
-                dopr_index(i) = 81
-                dopr_unit(i)  = 'K/s'
-                hom(:,2,81,:) = SPREAD( zu, 2, statistic_regions+1 )
-             ENDIF
-
-          CASE ( 'td_lsa_q' )
-             IF (  .NOT.  large_scale_forcing )  THEN
-                message_string = 'data_output_pr = ' //                        &
-                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
-                                 'lemented for large_scale_forcing = .FALSE.'
-                CALL message( 'check_parameters', 'PA0393', 1, 2, 0, 6, 0 )
-             ELSE
-                dopr_index(i) = 82
-                dopr_unit(i)  = 'kg/kgs'
-                hom(:,2,82,:) = SPREAD( zu, 2, statistic_regions+1 )
-             ENDIF
-
-          CASE ( 'td_sub_lpt' )
-             IF (  .NOT.  large_scale_forcing )  THEN
-                message_string = 'data_output_pr = ' //                        &
-                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
-                                 'lemented for large_scale_forcing = .FALSE.'
-                CALL message( 'check_parameters', 'PA0393', 1, 2, 0, 6, 0 )
-             ELSE
-                dopr_index(i) = 83
-                dopr_unit(i)  = 'K/s'
-                hom(:,2,83,:) = SPREAD( zu, 2, statistic_regions+1 )
-             ENDIF
-
-          CASE ( 'td_sub_q' )
-             IF (  .NOT.  large_scale_forcing )  THEN
-                message_string = 'data_output_pr = ' //                        &
-                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
-                                 'lemented for large_scale_forcing = .FALSE.'
-                CALL message( 'check_parameters', 'PA0393', 1, 2, 0, 6, 0 )
-             ELSE
-                dopr_index(i) = 84
-                dopr_unit(i)  = 'kg/kgs'
-                hom(:,2,84,:) = SPREAD( zu, 2, statistic_regions+1 )
-             ENDIF
-
-          CASE ( 'td_nud_lpt' )
-             IF (  .NOT.  nudging )  THEN
-                message_string = 'data_output_pr = ' //                        &
-                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
-                                 'lemented for nudging = .FALSE.'
-                CALL message( 'check_parameters', 'PA0394', 1, 2, 0, 6, 0 )
-             ELSE
-                dopr_index(i) = 85
-                dopr_unit(i)  = 'K/s'
-                hom(:,2,85,:) = SPREAD( zu, 2, statistic_regions+1 )
-             ENDIF
-
-          CASE ( 'td_nud_q' )
-             IF (  .NOT.  nudging )  THEN
-                message_string = 'data_output_pr = ' //                        &
-                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
-                                 'lemented for nudging = .FALSE.'
-                CALL message( 'check_parameters', 'PA0394', 1, 2, 0, 6, 0 )
-             ELSE
-                dopr_index(i) = 86
-                dopr_unit(i)  = 'kg/kgs'
-                hom(:,2,86,:) = SPREAD( zu, 2, statistic_regions+1 )
-             ENDIF
-
-          CASE ( 'td_nud_u' )
-             IF (  .NOT.  nudging )  THEN
-                message_string = 'data_output_pr = ' //                        &
-                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
-                                 'lemented for nudging = .FALSE.'
-                CALL message( 'check_parameters', 'PA0394', 1, 2, 0, 6, 0 )
-             ELSE
-                dopr_index(i) = 87
-                dopr_unit(i)  = 'm/s2'
-                hom(:,2,87,:) = SPREAD( zu, 2, statistic_regions+1 )
-             ENDIF
-
-          CASE ( 'td_nud_v' )
-             IF (  .NOT.  nudging )  THEN
-                message_string = 'data_output_pr = ' //                        &
-                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
-                                 'lemented for nudging = .FALSE.'
-                CALL message( 'check_parameters', 'PA0394', 1, 2, 0, 6, 0 )
-             ELSE
-                dopr_index(i) = 88
-                dopr_unit(i)  = 'm/s2'
-                hom(:,2,88,:) = SPREAD( zu, 2, statistic_regions+1 )
-             ENDIF
-
           CASE ( 's*2' )
              IF (  .NOT.  passive_scalar )  THEN
@@ -2962,4 +2877,7 @@
              CALL lsm_check_data_output_pr( data_output_pr(i), i, unit,        &
                                             dopr_unit(i) )
+
+             CALL lsf_nudging_check_data_output_pr( data_output_pr(i), i,      &
+                                                    unit, dopr_unit(i) )
 
              IF ( unit == 'illegal' )  THEN
@@ -3885,45 +3803,4 @@
 
 !
-!-- 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 ( large_scale_forcing  .AND.  ( bc_lr /= 'cyclic'  .OR.                 &
-                                    bc_ns /= 'cyclic' ) )  THEN
-       message_string = 'Non-cyclic lateral boundaries do not allow for &' //  &
-                        'the usage of large scale forcing from external file.'
-       CALL message( 'check_parameters', 'PA0375', 1, 2, 0, 6, 0 )
-    ENDIF
-
-    IF ( large_scale_forcing  .AND.  (  .NOT.  humidity ) )  THEN
-       message_string = 'The usage of large scale forcing from external &'//   &
-                        'file LSF_DATA requires humidity = .T..'
-       CALL message( 'check_parameters', 'PA0376', 1, 2, 0, 6, 0 )
-    ENDIF
-
-    IF ( large_scale_forcing  .AND.  passive_scalar )  THEN
-       message_string = 'The usage of large scale forcing from external &'//   &
-                        'file LSF_DATA is not implemented for passive scalars'
-       CALL message( 'check_parameters', 'PA0440', 1, 2, 0, 6, 0 )
-    ENDIF
-
-    IF ( large_scale_forcing  .AND.  topography /= 'flat'                      &
-                              .AND.  .NOT.  lsf_exception )  THEN
-       message_string = 'The usage of large scale forcing from external &'//   &
-                        'file LSF_DATA is not implemented for non-flat topography'
-       CALL message( 'check_parameters', 'PA0377', 1, 2, 0, 6, 0 )
-    ENDIF
-
-    IF ( large_scale_forcing  .AND.  ocean )  THEN
-       message_string = 'The usage of large scale forcing from external &'//   &
-                        'file LSF_DATA is not implemented for ocean runs'
-       CALL message( 'check_parameters', 'PA0378', 1, 2, 0, 6, 0 )
-    ENDIF
-
-!
 !-- Prevent empty time records in volume, cross-section and masked data in case
 !-- of non-parallel netcdf-output in restart runs
Index: /palm/trunk/SOURCE/flow_statistics.f90
===================================================================
--- /palm/trunk/SOURCE/flow_statistics.f90	(revision 2319)
+++ /palm/trunk/SOURCE/flow_statistics.f90	(revision 2320)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Modularize large-scale forcing and nudging
+! 
+! 2296 2017-06-28 07:53:56Z maronga
 ! Enabled output of radiation quantities for radiation_scheme = 'constant'
 ! 
@@ -243,6 +246,5 @@
                momentumflux_output_conversion, nc, nr, p, prho, prr, pt, q,    &
                qc, ql, qr, rho_air, rho_air_zw, rho_ocean, s,                  &
-               sa, td_lsa_lpt, td_lsa_q, td_sub_lpt, td_sub_q, time_vert, u,   &
-               ug, v, vg, vpt, w, w_subs, waterflux_output_conversion, zw
+               sa, u, ug, v, vg, vpt, w, w_subs, waterflux_output_conversion, zw
         
     USE cloud_parameters,                                                      &
@@ -272,4 +274,7 @@
     USE land_surface_model_mod,                                                &
         ONLY:   m_soil_h, nzb_soil, nzt_soil, t_soil_h
+
+    USE lsf_nudging_mod,                                                       &
+        ONLY:   td_lsa_lpt, td_lsa_q, td_sub_lpt, td_sub_q, time_vert
 
     USE netcdf_interface,                                                      &
Index: /palm/trunk/SOURCE/header.f90
===================================================================
--- /palm/trunk/SOURCE/header.f90	(revision 2319)
+++ /palm/trunk/SOURCE/header.f90	(revision 2320)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Modularize large-scale forcing and nudging
+! 
+! 2300 2017-06-29 13:31:14Z raasch
 ! host-specific code removed
 ! 
@@ -342,4 +345,7 @@
     USE land_surface_model_mod,                                                &
         ONLY: lsm_header
+
+    USE lsf_nudging_mod,                                                       &
+        ONLY:  lsf_nudging_header
 
     USE microphysics_mod,                                                      &
@@ -746,41 +752,5 @@
 !
 !-- Large scale forcing and nudging
-    WRITE ( io, 160 )
-    IF ( large_scale_forcing )  THEN
-       WRITE ( io, 162 )
-       WRITE ( io, 163 )
-
-       IF ( large_scale_subsidence )  THEN
-          IF ( .NOT. use_subsidence_tendencies )  THEN
-             WRITE ( io, 164 )
-          ELSE
-             WRITE ( io, 165 )
-          ENDIF 
-       ENDIF
-
-       IF ( bc_pt_b == 'dirichlet' )  THEN
-          WRITE ( io, 180 )
-       ELSEIF ( bc_pt_b == 'neumann' )  THEN
-          WRITE ( io, 181 )
-       ENDIF
-
-       IF ( bc_q_b == 'dirichlet' )  THEN
-          WRITE ( io, 182 )
-       ELSEIF ( bc_q_b == 'neumann' )  THEN
-          WRITE ( io, 183 )
-       ENDIF
-
-       WRITE ( io, 167 )
-       IF ( nudging )  THEN
-          WRITE ( io, 170 )
-       ENDIF
-    ELSE
-       WRITE ( io, 161 )
-       WRITE ( io, 171 )
-    ENDIF
-    IF ( large_scale_subsidence )  THEN
-       WRITE ( io, 168 )
-       WRITE ( io, 169 )
-    ENDIF
+    IF ( large_scale_forcing )  CALL lsf_nudging_header( io )
 
 !
@@ -2012,21 +1982,4 @@
            /'     ',2(1X,E12.5),'Pa/m in x/y direction', &
            /'     starting from dp_level_b =', F8.3, 'm', A /)
-160 FORMAT (//' Large scale forcing and nudging:'/ &
-              ' -------------------------------'/)
-161 FORMAT (' --> No large scale forcing from external is used (default) ')
-162 FORMAT (' --> Large scale forcing from external file LSF_DATA is used: ')
-163 FORMAT ('     - large scale advection tendencies ')
-164 FORMAT ('     - large scale subsidence velocity w_subs ')
-165 FORMAT ('     - large scale subsidence tendencies ')
-167 FORMAT ('     - and geostrophic wind components ug and vg')
-168 FORMAT (' --> Large-scale vertical motion is used in the ', &
-                  'prognostic equation(s) for')
-169 FORMAT ('     the scalar(s) only')
-170 FORMAT (' --> Nudging is used')
-171 FORMAT (' --> No nudging is used (default) ')
-180 FORMAT ('     - prescribed surface values for temperature')
-181 FORMAT ('     - prescribed surface fluxes for temperature')
-182 FORMAT ('     - prescribed surface values for humidity')
-183 FORMAT ('     - prescribed surface fluxes for humidity')
 200 FORMAT (//' Run time and time step information:'/ &
              ' ----------------------------------'/)
Index: /palm/trunk/SOURCE/init_3d_model.f90
===================================================================
--- /palm/trunk/SOURCE/init_3d_model.f90	(revision 2319)
+++ /palm/trunk/SOURCE/init_3d_model.f90	(revision 2320)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Modularize large-scale forcing and nudging
+! 
+! 2292 2017-06-20 09:51:42Z schwenkel
 ! Implementation of new microphysic scheme: cloud_scheme = 'morrison' 
 ! includes two more prognostic equations for cloud drop concentration (nc)  
@@ -372,5 +375,6 @@
         ONLY:  lsm_init, lsm_init_arrays
   
-    USE ls_forcing_mod
+    USE lsf_nudging_mod,                                                       &
+        ONLY:  lsf_init, ls_forcing_surf, nudge_init 
 
     USE microphysics_mod,                                                      &
@@ -956,4 +960,16 @@
 
 !
+!-- Initialize nudging if required
+    IF ( nudging )  THEN
+       CALL nudge_init
+    ENDIF
+
+!
+!-- Initialize reading of large scale forcing from external file - if required
+    IF ( large_scale_forcing )  THEN
+       CALL lsf_init
+    ENDIF
+
+!
 !-- Allocate arrays containing the RK coefficient for calculation of 
 !-- perturbation pressure and turbulent fluxes. At this point values are
@@ -1092,19 +1108,4 @@
 
           CALL location_message( 'initializing with constant profiles', .FALSE. )
-!
-!--       Overwrite initial profiles in case of nudging
-          IF ( nudging )  THEN
-             pt_init = ptnudge(:,1)
-             u_init  = unudge(:,1)
-             v_init  = vnudge(:,1)
-             IF ( humidity  )  THEN ! is passive_scalar correct???
-                q_init = qnudge(:,1)
-             ENDIF
-
-             WRITE( message_string, * ) 'Initial profiles of u, v and ',       &
-                 'scalars from NUDGING_DATA are used.'
-             CALL message( 'init_3d_model', 'PA0370', 0, 0, 0, 6, 0 )
-          ENDIF
-
 !
 !--       Overwrite initial profiles in case of synthetic turbulence generator
Index: /palm/trunk/SOURCE/large_scale_forcing_nudging_mod.f90
===================================================================
--- /palm/trunk/SOURCE/large_scale_forcing_nudging_mod.f90	(revision 2320)
+++ /palm/trunk/SOURCE/large_scale_forcing_nudging_mod.f90	(revision 2320)
@@ -0,0 +1,1394 @@
+!> @file lsf_nudging_mod.f90
+!------------------------------------------------------------------------------!
+! This file is part of PALM.
+!
+! PALM is free software: you can redistribute it and/or modify it under the
+! terms of the GNU General Public License as published by the Free Software
+! Foundation, either version 3 of the License, or (at your option) any later
+! version.
+!
+! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
+! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
+! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
+!
+! You should have received a copy of the GNU General Public License along with
+! PALM. If not, see <http://www.gnu.org/licenses/>.
+!
+! Copyright 1997-2017 Leibniz Universitaet Hannover
+!------------------------------------------------------------------------------!
+!
+! Current revisions:
+! ------------------
+! 
+! 
+! Former revisions:
+! -----------------
+! $Id$
+!
+! Description:
+! ------------
+!> Calculates large scale forcings (geostrophic wind and subsidence velocity) as 
+!> well as surfaces fluxes dependent on time given in an external file (LSF_DATA).
+!> Moreover, module contains nudging routines, where u, v, pt and q are nudged 
+!> to given profiles on a relaxation timescale tnudge. 
+!> Profiles are read in from NUDGIN_DATA. 
+!> Code is based on Neggers et al. (2012) and also in parts on DALES and UCLA-LES.
+!> @todo: Revise reading of ASCII-files
+!> @todo: Remove unused variables and control flags
+!> @todo: Revise large-scale facing of surface variables
+!> @todo: Revise control flags lsf_exception, lsf_surf, lsf_vert, etc. 
+!--------------------------------------------------------------------------------!
+ MODULE lsf_nudging_mod
+
+    USE arrays_3d,                                                             &
+        ONLY:  heatflux_input_conversion, pt, pt_init, q, q_init, tend,        &
+               u, u_init, ug, v, v_init, vg, w_subs,                           &
+               waterflux_input_conversion, zu                  
+
+    USE control_parameters,                                                    &
+        ONLY:  bc_lr, bc_ns, bc_pt_b, bc_q_b, data_output_pr, dt_3d, end_time, &  
+               humidity, intermediate_timestep_count, ibc_pt_b, ibc_q_b,       &
+               large_scale_forcing, large_scale_subsidence, lsf_surf, lsf_vert,&
+               lsf_exception, message_string, nudging, passive_scalar,         &
+               pt_surface, ocean, q_surface, surface_pressure, topography,     &
+               use_subsidence_tendencies
+
+
+    USE indices,                                                               &
+        ONLY:  ngp_sums_ls, nxl, nxr, nys, nyn, nzb, nz, nzt, wall_flags_0
+
+    USE kinds
+
+    USE surface_mod,                                                           &
+        ONLY:  surf_def_h, surf_lsm_h, surf_usm_h
+
+    USE statistics,                                                            &
+        ONLY:  hom, statistic_regions, sums_ls_l, weight_substep
+
+
+    INTEGER(iwp) ::  nlsf = 1000                       !< maximum number of profiles in LSF_DATA (large scale forcing)
+    INTEGER(iwp) ::  ntnudge = 1000                    !< maximum number of profiles in NUDGING_DATA (nudging)
+
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  ptnudge     !< vertical profile of pot. temperature interpolated to vertical grid (nudging)
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  qnudge      !< vertical profile of specific humidity interpolated to vertical grid (nudging) 
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  tnudge      !< vertical profile of nudging time scale interpolated to vertical grid (nudging)  
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  td_lsa_lpt  !< temperature tendency due to large scale advection (large scale forcing)
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  td_lsa_q    !< specific humidity tendency due to large scale advection (large scale forcing)
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  td_sub_lpt  !< temperature tendency due to subsidence/ascent (large scale forcing)
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  td_sub_q    !< specific humidity tendency due to subsidence/ascent (large scale forcing)
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  ug_vert     !< vertical profile of geostrophic wind component in x-direction interpolated to vertical grid (large scale forcing)
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  unudge      !< vertical profile of wind component in x-direction interpolated to vertical grid (nudging) 
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  vnudge      !< vertical profile of wind component in y-direction interpolated to vertical grid (nudging) 
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  vg_vert     !< vertical profile of geostrophic wind component in y-direction interpolated to vertical grid (large scale forcing)
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  wnudge      !< vertical profile of subsidence/ascent velocity interpolated to vertical grid (nudging) ???
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  wsubs_vert  !< vertical profile of wind component in z-direction interpolated to vertical grid (nudging) ???
+
+    REAL(wp), DIMENSION(:), ALLOCATABLE ::  shf_surf      !< time-dependent surface sensible heat flux (large scale forcing)
+    REAL(wp), DIMENSION(:), ALLOCATABLE ::  timenudge     !< times at which vertical profiles are defined in NUDGING_DATA (nudging)
+    REAL(wp), DIMENSION(:), ALLOCATABLE ::  time_surf     !< times at which surface values/fluxes are defined in LSF_DATA (large scale forcing)
+    REAL(wp), DIMENSION(:), ALLOCATABLE ::  time_vert     !< times at which vertical profiles are defined in LSF_DATA (large scale forcing)
+
+    REAL(wp), DIMENSION(:), ALLOCATABLE ::  tmp_tnudge    !< current nudging time scale
+
+    REAL(wp), DIMENSION(:), ALLOCATABLE ::  p_surf        !< time-dependent surface pressure (large scale forcing)
+    REAL(wp), DIMENSION(:), ALLOCATABLE ::  pt_surf       !< time-dependent surface temperature (large scale forcing)
+    REAL(wp), DIMENSION(:), ALLOCATABLE ::  qsws_surf     !< time-dependent surface latent heat flux (large scale forcing)
+    REAL(wp), DIMENSION(:), ALLOCATABLE ::  q_surf        !< time-dependent surface specific humidity (large scale forcing)
+
+    SAVE
+    PRIVATE
+!
+!-- Public subroutines
+    PUBLIC ls_forcing_surf, ls_forcing_vert, ls_advec, lsf_init,               &
+           lsf_nudging_check_parameters, nudge_init,                           &
+           lsf_nudging_check_data_output_pr, lsf_nudging_header,               &
+           calc_tnudge, nudge, nudge_ref
+!
+!-- Public variables
+    PUBLIC qsws_surf, shf_surf, td_lsa_lpt, td_lsa_q, td_sub_lpt,              &
+           td_sub_q, time_vert
+
+    INTERFACE ls_advec
+       MODULE PROCEDURE ls_advec
+       MODULE PROCEDURE ls_advec_ij
+    END INTERFACE ls_advec
+
+    INTERFACE nudge
+       MODULE PROCEDURE nudge
+       MODULE PROCEDURE nudge_ij
+    END INTERFACE nudge
+
+ CONTAINS
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> @todo Missing subroutine description.
+!------------------------------------------------------------------------------!
+    SUBROUTINE lsf_nudging_check_parameters
+
+       IMPLICIT NONE
+!
+!--    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 ( large_scale_forcing  .AND.  ( bc_lr /= 'cyclic'  .OR.              &
+                                       bc_ns /= 'cyclic' ) )  THEN
+          message_string = 'Non-cyclic lateral boundaries do not allow for &'//&
+                        'the usage of large scale forcing from external file.'
+          CALL message( 'check_parameters', 'PA0375', 1, 2, 0, 6, 0 )
+       ENDIF
+
+       IF ( large_scale_forcing  .AND.  (  .NOT.  humidity ) )  THEN
+          message_string = 'The usage of large scale forcing from external &'//&
+                        'file LSF_DATA requires humidity = .T..'
+          CALL message( 'check_parameters', 'PA0376', 1, 2, 0, 6, 0 )
+       ENDIF
+
+       IF ( large_scale_forcing  .AND.  passive_scalar )  THEN
+          message_string = 'The usage of large scale forcing from external &'//&
+                        'file LSF_DATA is not implemented for passive scalars'
+          CALL message( 'check_parameters', 'PA0440', 1, 2, 0, 6, 0 )
+       ENDIF
+
+       IF ( large_scale_forcing  .AND.  topography /= 'flat'                   &
+                              .AND.  .NOT.  lsf_exception )  THEN
+          message_string = 'The usage of large scale forcing from external &'//&
+                        'file LSF_DATA is not implemented for non-flat topography'
+          CALL message( 'check_parameters', 'PA0377', 1, 2, 0, 6, 0 )
+       ENDIF
+
+       IF ( large_scale_forcing  .AND.  ocean )  THEN
+          message_string = 'The usage of large scale forcing from external &'//&
+                        'file LSF_DATA is not implemented for ocean runs'
+          CALL message( 'check_parameters', 'PA0378', 1, 2, 0, 6, 0 )
+       ENDIF
+
+    END SUBROUTINE lsf_nudging_check_parameters
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Check data output of profiles for land surface model
+!------------------------------------------------------------------------------!
+    SUBROUTINE lsf_nudging_check_data_output_pr( variable, var_count, unit,    &
+                                                 dopr_unit )
+ 
+       USE profil_parameter
+
+       IMPLICIT NONE
+   
+       CHARACTER (LEN=*) ::  unit      !< 
+       CHARACTER (LEN=*) ::  variable  !< 
+       CHARACTER (LEN=*) ::  dopr_unit !< local value of dopr_unit
+ 
+       INTEGER(iwp) ::  user_pr_index !<    
+       INTEGER(iwp) ::  var_count     !< 
+
+       SELECT CASE ( TRIM( variable ) )
+       
+
+          CASE ( 'td_lsa_lpt' )
+             IF (  .NOT.  large_scale_forcing )  THEN
+                message_string = 'data_output_pr = ' //                        &
+                                 TRIM( data_output_pr(var_count) ) //          &
+                                 ' is not implemented for ' //                 &
+                                 'large_scale_forcing = .FALSE.'
+                CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393',    &
+                               1, 2, 0, 6, 0 )
+             ELSE
+                dopr_index(var_count) = 81
+                dopr_unit             = 'K/s'
+                unit                  = 'K/s'
+                hom(:,2,81,:) = SPREAD( zu, 2, statistic_regions+1 )
+             ENDIF
+
+          CASE ( 'td_lsa_q' )
+             IF (  .NOT.  large_scale_forcing )  THEN
+                message_string = 'data_output_pr = ' //                        &
+                                 TRIM( data_output_pr(var_count) ) //          &
+                                 ' is not implemented for ' //                 &
+                                 'large_scale_forcing = .FALSE.'
+                CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393',    &
+                               1, 2, 0, 6, 0 )
+             ELSE
+                dopr_index(var_count) = 82
+                dopr_unit             = 'kg/kgs'
+                unit                  = 'kg/kgs'
+                hom(:,2,82,:) = SPREAD( zu, 2, statistic_regions+1 )
+             ENDIF
+          CASE ( 'td_sub_lpt' )
+             IF (  .NOT.  large_scale_forcing )  THEN
+                message_string = 'data_output_pr = ' //                        &
+                                 TRIM( data_output_pr(var_count) ) //          &
+                                 ' is not implemented for ' //                 &
+                                 'large_scale_forcing = .FALSE.'
+                CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393',    &
+                               1, 2, 0, 6, 0 )
+             ELSE
+                dopr_index(var_count) = 83
+                dopr_unit             = 'K/s'
+                unit                  = 'K/s'
+                hom(:,2,83,:) = SPREAD( zu, 2, statistic_regions+1 )
+             ENDIF
+
+          CASE ( 'td_sub_q' )
+             IF (  .NOT.  large_scale_forcing )  THEN
+                message_string = 'data_output_pr = ' //                        &
+                                 TRIM( data_output_pr(var_count) ) //          &
+                                 ' is not implemented for ' //                 &
+                                 'large_scale_forcing = .FALSE.'
+                CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393',    &
+                               1, 2, 0, 6, 0 )
+             ELSE
+                dopr_index(var_count) = 84
+                dopr_unit             = 'kg/kgs'
+                unit                  = 'kg/kgs'
+                hom(:,2,84,:) = SPREAD( zu, 2, statistic_regions+1 )
+             ENDIF
+
+          CASE ( 'td_nud_lpt' )
+             IF (  .NOT.  nudging )  THEN
+                message_string = 'data_output_pr = ' //                        &
+                                 TRIM( data_output_pr(var_count) ) //          &
+                                 ' is not implemented for ' //                 &
+                                 'large_scale_forcing = .FALSE.'
+                CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393',    &
+                               1, 2, 0, 6, 0 )
+             ELSE
+                dopr_index(var_count) = 85
+                dopr_unit             = 'K/s'
+                unit                  = 'K/s'
+                hom(:,2,85,:) = SPREAD( zu, 2, statistic_regions+1 )
+             ENDIF
+
+          CASE ( 'td_nud_q' )
+             IF (  .NOT.  nudging )  THEN
+                message_string = 'data_output_pr = ' //                        &
+                                 TRIM( data_output_pr(var_count) ) //          &
+                                 ' is not implemented for ' //                 &
+                                 'large_scale_forcing = .FALSE.'
+                CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393',    &
+                               1, 2, 0, 6, 0 )
+             ELSE
+                dopr_index(var_count) = 86
+                dopr_unit             = 'kg/kgs'
+                unit                  = 'kg/kgs'
+                hom(:,2,86,:) = SPREAD( zu, 2, statistic_regions+1 )
+             ENDIF
+
+          CASE ( 'td_nud_u' )
+             IF (  .NOT.  nudging )  THEN
+                message_string = 'data_output_pr = ' //                        &
+                                 TRIM( data_output_pr(var_count) ) //          &
+                                 ' is not implemented for ' //                 &
+                                 'large_scale_forcing = .FALSE.'
+                CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393',    &
+                               1, 2, 0, 6, 0 )
+             ELSE
+                dopr_index(var_count) = 87
+                dopr_unit             = 'm/s2'
+                unit                  = 'm/s2'
+                hom(:,2,87,:) = SPREAD( zu, 2, statistic_regions+1 )
+             ENDIF
+
+          CASE ( 'td_nud_v' )
+             IF (  .NOT.  nudging )  THEN
+                message_string = 'data_output_pr = ' //                        &
+                                 TRIM( data_output_pr(var_count) ) //          &
+                                 ' is not implemented for ' //                 &
+                                 'large_scale_forcing = .FALSE.'
+                CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393',    &
+                               1, 2, 0, 6, 0 )
+             ELSE
+                dopr_index(var_count) = 88
+                dopr_unit             = 'm/s2'
+                unit                  = 'm/s2'
+                hom(:,2,88,:) = SPREAD( zu, 2, statistic_regions+1 )
+             ENDIF
+
+
+          CASE DEFAULT
+             unit = 'illegal'
+   
+       END SELECT
+
+    END SUBROUTINE lsf_nudging_check_data_output_pr
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> @todo Missing subroutine description.
+!------------------------------------------------------------------------------!
+    SUBROUTINE lsf_nudging_header ( io )
+
+       IMPLICIT NONE
+
+       INTEGER(iwp), INTENT(IN) ::  io !< Unit of the output file
+
+       WRITE ( io, 1 )
+       IF ( large_scale_forcing )  THEN
+          WRITE ( io, 3 )
+          WRITE ( io, 4 )
+
+          IF ( large_scale_subsidence )  THEN
+             IF ( .NOT. use_subsidence_tendencies )  THEN
+                WRITE ( io, 5 )
+             ELSE
+                WRITE ( io, 6 )
+             ENDIF 
+          ENDIF
+
+          IF ( bc_pt_b == 'dirichlet' )  THEN
+             WRITE ( io, 12 )
+          ELSEIF ( bc_pt_b == 'neumann' )  THEN
+             WRITE ( io, 13 )
+          ENDIF
+
+          IF ( bc_q_b == 'dirichlet' )  THEN
+             WRITE ( io, 14 )
+          ELSEIF ( bc_q_b == 'neumann' )  THEN
+             WRITE ( io, 15 )
+          ENDIF
+
+          WRITE ( io, 7 )
+          IF ( nudging )  THEN
+             WRITE ( io, 10 )
+          ENDIF
+       ELSE
+          WRITE ( io, 2 )
+          WRITE ( io, 11 )
+       ENDIF
+       IF ( large_scale_subsidence )  THEN
+          WRITE ( io, 8 )
+          WRITE ( io, 9 )
+       ENDIF
+
+
+1 FORMAT (//' Large scale forcing and nudging:'/ &
+              ' -------------------------------'/)
+2 FORMAT (' --> No large scale forcing from external is used (default) ')
+3 FORMAT (' --> Large scale forcing from external file LSF_DATA is used: ')
+4 FORMAT ('     - large scale advection tendencies ')
+5 FORMAT ('     - large scale subsidence velocity w_subs ')
+6 FORMAT ('     - large scale subsidence tendencies ')
+7 FORMAT ('     - and geostrophic wind components ug and vg')
+8 FORMAT (' --> Large-scale vertical motion is used in the ', &
+                  'prognostic equation(s) for')
+9 FORMAT ('     the scalar(s) only')
+10 FORMAT (' --> Nudging is used')
+11 FORMAT (' --> No nudging is used (default) ')
+12 FORMAT ('     - prescribed surface values for temperature')
+13 FORMAT ('     - prescribed surface fluxes for temperature')
+14 FORMAT ('     - prescribed surface values for humidity')
+15 FORMAT ('     - prescribed surface fluxes for humidity')
+
+    END SUBROUTINE lsf_nudging_header 
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> @todo Missing subroutine description.
+!------------------------------------------------------------------------------!
+    SUBROUTINE lsf_init
+
+       IMPLICIT NONE
+
+       CHARACTER(100) ::  chmess      !<
+       CHARACTER(1)   ::  hash        !<
+
+       INTEGER(iwp) ::  ierrn         !<
+       INTEGER(iwp) ::  finput = 90   !<
+       INTEGER(iwp) ::  k             !<
+       INTEGER(iwp) ::  nt             !<
+
+       REAL(wp) ::  fac               !<
+       REAL(wp) ::  highheight        !<
+       REAL(wp) ::  highug_vert       !<
+       REAL(wp) ::  highvg_vert       !<
+       REAL(wp) ::  highwsubs_vert    !<
+       REAL(wp) ::  lowheight         !<
+       REAL(wp) ::  lowug_vert        !<
+       REAL(wp) ::  lowvg_vert        !<
+       REAL(wp) ::  lowwsubs_vert     !<
+       REAL(wp) ::  high_td_lsa_lpt   !<
+       REAL(wp) ::  low_td_lsa_lpt    !<
+       REAL(wp) ::  high_td_lsa_q     !<
+       REAL(wp) ::  low_td_lsa_q      !<
+       REAL(wp) ::  high_td_sub_lpt   !<
+       REAL(wp) ::  low_td_sub_lpt    !<
+       REAL(wp) ::  high_td_sub_q     !<
+       REAL(wp) ::  low_td_sub_q      !<
+       REAL(wp) ::  r_dummy           !<
+
+       ALLOCATE( p_surf(0:nlsf), pt_surf(0:nlsf), q_surf(0:nlsf),              &
+                 qsws_surf(0:nlsf), shf_surf(0:nlsf),                          &
+                 td_lsa_lpt(nzb:nzt+1,0:nlsf), td_lsa_q(nzb:nzt+1,0:nlsf),     &
+                 td_sub_lpt(nzb:nzt+1,0:nlsf), td_sub_q(nzb:nzt+1,0:nlsf),     &
+                 time_vert(0:nlsf), time_surf(0:nlsf),                         &
+                 ug_vert(nzb:nzt+1,0:nlsf), vg_vert(nzb:nzt+1,0:nlsf),         &
+                 wsubs_vert(nzb:nzt+1,0:nlsf) )
+
+       p_surf = 0.0_wp; pt_surf = 0.0_wp; q_surf = 0.0_wp; qsws_surf = 0.0_wp
+       shf_surf = 0.0_wp; time_vert = 0.0_wp; td_lsa_lpt = 0.0_wp
+       td_lsa_q = 0.0_wp; td_sub_lpt = 0.0_wp; td_sub_q = 0.0_wp
+       time_surf = 0.0_wp; ug_vert = 0.0_wp; vg_vert = 0.0_wp
+       wsubs_vert = 0.0_wp
+
+!
+!--    Array for storing large scale forcing and nudging tendencies at each 
+!--    timestep for data output
+       ALLOCATE( sums_ls_l(nzb:nzt+1,0:7) )
+       sums_ls_l = 0.0_wp
+
+       ngp_sums_ls = (nz+2)*6
+
+       OPEN ( finput, FILE='LSF_DATA', STATUS='OLD', &
+              FORM='FORMATTED', IOSTAT=ierrn )
+
+       IF ( ierrn /= 0 )  THEN
+          message_string = 'file LSF_DATA does not exist'
+          CALL message( 'ls_forcing', 'PA0368', 1, 2, 0, 6, 0 )
+       ENDIF
+
+       ierrn = 0
+!
+!--    First three lines of LSF_DATA contain header
+       READ ( finput, FMT='(a100)', IOSTAT=ierrn ) chmess
+       READ ( finput, FMT='(a100)', IOSTAT=ierrn ) chmess
+       READ ( finput, FMT='(a100)', IOSTAT=ierrn ) chmess
+
+       IF ( ierrn /= 0 )  THEN
+          message_string = 'errors in file LSF_DATA'
+          CALL message( 'ls_forcing', 'PA0369', 1, 2, 0, 6, 0 )
+       ENDIF
+
+!
+!--    Surface values are read in
+       nt     = 0
+       ierrn = 0
+
+       DO WHILE ( time_surf(nt) < end_time )
+          nt = nt + 1
+          READ ( finput, *, IOSTAT = ierrn ) time_surf(nt), shf_surf(nt),      &
+                                             qsws_surf(nt), pt_surf(nt),       &
+                                             q_surf(nt), p_surf(nt)
+
+          IF ( ierrn < 0 )  THEN
+            WRITE ( message_string, * ) 'No time dependent surface variables ',&
+                              'in&LSF_DATA for end of run found'
+
+             CALL message( 'ls_forcing', 'PA0363', 1, 2, 0, 6, 0 )
+          ENDIF
+       ENDDO
+
+       IF ( time_surf(1) > end_time )  THEN
+          WRITE ( message_string, * ) 'Time dependent surface variables in ',  &
+                                      '&LSF_DATA set in after end of ' ,      &
+                                      'simulation - lsf_surf is set to FALSE'
+          CALL message( 'ls_forcing', 'PA0371', 0, 0, 0, 6, 0 )
+          lsf_surf = .FALSE.
+       ENDIF
+
+!
+!--    Go to the end of the list with surface variables
+       DO WHILE ( ierrn == 0 )
+          READ ( finput, *, IOSTAT = ierrn ) r_dummy
+       ENDDO
+
+!
+!--    Profiles of ug, vg and w_subs are read in (large scale forcing)
+
+       nt = 0
+       DO WHILE ( time_vert(nt) < end_time )
+          nt = nt + 1
+          hash = "#"
+          ierrn = 1 ! not zero
+!
+!--       Search for the next line consisting of "# time", 
+!--       from there onwards the profiles will be read
+          DO WHILE ( .NOT. ( hash == "#" .AND. ierrn == 0 ) ) 
+             READ ( finput, *, IOSTAT=ierrn ) hash, time_vert(nt)
+             IF ( ierrn < 0 )  THEN 
+                WRITE( message_string, * ) 'No time dependent vertical profiles',&
+                                 ' in&LSF_DATA for end of run found'
+                CALL message( 'ls_forcing', 'PA0372', 1, 2, 0, 6, 0 )
+             ENDIF
+          ENDDO
+
+          IF ( nt == 1 .AND. time_vert(nt) > end_time ) EXIT
+
+          READ ( finput, *, IOSTAT=ierrn ) lowheight, lowug_vert, lowvg_vert,  &
+                                           lowwsubs_vert, low_td_lsa_lpt,      &
+                                           low_td_lsa_q, low_td_sub_lpt,       &
+                                           low_td_sub_q
+          IF ( ierrn /= 0 )  THEN
+             message_string = 'errors in file LSF_DATA'
+             CALL message( 'ls_forcing', 'PA0369', 1, 2, 0, 6, 0 )
+          ENDIF
+
+          READ ( finput, *, IOSTAT=ierrn ) highheight, highug_vert,            &
+                                           highvg_vert, highwsubs_vert,        &
+                                           high_td_lsa_lpt, high_td_lsa_q,     &
+                                           high_td_sub_lpt, high_td_sub_q
+      
+          IF ( ierrn /= 0 )  THEN
+             message_string = 'errors in file LSF_DATA'
+             CALL message( 'ls_forcing', 'PA0369', 1, 2, 0, 6, 0 )
+          ENDIF
+
+
+          DO  k = nzb, nzt+1
+             IF ( highheight < zu(k) )  THEN
+                lowheight      = highheight
+                lowug_vert     = highug_vert
+                lowvg_vert     = highvg_vert
+                lowwsubs_vert  = highwsubs_vert
+                low_td_lsa_lpt = high_td_lsa_lpt
+                low_td_lsa_q   = high_td_lsa_q
+                low_td_sub_lpt = high_td_sub_lpt
+                low_td_sub_q   = high_td_sub_q
+
+                ierrn = 0
+                READ ( finput, *, IOSTAT=ierrn ) highheight, highug_vert,      &
+                                                 highvg_vert, highwsubs_vert,  &
+                                                 high_td_lsa_lpt,              &
+                                                 high_td_lsa_q,                &
+                                                 high_td_sub_lpt, high_td_sub_q
+
+                IF ( ierrn /= 0 )  THEN
+                   WRITE( message_string, * ) 'zu(',k,') = ', zu(k), 'm ', &
+                        'is higher than the maximum height in LSF_DATA which ',&
+                        'is ', lowheight, 'm. Interpolation on PALM ',         &
+                        'grid is not possible.'
+                   CALL message( 'ls_forcing', 'PA0395', 1, 2, 0, 6, 0 )
+                ENDIF
+
+             ENDIF
+
+!
+!--          Interpolation of prescribed profiles in space 
+             fac = (highheight-zu(k))/(highheight - lowheight)
+
+             ug_vert(k,nt)    = fac * lowug_vert                               &
+                                + ( 1.0_wp - fac ) * highug_vert
+             vg_vert(k,nt)    = fac * lowvg_vert                               &
+                                + ( 1.0_wp - fac ) * highvg_vert
+             wsubs_vert(k,nt) = fac * lowwsubs_vert                            &
+                                + ( 1.0_wp - fac ) * highwsubs_vert
+
+             td_lsa_lpt(k,nt) = fac * low_td_lsa_lpt                           &
+                                + ( 1.0_wp - fac ) * high_td_lsa_lpt
+             td_lsa_q(k,nt)   = fac * low_td_lsa_q                             &
+                                + ( 1.0_wp - fac ) * high_td_lsa_q
+             td_sub_lpt(k,nt) = fac * low_td_sub_lpt                           &
+                                + ( 1.0_wp - fac ) * high_td_sub_lpt
+             td_sub_q(k,nt)   = fac * low_td_sub_q                             &
+                                + ( 1.0_wp - fac ) * high_td_sub_q
+
+          ENDDO
+
+       ENDDO 
+
+!
+!--    Large scale vertical velocity has to be zero at the surface
+       wsubs_vert(nzb,:) = 0.0_wp
+ 
+       IF ( time_vert(1) > end_time )  THEN
+          WRITE ( message_string, * ) 'Time dependent large scale profile ',   &
+                             'forcing from&LSF_DATA sets in after end of ' ,   &
+                             'simulation - lsf_vert is set to FALSE'
+          CALL message( 'ls_forcing', 'PA0373', 0, 0, 0, 6, 0 )
+          lsf_vert = .FALSE.
+       ENDIF
+
+       CLOSE( finput )
+
+
+    END SUBROUTINE lsf_init
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> @todo Missing subroutine description.
+!------------------------------------------------------------------------------!
+    SUBROUTINE ls_forcing_surf ( time )
+
+       IMPLICIT NONE
+
+       INTEGER(iwp) ::  nt                     !<
+
+       REAL(wp)             :: dum_surf_flux  !<
+       REAL(wp)             :: fac            !<
+       REAL(wp), INTENT(in) :: time           !<
+
+!
+!--    Interpolation in time of LSF_DATA at the surface
+       nt = 1
+       DO WHILE ( time > time_surf(nt) )
+          nt = nt + 1
+       ENDDO
+       IF ( time /= time_surf(nt) )  THEN
+         nt = nt - 1
+       ENDIF
+
+       fac = ( time -time_surf(nt) ) / ( time_surf(nt+1) - time_surf(nt) )
+
+       IF ( ibc_pt_b == 0 )  THEN
+!
+!--       In case of Dirichlet boundary condition shf must not 
+!--       be set - it is calculated via MOST in prandtl_fluxes
+          pt_surface = pt_surf(nt) + fac * ( pt_surf(nt+1) - pt_surf(nt) )
+
+       ELSEIF ( ibc_pt_b == 1 )  THEN
+!
+!--       In case of Neumann boundary condition pt_surface is needed for 
+!--       calculation of reference density
+          dum_surf_flux = ( shf_surf(nt) + fac *                               &
+                            ( shf_surf(nt+1) - shf_surf(nt) )                  &
+                          ) * heatflux_input_conversion(nzb)
+!
+!--       Save surface sensible heat flux on default, natural and urban surface
+!--       type, if required 
+          IF ( surf_def_h(0)%ns >= 1 )  surf_def_h(0)%shf(:) = dum_surf_flux
+          IF ( surf_lsm_h%ns    >= 1 )  surf_lsm_h%shf(:)    = dum_surf_flux
+          IF ( surf_usm_h%ns    >= 1 )  surf_usm_h%shf(:)    = dum_surf_flux
+
+          pt_surface    = pt_surf(nt) + fac * ( pt_surf(nt+1) - pt_surf(nt) )
+
+       ENDIF
+
+       IF ( ibc_q_b == 0 )  THEN
+!
+!--       In case of Dirichlet boundary condition qsws must not 
+!--       be set - it is calculated via MOST in prandtl_fluxes
+          q_surface = q_surf(nt) + fac * ( q_surf(nt+1) - q_surf(nt) )
+
+       ELSEIF ( ibc_q_b == 1 )  THEN
+          dum_surf_flux = ( qsws_surf(nt) + fac *                              &
+                             ( qsws_surf(nt+1) - qsws_surf(nt) )               &
+                             ) * waterflux_input_conversion(nzb)
+!
+!--       Save surface sensible heat flux on default, natural and urban surface
+!--       type, if required 
+          IF ( surf_def_h(0)%ns >= 1 )  surf_def_h(0)%qsws(:) = dum_surf_flux
+          IF ( surf_lsm_h%ns    >= 1 )  surf_lsm_h%qsws(:)    = dum_surf_flux
+          IF ( surf_usm_h%ns    >= 1 )  surf_usm_h%qsws(:)    = dum_surf_flux
+
+       ENDIF
+
+       surface_pressure = p_surf(nt) + fac * ( p_surf(nt+1) - p_surf(nt) )
+
+    END SUBROUTINE ls_forcing_surf 
+
+
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> @todo Missing subroutine description.
+!------------------------------------------------------------------------------!
+    SUBROUTINE ls_forcing_vert ( time )
+
+
+       IMPLICIT NONE
+
+       INTEGER(iwp) ::  nt                     !<
+
+       REAL(wp)             ::  fac           !<
+       REAL(wp), INTENT(in) ::  time          !<
+
+!
+!--    Interpolation in time of LSF_DATA for ug, vg and w_subs
+       nt = 1
+       DO WHILE ( time > time_vert(nt) )
+          nt = nt + 1
+       ENDDO
+       IF ( time /= time_vert(nt) )  THEN
+         nt = nt - 1
+       ENDIF
+
+       fac = ( time-time_vert(nt) ) / ( time_vert(nt+1)-time_vert(nt) )
+
+       ug     = ug_vert(:,nt) + fac * ( ug_vert(:,nt+1) - ug_vert(:,nt) )
+       vg     = vg_vert(:,nt) + fac * ( vg_vert(:,nt+1) - vg_vert(:,nt) )
+
+       IF ( large_scale_subsidence )  THEN
+          w_subs = wsubs_vert(:,nt)                                            &
+                   + fac * ( wsubs_vert(:,nt+1) - wsubs_vert(:,nt) )
+       ENDIF
+
+    END SUBROUTINE ls_forcing_vert
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Call for all grid points
+!------------------------------------------------------------------------------!
+    SUBROUTINE ls_advec ( time, prog_var )
+      
+
+       IMPLICIT NONE
+
+       CHARACTER (LEN=*) ::  prog_var   !< 
+
+       REAL(wp), INTENT(in)  :: time    !< 
+       REAL(wp) :: fac                  !<  
+
+       INTEGER(iwp) ::  i               !< 
+       INTEGER(iwp) ::  j               !< 
+       INTEGER(iwp) ::  k               !< 
+       INTEGER(iwp) ::  nt               !< 
+
+!
+!--    Interpolation in time of LSF_DATA 
+       nt = 1
+       DO WHILE ( time > time_vert(nt) )
+          nt = nt + 1
+       ENDDO
+       IF ( time /= time_vert(nt) )  THEN
+         nt = nt - 1
+       ENDIF
+
+       fac = ( time-time_vert(nt) ) / ( time_vert(nt+1)-time_vert(nt) )
+
+!
+!--    Add horizontal large scale advection tendencies of pt and q 
+       SELECT CASE ( prog_var )
+
+          CASE ( 'pt' )
+
+             DO  i = nxl, nxr
+                DO  j = nys, nyn
+                   DO  k = nzb+1, nzt
+                      tend(k,j,i) = tend(k,j,i) + td_lsa_lpt(k,nt) + fac *     &
+                                    ( td_lsa_lpt(k,nt+1) - td_lsa_lpt(k,nt) ) *&
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 0 ) )
+                   ENDDO
+                ENDDO
+             ENDDO
+
+          CASE ( 'q' )
+
+             DO  i = nxl, nxr
+                DO  j = nys, nyn
+                   DO  k = nzb+1, nzt
+                      tend(k,j,i) = tend(k,j,i) + td_lsa_q(k,nt) + fac *       &
+                                    ( td_lsa_q(k,nt+1) - td_lsa_q(k,nt) ) *    &
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 0 ) )
+                   ENDDO
+                ENDDO
+             ENDDO
+
+       END SELECT
+
+!
+!--    Subsidence of pt and q with prescribed subsidence tendencies
+       IF ( large_scale_subsidence .AND. use_subsidence_tendencies )  THEN
+
+          SELECT CASE ( prog_var )
+
+             CASE ( 'pt' )
+
+                DO  i = nxl, nxr
+                   DO  j = nys, nyn
+                      DO  k = nzb+1, nzt
+                         tend(k,j,i) = tend(k,j,i) + td_sub_lpt(k,nt) + fac *  &
+                                     ( td_sub_lpt(k,nt+1) - td_sub_lpt(k,nt) )*&
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 0 ) )
+                      ENDDO
+                   ENDDO
+                ENDDO
+  
+             CASE ( 'q' )
+
+                DO  i = nxl, nxr
+                   DO  j = nys, nyn
+                      DO  k = nzb+1, nzt
+                         tend(k,j,i) = tend(k,j,i) + td_sub_q(k,nt) + fac *    &
+                                       ( td_sub_q(k,nt+1) - td_sub_q(k,nt) ) * &
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 0 ) )
+                      ENDDO
+                   ENDDO
+                ENDDO
+
+          END SELECT
+
+       ENDIF
+
+    END SUBROUTINE ls_advec
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Call for grid point i,j
+!------------------------------------------------------------------------------!
+    SUBROUTINE ls_advec_ij ( i, j, time, prog_var )
+
+       IMPLICIT NONE
+
+       CHARACTER (LEN=*) ::  prog_var   !< 
+
+       REAL(wp), INTENT(in)  :: time    !< 
+       REAL(wp) :: fac                  !< 
+
+       INTEGER(iwp) ::  i               !< 
+       INTEGER(iwp) ::  j               !< 
+       INTEGER(iwp) ::  k               !< 
+       INTEGER(iwp) ::  nt               !< 
+
+!
+!--    Interpolation in time of LSF_DATA 
+       nt = 1
+       DO WHILE ( time > time_vert(nt) )
+          nt = nt + 1
+       ENDDO
+       IF ( time /= time_vert(nt) )  THEN
+         nt = nt - 1
+       ENDIF
+
+       fac = ( time-time_vert(nt) ) / ( time_vert(nt+1)-time_vert(nt) )
+
+!
+!--    Add horizontal large scale advection tendencies of pt and q 
+       SELECT CASE ( prog_var )
+
+          CASE ( 'pt' )
+
+             DO  k = nzb+1, nzt
+                tend(k,j,i) = tend(k,j,i) + td_lsa_lpt(k,nt)                   &
+                             + fac * ( td_lsa_lpt(k,nt+1) - td_lsa_lpt(k,nt) )*&
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 0 ) )
+             ENDDO
+
+          CASE ( 'q' )
+
+             DO  k = nzb+1, nzt
+                tend(k,j,i) = tend(k,j,i) + td_lsa_q(k,nt)                     &
+                              + fac * ( td_lsa_q(k,nt+1) - td_lsa_q(k,nt) ) *  &
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 0 ) )
+             ENDDO
+
+       END SELECT
+
+!
+!--    Subsidence of pt and q with prescribed profiles
+       IF ( large_scale_subsidence .AND. use_subsidence_tendencies )  THEN
+
+          SELECT CASE ( prog_var )
+
+             CASE ( 'pt' )
+
+                DO  k = nzb+1, nzt
+                   tend(k,j,i) = tend(k,j,i) + td_sub_lpt(k,nt) + fac *        &
+                                 ( td_sub_lpt(k,nt+1) - td_sub_lpt(k,nt) ) *   &
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 0 ) )
+                ENDDO
+  
+             CASE ( 'q' )
+
+                DO  k = nzb+1, nzt
+                   tend(k,j,i) = tend(k,j,i) + td_sub_q(k,nt) + fac *          &
+                                 ( td_sub_q(k,nt+1) - td_sub_q(k,nt) ) *       &
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 0 ) )
+                ENDDO
+
+          END SELECT
+
+       ENDIF
+
+    END SUBROUTINE ls_advec_ij
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> @todo Missing subroutine description.
+!------------------------------------------------------------------------------!
+    SUBROUTINE nudge_init
+
+       IMPLICIT NONE
+
+
+       INTEGER(iwp) ::  finput = 90  !<
+       INTEGER(iwp) ::  ierrn        !<
+       INTEGER(iwp) ::  k            !<
+       INTEGER(iwp) ::  nt            !<
+
+       CHARACTER(1) ::  hash     !<
+
+       REAL(wp) ::  highheight   !<
+       REAL(wp) ::  highqnudge   !<
+       REAL(wp) ::  highptnudge  !<
+       REAL(wp) ::  highunudge   !<
+       REAL(wp) ::  highvnudge   !<
+       REAL(wp) ::  highwnudge   !<
+       REAL(wp) ::  hightnudge   !<
+
+       REAL(wp) ::  lowheight    !<
+       REAL(wp) ::  lowqnudge    !<
+       REAL(wp) ::  lowptnudge   !<
+       REAL(wp) ::  lowunudge    !<
+       REAL(wp) ::  lowvnudge    !<
+       REAL(wp) ::  lowwnudge    !<
+       REAL(wp) ::  lowtnudge    !<
+
+       REAL(wp) ::  fac          !<
+
+       ALLOCATE( ptnudge(nzb:nzt+1,1:ntnudge), qnudge(nzb:nzt+1,1:ntnudge), &
+                 tnudge(nzb:nzt+1,1:ntnudge), unudge(nzb:nzt+1,1:ntnudge),  &
+                 vnudge(nzb:nzt+1,1:ntnudge), wnudge(nzb:nzt+1,1:ntnudge)  )
+
+       ALLOCATE( tmp_tnudge(nzb:nzt) )
+
+       ALLOCATE( timenudge(0:ntnudge) )
+
+       ptnudge = 0.0_wp; qnudge = 0.0_wp; tnudge = 0.0_wp; unudge = 0.0_wp
+       vnudge = 0.0_wp; wnudge = 0.0_wp; timenudge = 0.0_wp
+!
+!--    Initialize array tmp_nudge with a current nudging time scale of 6 hours
+       tmp_tnudge = 21600.0_wp
+
+       nt = 0
+       OPEN ( finput, FILE='NUDGING_DATA', STATUS='OLD', &
+              FORM='FORMATTED', IOSTAT=ierrn )
+
+       IF ( ierrn /= 0 )  THEN
+          message_string = 'file NUDGING_DATA does not exist'
+          CALL message( 'nudging', 'PA0365', 1, 2, 0, 6, 0 )
+       ENDIF
+
+       ierrn = 0
+
+ rloop:DO
+          nt = nt + 1
+          hash = "#"
+          ierrn = 1 ! not zero
+!
+!--       Search for the next line consisting of "# time", 
+!--       from there onwards the profiles will be read
+          DO WHILE ( .NOT. ( hash == "#" .AND. ierrn == 0 ) ) 
+          
+            READ ( finput, *, IOSTAT=ierrn ) hash, timenudge(nt)
+            IF ( ierrn < 0 )  EXIT rloop
+
+          ENDDO
+
+          ierrn = 0
+          READ ( finput, *, IOSTAT=ierrn ) lowheight, lowtnudge, lowunudge,   &
+                                           lowvnudge, lowwnudge , lowptnudge, &
+                                           lowqnudge
+
+          IF ( ierrn /= 0 )  THEN
+             message_string = 'errors in file NUDGING_DATA'
+             CALL message( 'nudging', 'PA0366', 1, 2, 0, 6, 0 )
+          ENDIF
+
+          ierrn = 0
+          READ ( finput, *, IOSTAT=ierrn ) highheight, hightnudge, highunudge,   &
+                                           highvnudge, highwnudge , highptnudge, &
+                                           highqnudge
+
+          IF ( ierrn /= 0 )  THEN
+             message_string = 'errors in file NUDGING_DATA'
+             CALL message( 'nudging', 'PA0366', 1, 2, 0, 6, 0 )
+          ENDIF
+
+          DO  k = nzb, nzt+1
+             DO WHILE ( highheight < zu(k) )
+                lowheight  = highheight
+                lowtnudge  = hightnudge
+                lowunudge  = highunudge
+                lowvnudge  = highvnudge
+                lowwnudge  = highwnudge
+                lowptnudge = highptnudge
+                lowqnudge  = highqnudge
+ 
+                ierrn = 0
+                READ ( finput, *, IOSTAT=ierrn )  highheight , hightnudge , &
+                                                  highunudge , highvnudge , &
+                                                  highwnudge , highptnudge, &
+                                                  highqnudge
+                IF (ierrn /= 0 )  THEN
+                   WRITE( message_string, * ) 'zu(',k,') = ', zu(k), 'm is ',&
+                        'higher than the maximum height in NUDING_DATA which ',  &
+                        'is ', lowheight, 'm. Interpolation on PALM ',           &
+                        'grid is not possible.'
+                   CALL message( 'nudging', 'PA0364', 1, 2, 0, 6, 0 )
+                ENDIF
+             ENDDO
+
+!
+!--          Interpolation of prescribed profiles in space 
+
+             fac = ( highheight - zu(k) ) / ( highheight - lowheight )
+
+             tnudge(k,nt)  = fac * lowtnudge  + ( 1.0_wp - fac ) * hightnudge
+             unudge(k,nt)  = fac * lowunudge  + ( 1.0_wp - fac ) * highunudge
+             vnudge(k,nt)  = fac * lowvnudge  + ( 1.0_wp - fac ) * highvnudge
+             wnudge(k,nt)  = fac * lowwnudge  + ( 1.0_wp - fac ) * highwnudge
+             ptnudge(k,nt) = fac * lowptnudge + ( 1.0_wp - fac ) * highptnudge
+             qnudge(k,nt)  = fac * lowqnudge  + ( 1.0_wp - fac ) * highqnudge
+          ENDDO
+
+       ENDDO rloop
+
+       CLOSE ( finput )
+
+!
+!--    Overwrite initial profiles in case of nudging
+       IF ( nudging )  THEN
+          pt_init = ptnudge(:,1)
+          u_init  = unudge(:,1)
+          v_init  = vnudge(:,1)
+          IF ( humidity  )  THEN ! is passive_scalar correct???
+             q_init = qnudge(:,1)
+          ENDIF
+
+          WRITE( message_string, * ) 'Initial profiles of u, v and ',          &
+                                     'scalars from NUDGING_DATA are used.'
+          CALL message( 'large_scale_forcing_nudging', 'PA0370', 0, 0, 0, 6, 0 )
+       ENDIF
+
+
+    END SUBROUTINE nudge_init
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> @todo Missing subroutine description.
+!------------------------------------------------------------------------------!
+    SUBROUTINE calc_tnudge ( time )
+
+       IMPLICIT NONE
+
+
+       REAL(wp) ::  dtm         !<
+       REAL(wp) ::  dtp         !<
+       REAL(wp) ::  time        !<
+
+       INTEGER(iwp) ::  k   !<
+       INTEGER(iwp) ::  nt  !<
+
+       nt = 1
+       DO WHILE ( time > timenudge(nt) )
+         nt = nt+1
+       ENDDO
+       IF ( time /= timenudge(1) ) THEN
+         nt = nt-1
+       ENDIF
+
+       dtm = ( time - timenudge(nt) ) / ( timenudge(nt+1) - timenudge(nt) )
+       dtp = ( timenudge(nt+1) - time ) / ( timenudge(nt+1) - timenudge(nt) )
+
+       DO  k = nzb, nzt
+          tmp_tnudge(k) = MAX( dt_3d, tnudge(k,nt) * dtp + tnudge(k,nt+1) * dtm )
+       ENDDO
+
+    END SUBROUTINE calc_tnudge
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Call for all grid points
+!------------------------------------------------------------------------------!
+    SUBROUTINE nudge ( time, prog_var )
+
+       IMPLICIT NONE
+
+       CHARACTER (LEN=*) ::  prog_var  !<
+
+       REAL(wp) ::  tmp_tend    !<
+       REAL(wp) ::  dtm         !<
+       REAL(wp) ::  dtp         !<
+       REAL(wp) ::  time        !<
+
+       INTEGER(iwp) ::  i  !<
+       INTEGER(iwp) ::  j  !<
+       INTEGER(iwp) ::  k  !<
+       INTEGER(iwp) ::  nt  !<
+
+
+       nt = 1
+       DO WHILE ( time > timenudge(nt) )
+         nt = nt+1
+       ENDDO
+       IF ( time /= timenudge(1) ) THEN
+         nt = nt-1
+       ENDIF
+
+       dtm = ( time - timenudge(nt) ) / ( timenudge(nt+1) - timenudge(nt) )
+       dtp = ( timenudge(nt+1) - time ) / ( timenudge(nt+1) - timenudge(nt) )
+
+       SELECT CASE ( prog_var )
+
+          CASE ( 'u' )
+
+             DO  i = nxl, nxr
+                DO  j = nys, nyn
+
+                   DO  k = nzb+1, nzt
+
+                      tmp_tend = - ( hom(k,1,1,0) - ( unudge(k,nt) * dtp +     &
+                                     unudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
+
+                      tend(k,j,i) = tend(k,j,i) + tmp_tend *                   &
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 1 ) )
+
+                      sums_ls_l(k,6) = sums_ls_l(k,6) + tmp_tend *             &
+                                     weight_substep(intermediate_timestep_count)
+                   ENDDO
+                  
+                   sums_ls_l(nzt+1,6) = sums_ls_l(nzt,6)
+ 
+                ENDDO
+            ENDDO
+
+          CASE ( 'v' )
+
+             DO  i = nxl, nxr
+                DO  j = nys, nyn
+
+                   DO  k = nzb+1, nzt
+
+                      tmp_tend = - ( hom(k,1,2,0) - ( vnudge(k,nt) * dtp +     &
+                                     vnudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
+
+                      tend(k,j,i) = tend(k,j,i) + tmp_tend *                   &
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 2 ) )
+
+                      sums_ls_l(k,7) = sums_ls_l(k,7) + tmp_tend *             &
+                                     weight_substep(intermediate_timestep_count)
+                   ENDDO
+                  
+                   sums_ls_l(nzt+1,7) = sums_ls_l(nzt,7)
+
+                ENDDO
+            ENDDO
+
+          CASE ( 'pt' )
+
+             DO  i = nxl, nxr
+                DO  j = nys, nyn
+
+                   DO  k = nzb+1, nzt
+
+                      tmp_tend = - ( hom(k,1,4,0) - ( ptnudge(k,nt) * dtp +    &
+                                     ptnudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
+
+                      tend(k,j,i) = tend(k,j,i) + tmp_tend *                   &
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 0 ) )
+
+                      sums_ls_l(k,4) = sums_ls_l(k,4) + tmp_tend *             &
+                                     weight_substep(intermediate_timestep_count)
+                   ENDDO
+
+                   sums_ls_l(nzt+1,4) = sums_ls_l(nzt,4)
+
+                ENDDO
+            ENDDO
+
+          CASE ( 'q' )
+
+             DO  i = nxl, nxr
+                DO  j = nys, nyn
+
+                   DO  k = nzb+1, nzt
+
+                      tmp_tend = - ( hom(k,1,41,0) - ( qnudge(k,nt) * dtp +    &
+                                     qnudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
+
+                      tend(k,j,i) = tend(k,j,i) + tmp_tend *                   &
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 0 ) )
+
+                      sums_ls_l(k,5) = sums_ls_l(k,5) + tmp_tend *             &
+                                     weight_substep(intermediate_timestep_count)
+                   ENDDO
+                  
+                   sums_ls_l(nzt+1,5) = sums_ls_l(nzt,5)
+
+                ENDDO
+            ENDDO
+
+          CASE DEFAULT
+             message_string = 'unknown prognostic variable "' // prog_var // '"'
+             CALL message( 'nudge', 'PA0367', 1, 2, 0, 6, 0 )
+
+       END SELECT
+
+    END SUBROUTINE nudge
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Call for grid point i,j
+!------------------------------------------------------------------------------!
+
+    SUBROUTINE nudge_ij( i, j, time, prog_var )
+
+       IMPLICIT NONE
+
+
+       CHARACTER (LEN=*) ::  prog_var  !<
+
+       REAL(wp) ::  tmp_tend    !<
+       REAL(wp) ::  dtm         !<
+       REAL(wp) ::  dtp         !<
+       REAL(wp) ::  time        !<
+
+       INTEGER(iwp) ::  i  !<
+       INTEGER(iwp) ::  j  !<
+       INTEGER(iwp) ::  k  !<
+       INTEGER(iwp) ::  nt  !<
+
+
+       nt = 1
+       DO WHILE ( time > timenudge(nt) )
+         nt = nt+1
+       ENDDO
+       IF ( time /= timenudge(1) )  THEN
+         nt = nt-1
+       ENDIF
+
+       dtm = ( time - timenudge(nt) ) / ( timenudge(nt+1) - timenudge(nt) )
+       dtp = ( timenudge(nt+1) - time ) / ( timenudge(nt+1) - timenudge(nt) )
+
+       SELECT CASE ( prog_var )
+
+          CASE ( 'u' )
+
+             DO  k = nzb+1, nzt
+
+                tmp_tend = - ( hom(k,1,1,0) - ( unudge(k,nt) * dtp +           &
+                               unudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
+
+                tend(k,j,i) = tend(k,j,i) + tmp_tend *                         &
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 1 ) )
+
+                sums_ls_l(k,6) = sums_ls_l(k,6) + tmp_tend                     &
+                                 * weight_substep(intermediate_timestep_count)
+             ENDDO
+
+             sums_ls_l(nzt+1,6) = sums_ls_l(nzt,6)
+
+          CASE ( 'v' )
+
+             DO  k = nzb+1, nzt
+
+                tmp_tend = - ( hom(k,1,2,0) - ( vnudge(k,nt) * dtp +           &
+                               vnudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
+
+                tend(k,j,i) = tend(k,j,i) + tmp_tend *                         &
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 2 ) )
+
+                sums_ls_l(k,7) = sums_ls_l(k,7) + tmp_tend                     &
+                                 * weight_substep(intermediate_timestep_count)
+             ENDDO
+
+             sums_ls_l(nzt+1,7) = sums_ls_l(nzt,7)
+
+          CASE ( 'pt' )
+
+             DO  k = nzb+1, nzt
+
+                tmp_tend = - ( hom(k,1,4,0) - ( ptnudge(k,nt) * dtp +          &
+                               ptnudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
+
+                tend(k,j,i) = tend(k,j,i) + tmp_tend *                         &
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 0 ) )
+
+                sums_ls_l(k,4) = sums_ls_l(k,4) + tmp_tend                     &
+                                 * weight_substep(intermediate_timestep_count)
+             ENDDO
+
+             sums_ls_l(nzt+1,4) = sums_ls_l(nzt,4)
+
+
+          CASE ( 'q' )
+
+             DO  k = nzb+1, nzt
+
+                tmp_tend = - ( hom(k,1,41,0) - ( qnudge(k,nt) * dtp +          &
+                               qnudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
+
+                tend(k,j,i) = tend(k,j,i) + tmp_tend *                         &
+                                        MERGE( 1.0_wp, 0.0_wp,                 &
+                                               BTEST( wall_flags_0(k,j,i), 0 ) )
+
+                sums_ls_l(k,5) = sums_ls_l(k,5) + tmp_tend                     &
+                                 * weight_substep(intermediate_timestep_count)
+             ENDDO
+
+             sums_ls_l(nzt+1,5) = sums_ls_l(nzt,5)
+
+          CASE DEFAULT
+             message_string = 'unknown prognostic variable "' // prog_var // '"'
+             CALL message( 'nudge', 'PA0367', 1, 2, 0, 6, 0 )
+
+       END SELECT
+
+
+    END SUBROUTINE nudge_ij
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> @todo Missing subroutine description.
+!------------------------------------------------------------------------------!
+    SUBROUTINE nudge_ref ( time )
+
+       IMPLICIT NONE
+
+       INTEGER(iwp) ::  nt                    !<
+
+       REAL(wp)             ::  fac           !<
+       REAL(wp), INTENT(in) ::  time          !<
+
+!
+!--    Interpolation in time of NUDGING_DATA for pt_init and q_init. This is 
+!--    needed for correct upper boundary conditions for pt and q and in case that 
+!      large scale subsidence as well as scalar Rayleigh-damping are used
+       nt = 1
+       DO WHILE ( time > time_vert(nt) )
+          nt = nt + 1
+       ENDDO
+       IF ( time /= time_vert(nt) )  THEN
+        nt = nt - 1
+       ENDIF
+
+       fac = ( time-time_vert(nt) ) / ( time_vert(nt+1)-time_vert(nt) )
+
+       pt_init = ptnudge(:,nt) + fac * ( ptnudge(:,nt+1) - ptnudge(:,nt) )
+       q_init  = qnudge(:,nt) + fac * ( qnudge(:,nt+1) - qnudge(:,nt) )
+       u_init  = unudge(:,nt) + fac * ( unudge(:,nt+1) - unudge(:,nt) )
+       v_init  = vnudge(:,nt) + fac * ( vnudge(:,nt+1) - vnudge(:,nt) )
+
+    END SUBROUTINE nudge_ref
+
+ END MODULE lsf_nudging_mod
Index: lm/trunk/SOURCE/ls_forcing_mod.f90
===================================================================
--- /palm/trunk/SOURCE/ls_forcing_mod.f90	(revision 2319)
+++ 	(revision )
@@ -1,704 +1,0 @@
-!> @file ls_forcing_mod.f90
-!------------------------------------------------------------------------------!
-! This file is part of PALM.
-!
-! PALM is free software: you can redistribute it and/or modify it under the
-! terms of the GNU General Public License as published by the Free Software
-! Foundation, either version 3 of the License, or (at your option) any later
-! version.
-!
-! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
-! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
-! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
-!
-! You should have received a copy of the GNU General Public License along with
-! PALM. If not, see <http://www.gnu.org/licenses/>.
-!
-! Copyright 1997-2017 Leibniz Universitaet Hannover
-!------------------------------------------------------------------------------!
-!
-! Current revisions:
-! ------------------
-! 
-! 
-! Former revisions:
-! -----------------
-! $Id$
-! Error messages changed
-! 
-! 2233 2017-05-30 18:08:54Z suehring
-!
-! 2232 2017-05-30 17:47:52Z suehring
-! Adopt to new topography structure, even though no well-conceived topography 
-! concept concerning nudging and large-scale for exist so far. 
-!
-! Also adopt to new surface-structure, i.e. fluxes are obtained from data-types 
-! 
-! 2104 2017-01-06 16:01:15Z knoop
-! Bugfix for approximation related flux input conversion
-! 
-! 2071 2016-11-17 11:22:14Z maronga
-! Small bugfix (Resler)
-! 
-! 2037 2016-10-26 11:15:40Z knoop
-! Anelastic approximation implemented
-! 
-! 2000 2016-08-20 18:09:15Z knoop
-! Forced header and separation lines into 80 columns
-! 
-! 1850 2016-04-08 13:29:27Z maronga
-! Module renamed
-! 
-! 
-! 1682 2015-10-07 23:56:08Z knoop
-! Code annotations made doxygen readable
-! 
-! 1602 2015-06-22 07:50:56Z heinze 
-! PA0370 changed to PA0363
-!
-! 1382 2014-04-30 12:15:41Z boeske
-! Renamed variables which store large scale forcing tendencies
-! pt_lsa -> td_lsa_lpt, pt_subs -> td_sub_lpt, 
-! q_lsa  -> td_lsa_q,   q_subs  -> td_sub_q,
-! high|lowpt_lsa -> high|low_td_lsa_lpt, ...
-! 
-! 1365 2014-04-22 15:03:56Z boeske
-! Usage of large scale forcing for pt and q enabled:
-! Added new subroutine ls_advec for horizontal large scale advection and large 
-! scale subsidence, 
-! error message in init_ls_forcing specified,
-! variable t renamed nt
-! 
-! 1353 2014-04-08 15:21:23Z heinze
-! REAL constants provided with KIND-attribute
-!
-! 1320 2014-03-20 08:40:49Z raasch
-! ONLY-attribute added to USE-statements,
-! kind-parameters added to all INTEGER and REAL declaration statements,
-! kinds are defined in new module kinds,
-! comment fields (!:) to be used for variable explanations added to
-! all variable declaration statements
-!
-! 1318 2014-03-17 13:35:16Z raasch
-! module interfaces removed
-!
-! 1299 2014-03-06 13:15:21Z heinze
-! Ensure a zero large scale vertical velocity at the surface
-! Bugfix: typo in case of boundary condition in if-clause
-!
-! 1276 2014-01-15 13:40:41Z heinze
-! Use LSF_DATA also in case of Dirichlet bottom boundary condition for scalars
-!
-! 1249 2013-11-06 10:45:47Z heinze
-! remove call of user module
-! reformatting
-!
-! 1241 2013-10-30 11:36:58Z heinze
-! Initial revision
-!
-! Description:
-! ------------
-!> Calculates large scale forcings (geostrophic wind and subsidence velocity) as 
-!> well as surfaces fluxes dependent on time given in an external file (LSF_DATA).
-!> Code is based in parts on DALES and UCLA-LES.
-!--------------------------------------------------------------------------------!
- MODULE ls_forcing_mod
- 
-
-    PRIVATE
-    PUBLIC init_ls_forcing, ls_forcing_surf, ls_forcing_vert, ls_advec
-    SAVE
-
-    INTERFACE ls_advec
-       MODULE PROCEDURE ls_advec
-       MODULE PROCEDURE ls_advec_ij
-    END INTERFACE ls_advec
-
- CONTAINS
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> @todo Missing subroutine description.
-!------------------------------------------------------------------------------!
-    SUBROUTINE init_ls_forcing
-
-       USE arrays_3d,                                                          &
-           ONLY:  p_surf, pt_surf, q_surf, qsws_surf, shf_surf, td_lsa_lpt,    &
-                  td_lsa_q, td_sub_lpt, td_sub_q, time_surf, time_vert,        &
-                  ug_vert, vg_vert, wsubs_vert, zu
-
-       USE control_parameters,                                                 &
-           ONLY:  end_time, lsf_surf, lsf_vert, message_string, nlsf
-
-       USE indices,                                                            &
-           ONLY:  ngp_sums_ls, nzb, nz, nzt
-
-       USE kinds
-
-       USE statistics,                                                         &
-           ONLY:  sums_ls_l
-
-
-       IMPLICIT NONE
-
-       CHARACTER(100) ::  chmess      !<
-       CHARACTER(1)   ::  hash        !<
-
-       INTEGER(iwp) ::  ierrn         !<
-       INTEGER(iwp) ::  finput = 90   !<
-       INTEGER(iwp) ::  k             !<
-       INTEGER(iwp) ::  nt             !<
-
-       REAL(wp) ::  fac               !<
-       REAL(wp) ::  highheight        !<
-       REAL(wp) ::  highug_vert       !<
-       REAL(wp) ::  highvg_vert       !<
-       REAL(wp) ::  highwsubs_vert    !<
-       REAL(wp) ::  lowheight         !<
-       REAL(wp) ::  lowug_vert        !<
-       REAL(wp) ::  lowvg_vert        !<
-       REAL(wp) ::  lowwsubs_vert     !<
-       REAL(wp) ::  high_td_lsa_lpt   !<
-       REAL(wp) ::  low_td_lsa_lpt    !<
-       REAL(wp) ::  high_td_lsa_q     !<
-       REAL(wp) ::  low_td_lsa_q      !<
-       REAL(wp) ::  high_td_sub_lpt   !<
-       REAL(wp) ::  low_td_sub_lpt    !<
-       REAL(wp) ::  high_td_sub_q     !<
-       REAL(wp) ::  low_td_sub_q      !<
-       REAL(wp) ::  r_dummy           !<
-
-       ALLOCATE( p_surf(0:nlsf), pt_surf(0:nlsf), q_surf(0:nlsf),              &
-                 qsws_surf(0:nlsf), shf_surf(0:nlsf),                          &
-                 td_lsa_lpt(nzb:nzt+1,0:nlsf), td_lsa_q(nzb:nzt+1,0:nlsf),     &
-                 td_sub_lpt(nzb:nzt+1,0:nlsf), td_sub_q(nzb:nzt+1,0:nlsf),     &
-                 time_vert(0:nlsf), time_surf(0:nlsf),                         &
-                 ug_vert(nzb:nzt+1,0:nlsf), vg_vert(nzb:nzt+1,0:nlsf),         &
-                 wsubs_vert(nzb:nzt+1,0:nlsf) )
-
-       p_surf = 0.0_wp; pt_surf = 0.0_wp; q_surf = 0.0_wp; qsws_surf = 0.0_wp
-       shf_surf = 0.0_wp; time_vert = 0.0_wp; td_lsa_lpt = 0.0_wp
-       td_lsa_q = 0.0_wp; td_sub_lpt = 0.0_wp; td_sub_q = 0.0_wp
-       time_surf = 0.0_wp; ug_vert = 0.0_wp; vg_vert = 0.0_wp
-       wsubs_vert = 0.0_wp
-
-!
-!--    Array for storing large scale forcing and nudging tendencies at each 
-!--    timestep for data output
-       ALLOCATE( sums_ls_l(nzb:nzt+1,0:7) )
-       sums_ls_l = 0.0_wp
-
-       ngp_sums_ls = (nz+2)*6
-
-       OPEN ( finput, FILE='LSF_DATA', STATUS='OLD', &
-              FORM='FORMATTED', IOSTAT=ierrn )
-
-       IF ( ierrn /= 0 )  THEN
-          message_string = 'file LSF_DATA does not exist'
-          CALL message( 'ls_forcing', 'PA0368', 1, 2, 0, 6, 0 )
-       ENDIF
-
-       ierrn = 0
-!
-!--    First three lines of LSF_DATA contain header
-       READ ( finput, FMT='(a100)', IOSTAT=ierrn ) chmess
-       READ ( finput, FMT='(a100)', IOSTAT=ierrn ) chmess
-       READ ( finput, FMT='(a100)', IOSTAT=ierrn ) chmess
-
-       IF ( ierrn /= 0 )  THEN
-          message_string = 'errors in file LSF_DATA'
-          CALL message( 'ls_forcing', 'PA0369', 1, 2, 0, 6, 0 )
-       ENDIF
-
-!
-!--    Surface values are read in
-       nt     = 0
-       ierrn = 0
-
-       DO WHILE ( time_surf(nt) < end_time )
-          nt = nt + 1
-          READ ( finput, *, IOSTAT = ierrn ) time_surf(nt), shf_surf(nt),      &
-                                             qsws_surf(nt), pt_surf(nt),       &
-                                             q_surf(nt), p_surf(nt)
-
-          IF ( ierrn < 0 )  THEN
-            WRITE ( message_string, * ) 'No time dependent surface variables ',&
-                              'in&LSF_DATA for end of run found'
-
-             CALL message( 'ls_forcing', 'PA0363', 1, 2, 0, 6, 0 )
-          ENDIF
-       ENDDO
-
-       IF ( time_surf(1) > end_time )  THEN
-          WRITE ( message_string, * ) 'Time dependent surface variables in ',  &
-                                      '&LSF_DATA set in after end of ' ,      &
-                                      'simulation - lsf_surf is set to FALSE'
-          CALL message( 'ls_forcing', 'PA0371', 0, 0, 0, 6, 0 )
-          lsf_surf = .FALSE.
-       ENDIF
-
-!
-!--    Go to the end of the list with surface variables
-       DO WHILE ( ierrn == 0 )
-          READ ( finput, *, IOSTAT = ierrn ) r_dummy
-       ENDDO
-
-!
-!--    Profiles of ug, vg and w_subs are read in (large scale forcing)
-
-       nt = 0
-       DO WHILE ( time_vert(nt) < end_time )
-          nt = nt + 1
-          hash = "#"
-          ierrn = 1 ! not zero
-!
-!--       Search for the next line consisting of "# time", 
-!--       from there onwards the profiles will be read
-          DO WHILE ( .NOT. ( hash == "#" .AND. ierrn == 0 ) ) 
-             READ ( finput, *, IOSTAT=ierrn ) hash, time_vert(nt)
-             IF ( ierrn < 0 )  THEN 
-                WRITE( message_string, * ) 'No time dependent vertical profiles',&
-                                 ' in&LSF_DATA for end of run found'
-                CALL message( 'ls_forcing', 'PA0372', 1, 2, 0, 6, 0 )
-             ENDIF
-          ENDDO
-
-          IF ( nt == 1 .AND. time_vert(nt) > end_time ) EXIT
-
-          READ ( finput, *, IOSTAT=ierrn ) lowheight, lowug_vert, lowvg_vert,  &
-                                           lowwsubs_vert, low_td_lsa_lpt,      &
-                                           low_td_lsa_q, low_td_sub_lpt,       &
-                                           low_td_sub_q
-          IF ( ierrn /= 0 )  THEN
-             message_string = 'errors in file LSF_DATA'
-             CALL message( 'ls_forcing', 'PA0369', 1, 2, 0, 6, 0 )
-          ENDIF
-
-          READ ( finput, *, IOSTAT=ierrn ) highheight, highug_vert,            &
-                                           highvg_vert, highwsubs_vert,        &
-                                           high_td_lsa_lpt, high_td_lsa_q,     &
-                                           high_td_sub_lpt, high_td_sub_q
-      
-          IF ( ierrn /= 0 )  THEN
-             message_string = 'errors in file LSF_DATA'
-             CALL message( 'ls_forcing', 'PA0369', 1, 2, 0, 6, 0 )
-          ENDIF
-
-
-          DO  k = nzb, nzt+1
-             IF ( highheight < zu(k) )  THEN
-                lowheight      = highheight
-                lowug_vert     = highug_vert
-                lowvg_vert     = highvg_vert
-                lowwsubs_vert  = highwsubs_vert
-                low_td_lsa_lpt = high_td_lsa_lpt
-                low_td_lsa_q   = high_td_lsa_q
-                low_td_sub_lpt = high_td_sub_lpt
-                low_td_sub_q   = high_td_sub_q
-
-                ierrn = 0
-                READ ( finput, *, IOSTAT=ierrn ) highheight, highug_vert,      &
-                                                 highvg_vert, highwsubs_vert,  &
-                                                 high_td_lsa_lpt,              &
-                                                 high_td_lsa_q,                &
-                                                 high_td_sub_lpt, high_td_sub_q
-
-                IF ( ierrn /= 0 )  THEN
-                   WRITE( message_string, * ) 'zu(',k,') = ', zu(k), 'm ', &
-                        'is higher than the maximum height in LSF_DATA which ',&
-                        'is ', lowheight, 'm. Interpolation on PALM ',         &
-                        'grid is not possible.'
-                   CALL message( 'ls_forcing', 'PA0395', 1, 2, 0, 6, 0 )
-                ENDIF
-
-             ENDIF
-
-!
-!--          Interpolation of prescribed profiles in space 
-             fac = (highheight-zu(k))/(highheight - lowheight)
-
-             ug_vert(k,nt)    = fac * lowug_vert                               &
-                                + ( 1.0_wp - fac ) * highug_vert
-             vg_vert(k,nt)    = fac * lowvg_vert                               &
-                                + ( 1.0_wp - fac ) * highvg_vert
-             wsubs_vert(k,nt) = fac * lowwsubs_vert                            &
-                                + ( 1.0_wp - fac ) * highwsubs_vert
-
-             td_lsa_lpt(k,nt) = fac * low_td_lsa_lpt                           &
-                                + ( 1.0_wp - fac ) * high_td_lsa_lpt
-             td_lsa_q(k,nt)   = fac * low_td_lsa_q                             &
-                                + ( 1.0_wp - fac ) * high_td_lsa_q
-             td_sub_lpt(k,nt) = fac * low_td_sub_lpt                           &
-                                + ( 1.0_wp - fac ) * high_td_sub_lpt
-             td_sub_q(k,nt)   = fac * low_td_sub_q                             &
-                                + ( 1.0_wp - fac ) * high_td_sub_q
-
-          ENDDO
-
-       ENDDO 
-
-!
-!--    Large scale vertical velocity has to be zero at the surface
-       wsubs_vert(nzb,:) = 0.0_wp
- 
-       IF ( time_vert(1) > end_time )  THEN
-          WRITE ( message_string, * ) 'Time dependent large scale profile ',   &
-                             'forcing from&LSF_DATA sets in after end of ' ,   &
-                             'simulation - lsf_vert is set to FALSE'
-          CALL message( 'ls_forcing', 'PA0373', 0, 0, 0, 6, 0 )
-          lsf_vert = .FALSE.
-       ENDIF
-
-       CLOSE( finput )
-
-
-    END SUBROUTINE init_ls_forcing 
-
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> @todo Missing subroutine description.
-!------------------------------------------------------------------------------!
-    SUBROUTINE ls_forcing_surf ( time )
-
-       USE arrays_3d,                                                          &
-           ONLY:  p_surf, pt_surf, q_surf, qsws_surf, shf_surf,                &
-                  heatflux_input_conversion, waterflux_input_conversion,       &
-                  time_surf, time_vert, ug, ug_vert, vg, vg_vert
-
-       USE control_parameters,                                                 &
-           ONLY:  bc_q_b, ibc_pt_b, ibc_q_b, pt_surface, q_surface,            &
-                  surface_pressure
-
-       USE indices,                                                            &
-           ONLY:  nzb
-
-       USE kinds
-
-       USE surface_mod,                                                        &
-           ONLY:  surf_def_h, surf_lsm_h, surf_usm_h
-
-       IMPLICIT NONE
-
-       INTEGER(iwp) ::  nt                     !<
-
-       REAL(wp)             :: dum_surf_flux  !<
-       REAL(wp)             :: fac            !<
-       REAL(wp), INTENT(in) :: time           !<
-
-!
-!--    Interpolation in time of LSF_DATA at the surface
-       nt = 1
-       DO WHILE ( time > time_surf(nt) )
-          nt = nt + 1
-       ENDDO
-       IF ( time /= time_surf(nt) )  THEN
-         nt = nt - 1
-       ENDIF
-
-       fac = ( time -time_surf(nt) ) / ( time_surf(nt+1) - time_surf(nt) )
-
-       IF ( ibc_pt_b == 0 )  THEN
-!
-!--       In case of Dirichlet boundary condition shf must not 
-!--       be set - it is calculated via MOST in prandtl_fluxes
-          pt_surface = pt_surf(nt) + fac * ( pt_surf(nt+1) - pt_surf(nt) )
-
-       ELSEIF ( ibc_pt_b == 1 )  THEN
-!
-!--       In case of Neumann boundary condition pt_surface is needed for 
-!--       calculation of reference density
-          dum_surf_flux = ( shf_surf(nt) + fac *                               &
-                            ( shf_surf(nt+1) - shf_surf(nt) )                  &
-                          ) * heatflux_input_conversion(nzb)
-!
-!--       Save surface sensible heat flux on default, natural and urban surface
-!--       type, if required 
-          IF ( surf_def_h(0)%ns >= 1 )  surf_def_h(0)%shf(:) = dum_surf_flux
-          IF ( surf_lsm_h%ns    >= 1 )  surf_lsm_h%shf(:)    = dum_surf_flux
-          IF ( surf_usm_h%ns    >= 1 )  surf_usm_h%shf(:)    = dum_surf_flux
-
-          pt_surface    = pt_surf(nt) + fac * ( pt_surf(nt+1) - pt_surf(nt) )
-
-       ENDIF
-
-       IF ( ibc_q_b == 0 )  THEN
-!
-!--       In case of Dirichlet boundary condition qsws must not 
-!--       be set - it is calculated via MOST in prandtl_fluxes
-          q_surface = q_surf(nt) + fac * ( q_surf(nt+1) - q_surf(nt) )
-
-       ELSEIF ( ibc_q_b == 1 )  THEN
-          dum_surf_flux = ( qsws_surf(nt) + fac *                              &
-                             ( qsws_surf(nt+1) - qsws_surf(nt) )               &
-                             ) * waterflux_input_conversion(nzb)
-!
-!--       Save surface sensible heat flux on default, natural and urban surface
-!--       type, if required 
-          IF ( surf_def_h(0)%ns >= 1 )  surf_def_h(0)%qsws(:) = dum_surf_flux
-          IF ( surf_lsm_h%ns    >= 1 )  surf_lsm_h%qsws(:)    = dum_surf_flux
-          IF ( surf_usm_h%ns    >= 1 )  surf_usm_h%qsws(:)    = dum_surf_flux
-
-       ENDIF
-
-       surface_pressure = p_surf(nt) + fac * ( p_surf(nt+1) - p_surf(nt) )
-
-    END SUBROUTINE ls_forcing_surf 
-
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> @todo Missing subroutine description.
-!------------------------------------------------------------------------------!
-    SUBROUTINE ls_forcing_vert ( time )
-
-       USE arrays_3d,                                                          &
-           ONLY:  time_vert, ug, ug_vert, vg, vg_vert, w_subs, wsubs_vert
-
-       USE control_parameters,                                                 &
-           ONLY:  large_scale_subsidence
-
-       USE kinds
-
-       IMPLICIT NONE
-
-       INTEGER(iwp) ::  nt                     !<
-
-       REAL(wp)             ::  fac           !<
-       REAL(wp), INTENT(in) ::  time          !<
-
-!
-!--    Interpolation in time of LSF_DATA for ug, vg and w_subs
-       nt = 1
-       DO WHILE ( time > time_vert(nt) )
-          nt = nt + 1
-       ENDDO
-       IF ( time /= time_vert(nt) )  THEN
-         nt = nt - 1
-       ENDIF
-
-       fac = ( time-time_vert(nt) ) / ( time_vert(nt+1)-time_vert(nt) )
-
-       ug     = ug_vert(:,nt) + fac * ( ug_vert(:,nt+1) - ug_vert(:,nt) )
-       vg     = vg_vert(:,nt) + fac * ( vg_vert(:,nt+1) - vg_vert(:,nt) )
-
-       IF ( large_scale_subsidence )  THEN
-          w_subs = wsubs_vert(:,nt)                                            &
-                   + fac * ( wsubs_vert(:,nt+1) - wsubs_vert(:,nt) )
-       ENDIF
-
-    END SUBROUTINE ls_forcing_vert
-
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> Call for all grid points
-!------------------------------------------------------------------------------!
-    SUBROUTINE ls_advec ( time, prog_var )
-
-       USE arrays_3d,                                                          &
-           ONLY:  td_lsa_lpt, td_lsa_q, td_sub_lpt, td_sub_q, tend, time_vert        
-       
-       USE control_parameters,                                                 &
-           ONLY:  large_scale_subsidence, use_subsidence_tendencies
-       
-       USE indices
-       
-       USE kinds
-
-       IMPLICIT NONE
-
-       CHARACTER (LEN=*) ::  prog_var   !< 
-
-       REAL(wp), INTENT(in)  :: time    !< 
-       REAL(wp) :: fac                  !<  
-
-       INTEGER(iwp) ::  i               !< 
-       INTEGER(iwp) ::  j               !< 
-       INTEGER(iwp) ::  k               !< 
-       INTEGER(iwp) ::  nt               !< 
-
-!
-!--    Interpolation in time of LSF_DATA 
-       nt = 1
-       DO WHILE ( time > time_vert(nt) )
-          nt = nt + 1
-       ENDDO
-       IF ( time /= time_vert(nt) )  THEN
-         nt = nt - 1
-       ENDIF
-
-       fac = ( time-time_vert(nt) ) / ( time_vert(nt+1)-time_vert(nt) )
-
-!
-!--    Add horizontal large scale advection tendencies of pt and q 
-       SELECT CASE ( prog_var )
-
-          CASE ( 'pt' )
-
-             DO  i = nxl, nxr
-                DO  j = nys, nyn
-                   DO  k = nzb+1, nzt
-                      tend(k,j,i) = tend(k,j,i) + td_lsa_lpt(k,nt) + fac *     &
-                                    ( td_lsa_lpt(k,nt+1) - td_lsa_lpt(k,nt) ) *&
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 0 ) )
-                   ENDDO
-                ENDDO
-             ENDDO
-
-          CASE ( 'q' )
-
-             DO  i = nxl, nxr
-                DO  j = nys, nyn
-                   DO  k = nzb+1, nzt
-                      tend(k,j,i) = tend(k,j,i) + td_lsa_q(k,nt) + fac *       &
-                                    ( td_lsa_q(k,nt+1) - td_lsa_q(k,nt) ) *    &
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 0 ) )
-                   ENDDO
-                ENDDO
-             ENDDO
-
-       END SELECT
-
-!
-!--    Subsidence of pt and q with prescribed subsidence tendencies
-       IF ( large_scale_subsidence .AND. use_subsidence_tendencies )  THEN
-
-          SELECT CASE ( prog_var )
-
-             CASE ( 'pt' )
-
-                DO  i = nxl, nxr
-                   DO  j = nys, nyn
-                      DO  k = nzb+1, nzt
-                         tend(k,j,i) = tend(k,j,i) + td_sub_lpt(k,nt) + fac *  &
-                                     ( td_sub_lpt(k,nt+1) - td_sub_lpt(k,nt) )*&
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 0 ) )
-                      ENDDO
-                   ENDDO
-                ENDDO
-  
-             CASE ( 'q' )
-
-                DO  i = nxl, nxr
-                   DO  j = nys, nyn
-                      DO  k = nzb+1, nzt
-                         tend(k,j,i) = tend(k,j,i) + td_sub_q(k,nt) + fac *    &
-                                       ( td_sub_q(k,nt+1) - td_sub_q(k,nt) ) * &
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 0 ) )
-                      ENDDO
-                   ENDDO
-                ENDDO
-
-          END SELECT
-
-       ENDIF
-
-    END SUBROUTINE ls_advec
-
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> Call for grid point i,j
-!------------------------------------------------------------------------------!
-    SUBROUTINE ls_advec_ij ( i, j, time, prog_var )
-
-       USE arrays_3d,                                                          &
-           ONLY:  td_lsa_lpt, td_lsa_q, td_sub_lpt, td_sub_q, tend, time_vert        
-       
-       USE control_parameters,                                                 &
-           ONLY:  large_scale_subsidence, use_subsidence_tendencies
-       
-       USE indices
-       
-       USE kinds
-
-       IMPLICIT NONE
-
-       CHARACTER (LEN=*) ::  prog_var   !< 
-
-       REAL(wp), INTENT(in)  :: time    !< 
-       REAL(wp) :: fac                  !< 
-
-       INTEGER(iwp) ::  i               !< 
-       INTEGER(iwp) ::  j               !< 
-       INTEGER(iwp) ::  k               !< 
-       INTEGER(iwp) ::  nt               !< 
-
-!
-!--    Interpolation in time of LSF_DATA 
-       nt = 1
-       DO WHILE ( time > time_vert(nt) )
-          nt = nt + 1
-       ENDDO
-       IF ( time /= time_vert(nt) )  THEN
-         nt = nt - 1
-       ENDIF
-
-       fac = ( time-time_vert(nt) ) / ( time_vert(nt+1)-time_vert(nt) )
-
-!
-!--    Add horizontal large scale advection tendencies of pt and q 
-       SELECT CASE ( prog_var )
-
-          CASE ( 'pt' )
-
-             DO  k = nzb+1, nzt
-                tend(k,j,i) = tend(k,j,i) + td_lsa_lpt(k,nt)                   &
-                             + fac * ( td_lsa_lpt(k,nt+1) - td_lsa_lpt(k,nt) )*&
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 0 ) )
-             ENDDO
-
-          CASE ( 'q' )
-
-             DO  k = nzb+1, nzt
-                tend(k,j,i) = tend(k,j,i) + td_lsa_q(k,nt)                     &
-                              + fac * ( td_lsa_q(k,nt+1) - td_lsa_q(k,nt) ) *  &
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 0 ) )
-             ENDDO
-
-       END SELECT
-
-!
-!--    Subsidence of pt and q with prescribed profiles
-       IF ( large_scale_subsidence .AND. use_subsidence_tendencies )  THEN
-
-          SELECT CASE ( prog_var )
-
-             CASE ( 'pt' )
-
-                DO  k = nzb+1, nzt
-                   tend(k,j,i) = tend(k,j,i) + td_sub_lpt(k,nt) + fac *        &
-                                 ( td_sub_lpt(k,nt+1) - td_sub_lpt(k,nt) ) *   &
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 0 ) )
-                ENDDO
-  
-             CASE ( 'q' )
-
-                DO  k = nzb+1, nzt
-                   tend(k,j,i) = tend(k,j,i) + td_sub_q(k,nt) + fac *          &
-                                 ( td_sub_q(k,nt+1) - td_sub_q(k,nt) ) *       &
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 0 ) )
-                ENDDO
-
-          END SELECT
-
-       ENDIF
-
-    END SUBROUTINE ls_advec_ij
-
-
- END MODULE ls_forcing_mod
Index: /palm/trunk/SOURCE/modules.f90
===================================================================
--- /palm/trunk/SOURCE/modules.f90	(revision 2319)
+++ /palm/trunk/SOURCE/modules.f90	(revision 2320)
@@ -25,4 +25,9 @@
 ! -----------------
 ! $Id$
+! -ptnudge, qnudge, tnudge, td_lsa_lpt, td_lsa_q, td_sub_lpt, td_sub_q, ug_vert,
+!  vg_vert, unudge, vnudge, wsubs_vert, shf_surf, p_surf, pt_surf, q_surt,
+!  qsws_surf, tmp_tnudge, timenudge, time_surf, time_vert
+! 
+! 2300 2017-06-29 13:31:14Z raasch
 ! default value for host changed to '????', default value for loop_optimization
 ! changed to 'cache', default value for termination_time_needed set to 35.0
@@ -532,11 +537,7 @@
     REAL(wp), DIMENSION(:), ALLOCATABLE ::  ptdf_x                 !< damping factor for potential temperature in x-direction 
     REAL(wp), DIMENSION(:), ALLOCATABLE ::  ptdf_y                 !< damping factor for potential temperature in y-direction 
-    REAL(wp), DIMENSION(:), ALLOCATABLE ::  p_surf                 !< time-dependent surface pressure (large scale forcing)
-    REAL(wp), DIMENSION(:), ALLOCATABLE ::  pt_surf                !< time-dependent surface temperature (large scale forcing)
     REAL(wp), DIMENSION(:), ALLOCATABLE ::  pt_init                !< initial profile of potential temperature 
-    REAL(wp), DIMENSION(:), ALLOCATABLE ::  qsws_surf              !< time-dependent surface latent heat flux (large scale forcing)
     REAL(wp), DIMENSION(:), ALLOCATABLE ::  q_init                 !< initial profile of specific humidity 
                                                                    !< (or total water content with active cloud physics) 
-    REAL(wp), DIMENSION(:), ALLOCATABLE ::  q_surf                 !< time-dependent surface specific humidity (large scale forcing)
     REAL(wp), DIMENSION(:), ALLOCATABLE ::  rdf                    !< rayleigh damping factor for velocity components 
     REAL(wp), DIMENSION(:), ALLOCATABLE ::  rdf_sc                 !< rayleigh damping factor for scalar quantities 
@@ -545,9 +546,4 @@
     REAL(wp), DIMENSION(:), ALLOCATABLE ::  s_init                 !< initial profile of passive scalar concentration
     REAL(wp), DIMENSION(:), ALLOCATABLE ::  sa_init                !< initial profile of salinity (ocean)
-    REAL(wp), DIMENSION(:), ALLOCATABLE ::  shf_surf               !< time-dependent surface sensible heat flux (large scale forcing)
-    REAL(wp), DIMENSION(:), ALLOCATABLE ::  timenudge              !< times at which vertical profiles are defined in NUDGING_DATA (nudging)
-    REAL(wp), DIMENSION(:), ALLOCATABLE ::  time_surf              !< times at which surface values/fluxes are defined in LSF_DATA (large scale forcing)
-    REAL(wp), DIMENSION(:), ALLOCATABLE ::  time_vert              !< times at which vertical profiles are defined in LSF_DATA (large scale forcing)
-    REAL(wp), DIMENSION(:), ALLOCATABLE ::  tmp_tnudge             !< current nudging time scale
     REAL(wp), DIMENSION(:), ALLOCATABLE ::  ug                     !< geostrophic wind component in x-direction
     REAL(wp), DIMENSION(:), ALLOCATABLE ::  u_init                 !< initial profile of horizontal velocity component u
@@ -593,21 +589,8 @@
     REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  precipitation_amount  !< precipitation amount due to gravitational settling (bulk microphysics)
     REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  precipitation_rate    !< precipitation rate (bulk microphysics)
-    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  ptnudge               !< vertical profile of pot. temperature interpolated to vertical grid (nudging)
     REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  pt_slope_ref          !< potential temperature in rotated coordinate system 
                                                                     !< (in case of sloped surface) 
-    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  qnudge                !< vertical profile of specific humidity interpolated to vertical grid (nudging) 
-    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  tnudge                !< vertical profile of nudging time scale interpolated to vertical grid (nudging)  
-    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  td_lsa_lpt            !< temperature tendency due to large scale advection (large scale forcing)
-    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  td_lsa_q              !< specific humidity tendency due to large scale advection (large scale forcing)
-    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  td_sub_lpt            !< temperature tendency due to subsidence/ascent (large scale forcing)
-    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  td_sub_q              !< specific humidity tendency due to subsidence/ascent (large scale forcing)
     REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  total_2d_a            !< horizontal array to store the total domain data, used for atmosphere-ocean coupling (atmosphere data)
     REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  total_2d_o            !< horizontal array to store the total domain data, used for atmosphere-ocean coupling (ocean data)
-    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  ug_vert               !< vertical profile of geostrophic wind component in x-direction interpolated to vertical grid (large scale forcing)
-    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  unudge                !< vertical profile of wind component in x-direction interpolated to vertical grid (nudging) 
-    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  vnudge                !< vertical profile of wind component in y-direction interpolated to vertical grid (nudging) 
-    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  vg_vert               !< vertical profile of geostrophic wind component in y-direction interpolated to vertical grid (large scale forcing)
-    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  wnudge                !< vertical profile of subsidence/ascent velocity interpolated to vertical grid (nudging) ???
-    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  wsubs_vert            !< vertical profile of wind component in z-direction interpolated to vertical grid (nudging) ???
     
     REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  d          !< divergence
@@ -1061,6 +1044,4 @@
     INTEGER(iwp) ::  mg_switch_to_pe0_level = -1       !< namelist parameter
     INTEGER(iwp) ::  mid                               !< masked output running index
-    INTEGER(iwp) ::  nlsf = 1000                       !< maximum number of profiles in LSF_DATA (large scale forcing)
-    INTEGER(iwp) ::  ntnudge = 1000                    !< maximum number of profiles in NUDGING_DATA (nudging)
     INTEGER(iwp) ::  ngsrb = 2                         !< namelist parameter
     INTEGER(iwp) ::  nr_timesteps_this_run = 0         !< number of timesteps (cpu time measurements)
Index: lm/trunk/SOURCE/nudging_mod.f90
===================================================================
--- /palm/trunk/SOURCE/nudging_mod.f90	(revision 2319)
+++ 	(revision )
@@ -1,643 +1,0 @@
-!> @file nudging_mod.f90
-!------------------------------------------------------------------------------!
-! This file is part of PALM.
-!
-! PALM is free software: you can redistribute it and/or modify it under the
-! terms of the GNU General Public License as published by the Free Software
-! Foundation, either version 3 of the License, or (at your option) any later
-! version.
-!
-! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
-! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
-! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
-!
-! You should have received a copy of the GNU General Public License along with
-! PALM. If not, see <http://www.gnu.org/licenses/>.
-!
-! Copyright 1997-2017 Leibniz Universitaet Hannover
-!------------------------------------------------------------------------------!
-!
-! Current revisions:
-! ------------------
-! 
-! 
-! Former revisions:
-! -----------------
-! $Id$
-! Removed unused variables lptnudge, lqnudge, lunudge, lvnudge, lwnudge
-! 
-! 2271 2017-06-09 12:34:55Z sward
-! Error message changed
-! 
-! 2233 2017-05-30 18:08:54Z suehring
-!
-! 2232 2017-05-30 17:47:52Z suehring
-! Adjustments to new topography concept
-! 
-! 2000 2016-08-20 18:09:15Z knoop
-! Forced header and separation lines into 80 columns
-! 
-! 1850 2016-04-08 13:29:27Z maronga
-! Module renamed
-! 
-! 
-! 1757 2016-02-22 15:49:32Z maronga
-! Bugfix: allow for using higher vertical resolution in nudging file than grid
-! spacing in the LES model
-! 
-! 1682 2015-10-07 23:56:08Z knoop
-! Code annotations made doxygen readable 
-! 
-! 1398 2014-05-07 11:15:00Z heinze
-! Subroutine nudge_ref is extended to set u_init and v_init to the current 
-! nudging profiles
-! 
-! 1382 2014-04-30 12:15:41Z boeske
-! Changed the weighting factor that is used in the summation of nudging 
-! tendencies for profile data output from weight_pres to weight_substep,
-! added Neumann boundary conditions for profile data output of nudging terms at
-! nzt+1
-! 
-! 1380 2014-04-28 12:40:45Z heinze
-! Subroutine nudge_ref added to account for proper upper scalar boundary 
-! conditions in case of nudging
-! 
-! 1365 2014-04-22 15:03:56Z boeske
-! Variable t renamed nt, variable currtnudge renamed tmp_tnudge,
-! summation of nudging tendencies for data output added
-! +sums_ls_l, tmp_tend
-! Added new subroutine calc_tnudge, which calculates the current nudging time 
-! scale at each time step
-! 
-! 1355 2014-04-10 10:21:29Z heinze
-! Error message specified.
-! 
-! 1353 2014-04-08 15:21:23Z heinze
-! REAL constants provided with KIND-attribute 
-!
-! 1320 2014-03-20 08:40:49Z raasch
-! ONLY-attribute added to USE-statements,
-! kind-parameters added to all INTEGER and REAL declaration statements,
-! kinds are defined in new module kinds,
-! old module precision_kind is removed,
-! revision history before 2012 removed,
-! comment fields (!:) to be used for variable explanations added to
-! all variable declaration statements
-!
-! 1318 2014-03-17 13:35:16Z raasch
-! module interfaces removed
-!
-! 1268 2013-12-12 09:47:53Z heinze
-! bugfix: argument of calc_mean_profile corrected
-!
-! 1251 2013-11-07 08:14:30Z heinze
-! bugfix: calculate dtm and dtp also in vector version
-!
-! 1249 2013-11-06 10:45:47Z heinze
-! remove call of user module
-! reformatting
-!
-! 1241 2013-10-30 11:36:58Z heinze
-! Initial revision
-!
-! Description:
-! ------------
-!> Nudges u, v, pt and q to given profiles on a relaxation timescale tnudge. 
-!> Profiles are read in from NUDGIN_DATA. Code is based on Neggers et al. (2012) 
-!> and also part of DALES and UCLA-LES.
-!--------------------------------------------------------------------------------!
- MODULE nudge_mod
- 
-
-    PRIVATE
-    PUBLIC init_nudge, calc_tnudge, nudge, nudge_ref
-    SAVE
-
-    INTERFACE nudge
-       MODULE PROCEDURE nudge
-       MODULE PROCEDURE nudge_ij
-    END INTERFACE nudge
-
- CONTAINS
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> @todo Missing subroutine description.
-!------------------------------------------------------------------------------!
-    SUBROUTINE init_nudge
-
-       USE arrays_3d,                                                          &
-           ONLY:  ptnudge, qnudge, timenudge, tmp_tnudge, tnudge, unudge,      &
-                  vnudge, wnudge, zu
-
-       USE control_parameters,                                                 &
-           ONLY:  dt_3d, message_string, ntnudge
-
-       USE indices,                                                            &
-           ONLY:  nzb, nzt
-
-       USE kinds
-
-       IMPLICIT NONE
-
-
-       INTEGER(iwp) ::  finput = 90  !<
-       INTEGER(iwp) ::  ierrn        !<
-       INTEGER(iwp) ::  k            !<
-       INTEGER(iwp) ::  nt            !<
-
-       CHARACTER(1) ::  hash     !<
-
-       REAL(wp) ::  highheight   !<
-       REAL(wp) ::  highqnudge   !<
-       REAL(wp) ::  highptnudge  !<
-       REAL(wp) ::  highunudge   !<
-       REAL(wp) ::  highvnudge   !<
-       REAL(wp) ::  highwnudge   !<
-       REAL(wp) ::  hightnudge   !<
-
-       REAL(wp) ::  lowheight    !<
-       REAL(wp) ::  lowqnudge    !<
-       REAL(wp) ::  lowptnudge   !<
-       REAL(wp) ::  lowunudge    !<
-       REAL(wp) ::  lowvnudge    !<
-       REAL(wp) ::  lowwnudge    !<
-       REAL(wp) ::  lowtnudge    !<
-
-       REAL(wp) ::  fac          !<
-
-       ALLOCATE( ptnudge(nzb:nzt+1,1:ntnudge), qnudge(nzb:nzt+1,1:ntnudge), &
-                 tnudge(nzb:nzt+1,1:ntnudge), unudge(nzb:nzt+1,1:ntnudge),  &
-                 vnudge(nzb:nzt+1,1:ntnudge), wnudge(nzb:nzt+1,1:ntnudge)  )
-
-       ALLOCATE( tmp_tnudge(nzb:nzt) )
-
-       ALLOCATE( timenudge(0:ntnudge) )
-
-       ptnudge = 0.0_wp; qnudge = 0.0_wp; tnudge = 0.0_wp; unudge = 0.0_wp
-       vnudge = 0.0_wp; wnudge = 0.0_wp; timenudge = 0.0_wp
-!
-!--    Initialize array tmp_nudge with a current nudging time scale of 6 hours
-       tmp_tnudge = 21600.0_wp
-
-       nt = 0
-       OPEN ( finput, FILE='NUDGING_DATA', STATUS='OLD', &
-              FORM='FORMATTED', IOSTAT=ierrn )
-
-       IF ( ierrn /= 0 )  THEN
-          message_string = 'file NUDGING_DATA does not exist'
-          CALL message( 'nudging', 'PA0365', 1, 2, 0, 6, 0 )
-       ENDIF
-
-       ierrn = 0
-
- rloop:DO
-          nt = nt + 1
-          hash = "#"
-          ierrn = 1 ! not zero
-!
-!--       Search for the next line consisting of "# time", 
-!--       from there onwards the profiles will be read
-          DO WHILE ( .NOT. ( hash == "#" .AND. ierrn == 0 ) ) 
-          
-            READ ( finput, *, IOSTAT=ierrn ) hash, timenudge(nt)
-            IF ( ierrn < 0 )  EXIT rloop
-
-          ENDDO
-
-          ierrn = 0
-          READ ( finput, *, IOSTAT=ierrn ) lowheight, lowtnudge, lowunudge,   &
-                                           lowvnudge, lowwnudge , lowptnudge, &
-                                           lowqnudge
-
-          IF ( ierrn /= 0 )  THEN
-             message_string = 'errors in file NUDGING_DATA'
-             CALL message( 'nudging', 'PA0366', 1, 2, 0, 6, 0 )
-          ENDIF
-
-          ierrn = 0
-          READ ( finput, *, IOSTAT=ierrn ) highheight, hightnudge, highunudge,   &
-                                           highvnudge, highwnudge , highptnudge, &
-                                           highqnudge
-
-          IF ( ierrn /= 0 )  THEN
-             message_string = 'errors in file NUDGING_DATA'
-             CALL message( 'nudging', 'PA0366', 1, 2, 0, 6, 0 )
-          ENDIF
-
-          DO  k = nzb, nzt+1
-             DO WHILE ( highheight < zu(k) )
-                lowheight  = highheight
-                lowtnudge  = hightnudge
-                lowunudge  = highunudge
-                lowvnudge  = highvnudge
-                lowwnudge  = highwnudge
-                lowptnudge = highptnudge
-                lowqnudge  = highqnudge
- 
-                ierrn = 0
-                READ ( finput, *, IOSTAT=ierrn )  highheight , hightnudge , &
-                                                  highunudge , highvnudge , &
-                                                  highwnudge , highptnudge, &
-                                                  highqnudge
-                IF (ierrn /= 0 )  THEN
-                   WRITE( message_string, * ) 'zu(',k,') = ', zu(k), 'm is ',&
-                        'higher than the maximum height in NUDING_DATA which ',  &
-                        'is ', lowheight, 'm. Interpolation on PALM ',           &
-                        'grid is not possible.'
-                   CALL message( 'nudging', 'PA0364', 1, 2, 0, 6, 0 )
-                ENDIF
-             ENDDO
-
-!
-!--          Interpolation of prescribed profiles in space 
-
-             fac = ( highheight - zu(k) ) / ( highheight - lowheight )
-
-             tnudge(k,nt)  = fac * lowtnudge  + ( 1.0_wp - fac ) * hightnudge
-             unudge(k,nt)  = fac * lowunudge  + ( 1.0_wp - fac ) * highunudge
-             vnudge(k,nt)  = fac * lowvnudge  + ( 1.0_wp - fac ) * highvnudge
-             wnudge(k,nt)  = fac * lowwnudge  + ( 1.0_wp - fac ) * highwnudge
-             ptnudge(k,nt) = fac * lowptnudge + ( 1.0_wp - fac ) * highptnudge
-             qnudge(k,nt)  = fac * lowqnudge  + ( 1.0_wp - fac ) * highqnudge
-          ENDDO
-
-       ENDDO rloop
-
-       CLOSE ( finput )
-
-
-    END SUBROUTINE init_nudge
-
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> @todo Missing subroutine description.
-!------------------------------------------------------------------------------!
-    SUBROUTINE calc_tnudge ( time )
-
-       USE arrays_3d,                                                          &
-           ONLY:  timenudge, tmp_tnudge, tnudge
-
-       USE control_parameters,                                                 &
-           ONLY:  dt_3d 
-
-       USE indices,                                                            &
-           ONLY:  nzb, nzt
-
-       USE kinds
-
-       IMPLICIT NONE
-
-
-       REAL(wp) ::  dtm         !<
-       REAL(wp) ::  dtp         !<
-       REAL(wp) ::  time        !<
-
-       INTEGER(iwp) ::  k   !<
-       INTEGER(iwp) ::  nt  !<
-
-       nt = 1
-       DO WHILE ( time > timenudge(nt) )
-         nt = nt+1
-       ENDDO
-       IF ( time /= timenudge(1) ) THEN
-         nt = nt-1
-       ENDIF
-
-       dtm = ( time - timenudge(nt) ) / ( timenudge(nt+1) - timenudge(nt) )
-       dtp = ( timenudge(nt+1) - time ) / ( timenudge(nt+1) - timenudge(nt) )
-
-       DO  k = nzb, nzt
-          tmp_tnudge(k) = MAX( dt_3d, tnudge(k,nt) * dtp + tnudge(k,nt+1) * dtm )
-       ENDDO
-
-    END SUBROUTINE calc_tnudge
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> Call for all grid points
-!------------------------------------------------------------------------------!
-    SUBROUTINE nudge ( time, prog_var )
-
-       USE arrays_3d,                                                          &
-           ONLY:  pt, ptnudge, q, qnudge, tend, timenudge, tmp_tnudge, tnudge, &
-                  u, unudge, v, vnudge
-
-       USE control_parameters,                                                 &
-           ONLY:  dt_3d, intermediate_timestep_count, message_string
-
-       USE indices,                                                            &
-           ONLY:  nxl, nxr, nys, nyn, nzb, nzt, wall_flags_0
-
-       USE kinds
-
-       USE statistics,                                                         &
-           ONLY:  hom, sums_ls_l, weight_substep
-
-       IMPLICIT NONE
-
-       CHARACTER (LEN=*) ::  prog_var  !<
-
-       REAL(wp) ::  tmp_tend    !<
-       REAL(wp) ::  dtm         !<
-       REAL(wp) ::  dtp         !<
-       REAL(wp) ::  time        !<
-
-       INTEGER(iwp) ::  i  !<
-       INTEGER(iwp) ::  j  !<
-       INTEGER(iwp) ::  k  !<
-       INTEGER(iwp) ::  nt  !<
-
-
-       nt = 1
-       DO WHILE ( time > timenudge(nt) )
-         nt = nt+1
-       ENDDO
-       IF ( time /= timenudge(1) ) THEN
-         nt = nt-1
-       ENDIF
-
-       dtm = ( time - timenudge(nt) ) / ( timenudge(nt+1) - timenudge(nt) )
-       dtp = ( timenudge(nt+1) - time ) / ( timenudge(nt+1) - timenudge(nt) )
-
-       SELECT CASE ( prog_var )
-
-          CASE ( 'u' )
-
-             DO  i = nxl, nxr
-                DO  j = nys, nyn
-
-                   DO  k = nzb+1, nzt
-
-                      tmp_tend = - ( hom(k,1,1,0) - ( unudge(k,nt) * dtp +     &
-                                     unudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
-
-                      tend(k,j,i) = tend(k,j,i) + tmp_tend *                   &
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 1 ) )
-
-                      sums_ls_l(k,6) = sums_ls_l(k,6) + tmp_tend *             &
-                                     weight_substep(intermediate_timestep_count)
-                   ENDDO
-                  
-                   sums_ls_l(nzt+1,6) = sums_ls_l(nzt,6)
- 
-                ENDDO
-            ENDDO
-
-          CASE ( 'v' )
-
-             DO  i = nxl, nxr
-                DO  j = nys, nyn
-
-                   DO  k = nzb+1, nzt
-
-                      tmp_tend = - ( hom(k,1,2,0) - ( vnudge(k,nt) * dtp +     &
-                                     vnudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
-
-                      tend(k,j,i) = tend(k,j,i) + tmp_tend *                   &
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 2 ) )
-
-                      sums_ls_l(k,7) = sums_ls_l(k,7) + tmp_tend *             &
-                                     weight_substep(intermediate_timestep_count)
-                   ENDDO
-                  
-                   sums_ls_l(nzt+1,7) = sums_ls_l(nzt,7)
-
-                ENDDO
-            ENDDO
-
-          CASE ( 'pt' )
-
-             DO  i = nxl, nxr
-                DO  j = nys, nyn
-
-                   DO  k = nzb+1, nzt
-
-                      tmp_tend = - ( hom(k,1,4,0) - ( ptnudge(k,nt) * dtp +    &
-                                     ptnudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
-
-                      tend(k,j,i) = tend(k,j,i) + tmp_tend *                   &
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 0 ) )
-
-                      sums_ls_l(k,4) = sums_ls_l(k,4) + tmp_tend *             &
-                                     weight_substep(intermediate_timestep_count)
-                   ENDDO
-
-                   sums_ls_l(nzt+1,4) = sums_ls_l(nzt,4)
-
-                ENDDO
-            ENDDO
-
-          CASE ( 'q' )
-
-             DO  i = nxl, nxr
-                DO  j = nys, nyn
-
-                   DO  k = nzb+1, nzt
-
-                      tmp_tend = - ( hom(k,1,41,0) - ( qnudge(k,nt) * dtp +    &
-                                     qnudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
-
-                      tend(k,j,i) = tend(k,j,i) + tmp_tend *                   &
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 0 ) )
-
-                      sums_ls_l(k,5) = sums_ls_l(k,5) + tmp_tend *             &
-                                     weight_substep(intermediate_timestep_count)
-                   ENDDO
-                  
-                   sums_ls_l(nzt+1,5) = sums_ls_l(nzt,5)
-
-                ENDDO
-            ENDDO
-
-          CASE DEFAULT
-             message_string = 'unknown prognostic variable "' // prog_var // '"'
-             CALL message( 'nudge', 'PA0367', 1, 2, 0, 6, 0 )
-
-       END SELECT
-
-    END SUBROUTINE nudge
-
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> Call for grid point i,j
-!------------------------------------------------------------------------------!
-
-    SUBROUTINE nudge_ij( i, j, time, prog_var )
-
-       USE arrays_3d,                                                          &
-           ONLY:  pt, ptnudge, q, qnudge, tend, timenudge, tmp_tnudge, tnudge, &
-                  u, unudge, v, vnudge
-
-       USE control_parameters,                                                 &
-           ONLY:  dt_3d, intermediate_timestep_count, message_string
-
-       USE indices,                                                            &
-           ONLY:  nxl, nxr, nys, nyn, nzb, nzt, wall_flags_0
-
-       USE kinds
-
-       USE statistics,                                                         &
-           ONLY:  hom, sums_ls_l, weight_substep
-
-       IMPLICIT NONE
-
-
-       CHARACTER (LEN=*) ::  prog_var  !<
-
-       REAL(wp) ::  tmp_tend    !<
-       REAL(wp) ::  dtm         !<
-       REAL(wp) ::  dtp         !<
-       REAL(wp) ::  time        !<
-
-       INTEGER(iwp) ::  i  !<
-       INTEGER(iwp) ::  j  !<
-       INTEGER(iwp) ::  k  !<
-       INTEGER(iwp) ::  nt  !<
-
-
-       nt = 1
-       DO WHILE ( time > timenudge(nt) )
-         nt = nt+1
-       ENDDO
-       IF ( time /= timenudge(1) )  THEN
-         nt = nt-1
-       ENDIF
-
-       dtm = ( time - timenudge(nt) ) / ( timenudge(nt+1) - timenudge(nt) )
-       dtp = ( timenudge(nt+1) - time ) / ( timenudge(nt+1) - timenudge(nt) )
-
-       SELECT CASE ( prog_var )
-
-          CASE ( 'u' )
-
-             DO  k = nzb+1, nzt
-
-                tmp_tend = - ( hom(k,1,1,0) - ( unudge(k,nt) * dtp +           &
-                               unudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
-
-                tend(k,j,i) = tend(k,j,i) + tmp_tend *                         &
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 1 ) )
-
-                sums_ls_l(k,6) = sums_ls_l(k,6) + tmp_tend                     &
-                                 * weight_substep(intermediate_timestep_count)
-             ENDDO
-
-             sums_ls_l(nzt+1,6) = sums_ls_l(nzt,6)
-
-          CASE ( 'v' )
-
-             DO  k = nzb+1, nzt
-
-                tmp_tend = - ( hom(k,1,2,0) - ( vnudge(k,nt) * dtp +           &
-                               vnudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
-
-                tend(k,j,i) = tend(k,j,i) + tmp_tend *                         &
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 2 ) )
-
-                sums_ls_l(k,7) = sums_ls_l(k,7) + tmp_tend                     &
-                                 * weight_substep(intermediate_timestep_count)
-             ENDDO
-
-             sums_ls_l(nzt+1,7) = sums_ls_l(nzt,7)
-
-          CASE ( 'pt' )
-
-             DO  k = nzb+1, nzt
-
-                tmp_tend = - ( hom(k,1,4,0) - ( ptnudge(k,nt) * dtp +          &
-                               ptnudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
-
-                tend(k,j,i) = tend(k,j,i) + tmp_tend *                         &
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 0 ) )
-
-                sums_ls_l(k,4) = sums_ls_l(k,4) + tmp_tend                     &
-                                 * weight_substep(intermediate_timestep_count)
-             ENDDO
-
-             sums_ls_l(nzt+1,4) = sums_ls_l(nzt,4)
-
-
-          CASE ( 'q' )
-
-             DO  k = nzb+1, nzt
-
-                tmp_tend = - ( hom(k,1,41,0) - ( qnudge(k,nt) * dtp +          &
-                               qnudge(k,nt+1) * dtm ) ) / tmp_tnudge(k)
-
-                tend(k,j,i) = tend(k,j,i) + tmp_tend *                         &
-                                        MERGE( 1.0_wp, 0.0_wp,                 &
-                                               BTEST( wall_flags_0(k,j,i), 0 ) )
-
-                sums_ls_l(k,5) = sums_ls_l(k,5) + tmp_tend                     &
-                                 * weight_substep(intermediate_timestep_count)
-             ENDDO
-
-             sums_ls_l(nzt+1,5) = sums_ls_l(nzt,5)
-
-          CASE DEFAULT
-             message_string = 'unknown prognostic variable "' // prog_var // '"'
-             CALL message( 'nudge', 'PA0367', 1, 2, 0, 6, 0 )
-
-       END SELECT
-
-
-    END SUBROUTINE nudge_ij
-
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> @todo Missing subroutine description.
-!------------------------------------------------------------------------------!
-    SUBROUTINE nudge_ref ( time )
-
-       USE arrays_3d,                                                          &
-           ONLY:  time_vert, ptnudge, pt_init, qnudge, q_init, unudge, u_init, &
-                  vnudge, v_init
-
-       USE kinds
-
-
-       IMPLICIT NONE
-
-       INTEGER(iwp) ::  nt                    !<
-
-       REAL(wp)             ::  fac           !<
-       REAL(wp), INTENT(in) ::  time          !<
-
-!
-!--    Interpolation in time of NUDGING_DATA for pt_init and q_init. This is 
-!--    needed for correct upper boundary conditions for pt and q and in case that 
-!      large scale subsidence as well as scalar Rayleigh-damping are used
-       nt = 1
-       DO WHILE ( time > time_vert(nt) )
-          nt = nt + 1
-       ENDDO
-       IF ( time /= time_vert(nt) )  THEN
-        nt = nt - 1
-       ENDIF
-
-       fac = ( time-time_vert(nt) ) / ( time_vert(nt+1)-time_vert(nt) )
-
-       pt_init = ptnudge(:,nt) + fac * ( ptnudge(:,nt+1) - ptnudge(:,nt) )
-       q_init  = qnudge(:,nt) + fac * ( qnudge(:,nt+1) - qnudge(:,nt) )
-       u_init  = unudge(:,nt) + fac * ( unudge(:,nt+1) - unudge(:,nt) )
-       v_init  = vnudge(:,nt) + fac * ( vnudge(:,nt+1) - vnudge(:,nt) )
-
-    END SUBROUTINE nudge_ref
-
- END MODULE nudge_mod
Index: /palm/trunk/SOURCE/palm.f90
===================================================================
--- /palm/trunk/SOURCE/palm.f90	(revision 2319)
+++ /palm/trunk/SOURCE/palm.f90	(revision 2320)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Modularize large-scale forcing and nudging
+! 
+! 2298 2017-06-29 09:28:18Z raasch
 ! type of write_binary changed from CHARACTER to LOGICAL,
 ! user interface required revision updated, MPI2 related part removed
@@ -207,10 +210,4 @@
         ONLY:  lsm_last_actions
 
-    USE ls_forcing_mod,                                                        &
-        ONLY:  init_ls_forcing
-
-    USE nudge_mod,                                                             &
-        ONLY:  init_nudge
-
     USE particle_attributes,                                                   &
         ONLY:  particle_advection
@@ -330,16 +327,4 @@
 !-- Generate grid parameters
     CALL init_grid
-
-!
-!-- Initialize nudging if required
-    IF ( nudging )  THEN
-       CALL init_nudge
-    ENDIF
-
-!
-!-- Initialize reading of large scale forcing from external file - if required
-    IF ( large_scale_forcing )  THEN
-       CALL init_ls_forcing
-    ENDIF
 
 !
Index: /palm/trunk/SOURCE/prognostic_equations.f90
===================================================================
--- /palm/trunk/SOURCE/prognostic_equations.f90	(revision 2319)
+++ /palm/trunk/SOURCE/prognostic_equations.f90	(revision 2320)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Modularize large-scale forcing and nudging
+! 
+! 2292 2017-06-20 09:51:42Z schwenkel
 ! Implementation of new microphysic scheme: cloud_scheme = 'morrison' 
 ! includes two more prognostic equations for cloud drop concentration (nc)  
@@ -325,12 +328,9 @@
     USE kinds
 
-    USE ls_forcing_mod,                                                        &
-        ONLY:  ls_advec
+    USE lsf_nudging_mod,                                                       &
+        ONLY:  ls_advec, nudge
 
     USE microphysics_mod,                                                      &
         ONLY:  microphysics_control
-
-    USE nudge_mod,                                                             &
-        ONLY:  nudge
 
     USE plant_canopy_model_mod,                                                &
Index: /palm/trunk/SOURCE/read_var_list.f90
===================================================================
--- /palm/trunk/SOURCE/read_var_list.f90	(revision 2319)
+++ /palm/trunk/SOURCE/read_var_list.f90	(revision 2320)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Formatting adjustment
+! 
+! 2265 2017-06-08 16:58:28Z schwenkel
 ! Unused variables removed.
 ! 
@@ -504,5 +507,5 @@
           CASE ( 'large_scale_forcing' )
              READ ( 13 )  large_scale_forcing
-           CASE ( 'large_scale_subsidence' )
+          CASE ( 'large_scale_subsidence' )
              READ ( 13 )  large_scale_subsidence
           CASE ( 'limiter_sedimentation' )
Index: /palm/trunk/SOURCE/time_integration.f90
===================================================================
--- /palm/trunk/SOURCE/time_integration.f90	(revision 2319)
+++ /palm/trunk/SOURCE/time_integration.f90	(revision 2320)
@@ -343,12 +343,9 @@
                skip_time_do_lsm
 
-    USE ls_forcing_mod,                                                        &
-        ONLY:  ls_forcing_surf, ls_forcing_vert
+    USE lsf_nudging_mod,                                                       &
+        ONLY:  calc_tnudge, ls_forcing_surf, ls_forcing_vert, nudge_ref
 
     USE microphysics_mod,                                                      &
         ONLY: collision_turbulence
-
-    USE nudge_mod,                                                             &
-        ONLY:  calc_tnudge, nudge_ref
 
     USE particle_attributes,                                                   &
@@ -764,5 +761,5 @@
                  intermediate_timestep_count == intermediate_timestep_count_max )&
              THEN
-                CALL ls_forcing_surf ( simulated_time )
+                CALL ls_forcing_surf( simulated_time )
              ENDIF
 
