Index: /palm/trunk/SOURCE/Makefile
===================================================================
--- /palm/trunk/SOURCE/Makefile	(revision 2258)
+++ /palm/trunk/SOURCE/Makefile	(revision 2259)
@@ -25,4 +25,7 @@
 # -----------------
 # $Id$
+# Implemented synthetic turbulence generator
+#
+# 2256 2017-06-07 13:58:08Z suehring
 # Remove ring dependency in init_pegrid
 # 
@@ -366,5 +369,6 @@
 	set_slicer_attributes_dvrp.f90 singleton_mod.f90 sor.f90 spectra_mod.f90 \
 	subsidence_mod.f90 sum_up_3d_data.f90 surface_mod.f90 \
-	surface_coupler.f90 surface_layer_fluxes_mod.f90 swap_timelevel.f90 temperton_fft_mod.f90 \
+	surface_coupler.f90 surface_layer_fluxes_mod.f90 swap_timelevel.f90 \
+	synthetic_turbulence_generator_mod.f90 temperton_fft_mod.f90 \
 	time_integration.f90 time_to_string.f90 timestep.f90 \
 	timestep_scheme_steering.f90 transpose.f90 tridia_solver_mod.f90 \
@@ -438,7 +442,7 @@
 	netcdf_interface_mod.o posix_calls_from_fortran.o
 check_parameters.o: modules.o mod_kinds.o land_surface_model_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 microphysics_mod.o wind_turbine_model_mod.o \
-	urban_surface_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 \
+	wind_turbine_model_mod.o urban_surface_mod.o
 close_file.o: modules.o mod_kinds.o netcdf_interface_mod.o
 compute_vpt.o: modules.o mod_kinds.o
@@ -496,5 +500,5 @@
 init_masks.o: modules.o mod_kinds.o netcdf_interface_mod.o
 init_ocean.o: modules.o eqn_state_seawater.o mod_kinds.o
-init_pegrid.o: modules.o mod_kinds.o spectra_mod.o
+init_pegrid.o: modules.o mod_kinds.o spectra_mod.o synthetic_turbulence_generator_mod.o
 init_pt_anomaly.o: modules.o mod_kinds.o
 init_rankine.o: modules.o mod_kinds.o
@@ -547,5 +551,6 @@
 parin.o: modules.o cpulog_mod.o land_surface_model_mod.o mod_kinds.o netcdf_interface_mod.o \
    plant_canopy_model_mod.o pmc_interface_mod.o progress_bar_mod.o spectra_mod.o \
-   radiation_model_mod.o microphysics_mod.o virtual_flight_mod.o wind_turbine_model_mod.o
+   synthetic_turbulence_generator_mod.o radiation_model_mod.o microphysics_mod.o \
+   virtual_flight_mod.o wind_turbine_model_mod.o
 plant_canopy_model_mod.o: modules.o mod_kinds.o surface_mod.o
 pmc_interface_mod.o: modules.o mod_kinds.o pmc_child_mod.o pmc_general_mod.o \
@@ -596,4 +601,5 @@
 swap_timelevel.o: modules.o cpulog_mod.o mod_kinds.o land_surface_model_mod.o \
    pmc_interface_mod.o urban_surface_mod.o
+synthetic_turbulence_generator_mod.o: modules.o mod_kinds.o
 temperton_fft_mod.o: modules.o mod_kinds.o
 time_integration.o: modules.o advec_ws.o buoyancy.o calc_mean_profile.o \
@@ -601,5 +607,6 @@
 	ls_forcing_mod.o mod_kinds.o nudging_mod.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 surface_layer_fluxes_mod.o surface_mod.o \
+	spectra_mod.o user_actions.o microphysics_mod.o synthetic_turbulence_generator_mod.o \
+        surface_layer_fluxes_mod.o surface_mod.o \
 	urban_surface_mod.o virtual_flight_mod.o wind_turbine_model_mod.o
 time_to_string.o: mod_kinds.o
Index: /palm/trunk/SOURCE/check_parameters.f90
===================================================================
--- /palm/trunk/SOURCE/check_parameters.f90	(revision 2258)
+++ /palm/trunk/SOURCE/check_parameters.f90	(revision 2259)
@@ -20,9 +20,12 @@
 ! Current revisions:
 ! -----------------
-!  
+! 
 ! 
 ! Former revisions:
 ! -----------------
 ! $Id$
+! Implemented synthetic turbulence generator
+!
+! 2251 2017-06-06 15:10:46Z Giersch
 !
 ! 2250 2017-06-06 15:07:12Z Giersch
@@ -499,4 +502,6 @@
     USE subsidence_mod
     USE statistics
+    USE synthetic_turbulence_generator_mod,                                    &
+        ONLY:  stg_check_parameters
     USE transpose_indices
     USE urban_surface_mod,                                                     &
@@ -1207,4 +1212,7 @@
     ENDIF
 
+!
+!-- When synthetic turbulence generator is used, peform addtional checks
+    IF ( synthetic_turbulence_generator )  CALL stg_check_parameters
 !
 !-- When plant canopy model is used, peform addtional checks
Index: /palm/trunk/SOURCE/header.f90
===================================================================
--- /palm/trunk/SOURCE/header.f90	(revision 2258)
+++ /palm/trunk/SOURCE/header.f90	(revision 2259)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Implemented synthetic turbulence generator
+!
+! 2258 2017-06-08 07:55:13Z suehring
 ! Bugfix, add pre-preprocessor directives to enable non-parrallel mode
 ! 
@@ -356,9 +359,7 @@
     USE plant_canopy_model_mod,                                                &
         ONLY:  pcm_header, plant_canopy
-        
-#if defined( __parallel )
+
     USE pmc_handle_communicator,                                               &
         ONLY:  pmc_get_model_info
-#endif
 
     USE pmc_interface,                                                         &
@@ -373,4 +374,7 @@
     USE surface_mod,                                                           &
         ONLY:  surf_def_h
+
+    USE synthetic_turbulence_generator_mod,                                    &
+        ONLY:  stg_header
 
     IMPLICIT NONE
@@ -552,5 +556,4 @@
        WRITE ( io, 600 )  TRIM( nesting_mode ),                                &
                           TRIM( nesting_datatransfer_mode )
-#if defined( __parallel )
        CALL pmc_get_model_info( ncpl = ncpl, cpl_id = my_cpl_id )
 
@@ -570,5 +573,4 @@
                              TRIM( cpl_name )
        ENDDO
-#endif
     ENDIF
     WRITE ( io, 99 )
@@ -988,4 +990,6 @@
        ENDIF
     ENDIF
+
+    IF ( synthetic_turbulence_generator )  CALL stg_header ( io )
 
     IF ( plant_canopy )  CALL pcm_header ( io )
Index: /palm/trunk/SOURCE/init_3d_model.f90
===================================================================
--- /palm/trunk/SOURCE/init_3d_model.f90	(revision 2258)
+++ /palm/trunk/SOURCE/init_3d_model.f90	(revision 2259)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Implemented synthetic turbulence generator
+!
+! 2252 2017-06-07 09:35:37Z knoop
 ! rho_air now depending on surface_pressure even in Boussinesq mode
 ! 
@@ -390,5 +393,8 @@
                sums_l_l, sums_up_fraction_l, sums_wsts_bc_l, ts_value,         &
                weight_pres, weight_substep
- 
+
+    USE synthetic_turbulence_generator_mod,                                    &
+        ONLY:  stg_init, use_synthetic_turbulence_generator
+
     USE surface_layer_fluxes_mod,                                              &
         ONLY:  init_surface_layer_fluxes
@@ -1057,4 +1063,10 @@
 
 !
+!--       Overwrite initial profiles in case of synthetic turbulence generator
+          IF( use_synthetic_turbulence_generator ) THEN
+             CALL stg_init
+          ENDIF
+
+!
 !--       Use constructed initial profiles (velocity constant with height,
 !--       temperature profile with constant gradient)
@@ -1538,4 +1550,8 @@
        IF ( passive_scalar )  ts_m  = 0.0_wp
        IF ( ocean          )  tsa_m = 0.0_wp
+!
+!--    Initialize synthetic turbulence generator in case of restart.
+       IF ( TRIM( initializing_actions ) == 'read_restart_data'  .AND.         &
+            use_synthetic_turbulence_generator )  CALL stg_init
 
        CALL location_message( 'finished', .TRUE. )
Index: /palm/trunk/SOURCE/init_pegrid.f90
===================================================================
--- /palm/trunk/SOURCE/init_pegrid.f90	(revision 2258)
+++ /palm/trunk/SOURCE/init_pegrid.f90	(revision 2259)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Implemented synthetic turbulence generator
+!
+! 2238 2017-05-31 16:49:16Z suehring
 ! Remove unnecessary module load of pmc_interface
 ! 
@@ -209,4 +212,7 @@
     USE spectra_mod,                                                           &
         ONLY:  calculate_spectra, dt_dosp
+
+    USE synthetic_turbulence_generator_mod,                                    &
+        ONLY:  use_synthetic_turbulence_generator
 
     USE transpose_indices,                                                     &
@@ -423,5 +429,6 @@
 !-- 1. transposition  z --> x
 !-- This transposition is not neccessary in case of a 1d-decomposition along x
-    IF ( psolver == 'poisfft'  .OR.  calculate_spectra )  THEN
+    IF ( psolver == 'poisfft'  .OR.  calculate_spectra  .OR.                   &
+         use_synthetic_turbulence_generator )  THEN
 
        IF ( pdims(2) /= 1 )  THEN
Index: /palm/trunk/SOURCE/modules.f90
===================================================================
--- /palm/trunk/SOURCE/modules.f90	(revision 2258)
+++ /palm/trunk/SOURCE/modules.f90	(revision 2259)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Implemented synthetic turbulence generator
+!
+! 2256 2017-06-07 13:58:08Z suehring
 ! Change default value of zeta_min to -20
 ! Increase dimension for wall_heatflux, etc. 
@@ -1165,4 +1168,5 @@
     LOGICAL ::  stop_dt = .FALSE.                            !< 
     LOGICAL ::  synchronous_exchange = .FALSE.               !< 
+    LOGICAL ::  synthetic_turbulence_generator = .FALSE.     !< flag for synthetic turbulence generator module
     LOGICAL ::  terminate_run = .FALSE.                      !< 
     LOGICAL ::  transpose_compute_overlap = .FALSE.          !< 
Index: /palm/trunk/SOURCE/parin.f90
===================================================================
--- /palm/trunk/SOURCE/parin.f90	(revision 2258)
+++ /palm/trunk/SOURCE/parin.f90	(revision 2259)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Implemented synthetic turbulence generator
+!
+! 2233 2017-05-30 18:08:54Z suehring
 ! 
 ! 2232 2017-05-30 17:47:52Z suehring
@@ -315,4 +318,7 @@
     USE statistics,                                                            &
         ONLY:  hom, hom_sum, pr_palm, region, statistic_regions
+
+    USE synthetic_turbulence_generator_mod,                                    &
+        ONLY:  stg_parin
 
     USE urban_surface_mod,                                                     &
@@ -605,5 +611,8 @@
 !--       if required
           CALL flight_parin
-
+!
+!--       Check if synthetic turbulence generator is used and read stg_par if
+!--       required
+          CALL stg_parin
 !
 !--       Read user-defined variables
Index: /palm/trunk/SOURCE/read_var_list.f90
===================================================================
--- /palm/trunk/SOURCE/read_var_list.f90	(revision 2258)
+++ /palm/trunk/SOURCE/read_var_list.f90	(revision 2259)
@@ -25,20 +25,23 @@
 ! -----------------
 ! $Id$
+! Implemented synthetic turbulence generator
+!
+! 2233 2017-05-30 18:08:54Z suehring
 !
 ! 2232 2017-05-30 17:47:52Z suehring
 ! Replace wall_qflux, wall_sflux by wall_humidityflux and wall_scalarflux; add
 ! wall_salinityflux
-! +tunnel_height, tunnel_length, tunnel_width_x, tunnel_width_y, 
+! +tunnel_height, tunnel_length, tunnel_width_x, tunnel_width_y,
 !  tunnel_wall_depth
 !
 ! 2042 2016-11-02 13:47:31Z suehring
 ! Bugfix, read restart data for wall_heatflux, wall_qflux and wall_sflux
-! 
+!
 ! 2000 2016-08-20 18:09:15Z knoop
 ! Forced header and separation lines into 80 columns
-! 
+!
 ! 1992 2016-08-12 15:14:59Z suehring
 ! top_scalarflux added
-! 
+!
 ! 1960 2016-07-12 16:34:24Z suehring
 ! Separate humidity and passive scalar
@@ -66,12 +69,12 @@
 ! 1699 2015-10-29 08:02:35Z maronga
 ! Bugfix: update of binary version from 3.9b to 4.0 was missing
-! 
+!
 ! 1691 2015-10-26 16:17:44Z maronga
-! Added output of most_method, constant_flux_layer, zeta_min, zeta_max. Removed 
+! Added output of most_method, constant_flux_layer, zeta_min, zeta_max. Removed
 ! output of prandtl_layer and rif_min, rif_max.
-! 
+!
 ! 1682 2015-10-07 23:56:08Z knoop
 ! Code annotations made doxygen readable
-! 
+!
 ! 1615 2015-07-08 18:49:19Z suehring
 ! Enable turbulent inflow for passive_scalar and humidity
@@ -79,19 +82,19 @@
 ! 1585 2015-04-30 07:05:52Z maronga
 ! Adapted for RRTMG
-! 
+!
 ! 1560 2015-03-06 10:48:54Z keck
 ! +recycling_yshift
-! 
+!
 ! 1522 2015-01-14 10:53:12Z keck
 ! added query for checking if the advection scheme in the restart run is the
 ! same as the advection scheme in the corresponding initial run
-! 
+!
 ! 1502 2014-12-03 18:22:31Z kanani
-! Canopy module and parameters removed (parameters are always read from 
+! Canopy module and parameters removed (parameters are always read from
 ! canopy_par NAMELIST for initial and restart runs)
-! 
+!
 ! 1496 2014-12-02 17:25:50Z maronga
 ! Renamed "radiation" -> "cloud_top_radiation"
-! 
+!
 ! 1484 2014-10-21 10:53:05Z kanani
 ! Changes in the course of the canopy-model modularization:
@@ -118,5 +121,5 @@
 !
 ! 1253 2013-11-07 10:48:12Z fricke
-! Bugfix: add ref_state to read_parts_of_var_list, otherwise ref_state 
+! Bugfix: add ref_state to read_parts_of_var_list, otherwise ref_state
 ! is zero for initializing_actions = 'cyclic_fill'
 !
@@ -142,5 +145,5 @@
 !
 ! 1053 2012-11-13 17:11:03Z hoffmann
-! necessary expansions according to the two new prognostic equations (nr, qr) 
+! necessary expansions according to the two new prognostic equations (nr, qr)
 ! of the two-moment cloud physics scheme:
 ! +bc_*_b, +bc_*_t, +bc_*_t_val, *_init, *_surface, *_surface_initial_change,
@@ -148,6 +151,6 @@
 ! +surface_waterflux_*
 !
-! in addition, steering parameters parameters of the two-moment cloud physics 
-! scheme:   
+! in addition, steering parameters parameters of the two-moment cloud physics
+! scheme:
 ! +cloud_scheme, +drizzle, +mu_constant, +mu_constant_value, +ventilation_effect
 !
@@ -191,5 +194,5 @@
 !------------------------------------------------------------------------------!
  SUBROUTINE read_var_list
- 
+
 
     USE arrays_3d,                                                             &
@@ -201,5 +204,5 @@
     USE flight_mod,                                                            &
         ONLY:  flight_read_restart_data
-    
+
     USE grid_variables,                                                        &
         ONLY:  dx, dy
@@ -234,4 +237,6 @@
                v_max, v_max_ijk, w_max, w_max_ijk
 
+    USE synthetic_turbulence_generator_mod,                                    &
+        ONLY:  stg_read_restart_data
 
     IMPLICIT NONE
@@ -628,5 +633,5 @@
              READ ( 13 )  s_vertical_gradient_level
           CASE ( 's_vertical_gradient_level_ind' )
-             READ ( 13 )  s_vertical_gradient_level_ind 
+             READ ( 13 )  s_vertical_gradient_level_ind
           CASE ( 'sa_init' )
              READ ( 13 )  sa_init
@@ -646,7 +651,7 @@
              READ ( 13 )  surface_pressure
           CASE ( 'surface_scalarflux' )
-             READ ( 13 )  surface_scalarflux             
+             READ ( 13 )  surface_scalarflux
           CASE ( 'surface_waterflux' )
-             READ ( 13 )  surface_waterflux      
+             READ ( 13 )  surface_waterflux
           CASE ( 'time_coupling' )
              READ ( 13 )  time_coupling
@@ -770,5 +775,5 @@
              READ ( 13 )  vg_vertical_gradient_level_ind
           CASE ( 'virtual_flight' )
-             READ ( 13 ) virtual_flight  
+             READ ( 13 ) virtual_flight
           CASE ( 'volume_flow_area' )
              READ ( 13 )  volume_flow_area
@@ -805,5 +810,5 @@
              WRITE( message_string, * ) 'unknown variable named "',         &
                                         TRIM( variable_chr ), '" found in', &
-                                        ' data from prior run on PE ', myid 
+                                        ' data from prior run on PE ', myid
              CALL message( 'read_var_list', 'PA0302', 1, 2, 0, 6, 0 )
         END SELECT
@@ -813,7 +818,8 @@
 
     ENDDO
-    
+
     IF ( virtual_flight )  CALL flight_read_restart_data
 
+    IF ( synthetic_turbulence_generator )  CALL stg_read_restart_data
 
  END SUBROUTINE read_var_list
@@ -828,5 +834,5 @@
 !> run which used a smaller total domain or/and a different domain decomposition.
 !------------------------------------------------------------------------------!
- 
+
  SUBROUTINE read_parts_of_var_list
 
@@ -938,5 +944,5 @@
     READ ( 13 )  statistic_regions_on_file
     IF ( statistic_regions_on_file /= statistic_regions )  THEN
-       WRITE( message_string, * ) 'statistic regions on restart data file ',& 
+       WRITE( message_string, * ) 'statistic regions on restart data file ',&
                                   'differ from the current run',            &
                                   '&statistic regions on file    = "',      &
@@ -993,6 +999,6 @@
                 CALL message( 'read_parts_of_var_list', 'PA0100', &
                                                                  1, 2, 0, 6, 0 )
-             END IF                           
-             
+             END IF
+
           CASE ( 'nx' )
              READ ( 13 )  nx_on_file
@@ -1012,6 +1018,6 @@
                 CALL message( 'read_parts_of_var_list', 'PA0101', &
                                                                  1, 2, 0, 6, 0 )
-             END IF             
-             
+             END IF
+
           CASE DEFAULT
 
@@ -1040,13 +1046,13 @@
 !> Skipping the global control variables from restart-file (binary format)
 !------------------------------------------------------------------------------!
- 
+
  SUBROUTINE skip_var_list
 
     USE control_parameters,                                                    &
         ONLY:  virtual_flight
- 
+
     USE flight_mod,                                                            &
         ONLY:  flight_skip_var_list
- 
+
     IMPLICIT NONE
 
@@ -1068,8 +1074,8 @@
     ENDDO
 !
-!-- In case of virtual flights, skip also variables related to 
+!-- In case of virtual flights, skip also variables related to
 !-- this module.
     IF ( virtual_flight )  CALL flight_skip_var_list
-    
+
 
  END SUBROUTINE skip_var_list
Index: /palm/trunk/SOURCE/synthetic_turbulence_generator_mod.f90
===================================================================
--- /palm/trunk/SOURCE/synthetic_turbulence_generator_mod.f90	(revision 2259)
+++ /palm/trunk/SOURCE/synthetic_turbulence_generator_mod.f90	(revision 2259)
@@ -0,0 +1,995 @@
+!> @synthetic_turbulence_generator_mod.f90
+!------------------------------------------------------------------------------!
+! This file is part of PALM-4U.
+!
+! PALM-4U 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-4U 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 2017 Leibniz Universitaet Hannover
+!------------------------------------------------------------------------------!
+!
+! Current revisions:
+! -----------------
+! 
+! 
+! Former revisions:
+! -----------------
+! $Id$
+! Initial revision
+!
+! 
+!
+! Authors:
+! --------
+! @author Tobias Gronemeier, Atsushi Inagaki, Micha Gryschka, Christoph Knigge
+!
+!
+! Description:
+! ------------
+!> The module generates turbulence at the inflow boundary based on a method by
+!> Xie and Castro (2008) utilizing a Lund rotation (Lund, 1998) and a mass-flux
+!> correction by Kim et al. (2013).
+!> The turbulence is correlated based on length scales in y- and z-direction and
+!> a time scale for each velocity component. The profiles of length and time
+!> scales, mean u, v, w, e and pt, and all components of the Reynolds stress
+!> tensor are read from file STG_PROFILES.
+!>
+!> @todo test restart
+!>       enable cyclic_fill
+!>       implement turbulence generation for e and pt
+!> @note <Enter notes on the module>
+!> @bug  Height information from input file is not used. Profiles from input
+!>       must match with current PALM grid.
+!>       Transformation of length scales to number of gridpoints does not
+!>       consider grid stretching.
+!>       In case of restart, velocity seeds differ from precursor run if a11,
+!>       a22, or a33 are zero.
+!------------------------------------------------------------------------------!
+ MODULE synthetic_turbulence_generator_mod
+
+
+    USE arrays_3d,                                                             &
+        ONLY:  mean_inflow_profiles, u, v, w
+
+    USE constants,                                                             &
+        ONLY:  pi
+
+    USE control_parameters,                                                    &
+        ONLY:  initializing_actions, message_string,                           &
+               synthetic_turbulence_generator
+
+    USE cpulog,                                                                &
+        ONLY:  cpu_log, log_point, log_point_s
+
+    USE indices,                                                               &
+        ONLY:  nbgp, nzb, nzt, nyng, nysg
+
+    USE kinds
+
+    USE MPI
+
+    USE pegrid,                                                                &
+        ONLY:  comm1dx, comm1dy, comm2d, id_inflow, ierr, myidx, pdims
+
+    USE transpose_indices,                                                     &
+        ONLY: nzb_x, nzt_x
+
+
+    IMPLICIT NONE
+
+    LOGICAL :: velocity_seed_initialized = .FALSE.           !< true after first call of stg_main
+    LOGICAL :: use_synthetic_turbulence_generator = .FALSE.  !< switch to use synthetic turbulence generator
+
+    INTEGER(iwp) :: stg_type_yz        !< MPI type for full z range
+    INTEGER(iwp) :: stg_type_yz_small  !< MPI type for small z range
+    INTEGER(iwp) :: merg               !< maximum length scale (in gp)
+    INTEGER(iwp) :: mergp              !< merg + nbgp
+
+    REAL(wp) :: mc_factor    !< mass flux correction factor
+
+    INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: displs      !< displacement for MPI_GATHERV
+    INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: recv_count  !< receive count for MPI_GATHERV
+    INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: nuy         !< length scale of u in y direction (in gp)
+    INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: nuz         !< length scale of u in z direction (in gp)
+    INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: nvy         !< length scale of v in y direction (in gp)
+    INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: nvz         !< length scale of v in z direction (in gp)
+    INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: nwy         !< length scale of w in y direction (in gp)
+    INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: nwz         !< length scale of w in z direction (in gp)
+
+    INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: seed        !< seed of random number for rn-generator
+
+    REAL(wp), DIMENSION(:), ALLOCATABLE :: a11             !< coefficient for Lund rotation
+    REAL(wp), DIMENSION(:), ALLOCATABLE :: a21             !< coefficient for Lund rotation
+    REAL(wp), DIMENSION(:), ALLOCATABLE :: a22             !< coefficient for Lund rotation
+    REAL(wp), DIMENSION(:), ALLOCATABLE :: a31             !< coefficient for Lund rotation
+    REAL(wp), DIMENSION(:), ALLOCATABLE :: a32             !< coefficient for Lund rotation
+    REAL(wp), DIMENSION(:), ALLOCATABLE :: a33             !< coefficient for Lund rotation
+    REAL(wp), DIMENSION(:), ALLOCATABLE :: tu              !< Lagrangian time scale of u
+    REAL(wp), DIMENSION(:), ALLOCATABLE :: tv              !< Lagrangian time scale of v
+    REAL(wp), DIMENSION(:), ALLOCATABLE :: tw              !< Lagrangian time scale of w
+
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE :: buy           !< filter function for u in y direction
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE :: buz           !< filter function for u in z direction
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE :: bvy           !< filter function for v in y direction
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE :: bvz           !< filter function for v in z direction
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE :: bwy           !< filter function for w in y direction
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE :: bwz           !< filter function for w in z direction
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE :: fu            !< velocity seed for u
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE :: fuo           !< velocity seed for u with new random number
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE :: fv            !< velocity seed for v
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE :: fvo           !< velocity seed for v with new random number
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE :: fw            !< velocity seed for w
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE :: fwo           !< velocity seed for w with new random number
+
+
+!
+!-- PALM interfaces:
+!-- Input parameter checks to be done in check_parameters
+    INTERFACE stg_check_parameters
+       MODULE PROCEDURE stg_check_parameters
+    END INTERFACE stg_check_parameters
+
+!
+!-- Calculate filter functions
+    INTERFACE stg_filter_func
+       MODULE PROCEDURE stg_filter_func
+    END INTERFACE stg_filter_func
+
+!
+!-- Generate velocity seeds
+    INTERFACE stg_generate_seed_yz
+       MODULE PROCEDURE stg_generate_seed_yz
+    END INTERFACE stg_generate_seed_yz
+
+!
+!-- Output of information to the header file
+    INTERFACE stg_header
+       MODULE PROCEDURE stg_header
+    END INTERFACE stg_header
+
+!
+!-- Initialization actions
+    INTERFACE stg_init
+       MODULE PROCEDURE stg_init
+    END INTERFACE stg_init
+
+!
+!-- Main procedure of synth. turb. gen.
+    INTERFACE stg_main
+       MODULE PROCEDURE stg_main
+    END INTERFACE stg_main
+
+!
+!-- Reading of NAMELIST parameters
+    INTERFACE stg_parin
+       MODULE PROCEDURE stg_parin
+    END INTERFACE stg_parin
+
+!
+!-- Reading of parameters for restart runs
+    INTERFACE stg_read_restart_data
+       MODULE PROCEDURE stg_read_restart_data
+    END INTERFACE stg_read_restart_data
+
+!
+!-- Writing of binary output for restart runs
+    INTERFACE stg_write_restart_data
+       MODULE PROCEDURE stg_write_restart_data
+    END INTERFACE stg_write_restart_data
+
+    SAVE
+
+    PRIVATE
+
+!
+!-- Public interfaces
+    PUBLIC  stg_check_parameters, stg_header, stg_init, stg_main, stg_parin,   &
+            stg_read_restart_data, stg_write_restart_data
+
+!
+!-- Public variables
+    PUBLIC  use_synthetic_turbulence_generator
+
+
+ CONTAINS
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Check parameters routine for synthetic turbulence generator
+!------------------------------------------------------------------------------!
+ SUBROUTINE stg_check_parameters
+
+
+    USE control_parameters,                                                    &
+        ONLY:  bc_lr, bc_ns, turbulent_inflow
+
+
+    IMPLICIT NONE
+
+    IF ( use_synthetic_turbulence_generator )  THEN
+
+       IF ( INDEX( initializing_actions, 'set_constant_profiles' ) == 0  .AND. &
+            INDEX( initializing_actions, 'read_restart_data' ) == 0 )  THEN
+          message_string = 'Using synthetic turbulence generator ' //          &
+                           'requires &initializing_actions = '            //   &
+                           '"set_constant_profiles" or "read_restart_data"'
+          CALL message( 'stg_check_parameters', 'PA0015', 1, 2, 0, 6, 0 )
+       ENDIF
+
+       IF ( bc_lr /= 'dirichlet/radiation' )  THEN
+          message_string = 'Using synthetic turbulence generator ' //          &
+                           'requires &bc_lr = "dirichlet/radiation"'
+          CALL message( 'stg_check_parameters', 'PA0035', 1, 2, 0, 6, 0 )
+       ENDIF
+       IF ( bc_ns /= 'cyclic' )  THEN
+          message_string = 'Using synthetic turbulence generator ' //          &
+                           'requires &bc_ns = "cyclic"'
+          CALL message( 'stg_check_parameters', 'PA0037', 1, 2, 0, 6, 0 )
+       ENDIF
+       IF ( turbulent_inflow )  THEN
+          message_string = 'Using synthetic turbulence generator ' //          &
+                           'in combination &with turbulent_inflow = .T. ' //   &
+                           'is not allowed'
+          CALL message( 'stg_check_parameters', 'PA0039', 1, 2, 0, 6, 0 )
+       ENDIF
+
+    ENDIF
+
+ END SUBROUTINE stg_check_parameters
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Header output for synthetic turbulence generator
+!------------------------------------------------------------------------------!
+ SUBROUTINE stg_header ( io )
+
+
+    IMPLICIT NONE
+
+    INTEGER(iwp), INTENT(IN) ::  io   !< Unit of the output file
+
+!
+!-- Write synthetic turbulence generator Header
+    WRITE( io, 1 )
+    IF ( use_synthetic_turbulence_generator )  THEN
+       WRITE( io, 2 )
+    ELSE
+       WRITE( io, 3 )
+    ENDIF
+
+1   FORMAT (//' Synthetic turbulence generator information:'/                  &
+              ' ------------------------------------------'/)
+2   FORMAT ('    synthetic turbulence generator switched on')
+3   FORMAT ('    synthetic turbulence generator switched off')
+
+ END SUBROUTINE stg_header
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Initialization of the synthetic turbulence generator
+!------------------------------------------------------------------------------!
+ SUBROUTINE stg_init
+
+
+    USE arrays_3d,                                                             &
+        ONLY:  u_init, v_init
+
+    USE control_parameters,                                                    &
+        ONLY:  coupling_char, dz, e_init
+
+    USE grid_variables,                                                        &
+        ONLY:  ddy
+
+
+    IMPLICIT NONE
+
+    INTEGER(KIND=MPI_ADDRESS_KIND) :: extent !< extent of new MPI type
+    INTEGER(KIND=MPI_ADDRESS_KIND) :: tob=0  !< dummy variable
+
+    INTEGER(iwp) :: j                        !> loop index
+    INTEGER(iwp) :: k                        !< index
+    INTEGER(iwp) :: newtype                  !< dummy MPI type
+    INTEGER(iwp) :: realsize                 !< size of REAL variables
+    INTEGER(iwp) :: nseed                    !< dimension of random number seed
+    INTEGER(iwp) :: startseed = 1234567890   !< start seed for random number generator
+
+!
+!-- Dummy variables used for reading profiles
+    REAL(wp) :: d1      !< u profile
+    REAL(wp) :: d2      !< v profile
+    REAL(wp) :: d3      !< w profile
+    REAL(wp) :: d5      !< e profile
+    REAL(wp) :: d11     !< vertical interpolation for a11
+    REAL(wp) :: d21     !< vertical interpolation for a21
+    REAL(wp) :: d22     !< vertical interpolation for a22
+    REAL(wp) :: luy     !< length scale for u in y direction
+    REAL(wp) :: luz     !< length scale for u in z direction
+    REAL(wp) :: lvy     !< length scale for v in y direction
+    REAL(wp) :: lvz     !< length scale for v in z direction
+    REAL(wp) :: lwy     !< length scale for w in y direction
+    REAL(wp) :: lwz     !< length scale for w in z direction
+    REAL(wp) :: zz      !< height
+
+    REAL(wp),DIMENSION(nzb:nzt+1) :: r11  !< Reynolds parameter
+    REAL(wp),DIMENSION(nzb:nzt+1) :: r21  !< Reynolds parameter
+    REAL(wp),DIMENSION(nzb:nzt+1) :: r22  !< Reynolds parameter
+    REAL(wp),DIMENSION(nzb:nzt+1) :: r31  !< Reynolds parameter
+    REAL(wp),DIMENSION(nzb:nzt+1) :: r32  !< Reynolds parameter
+    REAL(wp),DIMENSION(nzb:nzt+1) :: r33  !< Reynolds parameter
+
+
+#if defined( __parallel )
+    CALL MPI_BARRIER( comm2d, ierr )
+#endif
+
+    CALL  cpu_log( log_point(57), 'synthetic_turbulence_gen', 'start' )
+
+    IF ( .NOT. ALLOCATED( mean_inflow_profiles ) )                             &
+       ALLOCATE( mean_inflow_profiles(nzb:nzt+1,5) )
+
+    ALLOCATE ( a11(nzb:nzt+1), a21(nzb:nzt+1), a22(nzb:nzt+1),                 &
+               a31(nzb:nzt+1), a32(nzb:nzt+1), a33(nzb:nzt+1),                 &
+               nuy(nzb:nzt+1), nuz(nzb:nzt+1), nvy(nzb:nzt+1),                 &
+               nvz(nzb:nzt+1), nwy(nzb:nzt+1), nwz(nzb:nzt+1),                 &
+               tu(nzb:nzt+1),  tv(nzb:nzt+1),  tw(nzb:nzt+1)   )
+
+!
+!-- Define MPI type used in stg_generate_seed_yz to gather vertical splitted
+!-- velocity seeds
+    CALL MPI_TYPE_SIZE( MPI_REAL, realsize, ierr )
+    extent = 1 * realsize
+
+    ! stg_type_yz: yz-slice with vertical bounds nzb:nzt+1
+    CALL MPI_TYPE_CREATE_SUBARRAY( 2, [nzt-nzb+2,nyng-nysg+1],                 &
+            [1,nyng-nysg+1], [0,0], MPI_ORDER_FORTRAN, MPI_REAL, newtype, ierr )
+    CALL MPI_TYPE_CREATE_RESIZED( newtype, tob, extent, stg_type_yz, ierr )
+    CALL MPI_TYPE_COMMIT( stg_type_yz, ierr )
+    CALL MPI_TYPE_FREE( newtype, ierr )
+
+    ! stg_type_yz_small: yz-slice with vertical bounds nzb_x:nzt_x+1
+    CALL MPI_TYPE_CREATE_SUBARRAY( 2, [nzt_x-nzb_x+2,nyng-nysg+1],             &
+            [1,nyng-nysg+1], [0,0], MPI_ORDER_FORTRAN, MPI_REAL, newtype, ierr )
+    CALL MPI_TYPE_CREATE_RESIZED( newtype, tob, extent, stg_type_yz_small, ierr )
+    CALL MPI_TYPE_COMMIT( stg_type_yz_small, ierr )
+    CALL MPI_TYPE_FREE( newtype, ierr )
+
+    ! receive count and displacement for MPI_GATHERV in stg_generate_seed_yz
+    ALLOCATE( recv_count(pdims(1)), displs(pdims(1)) )
+
+    recv_count           = nzt_x-nzb_x + 1
+    recv_count(pdims(1)) = recv_count(pdims(1)) + 1
+
+    DO  j = 1, pdims(1)
+       displs(j) = 0 + (nzt_x-nzb_x+1) * (j-1)
+    ENDDO
+
+!
+!-- Define seed of random number
+    CALL RANDOM_SEED()
+    CALL RANDOM_SEED( size=nseed )
+    ALLOCATE( seed(1:nseed) )
+    DO  j = 1, nseed
+       seed(j) = startseed + j
+    ENDDO
+    CALL RANDOM_SEED(put=seed)
+
+!
+!-- Read inflow profile
+!-- Height levels of profiles in input profile is as follows:
+!-- zu: luy, luz, tu, lvy, lvz, tv, r11, r21, r22, d1, d2, d5
+!-- zw: lwy, lwz, tw, r31, r32, r33, d3
+!-- WARNING: zz is not used at the moment
+    OPEN( 24, FILE='STG_PROFILES'//TRIM( coupling_char ), STATUS='OLD',        &
+                   FORM='FORMATTED')
+
+    ! Skip header
+    READ( 24, * )
+
+    DO  k = nzb, nzt+1
+       READ( 24, * ) zz, luy, luz, tu(k), lvy, lvz, tv(k), lwy, lwz, tw(k),    &
+                     r11(k), r21(k), r22(k), r31(k), r32(k), r33(k),           &
+                     d1, d2, d3, d5
+
+!
+!--    Convert length scales from meter to number of grid points
+       nuy(k) = INT( luy * ddy )
+       nuz(k) = INT( luz / dz  )
+       nvy(k) = INT( lvy * ddy )
+       nvz(k) = INT( lvz / dz  )
+       nwy(k) = INT( lwy * ddy )
+       nwz(k) = INT( lwz / dz  )
+
+!
+!--    Save Mean inflow profiles
+       IF ( TRIM( initializing_actions ) /= 'read_restart_data' ) THEN
+          mean_inflow_profiles(k,1) = d1
+          mean_inflow_profiles(k,2) = d2
+         !  mean_inflow_profiles(k,4) = d4
+          mean_inflow_profiles(k,5) = d5
+       ENDIF
+    ENDDO
+
+    CLOSE(24)
+
+!
+!--    Assign initial profiles
+    u_init = mean_inflow_profiles(:,1)
+    v_init = mean_inflow_profiles(:,2)
+   !  pt_init = mean_inflow_profiles(:,4)
+    e_init = MAXVAL( mean_inflow_profiles(:,5) )
+
+!
+!-- Calculate coefficient matrix from Reynolds stress tensor (Lund rotation)
+    DO  k = nzb, nzt+1
+       IF ( r11(k) > 0.0_wp )  THEN
+          a11(k) = SQRT( r11(k) )
+          a21(k) = r21(k) / a11(k)
+       ELSE
+          a11(k) = 0.0_wp
+          a21(k) = 0.0_wp
+       ENDIF
+
+       a22(k) = r22(k) - a21(k)**2
+       IF ( a22(k) > 0.0_wp )  THEN
+          a22(k) = SQRT( a22(k) )
+          a32(k) = ( r32(k) - a21(k) * a31(k) ) / a22(k)
+       ELSE
+          a22(k) = 0.0_wp
+          a32(k) = 0.0_wp
+       ENDIF
+
+!
+!--    a31, a32, a33 must be calculated with interpolated a11, a21, a22 (d11,
+!--    d21, d22) because of different vertical grid
+       IF ( k .le. nzt )  THEN
+          d11 = 0.5_wp * ( r11(k) + r11(k+1) )
+          IF ( d11 > 0.0_wp )  THEN
+             d11 = SQRT( d11 )
+             d21 = ( 0.5_wp * ( r21(k) + r21(k+1) ) ) / d11
+             a31(k) = r31(k) / d11
+          ELSE
+             d21 = 0.0_wp
+             a31(k) = 0.0_wp
+          ENDIF
+
+          d22 = 0.5_wp * ( r22(k) + r22(k+1) ) - d21 ** 2
+          IF ( d22 > 0.0_wp )  THEN
+             a32(k) = ( r32(k) - d21 * a31(k) ) / SQRT( d22 )
+          ELSE
+             a32(k) = 0.0_wp
+          ENDIF
+
+          a33(k) = r33(k) - a31(k) ** 2 - a32(k) ** 2
+          IF ( a33(k) > 0.0_wp )  THEN
+             a33(k) = SQRT( a33(k) )
+          ELSE
+             a33(k) = 0.0_wp
+          ENDIF
+       ELSE
+          a31(k) = a31(k-1)
+          a32(k) = a32(k-1)
+          a33(k) = a33(k-1)
+       ENDIF
+
+    ENDDO
+
+!
+!-- Define the size of the filter functions and allocate them.
+    merg = 0
+
+    ! arrays must be large enough to cover the largest length scale
+    DO  k = nzb, nzt+1
+       j = MAX( ABS(nuy(k)), ABS(nuz(k)), &
+                ABS(nvy(k)), ABS(nvz(k)), &
+                ABS(nwy(k)), ABS(nwz(k))  )
+       IF ( j > merg )  merg = j
+    ENDDO
+
+    merg  = 2 * merg
+    mergp = merg + nbgp
+
+    ALLOCATE ( buy(-merg:merg,nzb:nzt+1), buz(-merg:merg,nzb:nzt+1), &
+               bvy(-merg:merg,nzb:nzt+1), bvz(-merg:merg,nzb:nzt+1), &
+               bwy(-merg:merg,nzb:nzt+1), bwz(-merg:merg,nzb:nzt+1)  )
+
+!
+!-- Allocate velocity seeds
+    ALLOCATE ( fu( nzb:nzt+1,nysg:nyng), fuo(nzb:nzt+1,nysg:nyng), &
+               fv( nzb:nzt+1,nysg:nyng), fvo(nzb:nzt+1,nysg:nyng), &
+               fw( nzb:nzt+1,nysg:nyng), fwo(nzb:nzt+1,nysg:nyng)  )
+
+    fu  = 0._wp
+    fuo = 0._wp
+    fv  = 0._wp
+    fvo = 0._wp
+    fw  = 0._wp
+    fwo = 0._wp
+
+!
+!-- Create filter functions
+    CALL stg_filter_func( nuy, buy ) !filter uy
+    CALL stg_filter_func( nuz, buz ) !filter uz
+    CALL stg_filter_func( nvy, bvy ) !filter vy
+    CALL stg_filter_func( nvz, bvz ) !filter vz
+    CALL stg_filter_func( nwy, bwy ) !filter wy
+    CALL stg_filter_func( nwz, bwz ) !filter wz
+
+#if defined( __parallel )
+    CALL MPI_BARRIER( comm2d, ierr )
+#endif
+
+!
+!--    In case of restart, calculate velocity seeds fu, fv, fw from former
+!      time step.
+!      Bug: fu, fv, fw are different in those heights where a11, a22, a33
+!           are 0 compared to the prerun. This is mostly for k=nzt+1.
+    IF ( TRIM( initializing_actions ) == 'read_restart_data' )  THEN
+       IF ( myidx == id_inflow )  THEN
+          DO  j = nysg, nyng
+             DO  k = nzb, nzt+1
+
+                IF  ( a11(k) .NE. 0._wp ) THEN
+                fu(k,j) = ( u(k,j,-1) / mc_factor -                    &
+                            mean_inflow_profiles(k,1) ) / a11(k)
+                ELSE
+                   fu(k,j) = 0._wp
+                ENDIF
+
+                IF  ( a22(k) .NE. 0._wp ) THEN
+                fv(k,j) = ( v(k,j,-1) / mc_factor - a21(k) * fu(k,j) - &
+                            mean_inflow_profiles(k,2) ) / a22(k)
+                ELSE
+                   fv(k,j) = 0._wp
+                ENDIF
+
+                IF  ( a33(k) .NE. 0._wp ) THEN
+                fw(k,j) = ( w(k,j,-1) / mc_factor - a31(k) * fu(k,j) - &
+                            a32(k) * fv(k,j) ) / a33(k)
+                ELSE
+                   fw = 0._wp
+                ENDIF
+
+             ENDDO
+          ENDDO
+       ENDIF
+    ENDIF
+
+    CALL  cpu_log( log_point(57), 'synthetic_turbulence_gen', 'stop' )
+
+ END SUBROUTINE stg_init
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Calculate filter function bxx from length scale nxx following Eg.9 and 10
+!> (Xie and Castro, 2008)
+!------------------------------------------------------------------------------!
+ SUBROUTINE stg_filter_func( nxx, bxx )
+
+
+    IMPLICIT NONE
+
+    INTEGER(iwp) :: k         !< loop index
+    INTEGER(iwp) :: n_k       !< length scale nXX in height k
+    INTEGER(iwp) :: n_k2      !< n_k * 2
+    INTEGER(iwp) :: nf        !< index for length scales
+
+    REAL(wp) :: bdenom        !< denominator for filter functions bXX
+    REAL(wp) :: qsi = 1.0_wp  !< minimization factor
+
+    INTEGER(iwp), DIMENSION(:) :: nxx(nzb:nzt+1)           !< length scale (in gp)
+
+    REAL(wp), DIMENSION(:,:) :: bxx(-merg:merg,nzb:nzt+1)  !< filter function
+
+
+    bxx = 0.0_wp
+
+    DO  k = nzb, nzt+1
+       bdenom = 0.0_wp
+       n_k    = nxx(k)
+       IF ( n_k /= 0 )  THEN
+          n_k2 = n_k * 2
+
+!
+!--       ( Eq.10 )^2
+          DO  nf = -n_k2, n_k2
+             bdenom = bdenom + EXP( -qsi * pi * ABS(nf) / n_k )**2
+          ENDDO
+
+!
+!--       ( Eq.9 )
+          bdenom = SQRT( bdenom )
+          DO  nf = -n_k2, n_k2
+             bxx(nf,k) = EXP( -qsi * pi * ABS(nf) / n_k ) / bdenom
+          ENDDO
+       ENDIF
+    ENDDO
+
+END SUBROUTINE stg_filter_func
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Parin for &stg_par for synthetic turbulence generator
+!------------------------------------------------------------------------------!
+ SUBROUTINE stg_parin
+
+
+    IMPLICIT NONE
+
+    CHARACTER (LEN=80) ::  line   !< dummy string that contains the current line of the parameter file
+
+
+    NAMELIST /stg_par/   use_synthetic_turbulence_generator
+
+    line = ' '
+
+!
+!-- Try to find stg package
+    REWIND ( 11 )
+    line = ' '
+    DO WHILE ( INDEX( line, '&stg_par' ) == 0 )
+       READ ( 11, '(A)', END=10 )  line
+    ENDDO
+    BACKSPACE ( 11 )
+
+!
+!-- Read namelist
+    READ ( 11, stg_par )
+
+!
+!-- Set flag that indicates that the synthetic turbulence generator is switched
+!-- on
+    synthetic_turbulence_generator = .TRUE.
+
+ 10 CONTINUE
+
+ END SUBROUTINE stg_parin
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> This routine reads the respective restart data.
+!------------------------------------------------------------------------------!
+ SUBROUTINE stg_read_restart_data
+
+
+    IMPLICIT NONE
+
+    CHARACTER (LEN=30) ::  variable_chr  !< dummy variable to read string
+
+
+    READ ( 13 )  variable_chr
+    DO  WHILE ( TRIM( variable_chr ) /= '*** end stg module ***' )
+
+       SELECT CASE ( TRIM( variable_chr ) )
+
+          CASE ( 'use_synthetic_turbulence_generator ' )
+             READ ( 13 )  use_synthetic_turbulence_generator
+          CASE ( 'mc_factor' )
+             READ ( 13 )  mc_factor
+
+       END SELECT
+
+       READ ( 13 )  variable_chr
+
+    ENDDO
+
+ END SUBROUTINE stg_read_restart_data
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> This routine writes the respective restart data.
+!------------------------------------------------------------------------------!
+ SUBROUTINE stg_write_restart_data
+
+
+    IMPLICIT NONE
+
+    WRITE ( 14 )  'use_synthetic_turbulence_generator '
+    WRITE ( 14 )  use_synthetic_turbulence_generator
+    WRITE ( 14 )  'mc_factor           '
+    WRITE ( 14 )  mc_factor
+
+    WRITE ( 14 )  '*** end stg module ***            '
+
+END SUBROUTINE stg_write_restart_data
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Create turbulent inflow fields for u, v, w with prescribed length scales and
+!> Reynolds stress tensor after a method of Xie and Castro (2008), modified
+!> following suggestions of Kim et al. (2013), and using a Lund rotation
+!> (Lund, 1998).
+!------------------------------------------------------------------------------!
+ SUBROUTINE stg_main
+
+
+    USE arrays_3d,                                                             &
+        ONLY:  dzw
+
+    USE control_parameters,                                                    &
+        ONLY:  dt_3d, intermediate_timestep_count, simulated_time,             &
+               volume_flow_initial
+
+    USE indices,                                                               &
+        ONLY:  nyn, nys
+
+    USE statistics,                                                            &
+        ONLY:  weight_substep
+
+
+    IMPLICIT NONE
+
+    INTEGER(iwp) :: j           !< loop index in y-direction
+    INTEGER(iwp) :: k           !< loop index in z-direction
+
+    REAL(wp) :: dt_stg          !< wheighted subtimestep
+    REAL(wp) :: mc_factor_l     !< local mass flux correction factor
+
+    REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,5,nbgp) :: inflow_dist !< disturbance for inflow profiles
+
+
+    CALL  cpu_log( log_point(57), 'synthetic_turbulence_gen', 'start' )
+
+!
+!-- Calculate time step which is needed for filter functions
+    dt_stg = dt_3d * weight_substep(intermediate_timestep_count)
+
+!
+!-- Initial value of fu, fv, fw
+    IF ( simulated_time == 0.0_wp .AND. .NOT. velocity_seed_initialized )  THEN
+       CALL stg_generate_seed_yz( nuy, nuz, buy, buz, fu )
+       CALL stg_generate_seed_yz( nvy, nvz, bvy, bvz, fv )
+       CALL stg_generate_seed_yz( nwy, nwz, bwy, bwz, fw )
+       velocity_seed_initialized = .TRUE.
+    ENDIF
+!
+!-- New set of fu, fv, fw
+    CALL stg_generate_seed_yz( nuy, nuz, buy, buz, fuo )
+    CALL stg_generate_seed_yz( nvy, nvz, bvy, bvz, fvo )
+    CALL stg_generate_seed_yz( nwy, nwz, bwy, bwz, fwo )
+
+    IF ( myidx == id_inflow )  THEN
+
+       DO  j = nysg, nyng
+          DO  k = nzb, nzt + 1
+!
+!--          Update fu, fv, fw following Eq. 14 of Xie and Castro (2008)
+             IF ( tu(k) == 0.0_wp )  THEN
+                fu(k,j) = fuo(k,j)
+             ELSE
+                fu(k,j) = fu(k,j) * EXP( -pi * dt_stg * 0.5_wp / tu(k) ) +     &
+                         fuo(k,j) * SQRT( 1.0_wp - EXP( -pi * dt_stg / tu(k) ) )
+             ENDIF
+
+             IF ( tv(k) == 0.0_wp )  THEN
+                fv(k,j) = fvo(k,j)
+             ELSE
+                fv(k,j) = fv(k,j) * EXP( -pi * dt_stg * 0.5_wp / tv(k) ) +     &
+                         fvo(k,j) * SQRT( 1.0_wp - EXP( -pi * dt_stg / tv(k) ) )
+             ENDIF
+
+             IF ( tw(k) == 0.0_wp )  THEN
+                fw(k,j) = fwo(k,j)
+             ELSE
+                fw(k,j) = fw(k,j) * EXP( -pi * dt_stg * 0.5_wp / tw(k) ) +     &
+                         fwo(k,j) * SQRT( 1.0_wp - EXP( -pi * dt_stg / tw(k) ) )
+             ENDIF
+!
+!--          Lund rotation following Eq. 17 in Xie and Castro (2008).
+!--          Additional factors are added to improve the variance of v and w
+             IF( k == 0 )  THEN
+                inflow_dist(k,j,1,1) = 0.0_wp
+                inflow_dist(k,j,2,1) = 0.0_wp
+                inflow_dist(k,j,3,1) = 0.0_wp
+                inflow_dist(k,j,4,1) = 0.0_wp
+                inflow_dist(k,j,5,1) = 0.0_wp
+             ELSE
+                inflow_dist(k,j,1,1) = a11(k) * fu(k,j)
+                !experimental test of 1.2
+                inflow_dist(k,j,2,1) = ( SQRT( a22(k) / MAXVAL(a22) )          &
+                                         * 1.2_wp )                            &
+                                       * (   a21(k) * fu(k,j)                  &
+                                           + a22(k) * fv(k,j) )
+                inflow_dist(k,j,3,1) = ( SQRT(a33(k) / MAXVAL(a33) )           &
+                                         * 1.3_wp )                            &
+                                       * (   a31(k) * fu(k,j)                  &
+                                           + a32(k) * fv(k,j)                  &
+                                           + a33(k) * fw(k,j) )
+                ! Calculation for pt and e not yet implemented
+                inflow_dist(k,j,4,1) = 0.0_wp
+                inflow_dist(k,j,5,1) = 0.0_wp
+             ENDIF
+
+          ENDDO
+       ENDDO
+
+!
+!--    Mass flux correction following Kim et al. (2013)
+!--    This correction factor insures that the mass flux is preserved at the
+!--    inflow boundary
+       mc_factor_l = 0.0_wp
+       mc_factor   = 0.0_wp
+       DO  j = nys, nyn
+          DO  k = nzb+1, nzt
+             mc_factor_l = mc_factor_l + dzw(k)  *                             &
+                           ( mean_inflow_profiles(k,1) + inflow_dist(k,j,1,1) )
+          ENDDO
+       ENDDO
+
+#if defined( __parallel )
+       CALL MPI_ALLREDUCE( mc_factor_l, mc_factor,  &
+                           1, MPI_REAL, MPI_SUM, comm1dy, ierr )
+#else
+       mc_factor = mc_factor_l
+#endif
+
+       mc_factor = volume_flow_initial(1) / mc_factor
+
+!
+!--    Add disturbance at the inflow
+       DO  j = nysg, nyng
+          DO  k = nzb, nzt+1
+              u(k,j,-nbgp+1:0) = ( mean_inflow_profiles(k,1) +                 &
+                                   inflow_dist(k,j,1,1)      ) * mc_factor
+              v(k,j,-nbgp:-1)  = ( mean_inflow_profiles(k,2) +                 &
+                                   inflow_dist(k,j,2,1)      ) * mc_factor
+              w(k,j,-nbgp:-1)  =   inflow_dist(k,j,3,1)        * mc_factor
+          ENDDO
+       ENDDO
+
+    ENDIF
+
+    CALL  cpu_log( log_point(57), 'synthetic_turbulence_gen', 'stop' )
+
+ END SUBROUTINE stg_main
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Generate a set of random number rand_it wich is equal on each PE
+!> and calculate the velocity seed f_n.
+!> f_n is splitted in vertical direction by the number of PEs in x-direction and
+!> and each PE calculates a vertical subsection of f_n. At the the end, all
+!> parts are collected to form the full array.
+!------------------------------------------------------------------------------!
+ SUBROUTINE stg_generate_seed_yz( n_y, n_z, b_y, b_z, f_n )
+
+
+    USE indices,                                                               &
+        ONLY: ny
+
+
+    IMPLICIT NONE
+
+    INTEGER(iwp) :: j           !< loop index in y-direction
+    INTEGER(iwp) :: jj          !< loop index in y-direction
+    INTEGER(iwp) :: k           !< loop index in z-direction
+    INTEGER(iwp) :: kk          !< loop index in z-direction
+    INTEGER(iwp) :: send_count  !< send count for MPI_GATHERV
+
+    INTEGER(iwp), DIMENSION(nzb:nzt+1) :: n_y    !< length scale in y-direction
+    INTEGER(iwp), DIMENSION(nzb:nzt+1) :: n_z    !< length scale in z-direction
+    INTEGER(iwp), DIMENSION(nzb:nzt+1) :: n_y2   !< n_y*2
+    INTEGER(iwp), DIMENSION(nzb:nzt+1) :: n_z2   !< n_z*2
+
+    REAL(wp) :: nyz_inv         !< inverse of number of grid points in yz-slice
+    REAL(wp) :: rand_av         !< average of random number
+    REAL(wp) :: rand_sigma_inv  !< inverse of stdev of random number
+
+    REAL(wp), DIMENSION(-merg:merg,nzb:nzt+1)    :: b_y     !< filter func in y-dir
+    REAL(wp), DIMENSION(-merg:merg,nzb:nzt+1)    :: b_z     !< filter func in z-dir
+    REAL(wp), DIMENSION(nzb_x:nzt_x+1,nysg:nyng) :: f_n_l   !<  local velocity seed
+    REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng)     :: f_n     !<  velocity seed
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE        :: rand_it !<  random number
+
+
+!
+!-- Generate random numbers using a seed generated in stg_init.
+!-- The set of random numbers are modified to have an average of 0 and
+!-- unit variance.
+    ALLOCATE( rand_it(nzb-mergp:nzt+1+mergp,-mergp:ny+mergp) )
+
+    rand_av        = 0.0_wp
+    rand_sigma_inv = 0.0_wp
+    nyz_inv        = 1.0_wp / REAL( ( nzt+1 - nzb+1 ) * ( ny+1 ), KIND=wp )
+
+    DO  j = 0, ny
+       DO  k = nzb, nzt+1
+          CALL RANDOM_NUMBER( rand_it(k,j) )
+          rand_av = rand_av + rand_it(k,j)
+       ENDDO
+    ENDDO
+
+    rand_av = rand_av * nyz_inv
+
+    DO  j = 0, ny
+       DO  k = nzb, nzt+1
+          rand_it(k,j)   = rand_it(k,j) - rand_av
+          rand_sigma_inv = rand_sigma_inv + rand_it(k,j) ** 2
+       ENDDO
+    ENDDO
+
+    rand_sigma_inv = 1.0_wp / SQRT(rand_sigma_inv * nyz_inv)
+
+    DO  j = 0, ny
+       DO  k = nzb, nzt+1
+          rand_it(k,j) = rand_it(k,j) * rand_sigma_inv
+       ENDDO
+    ENDDO
+
+!
+!-- Periodic fill of random number in space
+    DO  j = 0, ny
+       DO  k = 1, mergp
+          rand_it(nzb  -k,j) = rand_it(nzt+2-k,j)    ! bottom margin
+          rand_it(nzt+1+k,j) = rand_it(nzb+k-1,j)    ! top margin
+       ENDDO
+    ENDDO
+    DO  j = 1, mergp
+       DO  k = nzb-mergp, nzt+1+mergp
+          rand_it(k,  -j) = rand_it(k,ny-j+1)        ! south margin
+          rand_it(k,ny+j) = rand_it(k,   j-1)        ! north margin
+       ENDDO
+    ENDDO
+
+!
+!-- Generate velocity seed following Eq.6 of Xie and Castro (2008)
+    n_y2 = n_y * 2
+    n_z2 = n_z * 2
+    f_n_l  = 0.0_wp
+
+    DO  j = nysg, nyng
+       DO  k = nzb_x, nzt_x+1
+          DO  jj = -n_y2(k), n_y2(k)
+             DO  kk = -n_z2(k), n_z2(k)
+                f_n_l(k,j) = f_n_l(k,j)                                        &
+                           + b_y(jj,k) * b_z(kk,k) * rand_it(k+kk,j+jj)
+             ENDDO
+          ENDDO
+       ENDDO
+    ENDDO
+
+    DEALLOCATE( rand_it )
+
+!
+!-- Gather velocity seeds of full subdomain
+    send_count = nzt_x - nzb_x + 1
+    IF ( nzt_x == nzt )  send_count = send_count + 1
+
+#if defined( __parallel )
+    CALL MPI_GATHERV( f_n_l(nzb_x,nysg), send_count, stg_type_yz_small,        &
+                      f_n(nzb+1,nysg), recv_count, displs, stg_type_yz,        &
+                      id_inflow, comm1dx, ierr )
+#else
+    f_n(nzb+1:nzt+1,nysg:nyng) = f_n_l(nzb_x:nzt_x+1,nysg:nyng)
+#endif
+
+
+ END SUBROUTINE stg_generate_seed_yz
+
+ END MODULE synthetic_turbulence_generator_mod
Index: /palm/trunk/SOURCE/time_integration.f90
===================================================================
--- /palm/trunk/SOURCE/time_integration.f90	(revision 2258)
+++ /palm/trunk/SOURCE/time_integration.f90	(revision 2259)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Implemented synthetic turbulence generator
+!
+! 2233 2017-05-30 18:08:54Z suehring
 !
 ! 2232 2017-05-30 17:47:52Z suehring
@@ -373,5 +376,8 @@
     USE urban_surface_mod,                                                     &
         ONLY:  usm_material_heat_model, usm_material_model,                    &
-               usm_radiation, usm_surface_energy_balance              
+               usm_radiation, usm_surface_energy_balance
+
+    USE synthetic_turbulence_generator_mod,                                    &
+        ONLY:  stg_main, use_synthetic_turbulence_generator
 
     USE user_actions_mod,                                                      &
@@ -665,4 +671,10 @@
 !--       Impose a turbulent inflow using the recycling method
           IF ( turbulent_inflow )  CALL  inflow_turbulence
+
+!
+!--       Impose a turbulent inflow using synthetic generated turbulence
+         IF ( use_synthetic_turbulence_generator ) THEN 
+            CALL  stg_main
+         ENDIF
 
 !
Index: /palm/trunk/SOURCE/write_var_list.f90
===================================================================
--- /palm/trunk/SOURCE/write_var_list.f90	(revision 2258)
+++ /palm/trunk/SOURCE/write_var_list.f90	(revision 2259)
@@ -25,23 +25,26 @@
 ! -----------------
 ! $Id$
+! Implemented synthetic turbulence generator
+!
+! 2233 2017-05-30 18:08:54Z suehring
 !
 ! 2232 2017-05-30 17:47:52Z suehring
 ! Replace wall_qflux, wall_sflux by wall_humidityflux and wall_scalarflux; add
 ! wall_salinityflux
-! +tunnel_height, tunnel_lenth, tunnel_width_x, tunnel_width_y, 
+! +tunnel_height, tunnel_lenth, tunnel_width_x, tunnel_width_y,
 ! tunnel_wall_depth
-! 
+!
 ! 2042 2016-11-02 13:47:31Z suehring
 ! Bugfix, write restart data for wall_heatflux, wall_qflux and wall_sflux
-! 
+!
 ! 2000 2016-08-20 18:09:15Z knoop
 ! Forced header and separation lines into 80 columns
-! 
+!
 ! 1992 2016-08-12 15:14:59Z suehring
-! top scalarflux added 
+! top scalarflux added
 !
 ! 1957 2016-07-07 10:43:48Z suehring
 ! flight module added
-! 
+!
 ! 1849 2016-04-08 11:33:18Z hoffmann
 ! Adapted for modularization of microphysics
@@ -50,6 +53,6 @@
 ! spectra_mod added
 !
-! 1831 2016-04-07 13:15:51Z hoffmann 
-! turbulence renamed collision_turbulence, drizzle renamed 
+! 1831 2016-04-07 13:15:51Z hoffmann
+! turbulence renamed collision_turbulence, drizzle renamed
 ! cloud_water_sedimentation
 !
@@ -59,27 +62,27 @@
 ! 1705 2015-11-02 14:28:56Z maronga
 ! Bugfix: two lines required swapping
-! 
+!
 ! 1691 2015-10-26 16:17:44Z maronga
-! Added output of most_method, constant_flux_layer, zeta_min, zeta_max. Removed 
+! Added output of most_method, constant_flux_layer, zeta_min, zeta_max. Removed
 ! output of prandtl_layer and rif_min, rif_max.
-! 
+!
 ! 1682 2015-10-07 23:56:08Z knoop
-! Code annotations made doxygen readable 
-! 
+! Code annotations made doxygen readable
+!
 ! 1585 2015-04-30 07:05:52Z maronga
 ! Adapted for RRTMG
-! 
+!
 ! 1551 2015-03-03 14:18:16Z maronga
 ! Typo removed
-! 
+!
 ! 1502 2014-12-03 18:22:31Z kanani
-! Canopy module and parameters removed (parameters are always read from 
+! Canopy module and parameters removed (parameters are always read from
 ! canopy_par NAMELIST for initial and restart runs),
-! Bugfix: added blanks in "cloud_top_radiation"-string to a total of 30 
+! Bugfix: added blanks in "cloud_top_radiation"-string to a total of 30
 ! characters
-! 
+!
 ! 1496 2014-12-02 17:25:50Z maronga
 ! Renamed "radiation" -> "cloud_top_radiation"
-! 
+!
 ! 1484 2014-10-21 10:53:05Z kanani
 ! Changes in the course of the canopy-model modularization:
@@ -90,10 +93,10 @@
 !
 ! 1324 2014-03-21 09:13:16Z suehring
-! Bugfix: ONLY statement for module netcdf_control removed 
+! Bugfix: ONLY statement for module netcdf_control removed
 !
 ! 1320 2014-03-20 08:40:49Z raasch
 ! revision history before 2012 removed,
 ! comment fields (!:) to be used for variable explanations added to
-! all variable declaration statements 
+! all variable declaration statements
 !
 ! 1308 2014-03-13 14:58:42Z fricke
@@ -116,5 +119,5 @@
 !
 ! 1053 2012-11-13 17:11:03Z hoffmann
-! necessary expansions according to the two new prognostic equations (nr, qr) 
+! necessary expansions according to the two new prognostic equations (nr, qr)
 ! of the two-moment cloud physics scheme:
 ! +bc_*_b, bc_*_t, bc_*_t_val, *_init, *_surface, *_surface_initial_change,
@@ -122,6 +125,6 @@
 ! +surface_waterflux_*
 !
-! in addition, steering parameters parameters of the two-moment cloud physics 
-! scheme:   
+! in addition, steering parameters parameters of the two-moment cloud physics
+! scheme:
 ! +cloud_scheme, +drizzle, +mu_constant, +mu_constant_value, +ventilation_effect
 !
@@ -166,5 +169,5 @@
 !------------------------------------------------------------------------------!
  SUBROUTINE write_var_list
- 
+
 
     USE arrays_3d,                                                             &
@@ -173,11 +176,11 @@
 
     USE control_parameters
-    
+
     USE flight_mod,                                                            &
         ONLY:  flight_write_restart_data
-    
+
     USE grid_variables,                                                        &
         ONLY:  dx, dy
-    
+
     USE indices,                                                               &
         ONLY:  nz, nx, ny
@@ -190,11 +193,11 @@
     USE model_1d,                                                              &
         ONLY:  damp_level_1d, dt_pr_1d, dt_run_control_1d, end_time_1d
-    
+
     USE netcdf_interface,                                                      &
         ONLY:  netcdf_precision, output_for_t0
-    
+
     USE particle_attributes,                                                   &
         ONLY:  curvature_solution_effects, time_sort_particles
-    
+
     USE pegrid
 
@@ -204,4 +207,7 @@
     USE spectra_mod,                                                           &
         ONLY:  average_count_sp
+
+    USE synthetic_turbulence_generator_mod,                                    &
+        ONLY:  stg_write_restart_data
 
     USE statistics,                                                            &
@@ -209,8 +215,8 @@
                v_max_ijk, w_max, w_max_ijk
 
-    
+
     IMPLICIT NONE
 
-    CHARACTER (LEN=10) ::  binary_version   !< 
+    CHARACTER (LEN=10) ::  binary_version   !<
 
 
@@ -531,5 +537,5 @@
     WRITE ( 14 )  s_vertical_gradient_level
     WRITE ( 14 )  's_vertical_gradient_level_ind '
-    WRITE ( 14 )  s_vertical_gradient_level_ind 
+    WRITE ( 14 )  s_vertical_gradient_level_ind
     WRITE ( 14 )  'sa_init                       '
     WRITE ( 14 )  sa_init
@@ -549,7 +555,7 @@
     WRITE ( 14 )  surface_pressure
     WRITE ( 14 )  'surface_scalarflux            '
-    WRITE ( 14 )  surface_scalarflux    
+    WRITE ( 14 )  surface_scalarflux
     WRITE ( 14 )  'surface_waterflux             '
-    WRITE ( 14 )  surface_waterflux    
+    WRITE ( 14 )  surface_waterflux
     WRITE ( 14 )  's_surface                     '
     WRITE ( 14 )  s_surface
@@ -711,5 +717,5 @@
     WRITE ( 14 )  zeta_min
     WRITE ( 14 )  'z0h_factor                    '
-    WRITE ( 14 )  z0h_factor  
+    WRITE ( 14 )  z0h_factor
 
 !
@@ -718,6 +724,10 @@
 !
 !-- If required, write restart data for virtual measurements.
-    IF ( virtual_flight )  CALL flight_write_restart_data  
-    
-    
+    IF ( virtual_flight )  CALL flight_write_restart_data
+
+!
+!-- If required, write restart data for virtual measurements.
+    IF ( synthetic_turbulence_generator )  CALL stg_write_restart_data
+
+
  END SUBROUTINE write_var_list
