source: palm/trunk/SOURCE/check_parameters.f90 @ 2007

Last change on this file since 2007 was 2007, checked in by kanani, 8 years ago

changes in the course of urban surface model implementation

  • Property svn:keywords set to Id
  • Property svn:mergeinfo set to False
    /palm/branches/forwind/SOURCE/check_parameters.f901564-1913
File size: 163.6 KB
Line 
1!> @file check_parameters.f90
2!------------------------------------------------------------------------------!
3! This file is part of PALM.
4!
5! PALM is free software: you can redistribute it and/or modify it under the
6! terms of the GNU General Public License as published by the Free Software
7! Foundation, either version 3 of the License, or (at your option) any later
8! version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 1997-2016 Leibniz Universitaet Hannover
18!------------------------------------------------------------------------------!
19!
20! Current revisions:
21! -----------------
22! Added checks for the urban surface model,
23! increased counter in DO WHILE loop over data_output (for urban surface output)
24!
25! Former revisions:
26! -----------------
27! $Id: check_parameters.f90 2007 2016-08-24 15:47:17Z kanani $
28!
29! 2000 2016-08-20 18:09:15Z knoop
30! Forced header and separation lines into 80 columns
31!
32! 1994 2016-08-15 09:52:21Z suehring
33! Add missing check for cloud_physics and cloud_droplets
34!
35! 1992 2016-08-12 15:14:59Z suehring
36! New checks for top_scalarflux
37! Bugfixes concerning data output of profiles in case of passive_scalar
38!
39! 1984 2016-08-01 14:48:05Z suehring
40! Bugfix: checking for bottom and top boundary condition for humidity and scalars
41!
42! 1972 2016-07-26 07:52:02Z maronga
43! Removed check of lai* (done in land surface module)
44!
45! 1970 2016-07-18 14:27:12Z suehring
46! Bugfix, check for ibc_q_b and constant_waterflux in case of land-surface scheme
47!
48! 1962 2016-07-13 09:16:44Z suehring
49! typo removed
50!
51! 1960 2016-07-12 16:34:24Z suehring
52! Separate checks and calculations for humidity and passive scalar. Therefore,
53! a subroutine to calculate vertical gradients is introduced, in order  to reduce
54! redundance.
55! Additional check - large-scale forcing is not implemented for passive scalars
56! so far.
57!
58! 1931 2016-06-10 12:06:59Z suehring
59! Rename multigrid into multigrid_noopt and multigrid_fast into multigrid
60!
61! 1929 2016-06-09 16:25:25Z suehring
62! Bugfix in check for use_upstream_for_tke
63!
64! 1918 2016-05-27 14:35:57Z raasch
65! setting of a fixed reference state ('single_value') for ocean runs removed
66!
67! 1914 2016-05-26 14:44:07Z witha
68! Added checks for the wind turbine model
69!
70! 1841 2016-04-07 19:14:06Z raasch
71! redundant location message removed
72!
73! 1833 2016-04-07 14:23:03Z raasch
74! check of spectra quantities moved to spectra_mod
75!
76! 1829 2016-04-07 12:16:29Z maronga
77! Bugfix: output of user defined data required reset of the string 'unit'
78!
79! 1826 2016-04-07 12:01:39Z maronga
80! Moved checks for radiation model to the respective module.
81! Moved checks for plant canopy model to the respective module.
82! Bugfix for check of too large roughness_length
83!
84! 1824 2016-04-07 09:55:29Z gronemeier
85! Check if roughness_length < dz/2. Moved location_message(finished) to the end.
86!
87! 1822 2016-04-07 07:49:42Z hoffmann
88! PALM collision kernel deleted. Collision algorithms are checked for correct
89! spelling.
90!
91! Tails removed.
92! !
93! Checks for use_sgs_for_particles adopted for the use of droplets with
94! use_sgs_for_particles adopted.
95!
96! Unused variables removed.
97!
98! 1817 2016-04-06 15:44:20Z maronga
99! Moved checks for land_surface model to the respective module
100!
101! 1806 2016-04-05 18:55:35Z gronemeier
102! Check for recycling_yshift
103!
104! 1804 2016-04-05 16:30:18Z maronga
105! Removed code for parameter file check (__check)
106!
107! 1795 2016-03-18 15:00:53Z raasch
108! Bugfix: setting of initial scalar profile in ocean
109!
110! 1788 2016-03-10 11:01:04Z maronga
111! Added check for use of most_method = 'lookup' in combination with water
112! surface presribed in the land surface model. Added output of z0q.
113! Syntax layout improved.
114!
115! 1786 2016-03-08 05:49:27Z raasch
116! cpp-direktives for spectra removed, check of spectra level removed
117!
118! 1783 2016-03-06 18:36:17Z raasch
119! netcdf variables and module name changed,
120! check of netcdf precision removed (is done in the netcdf module)
121!
122! 1764 2016-02-28 12:45:19Z raasch
123! output of nest id in run description header,
124! bugfix: check of validity of lateral boundary conditions moved to parin
125!
126! 1762 2016-02-25 12:31:13Z hellstea
127! Introduction of nested domain feature
128!
129! 1745 2016-02-05 13:06:51Z gronemeier
130! Bugfix: check data output intervals to be /= 0.0 in case of parallel NetCDF4
131!
132! 1701 2015-11-02 07:43:04Z maronga
133! Bugfix: definition of rad_net timeseries was missing
134!
135! 1691 2015-10-26 16:17:44Z maronga
136! Added output of Obukhov length (ol) and radiative heating rates for RRTMG.
137! Added checks for use of radiation / lsm with topography.
138!
139! 1682 2015-10-07 23:56:08Z knoop
140! Code annotations made doxygen readable
141!
142! 1606 2015-06-29 10:43:37Z maronga
143! Added check for use of RRTMG without netCDF.
144!
145! 1587 2015-05-04 14:19:01Z maronga
146! Added check for set of albedos when using RRTMG
147!
148! 1585 2015-04-30 07:05:52Z maronga
149! Added support for RRTMG
150!
151! 1575 2015-03-27 09:56:27Z raasch
152! multigrid_fast added as allowed pressure solver
153!
154! 1573 2015-03-23 17:53:37Z suehring
155! Bugfix: check for advection schemes in case of non-cyclic boundary conditions
156!
157! 1557 2015-03-05 16:43:04Z suehring
158! Added checks for monotonic limiter
159!
160! 1555 2015-03-04 17:44:27Z maronga
161! Added output of r_a and r_s. Renumbering of LSM PA-messages.
162!
163! 1553 2015-03-03 17:33:54Z maronga
164! Removed check for missing soil temperature profile as default values were added.
165!
166! 1551 2015-03-03 14:18:16Z maronga
167! Added various checks for land surface and radiation model. In the course of this
168! action, the length of the variable var has to be increased
169!
170! 1504 2014-12-04 09:23:49Z maronga
171! Bugfix: check for large_scale forcing got lost.
172!
173! 1500 2014-12-03 17:42:41Z maronga
174! Boundary conditions changed to dirichlet for land surface model
175!
176! 1496 2014-12-02 17:25:50Z maronga
177! Added checks for the land surface model
178!
179! 1484 2014-10-21 10:53:05Z kanani
180! Changes due to new module structure of the plant canopy model:
181!   module plant_canopy_model_mod added,
182!   checks regarding usage of new method for leaf area density profile
183!   construction added,
184!   lad-profile construction moved to new subroutine init_plant_canopy within
185!   the module plant_canopy_model_mod,
186!   drag_coefficient renamed to canopy_drag_coeff.
187! Missing KIND-attribute for REAL constant added
188!
189! 1455 2014-08-29 10:47:47Z heinze
190! empty time records in volume, cross-section and masked data output prevented 
191! in case of non-parallel netcdf-output in restart runs
192!
193! 1429 2014-07-15 12:53:45Z knoop
194! run_description_header exended to provide ensemble_member_nr if specified
195!
196! 1425 2014-07-05 10:57:53Z knoop
197! bugfix: perturbation domain modified for parallel random number generator
198!
199! 1402 2014-05-09 14:25:13Z raasch
200! location messages modified
201!
202! 1400 2014-05-09 14:03:54Z knoop
203! Check random generator extended by option random-parallel
204!
205! 1384 2014-05-02 14:31:06Z raasch
206! location messages added
207!
208! 1365 2014-04-22 15:03:56Z boeske
209! Usage of large scale forcing for pt and q enabled:
210! output for profiles of large scale advection (td_lsa_lpt, td_lsa_q),
211! large scale subsidence (td_sub_lpt, td_sub_q)
212! and nudging tendencies (td_nud_lpt, td_nud_q, td_nud_u and td_nud_v) added,
213! check if use_subsidence_tendencies is used correctly
214!
215! 1361 2014-04-16 15:17:48Z hoffmann
216! PA0363 removed
217! PA0362 changed
218!
219! 1359 2014-04-11 17:15:14Z hoffmann
220! Do not allow the execution of PALM with use_particle_tails, since particle
221! tails are currently not supported by our new particle structure.
222!
223! PA0084 not necessary for new particle structure
224!
225! 1353 2014-04-08 15:21:23Z heinze
226! REAL constants provided with KIND-attribute
227!
228! 1330 2014-03-24 17:29:32Z suehring
229! In case of SGS-particle velocity advection of TKE is also allowed with
230! dissipative 5th-order scheme.
231!
232! 1327 2014-03-21 11:00:16Z raasch
233! "baroclinicity" renamed "baroclinity", "ocean version" replaced by "ocean mode"
234! bugfix: duplicate error message 56 removed,
235! check of data_output_format and do3d_compress removed
236!
237! 1322 2014-03-20 16:38:49Z raasch
238! some REAL constants defined as wp-kind
239!
240! 1320 2014-03-20 08:40:49Z raasch
241! Kind-parameters added to all INTEGER and REAL declaration statements,
242! kinds are defined in new module kinds,
243! revision history before 2012 removed,
244! comment fields (!:) to be used for variable explanations added to
245! all variable declaration statements
246!
247! 1308 2014-03-13 14:58:42Z fricke
248! +netcdf_data_format_save
249! Calculate fixed number of output time levels for parallel netcdf output.
250! For masked data, parallel netcdf output is not tested so far, hence
251! netcdf_data_format is switched back to non-paralell output.
252!
253! 1299 2014-03-06 13:15:21Z heinze
254! enable usage of large_scale subsidence in combination with large_scale_forcing
255! output for profile of large scale vertical velocity w_subs added
256!
257! 1276 2014-01-15 13:40:41Z heinze
258! Use LSF_DATA also in case of Dirichlet bottom boundary condition for scalars
259!
260! 1241 2013-10-30 11:36:58Z heinze
261! output for profiles of ug and vg added
262! set surface_heatflux and surface_waterflux also in dependence on
263! large_scale_forcing
264! checks for nudging and large scale forcing from external file
265!
266! 1236 2013-09-27 07:21:13Z raasch
267! check number of spectra levels
268!
269! 1216 2013-08-26 09:31:42Z raasch
270! check for transpose_compute_overlap (temporary)
271!
272! 1214 2013-08-21 12:29:17Z kanani
273! additional check for simultaneous use of vertical grid stretching
274! and particle advection
275!
276! 1212 2013-08-15 08:46:27Z raasch
277! checks for poisfft_hybrid removed
278!
279! 1210 2013-08-14 10:58:20Z raasch
280! check for fftw
281!
282! 1179 2013-06-14 05:57:58Z raasch
283! checks and settings of buoyancy parameters and switches revised,
284! initial profile for rho added to hom (id=77)
285!
286! 1174 2013-05-31 10:28:08Z gryschka
287! Bugfix in computing initial profiles for ug, vg, lad, q in case of Atmosphere
288!
289! 1159 2013-05-21 11:58:22Z fricke
290! bc_lr/ns_dirneu/neudir removed
291!
292! 1115 2013-03-26 18:16:16Z hoffmann
293! unused variables removed
294! drizzle can be used without precipitation
295!
296! 1111 2013-03-08 23:54:10Z raasch
297! ibc_p_b = 2 removed
298!
299! 1103 2013-02-20 02:15:53Z raasch
300! Bugfix: turbulent inflow must not require cyclic fill in restart runs
301!
302! 1092 2013-02-02 11:24:22Z raasch
303! unused variables removed
304!
305! 1069 2012-11-28 16:18:43Z maronga
306! allow usage of topography in combination with cloud physics
307!
308! 1065 2012-11-22 17:42:36Z hoffmann
309! Bugfix: It is not allowed to use cloud_scheme = seifert_beheng without
310!         precipitation in order to save computational resources.
311!
312! 1060 2012-11-21 07:19:51Z raasch
313! additional check for parameter turbulent_inflow
314!
315! 1053 2012-11-13 17:11:03Z hoffmann
316! necessary changes for the new two-moment cloud physics scheme added:
317! - check cloud physics scheme (Kessler or Seifert and Beheng)
318! - plant_canopy is not allowed
319! - currently, only cache loop_optimization is allowed
320! - initial profiles of nr, qr
321! - boundary condition of nr, qr
322! - check output quantities (qr, nr, prr)
323!
324! 1036 2012-10-22 13:43:42Z raasch
325! code put under GPL (PALM 3.9)
326!
327! 1031/1034 2012-10-22 11:32:49Z raasch
328! check of netcdf4 parallel file support
329!
330! 1019 2012-09-28 06:46:45Z raasch
331! non-optimized version of prognostic_equations not allowed any more
332!
333! 1015 2012-09-27 09:23:24Z raasch
334! acc allowed for loop optimization,
335! checks for adjustment of mixing length to the Prandtl mixing length removed
336!
337! 1003 2012-09-14 14:35:53Z raasch
338! checks for cases with unequal subdomain sizes removed
339!
340! 1001 2012-09-13 14:08:46Z raasch
341! all actions concerning leapfrog- and upstream-spline-scheme removed
342!
343! 996 2012-09-07 10:41:47Z raasch
344! little reformatting
345
346! 978 2012-08-09 08:28:32Z fricke
347! setting of bc_lr/ns_dirneu/neudir
348! outflow damping layer removed
349! check for z0h*
350! check for pt_damping_width
351!
352! 964 2012-07-26 09:14:24Z raasch
353! check of old profil-parameters removed
354!
355! 940 2012-07-09 14:31:00Z raasch
356! checks for parameter neutral
357!
358! 924 2012-06-06 07:44:41Z maronga
359! Bugfix: preprocessor directives caused error during compilation
360!
361! 892 2012-05-02 13:51:44Z maronga
362! Bugfix for parameter file check ( excluding __netcdf4 )
363!
364! 866 2012-03-28 06:44:41Z raasch
365! use only 60% of the geostrophic wind as translation speed in case of Galilean
366! transformation and use_ug_for_galilei_tr = .T. in order to mimimize the
367! timestep
368!
369! 861 2012-03-26 14:18:34Z suehring
370! Check for topography and ws-scheme removed.
371! Check for loop_optimization = 'vector' and ws-scheme removed.
372!
373! 845 2012-03-07 10:23:05Z maronga
374! Bugfix: exclude __netcdf4 directive part from namelist file check compilation
375!
376! 828 2012-02-21 12:00:36Z raasch
377! check of collision_kernel extended
378!
379! 825 2012-02-19 03:03:44Z raasch
380! check for collision_kernel and curvature_solution_effects
381!
382! 809 2012-01-30 13:32:58Z maronga
383! Bugfix: replaced .AND. and .NOT. with && and ! in the preprocessor directives
384!
385! 807 2012-01-25 11:53:51Z maronga
386! New cpp directive "__check" implemented which is used by check_namelist_files
387!
388! Revision 1.1  1997/08/26 06:29:23  raasch
389! Initial revision
390!
391!
392! Description:
393! ------------
394!> Check control parameters and deduce further quantities.
395!------------------------------------------------------------------------------!
396 SUBROUTINE check_parameters
397 
398
399    USE arrays_3d
400    USE cloud_parameters
401    USE constants
402    USE control_parameters
403    USE dvrp_variables
404    USE grid_variables
405    USE indices
406    USE land_surface_model_mod,                                                &
407        ONLY: land_surface, lsm_check_data_output, lsm_check_data_output_pr,   &
408              lsm_check_parameters
409
410    USE kinds
411    USE model_1d
412    USE netcdf_interface,                                                      &
413        ONLY:  dopr_unit, do2d_unit, do3d_unit, netcdf_data_format,            &
414               netcdf_data_format_string
415
416    USE particle_attributes
417    USE pegrid
418    USE plant_canopy_model_mod,                                                &
419        ONLY:  pcm_check_parameters, plant_canopy
420       
421    USE pmc_interface,                                                         &
422        ONLY:  cpl_id, nested_run
423
424    USE profil_parameter
425    USE radiation_model_mod,                                                   &
426        ONLY:  radiation, radiation_check_data_output,                         &
427               radiation_check_data_output_pr, radiation_check_parameters
428    USE spectra_mod,                                                           &
429        ONLY:  calculate_spectra, spectra_check_parameters
430    USE statistics
431    USE subsidence_mod
432    USE statistics
433    USE transpose_indices
434    USE urban_surface_mod,                                                     &
435        ONLY:  urban_surface, usm_check_data_output, usm_check_parameters
436    USE wind_turbine_model_mod,                                                &
437        ONLY:  wtm_check_parameters, wind_turbine
438
439
440    IMPLICIT NONE
441
442    CHARACTER (LEN=1)   ::  sq                       !<
443    CHARACTER (LEN=15)  ::  var                      !<
444    CHARACTER (LEN=7)   ::  unit                     !<
445    CHARACTER (LEN=8)   ::  date                     !<
446    CHARACTER (LEN=10)  ::  time                     !<
447    CHARACTER (LEN=10)  ::  ensemble_string          !<
448    CHARACTER (LEN=15)  ::  nest_string              !<
449    CHARACTER (LEN=40)  ::  coupling_string          !<
450    CHARACTER (LEN=100) ::  action                   !<
451
452    INTEGER(iwp) ::  i                               !<
453    INTEGER(iwp) ::  ilen                            !<
454    INTEGER(iwp) ::  j                               !<
455    INTEGER(iwp) ::  k                               !<
456    INTEGER(iwp) ::  kk                              !<
457    INTEGER(iwp) ::  netcdf_data_format_save         !<
458    INTEGER(iwp) ::  position                        !<
459   
460    LOGICAL     ::  found                            !<
461   
462    REAL(wp)    ::  gradient                         !<
463    REAL(wp)    ::  remote = 0.0_wp                  !<
464    REAL(wp)    ::  simulation_time_since_reference  !<
465
466
467    CALL location_message( 'checking parameters', .FALSE. )
468
469!
470!-- Check for overlap combinations, which are not realized yet
471    IF ( transpose_compute_overlap )  THEN
472       IF ( numprocs == 1 )  STOP '+++ transpose-compute-overlap not implemented for single PE runs'
473#if defined( __openacc )
474       STOP '+++ transpose-compute-overlap not implemented for GPU usage'
475#endif
476    ENDIF
477
478!
479!-- Warning, if host is not set
480    IF ( host(1:1) == ' ' )  THEN
481       message_string = '"host" is not set. Please check that environment ' // &
482                        'variable "localhost" & is set before running PALM'
483       CALL message( 'check_parameters', 'PA0001', 0, 0, 0, 6, 0 )
484    ENDIF
485
486!
487!-- Check the coupling mode
488    IF ( coupling_mode /= 'uncoupled'            .AND.  &
489         coupling_mode /= 'atmosphere_to_ocean'  .AND.  &
490         coupling_mode /= 'ocean_to_atmosphere' )  THEN
491       message_string = 'illegal coupling mode: ' // TRIM( coupling_mode )
492       CALL message( 'check_parameters', 'PA0002', 1, 2, 0, 6, 0 )
493    ENDIF
494
495!
496!-- Check dt_coupling, restart_time, dt_restart, end_time, dx, dy, nx and ny
497    IF ( coupling_mode /= 'uncoupled')  THEN
498
499       IF ( dt_coupling == 9999999.9_wp )  THEN
500          message_string = 'dt_coupling is not set but required for coup' //   &
501                           'ling mode "' //  TRIM( coupling_mode ) // '"'
502          CALL message( 'check_parameters', 'PA0003', 1, 2, 0, 6, 0 )
503       ENDIF
504
505#if defined( __parallel )
506
507!
508!--    NOTE: coupled runs have not been implemented in the check_namelist_files
509!--    program.
510!--    check_namelist_files will need the following information of the other
511!--    model (atmosphere/ocean).
512!       dt_coupling = remote
513!       dt_max = remote
514!       restart_time = remote
515!       dt_restart= remote
516!       simulation_time_since_reference = remote
517!       dx = remote
518
519
520       IF ( myid == 0 ) THEN
521          CALL MPI_SEND( dt_coupling, 1, MPI_REAL, target_id, 11, comm_inter,  &
522                         ierr )
523          CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 11, comm_inter,       &
524                         status, ierr )
525       ENDIF
526       CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
527 
528       IF ( dt_coupling /= remote )  THEN
529          WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
530                 '": dt_coupling = ', dt_coupling, '& is not equal to ',       &
531                 'dt_coupling_remote = ', remote
532          CALL message( 'check_parameters', 'PA0004', 1, 2, 0, 6, 0 )
533       ENDIF
534       IF ( dt_coupling <= 0.0_wp )  THEN
535
536          IF ( myid == 0  ) THEN
537             CALL MPI_SEND( dt_max, 1, MPI_REAL, target_id, 19, comm_inter, ierr )
538             CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 19, comm_inter,    &
539                            status, ierr )
540          ENDIF   
541          CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
542     
543          dt_coupling = MAX( dt_max, remote )
544          WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
545                 '": dt_coupling <= 0.0 & is not allowed and is reset to ',    &
546                 'MAX(dt_max(A,O)) = ', dt_coupling
547          CALL message( 'check_parameters', 'PA0005', 0, 1, 0, 6, 0 )
548       ENDIF
549
550       IF ( myid == 0 ) THEN
551          CALL MPI_SEND( restart_time, 1, MPI_REAL, target_id, 12, comm_inter, &
552                         ierr )
553          CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 12, comm_inter,       &
554                         status, ierr )
555       ENDIF
556       CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
557 
558       IF ( restart_time /= remote )  THEN
559          WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
560                 '": restart_time = ', restart_time, '& is not equal to ',     &
561                 'restart_time_remote = ', remote
562          CALL message( 'check_parameters', 'PA0006', 1, 2, 0, 6, 0 )
563       ENDIF
564
565       IF ( myid == 0 ) THEN
566          CALL MPI_SEND( dt_restart, 1, MPI_REAL, target_id, 13, comm_inter,   &
567                         ierr )
568          CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 13, comm_inter,       &
569                         status, ierr )
570       ENDIF   
571       CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
572   
573       IF ( dt_restart /= remote )  THEN
574          WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
575                 '": dt_restart = ', dt_restart, '& is not equal to ',         &
576                 'dt_restart_remote = ', remote
577          CALL message( 'check_parameters', 'PA0007', 1, 2, 0, 6, 0 )
578       ENDIF
579
580       simulation_time_since_reference = end_time - coupling_start_time
581
582       IF  ( myid == 0 ) THEN
583          CALL MPI_SEND( simulation_time_since_reference, 1, MPI_REAL, target_id, &
584                         14, comm_inter, ierr )
585          CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 14, comm_inter,       &
586                         status, ierr )   
587       ENDIF
588       CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
589   
590       IF ( simulation_time_since_reference /= remote )  THEN
591          WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
592                 '": simulation_time_since_reference = ',                      &
593                 simulation_time_since_reference, '& is not equal to ',        &
594                 'simulation_time_since_reference_remote = ', remote
595          CALL message( 'check_parameters', 'PA0008', 1, 2, 0, 6, 0 )
596       ENDIF
597
598       IF ( myid == 0 ) THEN
599          CALL MPI_SEND( dx, 1, MPI_REAL, target_id, 15, comm_inter, ierr )
600          CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 15, comm_inter,       &
601                                                             status, ierr )
602       ENDIF
603       CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
604
605
606       IF ( coupling_mode == 'atmosphere_to_ocean') THEN
607
608          IF ( dx < remote ) THEN
609             WRITE( message_string, * ) 'coupling mode "',                     &
610                   TRIM( coupling_mode ),                                      &
611           '": dx in Atmosphere is not equal to or not larger then dx in ocean'
612             CALL message( 'check_parameters', 'PA0009', 1, 2, 0, 6, 0 )
613          ENDIF
614
615          IF ( (nx_a+1)*dx /= (nx_o+1)*remote )  THEN
616             WRITE( message_string, * ) 'coupling mode "',                     &
617                    TRIM( coupling_mode ),                                     &
618             '": Domain size in x-direction is not equal in ocean and atmosphere'
619             CALL message( 'check_parameters', 'PA0010', 1, 2, 0, 6, 0 )
620          ENDIF
621
622       ENDIF
623
624       IF ( myid == 0) THEN
625          CALL MPI_SEND( dy, 1, MPI_REAL, target_id, 16, comm_inter, ierr )
626          CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 16, comm_inter,       &
627                         status, ierr )
628       ENDIF
629       CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
630
631       IF ( coupling_mode == 'atmosphere_to_ocean') THEN
632
633          IF ( dy < remote )  THEN
634             WRITE( message_string, * ) 'coupling mode "',                     &
635                    TRIM( coupling_mode ),                                     &
636                 '": dy in Atmosphere is not equal to or not larger then dy in ocean'
637             CALL message( 'check_parameters', 'PA0011', 1, 2, 0, 6, 0 )
638          ENDIF
639
640          IF ( (ny_a+1)*dy /= (ny_o+1)*remote )  THEN
641             WRITE( message_string, * ) 'coupling mode "',                     &
642                   TRIM( coupling_mode ),                                      &
643             '": Domain size in y-direction is not equal in ocean and atmosphere'
644             CALL message( 'check_parameters', 'PA0012', 1, 2, 0, 6, 0 )
645          ENDIF
646
647          IF ( MOD(nx_o+1,nx_a+1) /= 0 )  THEN
648             WRITE( message_string, * ) 'coupling mode "',                     &
649                   TRIM( coupling_mode ),                                      &
650             '": nx+1 in ocean is not divisible without remainder with nx+1 in', & 
651             ' atmosphere'
652             CALL message( 'check_parameters', 'PA0339', 1, 2, 0, 6, 0 )
653          ENDIF
654
655          IF ( MOD(ny_o+1,ny_a+1) /= 0 )  THEN
656             WRITE( message_string, * ) 'coupling mode "',                     &
657                   TRIM( coupling_mode ),                                      &
658             '": ny+1 in ocean is not divisible without remainder with ny+1 in', & 
659             ' atmosphere'
660             CALL message( 'check_parameters', 'PA0340', 1, 2, 0, 6, 0 )
661          ENDIF
662
663       ENDIF
664#else
665       WRITE( message_string, * ) 'coupling requires PALM to be called with',  &
666            ' ''mrun -K parallel'''
667       CALL message( 'check_parameters', 'PA0141', 1, 2, 0, 6, 0 )
668#endif
669    ENDIF
670
671#if defined( __parallel )
672!
673!-- Exchange via intercommunicator
674    IF ( coupling_mode == 'atmosphere_to_ocean' .AND. myid == 0 )  THEN
675       CALL MPI_SEND( humidity, 1, MPI_LOGICAL, target_id, 19, comm_inter,     &
676                      ierr )
677    ELSEIF ( coupling_mode == 'ocean_to_atmosphere' .AND. myid == 0)  THEN
678       CALL MPI_RECV( humidity_remote, 1, MPI_LOGICAL, target_id, 19,          &
679                      comm_inter, status, ierr )
680    ENDIF
681    CALL MPI_BCAST( humidity_remote, 1, MPI_LOGICAL, 0, comm2d, ierr)
682   
683#endif
684
685
686!
687!-- Generate the file header which is used as a header for most of PALM's
688!-- output files
689    CALL DATE_AND_TIME( date, time )
690    run_date = date(7:8)//'-'//date(5:6)//'-'//date(3:4)
691    run_time = time(1:2)//':'//time(3:4)//':'//time(5:6)
692    IF ( coupling_mode == 'uncoupled' )  THEN
693       coupling_string = ''
694    ELSEIF ( coupling_mode == 'atmosphere_to_ocean' )  THEN
695       coupling_string = ' coupled (atmosphere)'
696    ELSEIF ( coupling_mode == 'ocean_to_atmosphere' )  THEN
697       coupling_string = ' coupled (ocean)'
698    ENDIF
699    IF ( ensemble_member_nr /= 0 )  THEN
700       WRITE( ensemble_string, '(2X,A,I2.2)' )  'en-no: ', ensemble_member_nr
701    ELSE
702       ensemble_string = ''
703    ENDIF
704    IF ( nested_run )  THEN
705       WRITE( nest_string, '(2X,A,I2.2)' )  'nest-id: ', cpl_id
706    ELSE
707       nest_string = ''
708    ENDIF
709
710    WRITE ( run_description_header,                                            &
711            '(A,2X,A,2X,A,A,A,I2.2,A,A,A,2X,A,A,2X,A,1X,A)' )                  &
712          TRIM( version ), TRIM( revision ), 'run: ',                          &
713          TRIM( run_identifier ), '.', runnr, TRIM( coupling_string ),         &
714          TRIM( nest_string ), TRIM( ensemble_string), 'host: ', TRIM( host ), &
715          run_date, run_time
716
717!
718!-- Check the general loop optimization method
719    IF ( loop_optimization == 'default' )  THEN
720       IF ( host(1:3) == 'nec' )  THEN
721          loop_optimization = 'vector'
722       ELSE
723          loop_optimization = 'cache'
724       ENDIF
725    ENDIF
726
727    SELECT CASE ( TRIM( loop_optimization ) )
728
729       CASE ( 'acc', 'cache', 'vector' )
730          CONTINUE
731
732       CASE DEFAULT
733          message_string = 'illegal value given for loop_optimization: "' //   &
734                           TRIM( loop_optimization ) // '"'
735          CALL message( 'check_parameters', 'PA0013', 1, 2, 0, 6, 0 )
736
737    END SELECT
738
739!
740!-- Check if vertical grid stretching is used together with particles
741    IF ( dz_stretch_level < 100000.0_wp .AND. particle_advection )  THEN
742       message_string = 'Vertical grid stretching is not allowed together ' // &
743                        'with particle advection.'
744       CALL message( 'check_parameters', 'PA0017', 1, 2, 0, 6, 0 )
745    ENDIF
746
747!
748!-- Check topography setting (check for illegal parameter combinations)
749    IF ( topography /= 'flat' )  THEN
750       action = ' '
751       IF ( scalar_advec /= 'pw-scheme' .AND. scalar_advec /= 'ws-scheme'      &
752      .AND. scalar_advec /= 'ws-scheme-mono' )  THEN
753          WRITE( action, '(A,A)' )  'scalar_advec = ', scalar_advec
754       ENDIF
755       IF ( momentum_advec /= 'pw-scheme' .AND. momentum_advec /= 'ws-scheme' )&
756       THEN
757          WRITE( action, '(A,A)' )  'momentum_advec = ', momentum_advec
758       ENDIF
759       IF ( psolver == 'sor' )  THEN
760          WRITE( action, '(A,A)' )  'psolver = ', psolver
761       ENDIF
762       IF ( sloping_surface )  THEN
763          WRITE( action, '(A)' )  'sloping surface = .TRUE.'
764       ENDIF
765       IF ( galilei_transformation )  THEN
766          WRITE( action, '(A)' )  'galilei_transformation = .TRUE.'
767       ENDIF
768       IF ( cloud_physics )  THEN
769          WRITE( action, '(A)' )  'cloud_physics = .TRUE.'
770       ENDIF
771       IF ( cloud_droplets )  THEN
772          WRITE( action, '(A)' )  'cloud_droplets = .TRUE.'
773       ENDIF
774       IF ( .NOT. constant_flux_layer )  THEN
775          WRITE( action, '(A)' )  'constant_flux_layer = .FALSE.'
776       ENDIF
777       IF ( action /= ' ' )  THEN
778          message_string = 'a non-flat topography does not allow ' //          &
779                           TRIM( action )
780          CALL message( 'check_parameters', 'PA0014', 1, 2, 0, 6, 0 )
781       ENDIF
782!
783!--    In case of non-flat topography, check whether the convention how to
784!--    define the topography grid has been set correctly, or whether the default
785!--    is applicable. If this is not possible, abort.
786       IF ( TRIM( topography_grid_convention ) == ' ' )  THEN
787          IF ( TRIM( topography ) /= 'single_building' .AND.                   &
788               TRIM( topography ) /= 'single_street_canyon' .AND.              &
789               TRIM( topography ) /= 'read_from_file' )  THEN
790!--          The default value is not applicable here, because it is only valid
791!--          for the two standard cases 'single_building' and 'read_from_file'
792!--          defined in init_grid.
793             WRITE( message_string, * )                                        &
794                  'The value for "topography_grid_convention" ',               &
795                  'is not set. Its default value is & only valid for ',        &
796                  '"topography" = ''single_building'', ',                      &
797                  '''single_street_canyon'' & or ''read_from_file''.',         &
798                  ' & Choose ''cell_edge'' or ''cell_center''.'
799             CALL message( 'check_parameters', 'PA0239', 1, 2, 0, 6, 0 )
800          ELSE
801!--          The default value is applicable here.
802!--          Set convention according to topography.
803             IF ( TRIM( topography ) == 'single_building' .OR.                 &
804                  TRIM( topography ) == 'single_street_canyon' )  THEN
805                topography_grid_convention = 'cell_edge'
806             ELSEIF ( TRIM( topography ) == 'read_from_file' )  THEN
807                topography_grid_convention = 'cell_center'
808             ENDIF
809          ENDIF
810       ELSEIF ( TRIM( topography_grid_convention ) /= 'cell_edge' .AND.        &
811                TRIM( topography_grid_convention ) /= 'cell_center' )  THEN
812          WRITE( message_string, * )                                           &
813               'The value for "topography_grid_convention" is ',               &
814               'not recognized. & Choose ''cell_edge'' or ''cell_center''.'
815          CALL message( 'check_parameters', 'PA0240', 1, 2, 0, 6, 0 )
816       ENDIF
817
818    ENDIF
819
820!
821!-- Check ocean setting
822    IF ( ocean )  THEN
823
824       action = ' '
825       IF ( action /= ' ' )  THEN
826          message_string = 'ocean = .T. does not allow ' // TRIM( action )
827          CALL message( 'check_parameters', 'PA0015', 1, 2, 0, 6, 0 )
828       ENDIF
829
830    ELSEIF ( TRIM( coupling_mode ) == 'uncoupled'  .AND.                       &
831             TRIM( coupling_char ) == '_O' )  THEN
832
833!
834!--    Check whether an (uncoupled) atmospheric run has been declared as an
835!--    ocean run (this setting is done via mrun-option -y)
836       message_string = 'ocean = .F. does not allow coupling_char = "' //      &
837                        TRIM( coupling_char ) // '" set by mrun-option "-y"'
838       CALL message( 'check_parameters', 'PA0317', 1, 2, 0, 6, 0 )
839
840    ENDIF
841!
842!-- Check cloud scheme
843    IF ( cloud_scheme == 'saturation_adjust' )  THEN
844       microphysics_sat_adjust = .TRUE.
845       microphysics_seifert    = .FALSE.
846       microphysics_kessler    = .FALSE.
847       precipitation           = .FALSE.
848    ELSEIF ( cloud_scheme == 'seifert_beheng' )  THEN
849       microphysics_sat_adjust = .FALSE.
850       microphysics_seifert    = .TRUE.
851       microphysics_kessler    = .FALSE.
852       precipitation           = .TRUE.
853    ELSEIF ( cloud_scheme == 'kessler' )  THEN
854       microphysics_sat_adjust = .FALSE.
855       microphysics_seifert    = .FALSE.
856       microphysics_kessler    = .TRUE.
857       precipitation           = .TRUE.
858    ELSE
859       message_string = 'unknown cloud microphysics scheme cloud_scheme ="' // &
860                        TRIM( cloud_scheme ) // '"'
861       CALL message( 'check_parameters', 'PA0357', 1, 2, 0, 6, 0 )
862    ENDIF
863!
864!-- Check whether there are any illegal values
865!-- Pressure solver:
866    IF ( psolver /= 'poisfft'  .AND.  psolver /= 'sor'  .AND.                  &
867         psolver /= 'multigrid'  .AND.  psolver /= 'multigrid_noopt' )  THEN
868       message_string = 'unknown solver for perturbation pressure: psolver' // &
869                        ' = "' // TRIM( psolver ) // '"'
870       CALL message( 'check_parameters', 'PA0016', 1, 2, 0, 6, 0 )
871    ENDIF
872
873    IF ( psolver(1:9) == 'multigrid' )  THEN
874       IF ( cycle_mg == 'w' )  THEN
875          gamma_mg = 2
876       ELSEIF ( cycle_mg == 'v' )  THEN
877          gamma_mg = 1
878       ELSE
879          message_string = 'unknown multigrid cycle: cycle_mg = "' //          &
880                           TRIM( cycle_mg ) // '"'
881          CALL message( 'check_parameters', 'PA0020', 1, 2, 0, 6, 0 )
882       ENDIF
883    ENDIF
884
885    IF ( fft_method /= 'singleton-algorithm'  .AND.                            &
886         fft_method /= 'temperton-algorithm'  .AND.                            &
887         fft_method /= 'fftw'                 .AND.                            &
888         fft_method /= 'system-specific' )  THEN
889       message_string = 'unknown fft-algorithm: fft_method = "' //             &
890                        TRIM( fft_method ) // '"'
891       CALL message( 'check_parameters', 'PA0021', 1, 2, 0, 6, 0 )
892    ENDIF
893   
894    IF( momentum_advec == 'ws-scheme' .AND.                                    & 
895        .NOT. call_psolver_at_all_substeps  ) THEN
896        message_string = 'psolver must be called at each RK3 substep when "'// &
897                      TRIM(momentum_advec) // ' "is used for momentum_advec'
898        CALL message( 'check_parameters', 'PA0344', 1, 2, 0, 6, 0 )
899    END IF
900!
901!-- Advection schemes:
902    IF ( momentum_advec /= 'pw-scheme'  .AND.  momentum_advec /= 'ws-scheme' ) &
903    THEN
904       message_string = 'unknown advection scheme: momentum_advec = "' //      &
905                        TRIM( momentum_advec ) // '"'
906       CALL message( 'check_parameters', 'PA0022', 1, 2, 0, 6, 0 )
907    ENDIF
908    IF ( ( momentum_advec == 'ws-scheme' .OR.  scalar_advec == 'ws-scheme'     &
909           .OR. scalar_advec == 'ws-scheme-mono' )                             &
910           .AND. ( timestep_scheme == 'euler' .OR.                             &
911                   timestep_scheme == 'runge-kutta-2' ) )                      &
912    THEN
913       message_string = 'momentum_advec or scalar_advec = "'                   &
914         // TRIM( momentum_advec ) // '" is not allowed with timestep_scheme = "' // &
915         TRIM( timestep_scheme ) // '"'
916       CALL message( 'check_parameters', 'PA0023', 1, 2, 0, 6, 0 )
917    ENDIF
918    IF ( scalar_advec /= 'pw-scheme'  .AND.  scalar_advec /= 'ws-scheme' .AND. &
919         scalar_advec /= 'ws-scheme-mono' .AND. scalar_advec /= 'bc-scheme' )  &
920    THEN
921       message_string = 'unknown advection scheme: scalar_advec = "' //        &
922                        TRIM( scalar_advec ) // '"'
923       CALL message( 'check_parameters', 'PA0024', 1, 2, 0, 6, 0 )
924    ENDIF
925    IF ( scalar_advec == 'bc-scheme'  .AND.  loop_optimization == 'cache' )    &
926    THEN
927       message_string = 'advection_scheme scalar_advec = "'                    &
928         // TRIM( scalar_advec ) // '" not implemented for & loop_optimization = "' // &
929         TRIM( loop_optimization ) // '"'
930       CALL message( 'check_parameters', 'PA0026', 1, 2, 0, 6, 0 )
931    ENDIF
932
933    IF ( use_sgs_for_particles  .AND.  .NOT. cloud_droplets  .AND.               &
934         .NOT. use_upstream_for_tke  .AND.                                       &
935         ( scalar_advec /= 'ws-scheme'  .AND.  scalar_advec /= 'ws-scheme-mono' )&
936       )  THEN
937       use_upstream_for_tke = .TRUE.
938       message_string = 'use_upstream_for_tke set .TRUE. because ' //          &
939                        'use_sgs_for_particles = .TRUE. '          //          &
940                        'and scalar_advec /= ws-scheme'
941       CALL message( 'check_parameters', 'PA0025', 0, 1, 0, 6, 0 )
942    ENDIF
943
944    IF ( use_sgs_for_particles  .AND.  curvature_solution_effects )  THEN
945       message_string = 'use_sgs_for_particles = .TRUE. not allowed with ' //  &
946                        'curvature_solution_effects = .TRUE.'
947       CALL message( 'check_parameters', 'PA0349', 1, 2, 0, 6, 0 )
948    ENDIF
949
950!
951!-- Set LOGICAL switches to enhance performance
952    IF ( momentum_advec == 'ws-scheme' )       ws_scheme_mom = .TRUE.
953    IF ( scalar_advec   == 'ws-scheme' .OR.                                    &
954         scalar_advec   == 'ws-scheme-mono' )  ws_scheme_sca = .TRUE.
955    IF ( scalar_advec   == 'ws-scheme-mono' )  monotonic_adjustment = .TRUE.
956
957
958!
959!-- Timestep schemes:
960    SELECT CASE ( TRIM( timestep_scheme ) )
961
962       CASE ( 'euler' )
963          intermediate_timestep_count_max = 1
964
965       CASE ( 'runge-kutta-2' )
966          intermediate_timestep_count_max = 2
967
968       CASE ( 'runge-kutta-3' )
969          intermediate_timestep_count_max = 3
970
971       CASE DEFAULT
972          message_string = 'unknown timestep scheme: timestep_scheme = "' //   &
973                           TRIM( timestep_scheme ) // '"'
974          CALL message( 'check_parameters', 'PA0027', 1, 2, 0, 6, 0 )
975
976    END SELECT
977
978    IF ( (momentum_advec /= 'pw-scheme' .AND. momentum_advec /= 'ws-scheme')   &
979         .AND. timestep_scheme(1:5) == 'runge' ) THEN
980       message_string = 'momentum advection scheme "' // &
981                        TRIM( momentum_advec ) // '" & does not work with ' // &
982                        'timestep_scheme "' // TRIM( timestep_scheme ) // '"'
983       CALL message( 'check_parameters', 'PA0029', 1, 2, 0, 6, 0 )
984    ENDIF
985!
986!-- Check for proper settings for microphysics
987    IF ( cloud_physics  .AND.  cloud_droplets )  THEN
988       message_string = 'cloud_physics = .TRUE. is not allowed with ' //  &
989                        'cloud_droplets = .TRUE.'
990       CALL message( 'check_parameters', 'PA0442', 1, 2, 0, 6, 0 )
991    ENDIF
992
993!
994!-- Collision kernels:
995    SELECT CASE ( TRIM( collision_kernel ) )
996
997       CASE ( 'hall', 'hall_fast' )
998          hall_kernel = .TRUE.
999
1000       CASE ( 'wang', 'wang_fast' )
1001          wang_kernel = .TRUE.
1002
1003       CASE ( 'none' )
1004
1005
1006       CASE DEFAULT
1007          message_string = 'unknown collision kernel: collision_kernel = "' // &
1008                           TRIM( collision_kernel ) // '"'
1009          CALL message( 'check_parameters', 'PA0350', 1, 2, 0, 6, 0 )
1010
1011    END SELECT
1012    IF ( collision_kernel(6:9) == 'fast' )  use_kernel_tables = .TRUE.
1013
1014!
1015!-- Collision algorithms:
1016    SELECT CASE ( TRIM( collision_algorithm ) )
1017
1018       CASE ( 'all_or_nothing' )
1019          all_or_nothing = .TRUE.
1020
1021       CASE ( 'average_impact' )
1022          average_impact = .TRUE.
1023
1024       CASE DEFAULT
1025          message_string = 'unknown collision algorithm: collision_algorithm = "' // &
1026                           TRIM( collision_algorithm ) // '"'
1027          CALL message( 'check_parameters', 'PA0420', 1, 2, 0, 6, 0 )
1028
1029       END SELECT
1030
1031    IF ( TRIM( initializing_actions ) /= 'read_restart_data'  .AND.            &
1032         TRIM( initializing_actions ) /= 'cyclic_fill' )  THEN
1033!
1034!--    No restart run: several initialising actions are possible
1035       action = initializing_actions
1036       DO  WHILE ( TRIM( action ) /= '' )
1037          position = INDEX( action, ' ' )
1038          SELECT CASE ( action(1:position-1) )
1039
1040             CASE ( 'set_constant_profiles', 'set_1d-model_profiles',          &
1041                    'by_user', 'initialize_vortex',     'initialize_ptanom' )
1042                action = action(position+1:)
1043
1044             CASE DEFAULT
1045                message_string = 'initializing_action = "' //                  &
1046                                 TRIM( action ) // '" unkown or not allowed'
1047                CALL message( 'check_parameters', 'PA0030', 1, 2, 0, 6, 0 )
1048
1049          END SELECT
1050       ENDDO
1051    ENDIF
1052
1053    IF ( TRIM( initializing_actions ) == 'initialize_vortex'  .AND.            &
1054         conserve_volume_flow ) THEN
1055         message_string = 'initializing_actions = "initialize_vortex"' //      &
1056                        ' ist not allowed with conserve_volume_flow = .T.'
1057       CALL message( 'check_parameters', 'PA0343', 1, 2, 0, 6, 0 )
1058    ENDIF       
1059
1060
1061    IF ( INDEX( initializing_actions, 'set_constant_profiles' ) /= 0  .AND.    &
1062         INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
1063       message_string = 'initializing_actions = "set_constant_profiles"' //    &
1064                        ' and "set_1d-model_profiles" are not allowed ' //     &
1065                        'simultaneously'
1066       CALL message( 'check_parameters', 'PA0031', 1, 2, 0, 6, 0 )
1067    ENDIF
1068
1069    IF ( INDEX( initializing_actions, 'set_constant_profiles' ) /= 0  .AND.    &
1070         INDEX( initializing_actions, 'by_user' ) /= 0 )  THEN
1071       message_string = 'initializing_actions = "set_constant_profiles"' //    &
1072                        ' and "by_user" are not allowed simultaneously'
1073       CALL message( 'check_parameters', 'PA0032', 1, 2, 0, 6, 0 )
1074    ENDIF
1075
1076    IF ( INDEX( initializing_actions, 'by_user' ) /= 0  .AND.                  &
1077         INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
1078       message_string = 'initializing_actions = "by_user" and ' //             &
1079                        '"set_1d-model_profiles" are not allowed simultaneously'
1080       CALL message( 'check_parameters', 'PA0033', 1, 2, 0, 6, 0 )
1081    ENDIF
1082
1083    IF ( cloud_physics  .AND.  .NOT.  humidity )  THEN
1084       WRITE( message_string, * ) 'cloud_physics = ', cloud_physics, ' is ',   &
1085              'not allowed with humidity = ', humidity
1086       CALL message( 'check_parameters', 'PA0034', 1, 2, 0, 6, 0 )
1087    ENDIF
1088
1089    IF ( humidity  .AND.  sloping_surface )  THEN
1090       message_string = 'humidity = .TRUE. and sloping_surface = .TRUE. ' //   &
1091                        'are not allowed simultaneously'
1092       CALL message( 'check_parameters', 'PA0036', 1, 2, 0, 6, 0 )
1093    ENDIF
1094
1095!
1096!-- When plant canopy model is used, peform addtional checks
1097    IF ( plant_canopy )  CALL pcm_check_parameters
1098   
1099!
1100!-- Additional checks for spectra
1101    IF ( calculate_spectra )  CALL spectra_check_parameters
1102!
1103!-- When land surface model is used, perform additional checks
1104    IF ( land_surface )  CALL lsm_check_parameters
1105
1106!
1107!-- When urban surface model is used, perform additional checks
1108    IF ( urban_surface )  CALL usm_check_parameters
1109
1110!
1111!-- If wind turbine model is used, perform additional checks
1112    IF ( wind_turbine )  CALL wtm_check_parameters
1113!
1114!-- When radiation model is used, peform addtional checks
1115    IF ( radiation )  CALL radiation_check_parameters
1116
1117
1118    IF ( .NOT. ( loop_optimization == 'cache'  .OR.                            &
1119                 loop_optimization == 'vector' )                               &
1120         .AND.  cloud_physics  .AND.  microphysics_seifert )  THEN
1121       message_string = 'cloud_scheme = seifert_beheng requires ' //           &
1122                        'loop_optimization = "cache" or "vector"'
1123       CALL message( 'check_parameters', 'PA0362', 1, 2, 0, 6, 0 )
1124    ENDIF 
1125
1126!
1127!-- In case of no model continuation run, check initialising parameters and
1128!-- deduce further quantities
1129    IF ( TRIM( initializing_actions ) /= 'read_restart_data' )  THEN
1130
1131!
1132!--    Initial profiles for 1D and 3D model, respectively (u,v further below)
1133       pt_init = pt_surface
1134       IF ( humidity       )  q_init  = q_surface
1135       IF ( ocean          )  sa_init = sa_surface
1136       IF ( passive_scalar )  s_init  = s_surface
1137
1138!
1139!--
1140!--    If required, compute initial profile of the geostrophic wind
1141!--    (component ug)
1142       i = 1
1143       gradient = 0.0_wp
1144
1145       IF (  .NOT.  ocean )  THEN
1146
1147          ug_vertical_gradient_level_ind(1) = 0
1148          ug(0) = ug_surface
1149          DO  k = 1, nzt+1
1150             IF ( i < 11 )  THEN
1151                IF ( ug_vertical_gradient_level(i) < zu(k)  .AND.              &
1152                     ug_vertical_gradient_level(i) >= 0.0_wp )  THEN
1153                   gradient = ug_vertical_gradient(i) / 100.0_wp
1154                   ug_vertical_gradient_level_ind(i) = k - 1
1155                   i = i + 1
1156                ENDIF
1157             ENDIF       
1158             IF ( gradient /= 0.0_wp )  THEN
1159                IF ( k /= 1 )  THEN
1160                   ug(k) = ug(k-1) + dzu(k) * gradient
1161                ELSE
1162                   ug(k) = ug_surface + dzu(k) * gradient
1163                ENDIF
1164             ELSE
1165                ug(k) = ug(k-1)
1166             ENDIF
1167          ENDDO
1168
1169       ELSE
1170
1171          ug_vertical_gradient_level_ind(1) = nzt+1
1172          ug(nzt+1) = ug_surface
1173          DO  k = nzt, nzb, -1
1174             IF ( i < 11 )  THEN
1175                IF ( ug_vertical_gradient_level(i) > zu(k)  .AND.              &
1176                     ug_vertical_gradient_level(i) <= 0.0_wp )  THEN
1177                   gradient = ug_vertical_gradient(i) / 100.0_wp
1178                   ug_vertical_gradient_level_ind(i) = k + 1
1179                   i = i + 1
1180                ENDIF
1181             ENDIF
1182             IF ( gradient /= 0.0_wp )  THEN
1183                IF ( k /= nzt )  THEN
1184                   ug(k) = ug(k+1) - dzu(k+1) * gradient
1185                ELSE
1186                   ug(k)   = ug_surface - 0.5_wp * dzu(k+1) * gradient
1187                   ug(k+1) = ug_surface + 0.5_wp * dzu(k+1) * gradient
1188                ENDIF
1189             ELSE
1190                ug(k) = ug(k+1)
1191             ENDIF
1192          ENDDO
1193
1194       ENDIF
1195
1196!
1197!--    In case of no given gradients for ug, choose a zero gradient
1198       IF ( ug_vertical_gradient_level(1) == -9999999.9_wp )  THEN
1199          ug_vertical_gradient_level(1) = 0.0_wp
1200       ENDIF 
1201
1202!
1203!--
1204!--    If required, compute initial profile of the geostrophic wind
1205!--    (component vg)
1206       i = 1
1207       gradient = 0.0_wp
1208
1209       IF (  .NOT.  ocean )  THEN
1210
1211          vg_vertical_gradient_level_ind(1) = 0
1212          vg(0) = vg_surface
1213          DO  k = 1, nzt+1
1214             IF ( i < 11 )  THEN
1215                IF ( vg_vertical_gradient_level(i) < zu(k)  .AND.              &
1216                     vg_vertical_gradient_level(i) >= 0.0_wp )  THEN
1217                   gradient = vg_vertical_gradient(i) / 100.0_wp
1218                   vg_vertical_gradient_level_ind(i) = k - 1
1219                   i = i + 1
1220                ENDIF
1221             ENDIF
1222             IF ( gradient /= 0.0_wp )  THEN
1223                IF ( k /= 1 )  THEN
1224                   vg(k) = vg(k-1) + dzu(k) * gradient
1225                ELSE
1226                   vg(k) = vg_surface + dzu(k) * gradient
1227                ENDIF
1228             ELSE
1229                vg(k) = vg(k-1)
1230             ENDIF
1231          ENDDO
1232
1233       ELSE
1234
1235          vg_vertical_gradient_level_ind(1) = nzt+1
1236          vg(nzt+1) = vg_surface
1237          DO  k = nzt, nzb, -1
1238             IF ( i < 11 )  THEN
1239                IF ( vg_vertical_gradient_level(i) > zu(k)  .AND.              &
1240                     vg_vertical_gradient_level(i) <= 0.0_wp )  THEN
1241                   gradient = vg_vertical_gradient(i) / 100.0_wp
1242                   vg_vertical_gradient_level_ind(i) = k + 1
1243                   i = i + 1
1244                ENDIF
1245             ENDIF
1246             IF ( gradient /= 0.0_wp )  THEN
1247                IF ( k /= nzt )  THEN
1248                   vg(k) = vg(k+1) - dzu(k+1) * gradient
1249                ELSE
1250                   vg(k)   = vg_surface - 0.5_wp * dzu(k+1) * gradient
1251                   vg(k+1) = vg_surface + 0.5_wp * dzu(k+1) * gradient
1252                ENDIF
1253             ELSE
1254                vg(k) = vg(k+1)
1255             ENDIF
1256          ENDDO
1257
1258       ENDIF
1259
1260!
1261!--    In case of no given gradients for vg, choose a zero gradient
1262       IF ( vg_vertical_gradient_level(1) == -9999999.9_wp )  THEN
1263          vg_vertical_gradient_level(1) = 0.0_wp
1264       ENDIF
1265
1266!
1267!--    Let the initial wind profiles be the calculated ug/vg profiles or
1268!--    interpolate them from wind profile data (if given)
1269       IF ( u_profile(1) == 9999999.9_wp  .AND.  v_profile(1) == 9999999.9_wp )  THEN
1270
1271          u_init = ug
1272          v_init = vg
1273
1274       ELSEIF ( u_profile(1) == 0.0_wp  .AND.  v_profile(1) == 0.0_wp )  THEN
1275
1276          IF ( uv_heights(1) /= 0.0_wp )  THEN
1277             message_string = 'uv_heights(1) must be 0.0'
1278             CALL message( 'check_parameters', 'PA0345', 1, 2, 0, 6, 0 )
1279          ENDIF
1280
1281          use_prescribed_profile_data = .TRUE.
1282
1283          kk = 1
1284          u_init(0) = 0.0_wp
1285          v_init(0) = 0.0_wp
1286
1287          DO  k = 1, nz+1
1288
1289             IF ( kk < 100 )  THEN
1290                DO  WHILE ( uv_heights(kk+1) <= zu(k) )
1291                   kk = kk + 1
1292                   IF ( kk == 100 )  EXIT
1293                ENDDO
1294             ENDIF
1295
1296             IF ( kk < 100  .AND.  uv_heights(kk+1) /= 9999999.9_wp )  THEN
1297                u_init(k) = u_profile(kk) + ( zu(k) - uv_heights(kk) ) /       &
1298                                       ( uv_heights(kk+1) - uv_heights(kk) ) * &
1299                                       ( u_profile(kk+1) - u_profile(kk) )
1300                v_init(k) = v_profile(kk) + ( zu(k) - uv_heights(kk) ) /       &
1301                                       ( uv_heights(kk+1) - uv_heights(kk) ) * &
1302                                       ( v_profile(kk+1) - v_profile(kk) )
1303             ELSE
1304                u_init(k) = u_profile(kk)
1305                v_init(k) = v_profile(kk)
1306             ENDIF
1307
1308          ENDDO
1309
1310       ELSE
1311
1312          message_string = 'u_profile(1) and v_profile(1) must be 0.0'
1313          CALL message( 'check_parameters', 'PA0346', 1, 2, 0, 6, 0 )
1314
1315       ENDIF
1316
1317!
1318!--    Compute initial temperature profile using the given temperature gradients
1319       IF (  .NOT.  neutral )  THEN
1320          CALL init_vertical_profiles( pt_vertical_gradient_level_ind,          &
1321                                       pt_vertical_gradient_level,              &
1322                                       pt_vertical_gradient, pt_init,           &
1323                                       pt_surface, bc_pt_t_val )
1324       ENDIF
1325!
1326!--    Compute initial humidity profile using the given humidity gradients
1327       IF ( humidity )  THEN
1328          CALL init_vertical_profiles( q_vertical_gradient_level_ind,          &
1329                                       q_vertical_gradient_level,              &
1330                                       q_vertical_gradient, q_init,            &
1331                                       q_surface, bc_q_t_val )
1332       ENDIF
1333!
1334!--    Compute initial scalar profile using the given scalar gradients
1335       IF ( passive_scalar )  THEN
1336          CALL init_vertical_profiles( s_vertical_gradient_level_ind,          &
1337                                       s_vertical_gradient_level,              &
1338                                       s_vertical_gradient, s_init,            &
1339                                       s_surface, bc_s_t_val )
1340       ENDIF
1341!
1342!--    If required, compute initial salinity profile using the given salinity
1343!--    gradients
1344       IF ( ocean )  THEN
1345          CALL init_vertical_profiles( sa_vertical_gradient_level_ind,          &
1346                                       sa_vertical_gradient_level,              &
1347                                       sa_vertical_gradient, s_init,            &
1348                                       sa_surface, -999.0_wp )
1349       ENDIF
1350
1351         
1352    ENDIF
1353
1354!
1355!-- Check if the control parameter use_subsidence_tendencies is used correctly
1356    IF ( use_subsidence_tendencies  .AND.  .NOT.  large_scale_subsidence )  THEN
1357       message_string = 'The usage of use_subsidence_tendencies ' //           &
1358                            'requires large_scale_subsidence = .T..'
1359       CALL message( 'check_parameters', 'PA0396', 1, 2, 0, 6, 0 )
1360    ELSEIF ( use_subsidence_tendencies  .AND.  .NOT. large_scale_forcing )  THEN
1361       message_string = 'The usage of use_subsidence_tendencies ' //           &
1362                            'requires large_scale_forcing = .T..'
1363       CALL message( 'check_parameters', 'PA0397', 1, 2, 0, 6, 0 )
1364    ENDIF
1365
1366!
1367!-- Initialize large scale subsidence if required
1368    If ( large_scale_subsidence )  THEN
1369       IF ( subs_vertical_gradient_level(1) /= -9999999.9_wp  .AND.            &
1370                                     .NOT.  large_scale_forcing )  THEN
1371          CALL init_w_subsidence
1372       ENDIF
1373!
1374!--    In case large_scale_forcing is used, profiles for subsidence velocity
1375!--    are read in from file LSF_DATA
1376
1377       IF ( subs_vertical_gradient_level(1) == -9999999.9_wp  .AND.            &
1378            .NOT.  large_scale_forcing )  THEN
1379          message_string = 'There is no default large scale vertical ' //      &
1380                           'velocity profile set. Specify the subsidence ' //  &
1381                           'velocity profile via subs_vertical_gradient ' //   &
1382                           'and subs_vertical_gradient_level.'
1383          CALL message( 'check_parameters', 'PA0380', 1, 2, 0, 6, 0 )
1384       ENDIF
1385    ELSE
1386        IF ( subs_vertical_gradient_level(1) /= -9999999.9_wp )  THEN
1387           message_string = 'Enable usage of large scale subsidence by ' //    &
1388                            'setting large_scale_subsidence = .T..'
1389          CALL message( 'check_parameters', 'PA0381', 1, 2, 0, 6, 0 )
1390        ENDIF
1391    ENDIF   
1392
1393!
1394!-- Compute Coriolis parameter
1395    f  = 2.0_wp * omega * SIN( phi / 180.0_wp * pi )
1396    fs = 2.0_wp * omega * COS( phi / 180.0_wp * pi )
1397
1398!
1399!-- Check and set buoyancy related parameters and switches
1400    IF ( reference_state == 'horizontal_average' )  THEN
1401       CONTINUE
1402    ELSEIF ( reference_state == 'initial_profile' )  THEN
1403       use_initial_profile_as_reference = .TRUE.
1404    ELSEIF ( reference_state == 'single_value' )  THEN
1405       use_single_reference_value = .TRUE.
1406       IF ( pt_reference == 9999999.9_wp )  pt_reference = pt_surface
1407       vpt_reference = pt_reference * ( 1.0_wp + 0.61_wp * q_surface )
1408    ELSE
1409       message_string = 'illegal value for reference_state: "' //              &
1410                        TRIM( reference_state ) // '"'
1411       CALL message( 'check_parameters', 'PA0056', 1, 2, 0, 6, 0 )
1412    ENDIF
1413
1414!
1415!-- Sign of buoyancy/stability terms
1416    IF ( ocean )  atmos_ocean_sign = -1.0_wp
1417
1418!
1419!-- Ocean version must use flux boundary conditions at the top
1420    IF ( ocean .AND. .NOT. use_top_fluxes )  THEN
1421       message_string = 'use_top_fluxes must be .TRUE. in ocean mode'
1422       CALL message( 'check_parameters', 'PA0042', 1, 2, 0, 6, 0 )
1423    ENDIF
1424
1425!
1426!-- In case of a given slope, compute the relevant quantities
1427    IF ( alpha_surface /= 0.0_wp )  THEN
1428       IF ( ABS( alpha_surface ) > 90.0_wp )  THEN
1429          WRITE( message_string, * ) 'ABS( alpha_surface = ', alpha_surface,   &
1430                                     ' ) must be < 90.0'
1431          CALL message( 'check_parameters', 'PA0043', 1, 2, 0, 6, 0 )
1432       ENDIF
1433       sloping_surface = .TRUE.
1434       cos_alpha_surface = COS( alpha_surface / 180.0_wp * pi )
1435       sin_alpha_surface = SIN( alpha_surface / 180.0_wp * pi )
1436    ENDIF
1437
1438!
1439!-- Check time step and cfl_factor
1440    IF ( dt /= -1.0_wp )  THEN
1441       IF ( dt <= 0.0_wp  .AND.  dt /= -1.0_wp )  THEN
1442          WRITE( message_string, * ) 'dt = ', dt , ' <= 0.0'
1443          CALL message( 'check_parameters', 'PA0044', 1, 2, 0, 6, 0 )
1444       ENDIF
1445       dt_3d = dt
1446       dt_fixed = .TRUE.
1447    ENDIF
1448
1449    IF ( cfl_factor <= 0.0_wp  .OR.  cfl_factor > 1.0_wp )  THEN
1450       IF ( cfl_factor == -1.0_wp )  THEN
1451          IF ( timestep_scheme == 'runge-kutta-2' )  THEN
1452             cfl_factor = 0.8_wp
1453          ELSEIF ( timestep_scheme == 'runge-kutta-3' )  THEN
1454             cfl_factor = 0.9_wp
1455          ELSE
1456             cfl_factor = 0.9_wp
1457          ENDIF
1458       ELSE
1459          WRITE( message_string, * ) 'cfl_factor = ', cfl_factor,              &
1460                 ' out of range & 0.0 < cfl_factor <= 1.0 is required'
1461          CALL message( 'check_parameters', 'PA0045', 1, 2, 0, 6, 0 )
1462       ENDIF
1463    ENDIF
1464
1465!
1466!-- Store simulated time at begin
1467    simulated_time_at_begin = simulated_time
1468
1469!
1470!-- Store reference time for coupled runs and change the coupling flag,
1471!-- if ...
1472    IF ( simulated_time == 0.0_wp )  THEN
1473       IF ( coupling_start_time == 0.0_wp )  THEN
1474          time_since_reference_point = 0.0_wp
1475       ELSEIF ( time_since_reference_point < 0.0_wp )  THEN
1476          run_coupled = .FALSE.
1477       ENDIF
1478    ENDIF
1479
1480!
1481!-- Set wind speed in the Galilei-transformed system
1482    IF ( galilei_transformation )  THEN
1483       IF ( use_ug_for_galilei_tr                    .AND.                     &
1484            ug_vertical_gradient_level(1) == 0.0_wp  .AND.                     &
1485            ug_vertical_gradient(1) == 0.0_wp        .AND.                     & 
1486            vg_vertical_gradient_level(1) == 0.0_wp  .AND.                     &
1487            vg_vertical_gradient(1) == 0.0_wp )  THEN
1488          u_gtrans = ug_surface * 0.6_wp
1489          v_gtrans = vg_surface * 0.6_wp
1490       ELSEIF ( use_ug_for_galilei_tr  .AND.                                   &
1491                ( ug_vertical_gradient_level(1) /= 0.0_wp  .OR.                &
1492                ug_vertical_gradient(1) /= 0.0_wp ) )  THEN
1493          message_string = 'baroclinity (ug) not allowed simultaneously' //    &
1494                           ' with galilei transformation'
1495          CALL message( 'check_parameters', 'PA0046', 1, 2, 0, 6, 0 )
1496       ELSEIF ( use_ug_for_galilei_tr  .AND.                                   &
1497                ( vg_vertical_gradient_level(1) /= 0.0_wp  .OR.                &
1498                vg_vertical_gradient(1) /= 0.0_wp ) )  THEN
1499          message_string = 'baroclinity (vg) not allowed simultaneously' //    &
1500                           ' with galilei transformation'
1501          CALL message( 'check_parameters', 'PA0047', 1, 2, 0, 6, 0 )
1502       ELSE
1503          message_string = 'variable translation speed used for galilei-' //   &
1504             'transformation, which may cause & instabilities in stably ' //   &
1505             'stratified regions'
1506          CALL message( 'check_parameters', 'PA0048', 0, 1, 0, 6, 0 )
1507       ENDIF
1508    ENDIF
1509
1510!
1511!-- In case of using a prandtl-layer, calculated (or prescribed) surface
1512!-- fluxes have to be used in the diffusion-terms
1513    IF ( constant_flux_layer )  use_surface_fluxes = .TRUE.
1514
1515!
1516!-- Check boundary conditions and set internal variables:
1517!-- Attention: the lateral boundary conditions have been already checked in
1518!-- parin
1519!
1520!-- Non-cyclic lateral boundaries require the multigrid method and Piascek-
1521!-- Willimas or Wicker - Skamarock advection scheme. Several schemes
1522!-- and tools do not work with non-cyclic boundary conditions.
1523    IF ( bc_lr /= 'cyclic'  .OR.  bc_ns /= 'cyclic' )  THEN
1524       IF ( psolver(1:9) /= 'multigrid' )  THEN
1525          message_string = 'non-cyclic lateral boundaries do not allow ' //    &
1526                           'psolver = "' // TRIM( psolver ) // '"'
1527          CALL message( 'check_parameters', 'PA0051', 1, 2, 0, 6, 0 )
1528       ENDIF
1529       IF ( momentum_advec /= 'pw-scheme'  .AND.                               &
1530          ( momentum_advec /= 'ws-scheme'  .AND.                               &
1531            momentum_advec /= 'ws-scheme-mono' )                               &
1532          )  THEN
1533
1534          message_string = 'non-cyclic lateral boundaries do not allow ' //    &
1535                           'momentum_advec = "' // TRIM( momentum_advec ) // '"'
1536          CALL message( 'check_parameters', 'PA0052', 1, 2, 0, 6, 0 )
1537       ENDIF
1538       IF ( scalar_advec /= 'pw-scheme'  .AND.                                 &
1539          ( scalar_advec /= 'ws-scheme'  .AND.                                 &
1540            scalar_advec /= 'ws-scheme-mono' )                                 &
1541          )  THEN
1542          message_string = 'non-cyclic lateral boundaries do not allow ' //    &
1543                           'scalar_advec = "' // TRIM( scalar_advec ) // '"'
1544          CALL message( 'check_parameters', 'PA0053', 1, 2, 0, 6, 0 )
1545       ENDIF
1546       IF ( galilei_transformation )  THEN
1547          message_string = 'non-cyclic lateral boundaries do not allow ' //    &
1548                           'galilei_transformation = .T.'
1549          CALL message( 'check_parameters', 'PA0054', 1, 2, 0, 6, 0 )
1550       ENDIF
1551    ENDIF
1552
1553!
1554!-- Bottom boundary condition for the turbulent Kinetic energy
1555    IF ( bc_e_b == 'neumann' )  THEN
1556       ibc_e_b = 1
1557    ELSEIF ( bc_e_b == '(u*)**2+neumann' )  THEN
1558       ibc_e_b = 2
1559       IF ( .NOT. constant_flux_layer )  THEN
1560          bc_e_b = 'neumann'
1561          ibc_e_b = 1
1562          message_string = 'boundary condition bc_e_b changed to "' //         &
1563                           TRIM( bc_e_b ) // '"'
1564          CALL message( 'check_parameters', 'PA0057', 0, 1, 0, 6, 0 )
1565       ENDIF
1566    ELSE
1567       message_string = 'unknown boundary condition: bc_e_b = "' //            &
1568                        TRIM( bc_e_b ) // '"'
1569       CALL message( 'check_parameters', 'PA0058', 1, 2, 0, 6, 0 )
1570    ENDIF
1571
1572!
1573!-- Boundary conditions for perturbation pressure
1574    IF ( bc_p_b == 'dirichlet' )  THEN
1575       ibc_p_b = 0
1576    ELSEIF ( bc_p_b == 'neumann' )  THEN
1577       ibc_p_b = 1
1578    ELSE
1579       message_string = 'unknown boundary condition: bc_p_b = "' //            &
1580                        TRIM( bc_p_b ) // '"'
1581       CALL message( 'check_parameters', 'PA0059', 1, 2, 0, 6, 0 )
1582    ENDIF
1583
1584    IF ( bc_p_t == 'dirichlet' )  THEN
1585       ibc_p_t = 0
1586!-- TO_DO: later set bc_p_t to neumann before, in case of nested domain
1587    ELSEIF ( bc_p_t == 'neumann' .OR. bc_p_t == 'nested' )  THEN
1588       ibc_p_t = 1
1589    ELSE
1590       message_string = 'unknown boundary condition: bc_p_t = "' //            &
1591                        TRIM( bc_p_t ) // '"'
1592       CALL message( 'check_parameters', 'PA0061', 1, 2, 0, 6, 0 )
1593    ENDIF
1594
1595!
1596!-- Boundary conditions for potential temperature
1597    IF ( coupling_mode == 'atmosphere_to_ocean' )  THEN
1598       ibc_pt_b = 2
1599    ELSE
1600       IF ( bc_pt_b == 'dirichlet' )  THEN
1601          ibc_pt_b = 0
1602       ELSEIF ( bc_pt_b == 'neumann' )  THEN
1603          ibc_pt_b = 1
1604       ELSE
1605          message_string = 'unknown boundary condition: bc_pt_b = "' //        &
1606                           TRIM( bc_pt_b ) // '"'
1607          CALL message( 'check_parameters', 'PA0062', 1, 2, 0, 6, 0 )
1608       ENDIF
1609    ENDIF
1610
1611    IF ( bc_pt_t == 'dirichlet' )  THEN
1612       ibc_pt_t = 0
1613    ELSEIF ( bc_pt_t == 'neumann' )  THEN
1614       ibc_pt_t = 1
1615    ELSEIF ( bc_pt_t == 'initial_gradient' )  THEN
1616       ibc_pt_t = 2
1617    ELSEIF ( bc_pt_t == 'nested' )  THEN
1618       ibc_pt_t = 3
1619    ELSE
1620       message_string = 'unknown boundary condition: bc_pt_t = "' //           &
1621                        TRIM( bc_pt_t ) // '"'
1622       CALL message( 'check_parameters', 'PA0063', 1, 2, 0, 6, 0 )
1623    ENDIF
1624
1625!
1626!   This IF clause needs revision, got too complex!!
1627    IF ( surface_heatflux == 9999999.9_wp  )  THEN
1628       constant_heatflux = .FALSE.
1629       IF ( large_scale_forcing  .OR.  land_surface  .OR.  urban_surface )  THEN
1630          IF ( ibc_pt_b == 0 )  THEN
1631             constant_heatflux = .FALSE.
1632          ELSEIF ( ibc_pt_b == 1 )  THEN
1633             constant_heatflux = .TRUE.
1634             IF ( TRIM( initializing_actions ) /= 'read_restart_data'  .AND.   &
1635                  .NOT.  land_surface )  THEN
1636                surface_heatflux = shf_surf(1)
1637             ELSE
1638                surface_heatflux = 0.0_wp
1639             ENDIF
1640          ENDIF
1641       ENDIF
1642    ELSE
1643        constant_heatflux = .TRUE.
1644        IF ( TRIM( initializing_actions ) /= 'read_restart_data'  .AND.        &
1645             large_scale_forcing  .AND.  .NOT.  land_surface )  THEN
1646           surface_heatflux = shf_surf(1)
1647        ENDIF
1648    ENDIF
1649
1650    IF ( top_heatflux     == 9999999.9_wp )  constant_top_heatflux = .FALSE.
1651
1652    IF ( neutral )  THEN
1653
1654       IF ( surface_heatflux /= 0.0_wp  .AND.                                  &
1655            surface_heatflux /= 9999999.9_wp )  THEN
1656          message_string = 'heatflux must not be set for pure neutral flow'
1657          CALL message( 'check_parameters', 'PA0351', 1, 2, 0, 6, 0 )
1658       ENDIF
1659
1660       IF ( top_heatflux /= 0.0_wp  .AND.  top_heatflux /= 9999999.9_wp )      &
1661       THEN
1662          message_string = 'heatflux must not be set for pure neutral flow'
1663          CALL message( 'check_parameters', 'PA0351', 1, 2, 0, 6, 0 )
1664       ENDIF
1665
1666    ENDIF
1667
1668    IF ( top_momentumflux_u /= 9999999.9_wp  .AND.                             &
1669         top_momentumflux_v /= 9999999.9_wp )  THEN
1670       constant_top_momentumflux = .TRUE.
1671    ELSEIF (  .NOT. ( top_momentumflux_u == 9999999.9_wp  .AND.                &
1672           top_momentumflux_v == 9999999.9_wp ) )  THEN
1673       message_string = 'both, top_momentumflux_u AND top_momentumflux_v ' //  &
1674                        'must be set'
1675       CALL message( 'check_parameters', 'PA0064', 1, 2, 0, 6, 0 )
1676    ENDIF
1677
1678!
1679!-- A given surface temperature implies Dirichlet boundary condition for
1680!-- temperature. In this case specification of a constant heat flux is
1681!-- forbidden.
1682    IF ( ibc_pt_b == 0  .AND.  constant_heatflux  .AND.                        &
1683         surface_heatflux /= 0.0_wp )  THEN
1684       message_string = 'boundary_condition: bc_pt_b = "' // TRIM( bc_pt_b ) //&
1685                        '& is not allowed with constant_heatflux = .TRUE.'
1686       CALL message( 'check_parameters', 'PA0065', 1, 2, 0, 6, 0 )
1687    ENDIF
1688    IF ( constant_heatflux  .AND.  pt_surface_initial_change /= 0.0_wp )  THEN
1689       WRITE ( message_string, * )  'constant_heatflux = .TRUE. is not allo',  &
1690               'wed with pt_surface_initial_change (/=0) = ',                  &
1691               pt_surface_initial_change
1692       CALL message( 'check_parameters', 'PA0066', 1, 2, 0, 6, 0 )
1693    ENDIF
1694
1695!
1696!-- A given temperature at the top implies Dirichlet boundary condition for
1697!-- temperature. In this case specification of a constant heat flux is
1698!-- forbidden.
1699    IF ( ibc_pt_t == 0  .AND.  constant_top_heatflux  .AND.                    &
1700         top_heatflux /= 0.0_wp )  THEN
1701       message_string = 'boundary_condition: bc_pt_t = "' // TRIM( bc_pt_t ) //&
1702                        '" is not allowed with constant_top_heatflux = .TRUE.'
1703       CALL message( 'check_parameters', 'PA0067', 1, 2, 0, 6, 0 )
1704    ENDIF
1705
1706!
1707!-- Boundary conditions for salinity
1708    IF ( ocean )  THEN
1709       IF ( bc_sa_t == 'dirichlet' )  THEN
1710          ibc_sa_t = 0
1711       ELSEIF ( bc_sa_t == 'neumann' )  THEN
1712          ibc_sa_t = 1
1713       ELSE
1714          message_string = 'unknown boundary condition: bc_sa_t = "' //        &
1715                           TRIM( bc_sa_t ) // '"'
1716          CALL message( 'check_parameters', 'PA0068', 1, 2, 0, 6, 0 )
1717       ENDIF
1718
1719       IF ( top_salinityflux == 9999999.9_wp )  constant_top_salinityflux = .FALSE.
1720       IF ( ibc_sa_t == 1  .AND.  top_salinityflux == 9999999.9_wp )  THEN
1721          message_string = 'boundary condition: bc_sa_t = "' //                &
1722                           TRIM( bc_sa_t ) // '" requires to set ' //          &
1723                           'top_salinityflux'
1724          CALL message( 'check_parameters', 'PA0069', 1, 2, 0, 6, 0 )
1725       ENDIF
1726
1727!
1728!--    A fixed salinity at the top implies Dirichlet boundary condition for
1729!--    salinity. In this case specification of a constant salinity flux is
1730!--    forbidden.
1731       IF ( ibc_sa_t == 0  .AND.  constant_top_salinityflux  .AND.             &
1732            top_salinityflux /= 0.0_wp )  THEN
1733          message_string = 'boundary condition: bc_sa_t = "' //                &
1734                           TRIM( bc_sa_t ) // '" is not allowed with ' //      &
1735                           'top_salinityflux /= 0.0'
1736          CALL message( 'check_parameters', 'PA0070', 1, 2, 0, 6, 0 )
1737       ENDIF
1738
1739    ENDIF
1740
1741!
1742!-- Set boundary conditions for total water content
1743    IF ( humidity )  THEN
1744       CALL set_bc_scalars( 'q', bc_q_b, bc_q_t, ibc_q_b, ibc_q_t,           &
1745                            'PA0071', 'PA0072' )
1746
1747       IF ( surface_waterflux == 9999999.9_wp  )  THEN
1748          constant_waterflux = .FALSE.
1749          IF ( large_scale_forcing .OR. land_surface )  THEN
1750             IF ( ibc_q_b == 0 )  THEN
1751                constant_waterflux = .FALSE.
1752             ELSEIF ( ibc_q_b == 1 )  THEN
1753                constant_waterflux = .TRUE.
1754                IF ( TRIM( initializing_actions ) /= 'read_restart_data' )  THEN
1755                   surface_waterflux = qsws_surf(1)
1756                ENDIF
1757             ENDIF
1758          ENDIF
1759       ELSE
1760          constant_waterflux = .TRUE.
1761          IF ( TRIM( initializing_actions ) /= 'read_restart_data'  .AND.      &
1762               large_scale_forcing )  THEN
1763             surface_waterflux = qsws_surf(1)
1764          ENDIF
1765       ENDIF
1766
1767       CALL check_bc_scalars( 'q', bc_q_b, ibc_q_b, 'PA0073', 'PA0074',        &
1768                              constant_waterflux, q_surface_initial_change )
1769
1770    ENDIF
1771   
1772    IF ( passive_scalar )  THEN
1773       CALL set_bc_scalars( 's', bc_s_b, bc_s_t, ibc_s_b, ibc_s_t,             &
1774                            'PA0071', 'PA0072' )
1775
1776       CALL check_bc_scalars( 's', bc_s_b, ibc_s_b, 'PA0073', 'PA0074',        &
1777                              constant_scalarflux, s_surface_initial_change )
1778
1779       IF ( top_scalarflux == 9999999.9_wp )  constant_top_scalarflux = .FALSE.
1780!
1781!--    A fixed scalar concentration at the top implies Dirichlet boundary
1782!--    condition for scalar. Hence, in this case specification of a constant
1783!--    scalar flux is forbidden.
1784       IF ( ( ibc_s_t == 0 .OR. ibc_s_t == 2 )  .AND.  constant_top_scalarflux &
1785               .AND.  top_scalarflux /= 0.0_wp )  THEN
1786          message_string = 'boundary condition: bc_s_t = "' //                 &
1787                           TRIM( bc_sa_t ) // '" is not allowed with ' //      &
1788                           'top_scalarflux /= 0.0'
1789          CALL message( 'check_parameters', 'PA0441', 1, 2, 0, 6, 0 )
1790       ENDIF
1791    ENDIF
1792!
1793!-- Boundary conditions for horizontal components of wind speed
1794    IF ( bc_uv_b == 'dirichlet' )  THEN
1795       ibc_uv_b = 0
1796    ELSEIF ( bc_uv_b == 'neumann' )  THEN
1797       ibc_uv_b = 1
1798       IF ( constant_flux_layer )  THEN
1799          message_string = 'boundary condition: bc_uv_b = "' //                &
1800               TRIM( bc_uv_b ) // '" is not allowed with constant_flux_layer'  &
1801               // ' = .TRUE.'
1802          CALL message( 'check_parameters', 'PA0075', 1, 2, 0, 6, 0 )
1803       ENDIF
1804    ELSE
1805       message_string = 'unknown boundary condition: bc_uv_b = "' //           &
1806                        TRIM( bc_uv_b ) // '"'
1807       CALL message( 'check_parameters', 'PA0076', 1, 2, 0, 6, 0 )
1808    ENDIF
1809!
1810!-- In case of coupled simulations u and v at the ground in atmosphere will be
1811!-- assigned with the u and v values of the ocean surface
1812    IF ( coupling_mode == 'atmosphere_to_ocean' )  THEN
1813       ibc_uv_b = 2
1814    ENDIF
1815
1816    IF ( coupling_mode == 'ocean_to_atmosphere' )  THEN
1817       bc_uv_t = 'neumann'
1818       ibc_uv_t = 1
1819    ELSE
1820       IF ( bc_uv_t == 'dirichlet' .OR. bc_uv_t == 'dirichlet_0' )  THEN
1821          ibc_uv_t = 0
1822          IF ( bc_uv_t == 'dirichlet_0' )  THEN
1823!
1824!--          Velocities for the initial u,v-profiles are set zero at the top
1825!--          in case of dirichlet_0 conditions
1826             u_init(nzt+1)    = 0.0_wp
1827             v_init(nzt+1)    = 0.0_wp
1828          ENDIF
1829       ELSEIF ( bc_uv_t == 'neumann' )  THEN
1830          ibc_uv_t = 1
1831       ELSEIF ( bc_uv_t == 'nested' )  THEN
1832          ibc_uv_t = 3
1833       ELSE
1834          message_string = 'unknown boundary condition: bc_uv_t = "' //        &
1835                           TRIM( bc_uv_t ) // '"'
1836          CALL message( 'check_parameters', 'PA0077', 1, 2, 0, 6, 0 )
1837       ENDIF
1838    ENDIF
1839
1840!
1841!-- Compute and check, respectively, the Rayleigh Damping parameter
1842    IF ( rayleigh_damping_factor == -1.0_wp )  THEN
1843       rayleigh_damping_factor = 0.0_wp
1844    ELSE
1845       IF ( rayleigh_damping_factor < 0.0_wp  .OR.                             &
1846            rayleigh_damping_factor > 1.0_wp )  THEN
1847          WRITE( message_string, * )  'rayleigh_damping_factor = ',            &
1848                              rayleigh_damping_factor, ' out of range [0.0,1.0]'
1849          CALL message( 'check_parameters', 'PA0078', 1, 2, 0, 6, 0 )
1850       ENDIF
1851    ENDIF
1852
1853    IF ( rayleigh_damping_height == -1.0_wp )  THEN
1854       IF (  .NOT.  ocean )  THEN
1855          rayleigh_damping_height = 0.66666666666_wp * zu(nzt)
1856       ELSE
1857          rayleigh_damping_height = 0.66666666666_wp * zu(nzb)
1858       ENDIF
1859    ELSE
1860       IF (  .NOT.  ocean )  THEN
1861          IF ( rayleigh_damping_height < 0.0_wp  .OR.                          &
1862               rayleigh_damping_height > zu(nzt) )  THEN
1863             WRITE( message_string, * )  'rayleigh_damping_height = ',         &
1864                   rayleigh_damping_height, ' out of range [0.0,', zu(nzt), ']'
1865             CALL message( 'check_parameters', 'PA0079', 1, 2, 0, 6, 0 )
1866          ENDIF
1867       ELSE
1868          IF ( rayleigh_damping_height > 0.0_wp  .OR.                          &
1869               rayleigh_damping_height < zu(nzb) )  THEN
1870             WRITE( message_string, * )  'rayleigh_damping_height = ',         &
1871                   rayleigh_damping_height, ' out of range [0.0,', zu(nzb), ']'
1872             CALL message( 'check_parameters', 'PA0079', 1, 2, 0, 6, 0 )
1873          ENDIF
1874       ENDIF
1875    ENDIF
1876
1877!
1878!-- Check number of chosen statistic regions. More than 10 regions are not
1879!-- allowed, because so far no more than 10 corresponding output files can
1880!-- be opened (cf. check_open)
1881    IF ( statistic_regions > 9  .OR.  statistic_regions < 0 )  THEN
1882       WRITE ( message_string, * ) 'number of statistic_regions = ',           &
1883                   statistic_regions+1, ' but only 10 regions are allowed'
1884       CALL message( 'check_parameters', 'PA0082', 1, 2, 0, 6, 0 )
1885    ENDIF
1886    IF ( normalizing_region > statistic_regions  .OR.                          &
1887         normalizing_region < 0)  THEN
1888       WRITE ( message_string, * ) 'normalizing_region = ',                    &
1889                normalizing_region, ' must be >= 0 and <= ',statistic_regions, &
1890                ' (value of statistic_regions)'
1891       CALL message( 'check_parameters', 'PA0083', 1, 2, 0, 6, 0 )
1892    ENDIF
1893
1894!
1895!-- Set the default intervals for data output, if necessary
1896!-- NOTE: dt_dosp has already been set in spectra_parin
1897    IF ( dt_data_output /= 9999999.9_wp )  THEN
1898       IF ( dt_dopr           == 9999999.9_wp )  dt_dopr           = dt_data_output
1899       IF ( dt_dopts          == 9999999.9_wp )  dt_dopts          = dt_data_output
1900       IF ( dt_do2d_xy        == 9999999.9_wp )  dt_do2d_xy        = dt_data_output
1901       IF ( dt_do2d_xz        == 9999999.9_wp )  dt_do2d_xz        = dt_data_output
1902       IF ( dt_do2d_yz        == 9999999.9_wp )  dt_do2d_yz        = dt_data_output
1903       IF ( dt_do3d           == 9999999.9_wp )  dt_do3d           = dt_data_output
1904       IF ( dt_data_output_av == 9999999.9_wp )  dt_data_output_av = dt_data_output
1905       DO  mid = 1, max_masks
1906          IF ( dt_domask(mid) == 9999999.9_wp )  dt_domask(mid)    = dt_data_output
1907       ENDDO
1908    ENDIF
1909
1910!
1911!-- Set the default skip time intervals for data output, if necessary
1912    IF ( skip_time_dopr    == 9999999.9_wp )                                   &
1913                                       skip_time_dopr    = skip_time_data_output
1914    IF ( skip_time_do2d_xy == 9999999.9_wp )                                   &
1915                                       skip_time_do2d_xy = skip_time_data_output
1916    IF ( skip_time_do2d_xz == 9999999.9_wp )                                   &
1917                                       skip_time_do2d_xz = skip_time_data_output
1918    IF ( skip_time_do2d_yz == 9999999.9_wp )                                   &
1919                                       skip_time_do2d_yz = skip_time_data_output
1920    IF ( skip_time_do3d    == 9999999.9_wp )                                   &
1921                                       skip_time_do3d    = skip_time_data_output
1922    IF ( skip_time_data_output_av == 9999999.9_wp )                            &
1923                                skip_time_data_output_av = skip_time_data_output
1924    DO  mid = 1, max_masks
1925       IF ( skip_time_domask(mid) == 9999999.9_wp )                            &
1926                                skip_time_domask(mid)    = skip_time_data_output
1927    ENDDO
1928
1929!
1930!-- Check the average intervals (first for 3d-data, then for profiles)
1931    IF ( averaging_interval > dt_data_output_av )  THEN
1932       WRITE( message_string, * )  'averaging_interval = ',                    &
1933             averaging_interval, ' must be <= dt_data_output = ', dt_data_output
1934       CALL message( 'check_parameters', 'PA0085', 1, 2, 0, 6, 0 )
1935    ENDIF
1936
1937    IF ( averaging_interval_pr == 9999999.9_wp )  THEN
1938       averaging_interval_pr = averaging_interval
1939    ENDIF
1940
1941    IF ( averaging_interval_pr > dt_dopr )  THEN
1942       WRITE( message_string, * )  'averaging_interval_pr = ',                 &
1943             averaging_interval_pr, ' must be <= dt_dopr = ', dt_dopr
1944       CALL message( 'check_parameters', 'PA0086', 1, 2, 0, 6, 0 )
1945    ENDIF
1946
1947!
1948!-- Set the default interval for profiles entering the temporal average
1949    IF ( dt_averaging_input_pr == 9999999.9_wp )  THEN
1950       dt_averaging_input_pr = dt_averaging_input
1951    ENDIF
1952
1953!
1954!-- Set the default interval for the output of timeseries to a reasonable
1955!-- value (tries to minimize the number of calls of flow_statistics)
1956    IF ( dt_dots == 9999999.9_wp )  THEN
1957       IF ( averaging_interval_pr == 0.0_wp )  THEN
1958          dt_dots = MIN( dt_run_control, dt_dopr )
1959       ELSE
1960          dt_dots = MIN( dt_run_control, dt_averaging_input_pr )
1961       ENDIF
1962    ENDIF
1963
1964!
1965!-- Check the sample rate for averaging (first for 3d-data, then for profiles)
1966    IF ( dt_averaging_input > averaging_interval )  THEN
1967       WRITE( message_string, * )  'dt_averaging_input = ',                    &
1968                dt_averaging_input, ' must be <= averaging_interval = ',       &
1969                averaging_interval
1970       CALL message( 'check_parameters', 'PA0088', 1, 2, 0, 6, 0 )
1971    ENDIF
1972
1973    IF ( dt_averaging_input_pr > averaging_interval_pr )  THEN
1974       WRITE( message_string, * )  'dt_averaging_input_pr = ',                 &
1975                dt_averaging_input_pr, ' must be <= averaging_interval_pr = ', &
1976                averaging_interval_pr
1977       CALL message( 'check_parameters', 'PA0089', 1, 2, 0, 6, 0 )
1978    ENDIF
1979
1980!
1981!-- Set the default value for the integration interval of precipitation amount
1982    IF ( microphysics_seifert  .OR.  microphysics_kessler )  THEN
1983       IF ( precipitation_amount_interval == 9999999.9_wp )  THEN
1984          precipitation_amount_interval = dt_do2d_xy
1985       ELSE
1986          IF ( precipitation_amount_interval > dt_do2d_xy )  THEN
1987             WRITE( message_string, * )  'precipitation_amount_interval = ',   &
1988                 precipitation_amount_interval, ' must not be larger than ',   &
1989                 'dt_do2d_xy = ', dt_do2d_xy
1990             CALL message( 'check_parameters', 'PA0090', 1, 2, 0, 6, 0 )
1991          ENDIF
1992       ENDIF
1993    ENDIF
1994
1995!
1996!-- Determine the number of output profiles and check whether they are
1997!-- permissible
1998    DO  WHILE ( data_output_pr(dopr_n+1) /= '          ' )
1999
2000       dopr_n = dopr_n + 1
2001       i = dopr_n
2002
2003!
2004!--    Determine internal profile number (for hom, homs)
2005!--    and store height levels
2006       SELECT CASE ( TRIM( data_output_pr(i) ) )
2007
2008          CASE ( 'u', '#u' )
2009             dopr_index(i) = 1
2010             dopr_unit(i)  = 'm/s'
2011             hom(:,2,1,:)  = SPREAD( zu, 2, statistic_regions+1 )
2012             IF ( data_output_pr(i)(1:1) == '#' )  THEN
2013                dopr_initial_index(i) = 5
2014                hom(:,2,5,:)          = SPREAD( zu, 2, statistic_regions+1 )
2015                data_output_pr(i)     = data_output_pr(i)(2:)
2016             ENDIF
2017
2018          CASE ( 'v', '#v' )
2019             dopr_index(i) = 2
2020             dopr_unit(i)  = 'm/s'
2021             hom(:,2,2,:)  = SPREAD( zu, 2, statistic_regions+1 )
2022             IF ( data_output_pr(i)(1:1) == '#' )  THEN
2023                dopr_initial_index(i) = 6
2024                hom(:,2,6,:)          = SPREAD( zu, 2, statistic_regions+1 )
2025                data_output_pr(i)     = data_output_pr(i)(2:)
2026             ENDIF
2027
2028          CASE ( 'w' )
2029             dopr_index(i) = 3
2030             dopr_unit(i)  = 'm/s'
2031             hom(:,2,3,:)  = SPREAD( zw, 2, statistic_regions+1 )
2032
2033          CASE ( 'pt', '#pt' )
2034             IF ( .NOT. cloud_physics ) THEN
2035                dopr_index(i) = 4
2036                dopr_unit(i)  = 'K'
2037                hom(:,2,4,:)  = SPREAD( zu, 2, statistic_regions+1 )
2038                IF ( data_output_pr(i)(1:1) == '#' )  THEN
2039                   dopr_initial_index(i) = 7
2040                   hom(:,2,7,:)          = SPREAD( zu, 2, statistic_regions+1 )
2041                   hom(nzb,2,7,:)        = 0.0_wp    ! because zu(nzb) is negative
2042                   data_output_pr(i)     = data_output_pr(i)(2:)
2043                ENDIF
2044             ELSE
2045                dopr_index(i) = 43
2046                dopr_unit(i)  = 'K'
2047                hom(:,2,43,:)  = SPREAD( zu, 2, statistic_regions+1 )
2048                IF ( data_output_pr(i)(1:1) == '#' )  THEN
2049                   dopr_initial_index(i) = 28
2050                   hom(:,2,28,:)         = SPREAD( zu, 2, statistic_regions+1 )
2051                   hom(nzb,2,28,:)       = 0.0_wp    ! because zu(nzb) is negative
2052                   data_output_pr(i)     = data_output_pr(i)(2:)
2053                ENDIF
2054             ENDIF
2055
2056          CASE ( 'e' )
2057             dopr_index(i)  = 8
2058             dopr_unit(i)   = 'm2/s2'
2059             hom(:,2,8,:)   = SPREAD( zu, 2, statistic_regions+1 )
2060             hom(nzb,2,8,:) = 0.0_wp
2061
2062          CASE ( 'km', '#km' )
2063             dopr_index(i)  = 9
2064             dopr_unit(i)   = 'm2/s'
2065             hom(:,2,9,:)   = SPREAD( zu, 2, statistic_regions+1 )
2066             hom(nzb,2,9,:) = 0.0_wp
2067             IF ( data_output_pr(i)(1:1) == '#' )  THEN
2068                dopr_initial_index(i) = 23
2069                hom(:,2,23,:)         = hom(:,2,9,:)
2070                data_output_pr(i)     = data_output_pr(i)(2:)
2071             ENDIF
2072
2073          CASE ( 'kh', '#kh' )
2074             dopr_index(i)   = 10
2075             dopr_unit(i)    = 'm2/s'
2076             hom(:,2,10,:)   = SPREAD( zu, 2, statistic_regions+1 )
2077             hom(nzb,2,10,:) = 0.0_wp
2078             IF ( data_output_pr(i)(1:1) == '#' )  THEN
2079                dopr_initial_index(i) = 24
2080                hom(:,2,24,:)         = hom(:,2,10,:)
2081                data_output_pr(i)     = data_output_pr(i)(2:)
2082             ENDIF
2083
2084          CASE ( 'l', '#l' )
2085             dopr_index(i)   = 11
2086             dopr_unit(i)    = 'm'
2087             hom(:,2,11,:)   = SPREAD( zu, 2, statistic_regions+1 )
2088             hom(nzb,2,11,:) = 0.0_wp
2089             IF ( data_output_pr(i)(1:1) == '#' )  THEN
2090                dopr_initial_index(i) = 25
2091                hom(:,2,25,:)         = hom(:,2,11,:)
2092                data_output_pr(i)     = data_output_pr(i)(2:)
2093             ENDIF
2094
2095          CASE ( 'w"u"' )
2096             dopr_index(i) = 12
2097             dopr_unit(i)  = 'm2/s2'
2098             hom(:,2,12,:) = SPREAD( zw, 2, statistic_regions+1 )
2099             IF ( constant_flux_layer )  hom(nzb,2,12,:) = zu(1)
2100
2101          CASE ( 'w*u*' )
2102             dopr_index(i) = 13
2103             dopr_unit(i)  = 'm2/s2'
2104             hom(:,2,13,:) = SPREAD( zw, 2, statistic_regions+1 )
2105
2106          CASE ( 'w"v"' )
2107             dopr_index(i) = 14
2108             dopr_unit(i)  = 'm2/s2'
2109             hom(:,2,14,:) = SPREAD( zw, 2, statistic_regions+1 )
2110             IF ( constant_flux_layer )  hom(nzb,2,14,:) = zu(1)
2111
2112          CASE ( 'w*v*' )
2113             dopr_index(i) = 15
2114             dopr_unit(i)  = 'm2/s2'
2115             hom(:,2,15,:) = SPREAD( zw, 2, statistic_regions+1 )
2116
2117          CASE ( 'w"pt"' )
2118             dopr_index(i) = 16
2119             dopr_unit(i)  = 'K m/s'
2120             hom(:,2,16,:) = SPREAD( zw, 2, statistic_regions+1 )
2121
2122          CASE ( 'w*pt*' )
2123             dopr_index(i) = 17
2124             dopr_unit(i)  = 'K m/s'
2125             hom(:,2,17,:) = SPREAD( zw, 2, statistic_regions+1 )
2126
2127          CASE ( 'wpt' )
2128             dopr_index(i) = 18
2129             dopr_unit(i)  = 'K m/s'
2130             hom(:,2,18,:) = SPREAD( zw, 2, statistic_regions+1 )
2131
2132          CASE ( 'wu' )
2133             dopr_index(i) = 19
2134             dopr_unit(i)  = 'm2/s2'
2135             hom(:,2,19,:) = SPREAD( zw, 2, statistic_regions+1 )
2136             IF ( constant_flux_layer )  hom(nzb,2,19,:) = zu(1)
2137
2138          CASE ( 'wv' )
2139             dopr_index(i) = 20
2140             dopr_unit(i)  = 'm2/s2'
2141             hom(:,2,20,:) = SPREAD( zw, 2, statistic_regions+1 )
2142             IF ( constant_flux_layer )  hom(nzb,2,20,:) = zu(1)
2143
2144          CASE ( 'w*pt*BC' )
2145             dopr_index(i) = 21
2146             dopr_unit(i)  = 'K m/s'
2147             hom(:,2,21,:) = SPREAD( zw, 2, statistic_regions+1 )
2148
2149          CASE ( 'wptBC' )
2150             dopr_index(i) = 22
2151             dopr_unit(i)  = 'K m/s'
2152             hom(:,2,22,:) = SPREAD( zw, 2, statistic_regions+1 )
2153
2154          CASE ( 'sa', '#sa' )
2155             IF ( .NOT. ocean )  THEN
2156                message_string = 'data_output_pr = ' // &
2157                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2158                                 'lemented for ocean = .FALSE.'
2159                CALL message( 'check_parameters', 'PA0091', 1, 2, 0, 6, 0 )
2160             ELSE
2161                dopr_index(i) = 23
2162                dopr_unit(i)  = 'psu'
2163                hom(:,2,23,:) = SPREAD( zu, 2, statistic_regions+1 )
2164                IF ( data_output_pr(i)(1:1) == '#' )  THEN
2165                   dopr_initial_index(i) = 26
2166                   hom(:,2,26,:)         = SPREAD( zu, 2, statistic_regions+1 )
2167                   hom(nzb,2,26,:)       = 0.0_wp    ! because zu(nzb) is negative
2168                   data_output_pr(i)     = data_output_pr(i)(2:)
2169                ENDIF
2170             ENDIF
2171
2172          CASE ( 'u*2' )
2173             dopr_index(i) = 30
2174             dopr_unit(i)  = 'm2/s2'
2175             hom(:,2,30,:) = SPREAD( zu, 2, statistic_regions+1 )
2176
2177          CASE ( 'v*2' )
2178             dopr_index(i) = 31
2179             dopr_unit(i)  = 'm2/s2'
2180             hom(:,2,31,:) = SPREAD( zu, 2, statistic_regions+1 )
2181
2182          CASE ( 'w*2' )
2183             dopr_index(i) = 32
2184             dopr_unit(i)  = 'm2/s2'
2185             hom(:,2,32,:) = SPREAD( zw, 2, statistic_regions+1 )
2186
2187          CASE ( 'pt*2' )
2188             dopr_index(i) = 33
2189             dopr_unit(i)  = 'K2'
2190             hom(:,2,33,:) = SPREAD( zu, 2, statistic_regions+1 )
2191
2192          CASE ( 'e*' )
2193             dopr_index(i) = 34
2194             dopr_unit(i)  = 'm2/s2'
2195             hom(:,2,34,:) = SPREAD( zu, 2, statistic_regions+1 )
2196
2197          CASE ( 'w*2pt*' )
2198             dopr_index(i) = 35
2199             dopr_unit(i)  = 'K m2/s2'
2200             hom(:,2,35,:) = SPREAD( zw, 2, statistic_regions+1 )
2201
2202          CASE ( 'w*pt*2' )
2203             dopr_index(i) = 36
2204             dopr_unit(i)  = 'K2 m/s'
2205             hom(:,2,36,:) = SPREAD( zw, 2, statistic_regions+1 )
2206
2207          CASE ( 'w*e*' )
2208             dopr_index(i) = 37
2209             dopr_unit(i)  = 'm3/s3'
2210             hom(:,2,37,:) = SPREAD( zw, 2, statistic_regions+1 )
2211
2212          CASE ( 'w*3' )
2213             dopr_index(i) = 38
2214             dopr_unit(i)  = 'm3/s3'
2215             hom(:,2,38,:) = SPREAD( zw, 2, statistic_regions+1 )
2216
2217          CASE ( 'Sw' )
2218             dopr_index(i) = 39
2219             dopr_unit(i)  = 'none'
2220             hom(:,2,39,:) = SPREAD( zw, 2, statistic_regions+1 )
2221
2222          CASE ( 'p' )
2223             dopr_index(i) = 40
2224             dopr_unit(i)  = 'Pa'
2225             hom(:,2,40,:) = SPREAD( zu, 2, statistic_regions+1 )
2226
2227          CASE ( 'q', '#q' )
2228             IF ( .NOT. humidity )  THEN
2229                message_string = 'data_output_pr = ' // &
2230                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2231                                 'lemented for humidity = .FALSE.'
2232                CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
2233             ELSE
2234                dopr_index(i) = 41
2235                dopr_unit(i)  = 'kg/kg'
2236                hom(:,2,41,:) = SPREAD( zu, 2, statistic_regions+1 )
2237                IF ( data_output_pr(i)(1:1) == '#' )  THEN
2238                   dopr_initial_index(i) = 26
2239                   hom(:,2,26,:)         = SPREAD( zu, 2, statistic_regions+1 )
2240                   hom(nzb,2,26,:)       = 0.0_wp    ! because zu(nzb) is negative
2241                   data_output_pr(i)     = data_output_pr(i)(2:)
2242                ENDIF
2243             ENDIF
2244
2245          CASE ( 's', '#s' )
2246             IF ( .NOT. passive_scalar )  THEN
2247                message_string = 'data_output_pr = ' //                        &
2248                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2249                                 'lemented for passive_scalar = .FALSE.'
2250                CALL message( 'check_parameters', 'PA0093', 1, 2, 0, 6, 0 )
2251             ELSE
2252                dopr_index(i) = 117
2253                dopr_unit(i)  = 'kg/m3'
2254                hom(:,2,117,:) = SPREAD( zu, 2, statistic_regions+1 )
2255                IF ( data_output_pr(i)(1:1) == '#' )  THEN
2256                   dopr_initial_index(i) = 117
2257                   hom(:,2,117,:)        = SPREAD( zu, 2, statistic_regions+1 )
2258                   hom(nzb,2,117,:)      = 0.0_wp    ! because zu(nzb) is negative
2259                   data_output_pr(i)     = data_output_pr(i)(2:)
2260                ENDIF
2261             ENDIF
2262
2263          CASE ( 'qv', '#qv' )
2264             IF ( .NOT. cloud_physics ) THEN
2265                dopr_index(i) = 41
2266                dopr_unit(i)  = 'kg/kg'
2267                hom(:,2,41,:) = SPREAD( zu, 2, statistic_regions+1 )
2268                IF ( data_output_pr(i)(1:1) == '#' )  THEN
2269                   dopr_initial_index(i) = 26
2270                   hom(:,2,26,:)         = SPREAD( zu, 2, statistic_regions+1 )
2271                   hom(nzb,2,26,:)       = 0.0_wp    ! because zu(nzb) is negative
2272                   data_output_pr(i)     = data_output_pr(i)(2:)
2273                ENDIF
2274             ELSE
2275                dopr_index(i) = 42
2276                dopr_unit(i)  = 'kg/kg'
2277                hom(:,2,42,:) = SPREAD( zu, 2, statistic_regions+1 )
2278                IF ( data_output_pr(i)(1:1) == '#' )  THEN
2279                   dopr_initial_index(i) = 27
2280                   hom(:,2,27,:)         = SPREAD( zu, 2, statistic_regions+1 )
2281                   hom(nzb,2,27,:)       = 0.0_wp   ! because zu(nzb) is negative
2282                   data_output_pr(i)     = data_output_pr(i)(2:)
2283                ENDIF
2284             ENDIF
2285
2286          CASE ( 'lpt', '#lpt' )
2287             IF ( .NOT. cloud_physics ) THEN
2288                message_string = 'data_output_pr = ' //                        &
2289                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2290                                 'lemented for cloud_physics = .FALSE.'
2291                CALL message( 'check_parameters', 'PA0094', 1, 2, 0, 6, 0 )
2292             ELSE
2293                dopr_index(i) = 4
2294                dopr_unit(i)  = 'K'
2295                hom(:,2,4,:)  = SPREAD( zu, 2, statistic_regions+1 )
2296                IF ( data_output_pr(i)(1:1) == '#' )  THEN
2297                   dopr_initial_index(i) = 7
2298                   hom(:,2,7,:)          = SPREAD( zu, 2, statistic_regions+1 )
2299                   hom(nzb,2,7,:)        = 0.0_wp    ! because zu(nzb) is negative
2300                   data_output_pr(i)     = data_output_pr(i)(2:)
2301                ENDIF
2302             ENDIF
2303
2304          CASE ( 'vpt', '#vpt' )
2305             dopr_index(i) = 44
2306             dopr_unit(i)  = 'K'
2307             hom(:,2,44,:) = SPREAD( zu, 2, statistic_regions+1 )
2308             IF ( data_output_pr(i)(1:1) == '#' )  THEN
2309                dopr_initial_index(i) = 29
2310                hom(:,2,29,:)         = SPREAD( zu, 2, statistic_regions+1 )
2311                hom(nzb,2,29,:)       = 0.0_wp    ! because zu(nzb) is negative
2312                data_output_pr(i)     = data_output_pr(i)(2:)
2313             ENDIF
2314
2315          CASE ( 'w"vpt"' )
2316             dopr_index(i) = 45
2317             dopr_unit(i)  = 'K m/s'
2318             hom(:,2,45,:) = SPREAD( zw, 2, statistic_regions+1 )
2319
2320          CASE ( 'w*vpt*' )
2321             dopr_index(i) = 46
2322             dopr_unit(i)  = 'K m/s'
2323             hom(:,2,46,:) = SPREAD( zw, 2, statistic_regions+1 )
2324
2325          CASE ( 'wvpt' )
2326             dopr_index(i) = 47
2327             dopr_unit(i)  = 'K m/s'
2328             hom(:,2,47,:) = SPREAD( zw, 2, statistic_regions+1 )
2329
2330          CASE ( 'w"q"' )
2331             IF (  .NOT.  humidity )  THEN
2332                message_string = 'data_output_pr = ' //                        &
2333                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2334                                 'lemented for humidity = .FALSE.'
2335                CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
2336             ELSE
2337                dopr_index(i) = 48
2338                dopr_unit(i)  = 'kg/kg m/s'
2339                hom(:,2,48,:) = SPREAD( zw, 2, statistic_regions+1 )
2340             ENDIF
2341
2342          CASE ( 'w*q*' )
2343             IF (  .NOT.  humidity )  THEN
2344                message_string = 'data_output_pr = ' //                        &
2345                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2346                                 'lemented for humidity = .FALSE.'
2347                CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
2348             ELSE
2349                dopr_index(i) = 49
2350                dopr_unit(i)  = 'kg/kg m/s'
2351                hom(:,2,49,:) = SPREAD( zw, 2, statistic_regions+1 )
2352             ENDIF
2353
2354          CASE ( 'wq' )
2355             IF (  .NOT.  humidity )  THEN
2356                message_string = 'data_output_pr = ' //                        &
2357                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2358                                 'lemented for humidity = .FALSE.'
2359                CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
2360             ELSE
2361                dopr_index(i) = 50
2362                dopr_unit(i)  = 'kg/kg m/s'
2363                hom(:,2,50,:) = SPREAD( zw, 2, statistic_regions+1 )
2364             ENDIF
2365
2366          CASE ( 'w"s"' )
2367             IF (  .NOT.  passive_scalar )  THEN
2368                message_string = 'data_output_pr = ' //                        &
2369                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2370                                 'lemented for passive_scalar = .FALSE.'
2371                CALL message( 'check_parameters', 'PA0093', 1, 2, 0, 6, 0 )
2372             ELSE
2373                dopr_index(i) = 119
2374                dopr_unit(i)  = 'kg/m3 m/s'
2375                hom(:,2,119,:) = SPREAD( zw, 2, statistic_regions+1 )
2376             ENDIF
2377
2378          CASE ( 'w*s*' )
2379             IF (  .NOT.  passive_scalar )  THEN
2380                message_string = 'data_output_pr = ' //                        &
2381                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2382                                 'lemented for passive_scalar = .FALSE.'
2383                CALL message( 'check_parameters', 'PA0093', 1, 2, 0, 6, 0 )
2384             ELSE
2385                dopr_index(i) = 116
2386                dopr_unit(i)  = 'kg/m3 m/s'
2387                hom(:,2,116,:) = SPREAD( zw, 2, statistic_regions+1 )
2388             ENDIF
2389
2390          CASE ( 'ws' )
2391             IF (  .NOT.  passive_scalar )  THEN
2392                message_string = 'data_output_pr = ' //                        &
2393                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2394                                 'lemented for passive_scalar = .FALSE.'
2395                CALL message( 'check_parameters', 'PA0093', 1, 2, 0, 6, 0 )
2396             ELSE
2397                dopr_index(i) = 120
2398                dopr_unit(i)  = 'kg/m3 m/s'
2399                hom(:,2,120,:) = SPREAD( zw, 2, statistic_regions+1 )
2400             ENDIF
2401
2402          CASE ( 'w"qv"' )
2403             IF ( humidity  .AND.  .NOT.  cloud_physics )  THEN
2404                dopr_index(i) = 48
2405                dopr_unit(i)  = 'kg/kg m/s'
2406                hom(:,2,48,:) = SPREAD( zw, 2, statistic_regions+1 )
2407             ELSEIF ( humidity  .AND.  cloud_physics )  THEN
2408                dopr_index(i) = 51
2409                dopr_unit(i)  = 'kg/kg m/s'
2410                hom(:,2,51,:) = SPREAD( zw, 2, statistic_regions+1 )
2411             ELSE
2412                message_string = 'data_output_pr = ' //                        &
2413                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2414                                 'lemented for cloud_physics = .FALSE. an&' // &
2415                                 'd humidity = .FALSE.'
2416                CALL message( 'check_parameters', 'PA0095', 1, 2, 0, 6, 0 )
2417             ENDIF
2418
2419          CASE ( 'w*qv*' )
2420             IF ( humidity  .AND.  .NOT. cloud_physics )                       &
2421             THEN
2422                dopr_index(i) = 49
2423                dopr_unit(i)  = 'kg/kg m/s'
2424                hom(:,2,49,:) = SPREAD( zw, 2, statistic_regions+1 )
2425             ELSEIF( humidity .AND. cloud_physics ) THEN
2426                dopr_index(i) = 52
2427                dopr_unit(i)  = 'kg/kg m/s'
2428                hom(:,2,52,:) = SPREAD( zw, 2, statistic_regions+1 )
2429             ELSE
2430                message_string = 'data_output_pr = ' //                        &
2431                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2432                                 'lemented for cloud_physics = .FALSE. an&' // &
2433                                 'd humidity = .FALSE.'
2434                CALL message( 'check_parameters', 'PA0095', 1, 2, 0, 6, 0 )
2435             ENDIF
2436
2437          CASE ( 'wqv' )
2438             IF ( humidity  .AND.  .NOT.  cloud_physics )  THEN
2439                dopr_index(i) = 50
2440                dopr_unit(i)  = 'kg/kg m/s'
2441                hom(:,2,50,:) = SPREAD( zw, 2, statistic_regions+1 )
2442             ELSEIF ( humidity  .AND.  cloud_physics )  THEN
2443                dopr_index(i) = 53
2444                dopr_unit(i)  = 'kg/kg m/s'
2445                hom(:,2,53,:) = SPREAD( zw, 2, statistic_regions+1 )
2446             ELSE
2447                message_string = 'data_output_pr = ' //                        &
2448                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2449                                 'lemented for cloud_physics = .FALSE. an&' // &
2450                                 'd humidity = .FALSE.'
2451                CALL message( 'check_parameters', 'PA0095', 1, 2, 0, 6, 0 )
2452             ENDIF
2453
2454          CASE ( 'ql' )
2455             IF (  .NOT.  cloud_physics  .AND.  .NOT.  cloud_droplets )  THEN
2456                message_string = 'data_output_pr = ' //                        &
2457                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2458                                 'lemented for cloud_physics = .FALSE. or'  // &
2459                                 '&cloud_droplets = .FALSE.'
2460                CALL message( 'check_parameters', 'PA0096', 1, 2, 0, 6, 0 )
2461             ELSE
2462                dopr_index(i) = 54
2463                dopr_unit(i)  = 'kg/kg'
2464                hom(:,2,54,:)  = SPREAD( zu, 2, statistic_regions+1 )
2465             ENDIF
2466
2467          CASE ( 'w*u*u*:dz' )
2468             dopr_index(i) = 55
2469             dopr_unit(i)  = 'm2/s3'
2470             hom(:,2,55,:) = SPREAD( zu, 2, statistic_regions+1 )
2471
2472          CASE ( 'w*p*:dz' )
2473             dopr_index(i) = 56
2474             dopr_unit(i)  = 'm2/s3'
2475             hom(:,2,56,:) = SPREAD( zw, 2, statistic_regions+1 )
2476
2477          CASE ( 'w"e:dz' )
2478             dopr_index(i) = 57
2479             dopr_unit(i)  = 'm2/s3'
2480             hom(:,2,57,:) = SPREAD( zu, 2, statistic_regions+1 )
2481
2482
2483          CASE ( 'u"pt"' )
2484             dopr_index(i) = 58
2485             dopr_unit(i)  = 'K m/s'
2486             hom(:,2,58,:) = SPREAD( zu, 2, statistic_regions+1 )
2487
2488          CASE ( 'u*pt*' )
2489             dopr_index(i) = 59
2490             dopr_unit(i)  = 'K m/s'
2491             hom(:,2,59,:) = SPREAD( zu, 2, statistic_regions+1 )
2492
2493          CASE ( 'upt_t' )
2494             dopr_index(i) = 60
2495             dopr_unit(i)  = 'K m/s'
2496             hom(:,2,60,:) = SPREAD( zu, 2, statistic_regions+1 )
2497
2498          CASE ( 'v"pt"' )
2499             dopr_index(i) = 61
2500             dopr_unit(i)  = 'K m/s'
2501             hom(:,2,61,:) = SPREAD( zu, 2, statistic_regions+1 )
2502             
2503          CASE ( 'v*pt*' )
2504             dopr_index(i) = 62
2505             dopr_unit(i)  = 'K m/s'
2506             hom(:,2,62,:) = SPREAD( zu, 2, statistic_regions+1 )
2507
2508          CASE ( 'vpt_t' )
2509             dopr_index(i) = 63
2510             dopr_unit(i)  = 'K m/s'
2511             hom(:,2,63,:) = SPREAD( zu, 2, statistic_regions+1 )
2512
2513          CASE ( 'rho' )
2514             IF (  .NOT.  ocean ) THEN
2515                message_string = 'data_output_pr = ' //                        &
2516                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2517                                 'lemented for ocean = .FALSE.'
2518                CALL message( 'check_parameters', 'PA0091', 1, 2, 0, 6, 0 )
2519             ELSE
2520                dopr_index(i) = 64
2521                dopr_unit(i)  = 'kg/m3'
2522                hom(:,2,64,:) = SPREAD( zu, 2, statistic_regions+1 )
2523                IF ( data_output_pr(i)(1:1) == '#' )  THEN
2524                   dopr_initial_index(i) = 77
2525                   hom(:,2,77,:)         = SPREAD( zu, 2, statistic_regions+1 )
2526                   hom(nzb,2,77,:)       = 0.0_wp    ! because zu(nzb) is negative
2527                   data_output_pr(i)     = data_output_pr(i)(2:)
2528                ENDIF
2529             ENDIF
2530
2531          CASE ( 'w"sa"' )
2532             IF (  .NOT.  ocean ) THEN
2533                message_string = 'data_output_pr = ' //                        &
2534                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2535                                 'lemented for ocean = .FALSE.'
2536                CALL message( 'check_parameters', 'PA0091', 1, 2, 0, 6, 0 )
2537             ELSE
2538                dopr_index(i) = 65
2539                dopr_unit(i)  = 'psu m/s'
2540                hom(:,2,65,:) = SPREAD( zw, 2, statistic_regions+1 )
2541             ENDIF
2542
2543          CASE ( 'w*sa*' )
2544             IF (  .NOT. ocean  ) THEN
2545                message_string = 'data_output_pr = ' //                        &
2546                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2547                                 'lemented for ocean = .FALSE.'
2548                CALL message( 'check_parameters', 'PA0091', 1, 2, 0, 6, 0 )
2549             ELSE
2550                dopr_index(i) = 66
2551                dopr_unit(i)  = 'psu m/s'
2552                hom(:,2,66,:) = SPREAD( zw, 2, statistic_regions+1 )
2553             ENDIF
2554
2555          CASE ( 'wsa' )
2556             IF (  .NOT.  ocean ) THEN
2557                message_string = 'data_output_pr = ' //                        &
2558                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2559                                 'lemented for ocean = .FALSE.'
2560                CALL message( 'check_parameters', 'PA0091', 1, 2, 0, 6, 0 )
2561             ELSE
2562                dopr_index(i) = 67
2563                dopr_unit(i)  = 'psu m/s'
2564                hom(:,2,67,:) = SPREAD( zw, 2, statistic_regions+1 )
2565             ENDIF
2566
2567          CASE ( 'w*p*' )
2568             dopr_index(i) = 68
2569             dopr_unit(i)  = 'm3/s3'
2570             hom(:,2,68,:) = SPREAD( zu, 2, statistic_regions+1 )
2571
2572          CASE ( 'w"e' )
2573             dopr_index(i) = 69
2574             dopr_unit(i)  = 'm3/s3'
2575             hom(:,2,69,:) = SPREAD( zu, 2, statistic_regions+1 )
2576
2577          CASE ( 'q*2' )
2578             IF (  .NOT.  humidity )  THEN
2579                message_string = 'data_output_pr = ' //                        &
2580                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2581                                 'lemented for humidity = .FALSE.'
2582                CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
2583             ELSE
2584                dopr_index(i) = 70
2585                dopr_unit(i)  = 'kg2/kg2'
2586                hom(:,2,70,:) = SPREAD( zu, 2, statistic_regions+1 )
2587             ENDIF
2588
2589          CASE ( 'prho' )
2590             IF (  .NOT.  ocean ) THEN
2591                message_string = 'data_output_pr = ' //                        &
2592                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2593                                 'lemented for ocean = .FALSE.'
2594                CALL message( 'check_parameters', 'PA0091', 1, 2, 0, 6, 0 )
2595             ELSE
2596                dopr_index(i) = 71
2597                dopr_unit(i)  = 'kg/m3'
2598                hom(:,2,71,:) = SPREAD( zu, 2, statistic_regions+1 )
2599             ENDIF
2600
2601          CASE ( 'hyp' )
2602             dopr_index(i) = 72
2603             dopr_unit(i)  = 'dbar'
2604             hom(:,2,72,:) = SPREAD( zu, 2, statistic_regions+1 )
2605
2606          CASE ( 'nr' )
2607             IF (  .NOT.  cloud_physics )  THEN
2608                message_string = 'data_output_pr = ' //                        &
2609                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2610                                 'lemented for cloud_physics = .FALSE.'
2611                CALL message( 'check_parameters', 'PA0094', 1, 2, 0, 6, 0 )
2612             ELSEIF ( .NOT.  microphysics_seifert )  THEN
2613                message_string = 'data_output_pr = ' //                        &
2614                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2615                                 'lemented for cloud_scheme /= seifert_beheng'
2616                CALL message( 'check_parameters', 'PA0358', 1, 2, 0, 6, 0 )
2617             ELSE
2618                dopr_index(i) = 73
2619                dopr_unit(i)  = '1/m3'
2620                hom(:,2,73,:)  = SPREAD( zu, 2, statistic_regions+1 )
2621             ENDIF
2622
2623          CASE ( 'qr' )
2624             IF (  .NOT.  cloud_physics )  THEN
2625                message_string = 'data_output_pr = ' //                        &
2626                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2627                                 'lemented for cloud_physics = .FALSE.'
2628                CALL message( 'check_parameters', 'PA0094', 1, 2, 0, 6, 0 )
2629             ELSEIF ( .NOT.  microphysics_seifert )  THEN
2630                message_string = 'data_output_pr = ' //                        &
2631                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2632                                 'lemented for cloud_scheme /= seifert_beheng'
2633                CALL message( 'check_parameters', 'PA0358', 1, 2, 0, 6, 0 )
2634             ELSE
2635                dopr_index(i) = 74
2636                dopr_unit(i)  = 'kg/kg'
2637                hom(:,2,74,:)  = SPREAD( zu, 2, statistic_regions+1 )
2638             ENDIF
2639
2640          CASE ( 'qc' )
2641             IF (  .NOT.  cloud_physics )  THEN
2642                message_string = 'data_output_pr = ' //                        &
2643                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2644                                 'lemented for cloud_physics = .FALSE.'
2645                CALL message( 'check_parameters', 'PA0094', 1, 2, 0, 6, 0 )
2646             ELSE
2647                dopr_index(i) = 75
2648                dopr_unit(i)  = 'kg/kg'
2649                hom(:,2,75,:)  = SPREAD( zu, 2, statistic_regions+1 )
2650             ENDIF
2651
2652          CASE ( 'prr' )
2653             IF (  .NOT.  cloud_physics )  THEN
2654                message_string = 'data_output_pr = ' //                        &
2655                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2656                                 'lemented for cloud_physics = .FALSE.'
2657                CALL message( 'check_parameters', 'PA0094', 1, 2, 0, 6, 0 )
2658             ELSEIF ( microphysics_sat_adjust )  THEN
2659                message_string = 'data_output_pr = ' //                        &
2660                                 TRIM( data_output_pr(i) ) // ' is not ava' // &
2661                                 'ilable for cloud_scheme = saturation_adjust'
2662                CALL message( 'check_parameters', 'PA0422', 1, 2, 0, 6, 0 )
2663             ELSE
2664                dopr_index(i) = 76
2665                dopr_unit(i)  = 'kg/kg m/s'
2666                hom(:,2,76,:)  = SPREAD( zu, 2, statistic_regions+1 )
2667             ENDIF
2668
2669          CASE ( 'ug' )
2670             dopr_index(i) = 78
2671             dopr_unit(i)  = 'm/s'
2672             hom(:,2,78,:) = SPREAD( zu, 2, statistic_regions+1 )
2673
2674          CASE ( 'vg' )
2675             dopr_index(i) = 79
2676             dopr_unit(i)  = 'm/s'
2677             hom(:,2,79,:) = SPREAD( zu, 2, statistic_regions+1 )
2678
2679          CASE ( 'w_subs' )
2680             IF (  .NOT.  large_scale_subsidence )  THEN
2681                message_string = 'data_output_pr = ' //                        &
2682                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2683                                 'lemented for large_scale_subsidence = .FALSE.'
2684                CALL message( 'check_parameters', 'PA0382', 1, 2, 0, 6, 0 )
2685             ELSE
2686                dopr_index(i) = 80
2687                dopr_unit(i)  = 'm/s'
2688                hom(:,2,80,:) = SPREAD( zu, 2, statistic_regions+1 )
2689             ENDIF
2690
2691          CASE ( 'td_lsa_lpt' )
2692             IF (  .NOT.  large_scale_forcing )  THEN
2693                message_string = 'data_output_pr = ' //                        &
2694                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2695                                 'lemented for large_scale_forcing = .FALSE.'
2696                CALL message( 'check_parameters', 'PA0393', 1, 2, 0, 6, 0 )
2697             ELSE
2698                dopr_index(i) = 81
2699                dopr_unit(i)  = 'K/s'
2700                hom(:,2,81,:) = SPREAD( zu, 2, statistic_regions+1 )
2701             ENDIF
2702
2703          CASE ( 'td_lsa_q' )
2704             IF (  .NOT.  large_scale_forcing )  THEN
2705                message_string = 'data_output_pr = ' //                        &
2706                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2707                                 'lemented for large_scale_forcing = .FALSE.'
2708                CALL message( 'check_parameters', 'PA0393', 1, 2, 0, 6, 0 )
2709             ELSE
2710                dopr_index(i) = 82
2711                dopr_unit(i)  = 'kg/kgs'
2712                hom(:,2,82,:) = SPREAD( zu, 2, statistic_regions+1 )
2713             ENDIF
2714
2715          CASE ( 'td_sub_lpt' )
2716             IF (  .NOT.  large_scale_forcing )  THEN
2717                message_string = 'data_output_pr = ' //                        &
2718                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2719                                 'lemented for large_scale_forcing = .FALSE.'
2720                CALL message( 'check_parameters', 'PA0393', 1, 2, 0, 6, 0 )
2721             ELSE
2722                dopr_index(i) = 83
2723                dopr_unit(i)  = 'K/s'
2724                hom(:,2,83,:) = SPREAD( zu, 2, statistic_regions+1 )
2725             ENDIF
2726
2727          CASE ( 'td_sub_q' )
2728             IF (  .NOT.  large_scale_forcing )  THEN
2729                message_string = 'data_output_pr = ' //                        &
2730                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2731                                 'lemented for large_scale_forcing = .FALSE.'
2732                CALL message( 'check_parameters', 'PA0393', 1, 2, 0, 6, 0 )
2733             ELSE
2734                dopr_index(i) = 84
2735                dopr_unit(i)  = 'kg/kgs'
2736                hom(:,2,84,:) = SPREAD( zu, 2, statistic_regions+1 )
2737             ENDIF
2738
2739          CASE ( 'td_nud_lpt' )
2740             IF (  .NOT.  nudging )  THEN
2741                message_string = 'data_output_pr = ' //                        &
2742                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2743                                 'lemented for nudging = .FALSE.'
2744                CALL message( 'check_parameters', 'PA0394', 1, 2, 0, 6, 0 )
2745             ELSE
2746                dopr_index(i) = 85
2747                dopr_unit(i)  = 'K/s'
2748                hom(:,2,85,:) = SPREAD( zu, 2, statistic_regions+1 )
2749             ENDIF
2750
2751          CASE ( 'td_nud_q' )
2752             IF (  .NOT.  nudging )  THEN
2753                message_string = 'data_output_pr = ' //                        &
2754                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2755                                 'lemented for nudging = .FALSE.'
2756                CALL message( 'check_parameters', 'PA0394', 1, 2, 0, 6, 0 )
2757             ELSE
2758                dopr_index(i) = 86
2759                dopr_unit(i)  = 'kg/kgs'
2760                hom(:,2,86,:) = SPREAD( zu, 2, statistic_regions+1 )
2761             ENDIF
2762
2763          CASE ( 'td_nud_u' )
2764             IF (  .NOT.  nudging )  THEN
2765                message_string = 'data_output_pr = ' //                        &
2766                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2767                                 'lemented for nudging = .FALSE.'
2768                CALL message( 'check_parameters', 'PA0394', 1, 2, 0, 6, 0 )
2769             ELSE
2770                dopr_index(i) = 87
2771                dopr_unit(i)  = 'm/s2'
2772                hom(:,2,87,:) = SPREAD( zu, 2, statistic_regions+1 )
2773             ENDIF
2774
2775          CASE ( 'td_nud_v' )
2776             IF (  .NOT.  nudging )  THEN
2777                message_string = 'data_output_pr = ' //                        &
2778                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
2779                                 'lemented for nudging = .FALSE.'
2780                CALL message( 'check_parameters', 'PA0394', 1, 2, 0, 6, 0 )
2781             ELSE
2782                dopr_index(i) = 88
2783                dopr_unit(i)  = 'm/s2'
2784                hom(:,2,88,:) = SPREAD( zu, 2, statistic_regions+1 )
2785             ENDIF
2786
2787 
2788
2789          CASE DEFAULT
2790
2791             CALL lsm_check_data_output_pr( data_output_pr(i), i, unit,        &
2792                                            dopr_unit(i) )
2793
2794             IF ( unit == 'illegal' )  THEN
2795                CALL radiation_check_data_output_pr( data_output_pr(i), i,     &
2796                                                     unit, dopr_unit(i) )
2797             ENDIF
2798             
2799             IF ( unit == 'illegal' )  THEN
2800                unit = ''
2801                CALL user_check_data_output_pr( data_output_pr(i), i, unit )
2802             ENDIF
2803
2804             IF ( unit == 'illegal' )  THEN
2805                IF ( data_output_pr_user(1) /= ' ' )  THEN
2806                   message_string = 'illegal value for data_output_pr or ' //  &
2807                                    'data_output_pr_user = "' //               &
2808                                    TRIM( data_output_pr(i) ) // '"'
2809                   CALL message( 'check_parameters', 'PA0097', 1, 2, 0, 6, 0 )
2810                ELSE
2811                   message_string = 'illegal value for data_output_pr = "' //  &
2812                                    TRIM( data_output_pr(i) ) // '"'
2813                   CALL message( 'check_parameters', 'PA0098', 1, 2, 0, 6, 0 )
2814                ENDIF
2815             ENDIF
2816
2817       END SELECT
2818
2819    ENDDO
2820
2821
2822!
2823!-- Append user-defined data output variables to the standard data output
2824    IF ( data_output_user(1) /= ' ' )  THEN
2825       i = 1
2826       DO  WHILE ( data_output(i) /= ' '  .AND.  i <= 100 )
2827          i = i + 1
2828       ENDDO
2829       j = 1
2830       DO  WHILE ( data_output_user(j) /= ' '  .AND.  j <= 100 )
2831          IF ( i > 100 )  THEN
2832             message_string = 'number of output quantitities given by data' // &
2833                '_output and data_output_user exceeds the limit of 100'
2834             CALL message( 'check_parameters', 'PA0102', 1, 2, 0, 6, 0 )
2835          ENDIF
2836          data_output(i) = data_output_user(j)
2837          i = i + 1
2838          j = j + 1
2839       ENDDO
2840    ENDIF
2841
2842!
2843!-- Check and set steering parameters for 2d/3d data output and averaging
2844    i   = 1
2845    DO  WHILE ( data_output(i) /= ' '  .AND.  i <= 500 )
2846!
2847!--    Check for data averaging
2848       ilen = LEN_TRIM( data_output(i) )
2849       j = 0                                                 ! no data averaging
2850       IF ( ilen > 3 )  THEN
2851          IF ( data_output(i)(ilen-2:ilen) == '_av' )  THEN
2852             j = 1                                           ! data averaging
2853             data_output(i) = data_output(i)(1:ilen-3)
2854          ENDIF
2855       ENDIF
2856!
2857!--    Check for cross section or volume data
2858       ilen = LEN_TRIM( data_output(i) )
2859       k = 0                                                   ! 3d data
2860       var = data_output(i)(1:ilen)
2861       IF ( ilen > 3 )  THEN
2862          IF ( data_output(i)(ilen-2:ilen) == '_xy'  .OR.                      &
2863               data_output(i)(ilen-2:ilen) == '_xz'  .OR.                      &
2864               data_output(i)(ilen-2:ilen) == '_yz' )  THEN
2865             k = 1                                             ! 2d data
2866             var = data_output(i)(1:ilen-3)
2867          ENDIF
2868       ENDIF
2869
2870!
2871!--    Check for allowed value and set units
2872       SELECT CASE ( TRIM( var ) )
2873
2874          CASE ( 'e' )
2875             IF ( constant_diffusion )  THEN
2876                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2877                                 'res constant_diffusion = .FALSE.'
2878                CALL message( 'check_parameters', 'PA0103', 1, 2, 0, 6, 0 )
2879             ENDIF
2880             unit = 'm2/s2'
2881
2882          CASE ( 'lpt' )
2883             IF (  .NOT.  cloud_physics )  THEN
2884                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2885                         'res cloud_physics = .TRUE.'
2886                CALL message( 'check_parameters', 'PA0108', 1, 2, 0, 6, 0 )
2887             ENDIF
2888             unit = 'K'
2889
2890          CASE ( 'nr' )
2891             IF (  .NOT.  cloud_physics )  THEN
2892                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2893                         'res cloud_physics = .TRUE.'
2894                CALL message( 'check_parameters', 'PA0108', 1, 2, 0, 6, 0 )
2895             ELSEIF ( .NOT.  microphysics_seifert )  THEN
2896                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2897                         'res cloud_scheme = seifert_beheng'
2898                CALL message( 'check_parameters', 'PA0359', 1, 2, 0, 6, 0 )
2899             ENDIF
2900             unit = '1/m3'
2901
2902          CASE ( 'pc', 'pr' )
2903             IF (  .NOT.  particle_advection )  THEN
2904                message_string = 'output of "' // TRIM( var ) // '" requir' // &
2905                   'es a "particles_par"-NAMELIST in the parameter file (PARIN)'
2906                CALL message( 'check_parameters', 'PA0104', 1, 2, 0, 6, 0 )
2907             ENDIF
2908             IF ( TRIM( var ) == 'pc' )  unit = 'number'
2909             IF ( TRIM( var ) == 'pr' )  unit = 'm'
2910
2911          CASE ( 'prr' )
2912             IF (  .NOT.  cloud_physics )  THEN
2913                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2914                         'res cloud_physics = .TRUE.'
2915                CALL message( 'check_parameters', 'PA0108', 1, 2, 0, 6, 0 )
2916             ELSEIF ( microphysics_sat_adjust )  THEN
2917                message_string = 'output of "' // TRIM( var ) // '" is ' //  &
2918                         'not available for cloud_scheme = saturation_adjust'
2919                CALL message( 'check_parameters', 'PA0423', 1, 2, 0, 6, 0 )
2920             ENDIF
2921             unit = 'kg/kg m/s'
2922
2923          CASE ( 'q', 'vpt' )
2924             IF (  .NOT.  humidity )  THEN
2925                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2926                                 'res humidity = .TRUE.'
2927                CALL message( 'check_parameters', 'PA0105', 1, 2, 0, 6, 0 )
2928             ENDIF
2929             IF ( TRIM( var ) == 'q'   )  unit = 'kg/kg'
2930             IF ( TRIM( var ) == 'vpt' )  unit = 'K'
2931
2932          CASE ( 'qc' )
2933             IF (  .NOT.  cloud_physics )  THEN
2934                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2935                         'res cloud_physics = .TRUE.'
2936                CALL message( 'check_parameters', 'PA0108', 1, 2, 0, 6, 0 )
2937             ENDIF
2938             unit = 'kg/kg'
2939
2940          CASE ( 'ql' )
2941             IF ( .NOT.  ( cloud_physics  .OR.  cloud_droplets ) )  THEN
2942                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2943                         'res cloud_physics = .TRUE. or cloud_droplets = .TRUE.'
2944                CALL message( 'check_parameters', 'PA0106', 1, 2, 0, 6, 0 )
2945             ENDIF
2946             unit = 'kg/kg'
2947
2948          CASE ( 'ql_c', 'ql_v', 'ql_vp' )
2949             IF (  .NOT.  cloud_droplets )  THEN
2950                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2951                                 'res cloud_droplets = .TRUE.'
2952                CALL message( 'check_parameters', 'PA0107', 1, 2, 0, 6, 0 )
2953             ENDIF
2954             IF ( TRIM( var ) == 'ql_c'  )  unit = 'kg/kg'
2955             IF ( TRIM( var ) == 'ql_v'  )  unit = 'm3'
2956             IF ( TRIM( var ) == 'ql_vp' )  unit = 'none'
2957
2958          CASE ( 'qr' )
2959             IF (  .NOT.  cloud_physics )  THEN
2960                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2961                         'res cloud_physics = .TRUE.'
2962                CALL message( 'check_parameters', 'PA0108', 1, 2, 0, 6, 0 )
2963             ELSEIF ( .NOT.  microphysics_seifert ) THEN
2964                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2965                         'res cloud_scheme = seifert_beheng'
2966                CALL message( 'check_parameters', 'PA0359', 1, 2, 0, 6, 0 )
2967             ENDIF
2968             unit = 'kg/kg'
2969
2970          CASE ( 'qv' )
2971             IF (  .NOT.  cloud_physics )  THEN
2972                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2973                                 'res cloud_physics = .TRUE.'
2974                CALL message( 'check_parameters', 'PA0108', 1, 2, 0, 6, 0 )
2975             ENDIF
2976             unit = 'kg/kg'
2977
2978          CASE ( 'rho' )
2979             IF (  .NOT.  ocean )  THEN
2980                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2981                                 'res ocean = .TRUE.'
2982                CALL message( 'check_parameters', 'PA0109', 1, 2, 0, 6, 0 )
2983             ENDIF
2984             unit = 'kg/m3'
2985
2986          CASE ( 's' )
2987             IF (  .NOT.  passive_scalar )  THEN
2988                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2989                                 'res passive_scalar = .TRUE.'
2990                CALL message( 'check_parameters', 'PA0110', 1, 2, 0, 6, 0 )
2991             ENDIF
2992             unit = 'conc'
2993
2994          CASE ( 'sa' )
2995             IF (  .NOT.  ocean )  THEN
2996                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
2997                                 'res ocean = .TRUE.'
2998                CALL message( 'check_parameters', 'PA0109', 1, 2, 0, 6, 0 )
2999             ENDIF
3000             unit = 'psu'
3001
3002          CASE ( 'p', 'pt', 'u', 'v', 'w' )
3003             IF ( TRIM( var ) == 'p'  )  unit = 'Pa'
3004             IF ( TRIM( var ) == 'pt' )  unit = 'K'
3005             IF ( TRIM( var ) == 'u'  )  unit = 'm/s'
3006             IF ( TRIM( var ) == 'v'  )  unit = 'm/s'
3007             IF ( TRIM( var ) == 'w'  )  unit = 'm/s'
3008             CONTINUE
3009
3010          CASE ( 'lwp*', 'ol*', 'pra*', 'prr*', 'qsws*', 'shf*', 't*', &
3011                 'u*', 'z0*', 'z0h*', 'z0q*' )
3012             IF ( k == 0  .OR.  data_output(i)(ilen-2:ilen) /= '_xy' )  THEN
3013                message_string = 'illegal value for data_output: "' //         &
3014                                 TRIM( var ) // '" & only 2d-horizontal ' //   &
3015                                 'cross sections are allowed for this value'
3016                CALL message( 'check_parameters', 'PA0111', 1, 2, 0, 6, 0 )
3017             ENDIF
3018
3019             IF ( TRIM( var ) == 'lwp*'  .AND.  .NOT. cloud_physics )  THEN
3020                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
3021                                 'res cloud_physics = .TRUE.'
3022                CALL message( 'check_parameters', 'PA0108', 1, 2, 0, 6, 0 )
3023             ENDIF
3024             IF ( TRIM( var ) == 'pra*'  .AND.                                 &
3025                  .NOT. ( microphysics_kessler .OR. microphysics_seifert ) )  THEN
3026                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
3027                                 'res cloud_scheme = kessler or seifert_beheng'
3028                CALL message( 'check_parameters', 'PA0112', 1, 2, 0, 6, 0 )
3029             ENDIF
3030             IF ( TRIM( var ) == 'pra*'  .AND.  j == 1 )  THEN
3031                message_string = 'temporal averaging of precipitation ' //     &
3032                          'amount "' // TRIM( var ) // '" is not possible'
3033                CALL message( 'check_parameters', 'PA0113', 1, 2, 0, 6, 0 )
3034             ENDIF
3035             IF ( TRIM( var ) == 'prr*'  .AND.                                 &
3036                  .NOT. ( microphysics_kessler .OR. microphysics_seifert ) )  THEN
3037                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
3038                                 'res cloud_scheme = kessler or seifert_beheng'
3039                CALL message( 'check_parameters', 'PA0112', 1, 2, 0, 6, 0 )
3040             ENDIF
3041             IF ( TRIM( var ) == 'qsws*'  .AND.  .NOT.  humidity )  THEN
3042                message_string = 'output of "' // TRIM( var ) // '" requi' //  &
3043                                 'res humidity = .TRUE.'
3044                CALL message( 'check_parameters', 'PA0322', 1, 2, 0, 6, 0 )
3045             ENDIF
3046
3047             IF ( TRIM( var ) == 'lwp*'   )  unit = 'kg/m2'
3048             IF ( TRIM( var ) == 'ol*'   )   unit = 'm'
3049             IF ( TRIM( var ) == 'pra*'   )  unit = 'mm'
3050             IF ( TRIM( var ) == 'prr*'   )  unit = 'mm/s'
3051             IF ( TRIM( var ) == 'qsws*'  )  unit = 'kgm/kgs'
3052             IF ( TRIM( var ) == 'shf*'   )  unit = 'K*m/s'
3053             IF ( TRIM( var ) == 't*'     )  unit = 'K'
3054             IF ( TRIM( var ) == 'u*'     )  unit = 'm/s'
3055             IF ( TRIM( var ) == 'z0*'    )  unit = 'm'
3056             IF ( TRIM( var ) == 'z0h*'   )  unit = 'm' 
3057             
3058          CASE DEFAULT
3059
3060             CALL lsm_check_data_output ( var, unit, i, ilen, k )
3061 
3062             IF ( unit == 'illegal' )  THEN
3063                CALL radiation_check_data_output( var, unit, i, ilen, k )
3064             ENDIF
3065
3066!
3067!--          Block of urban surface model outputs
3068             IF ( unit == 'illegal' .AND. urban_surface .AND. var(1:3) == 'us_' ) THEN
3069                 CALL usm_check_data_output( var, unit )
3070             ENDIF
3071
3072             IF ( unit == 'illegal' )  THEN
3073                unit = ''
3074                CALL user_check_data_output( var, unit )
3075             ENDIF
3076
3077             IF ( unit == 'illegal' )  THEN
3078                IF ( data_output_user(1) /= ' ' )  THEN
3079                   message_string = 'illegal value for data_output or ' //     &
3080                         'data_output_user = "' // TRIM( data_output(i) ) // '"'
3081                   CALL message( 'check_parameters', 'PA0114', 1, 2, 0, 6, 0 )
3082                ELSE
3083                   message_string = 'illegal value for data_output = "' //     &
3084                                    TRIM( data_output(i) ) // '"'
3085                   CALL message( 'check_parameters', 'PA0115', 1, 2, 0, 6, 0 )
3086                ENDIF
3087             ENDIF
3088
3089       END SELECT
3090!
3091!--    Set the internal steering parameters appropriately
3092       IF ( k == 0 )  THEN
3093          do3d_no(j)              = do3d_no(j) + 1
3094          do3d(j,do3d_no(j))      = data_output(i)
3095          do3d_unit(j,do3d_no(j)) = unit
3096       ELSE
3097          do2d_no(j)              = do2d_no(j) + 1
3098          do2d(j,do2d_no(j))      = data_output(i)
3099          do2d_unit(j,do2d_no(j)) = unit
3100          IF ( data_output(i)(ilen-2:ilen) == '_xy' )  THEN
3101             data_output_xy(j) = .TRUE.
3102          ENDIF
3103          IF ( data_output(i)(ilen-2:ilen) == '_xz' )  THEN
3104             data_output_xz(j) = .TRUE.
3105          ENDIF
3106          IF ( data_output(i)(ilen-2:ilen) == '_yz' )  THEN
3107             data_output_yz(j) = .TRUE.
3108          ENDIF
3109       ENDIF
3110
3111       IF ( j == 1 )  THEN
3112!
3113!--       Check, if variable is already subject to averaging
3114          found = .FALSE.
3115          DO  k = 1, doav_n
3116             IF ( TRIM( doav(k) ) == TRIM( var ) )  found = .TRUE.
3117          ENDDO
3118
3119          IF ( .NOT. found )  THEN
3120             doav_n = doav_n + 1
3121             doav(doav_n) = var
3122          ENDIF
3123       ENDIF
3124
3125       i = i + 1
3126    ENDDO
3127
3128!
3129!-- Averaged 2d or 3d output requires that an averaging interval has been set
3130    IF ( doav_n > 0  .AND.  averaging_interval == 0.0_wp )  THEN
3131       WRITE( message_string, * )  'output of averaged quantity "',            &
3132                                   TRIM( doav(1) ), '_av" requires to set a ', &
3133                                   'non-zero & averaging interval'
3134       CALL message( 'check_parameters', 'PA0323', 1, 2, 0, 6, 0 )
3135    ENDIF
3136
3137!
3138!-- Check sectional planes and store them in one shared array
3139    IF ( ANY( section_xy > nz + 1 ) )  THEN
3140       WRITE( message_string, * )  'section_xy must be <= nz + 1 = ', nz + 1
3141       CALL message( 'check_parameters', 'PA0319', 1, 2, 0, 6, 0 )
3142    ENDIF
3143    IF ( ANY( section_xz > ny + 1 ) )  THEN
3144       WRITE( message_string, * )  'section_xz must be <= ny + 1 = ', ny + 1
3145       CALL message( 'check_parameters', 'PA0320', 1, 2, 0, 6, 0 )
3146    ENDIF
3147    IF ( ANY( section_yz > nx + 1 ) )  THEN
3148       WRITE( message_string, * )  'section_yz must be <= nx + 1 = ', nx + 1
3149       CALL message( 'check_parameters', 'PA0321', 1, 2, 0, 6, 0 )
3150    ENDIF
3151    section(:,1) = section_xy
3152    section(:,2) = section_xz
3153    section(:,3) = section_yz
3154
3155!
3156!-- Upper plot limit for 2D vertical sections
3157    IF ( z_max_do2d == -1.0_wp )  z_max_do2d = zu(nzt)
3158    IF ( z_max_do2d < zu(nzb+1)  .OR.  z_max_do2d > zu(nzt) )  THEN
3159       WRITE( message_string, * )  'z_max_do2d = ', z_max_do2d,                &
3160                    ' must be >= ', zu(nzb+1), '(zu(nzb+1)) and <= ', zu(nzt), &
3161                    ' (zu(nzt))'
3162       CALL message( 'check_parameters', 'PA0116', 1, 2, 0, 6, 0 )
3163    ENDIF
3164
3165!
3166!-- Upper plot limit for 3D arrays
3167    IF ( nz_do3d == -9999 )  nz_do3d = nzt + 1
3168
3169!
3170!-- Set output format string (used in header)
3171    SELECT CASE ( netcdf_data_format )
3172       CASE ( 1 )
3173          netcdf_data_format_string = 'netCDF classic'
3174       CASE ( 2 )
3175          netcdf_data_format_string = 'netCDF 64bit offset'
3176       CASE ( 3 )
3177          netcdf_data_format_string = 'netCDF4/HDF5'
3178       CASE ( 4 )
3179          netcdf_data_format_string = 'netCDF4/HDF5 classic'
3180       CASE ( 5 )
3181          netcdf_data_format_string = 'parallel netCDF4/HDF5'
3182       CASE ( 6 )
3183          netcdf_data_format_string = 'parallel netCDF4/HDF5 classic'
3184
3185    END SELECT
3186
3187!
3188!-- Check mask conditions
3189    DO mid = 1, max_masks
3190       IF ( data_output_masks(mid,1) /= ' '  .OR.                              &
3191            data_output_masks_user(mid,1) /= ' ' ) THEN
3192          masks = masks + 1
3193       ENDIF
3194    ENDDO
3195   
3196    IF ( masks < 0  .OR.  masks > max_masks )  THEN
3197       WRITE( message_string, * )  'illegal value: masks must be >= 0 and ',   &
3198            '<= ', max_masks, ' (=max_masks)'
3199       CALL message( 'check_parameters', 'PA0325', 1, 2, 0, 6, 0 )
3200    ENDIF
3201    IF ( masks > 0 )  THEN
3202       mask_scale(1) = mask_scale_x
3203       mask_scale(2) = mask_scale_y
3204       mask_scale(3) = mask_scale_z
3205       IF ( ANY( mask_scale <= 0.0_wp ) )  THEN
3206          WRITE( message_string, * )                                           &
3207               'illegal value: mask_scale_x, mask_scale_y and mask_scale_z',   &
3208               'must be > 0.0'
3209          CALL message( 'check_parameters', 'PA0326', 1, 2, 0, 6, 0 )
3210       ENDIF
3211!
3212!--    Generate masks for masked data output
3213!--    Parallel netcdf output is not tested so far for masked data, hence
3214!--    netcdf_data_format is switched back to non-paralell output.
3215       netcdf_data_format_save = netcdf_data_format
3216       IF ( netcdf_data_format > 4 )  THEN
3217          IF ( netcdf_data_format == 5 ) netcdf_data_format = 3
3218          IF ( netcdf_data_format == 6 ) netcdf_data_format = 4
3219          message_string = 'netCDF file formats '//                            &
3220                           '5 (parallel netCDF 4) and ' //                     &
3221                           '6 (parallel netCDF 4 Classic model) '//            &
3222                           '&are currently not supported (not yet tested) ' // &
3223                           'for masked data.&Using respective non-parallel' // & 
3224                           ' output for masked data.'
3225          CALL message( 'check_parameters', 'PA0383', 0, 0, 0, 6, 0 )
3226       ENDIF
3227       CALL init_masks
3228       netcdf_data_format = netcdf_data_format_save
3229    ENDIF
3230
3231!
3232!-- Check the NetCDF data format
3233    IF ( netcdf_data_format > 2 )  THEN
3234#if defined( __netcdf4 )
3235       CONTINUE
3236#else
3237       message_string = 'netCDF: netCDF4 format requested but no ' //          &
3238                        'cpp-directive __netcdf4 given & switch '  //          &
3239                        'back to 64-bit offset format'
3240       CALL message( 'check_parameters', 'PA0171', 0, 1, 0, 6, 0 )
3241       netcdf_data_format = 2
3242#endif
3243    ENDIF
3244    IF ( netcdf_data_format > 4 )  THEN
3245#if defined( __netcdf4 ) && defined( __netcdf4_parallel )
3246       CONTINUE
3247#else
3248       message_string = 'netCDF: netCDF4 parallel output requested but no ' // &
3249                        'cpp-directive __netcdf4_parallel given & switch '  // &
3250                        'back to netCDF4 non-parallel output'
3251       CALL message( 'check_parameters', 'PA0099', 0, 1, 0, 6, 0 )
3252       netcdf_data_format = netcdf_data_format - 2
3253#endif
3254    ENDIF
3255
3256!
3257!-- Calculate fixed number of output time levels for parallel netcdf output.
3258!-- The time dimension has to be defined as limited for parallel output,
3259!-- because otherwise the I/O performance drops significantly.
3260    IF ( netcdf_data_format > 4 )  THEN
3261
3262!
3263!--    Check if any of the follwoing data output interval is 0.0s, which is
3264!--    not allowed for parallel output.
3265       CALL check_dt_do( dt_do3d,    'dt_do3d'    )
3266       CALL check_dt_do( dt_do2d_xy, 'dt_do2d_xy' )
3267       CALL check_dt_do( dt_do2d_xz, 'dt_do2d_xz' )
3268       CALL check_dt_do( dt_do2d_yz, 'dt_do2d_yz' )
3269
3270       ntdim_3d(0)    = INT( ( end_time - skip_time_do3d ) / dt_do3d )
3271       IF ( do3d_at_begin ) ntdim_3d(0) = ntdim_3d(0) + 1
3272       ntdim_3d(1)    = INT( ( end_time - skip_time_data_output_av )           &
3273                             / dt_data_output_av )
3274       ntdim_2d_xy(0) = INT( ( end_time - skip_time_do2d_xy ) / dt_do2d_xy )
3275       ntdim_2d_xz(0) = INT( ( end_time - skip_time_do2d_xz ) / dt_do2d_xz )
3276       ntdim_2d_yz(0) = INT( ( end_time - skip_time_do2d_yz ) / dt_do2d_yz )
3277       IF ( do2d_at_begin )  THEN
3278          ntdim_2d_xy(0) = ntdim_2d_xy(0) + 1
3279          ntdim_2d_xz(0) = ntdim_2d_xz(0) + 1
3280          ntdim_2d_yz(0) = ntdim_2d_yz(0) + 1
3281       ENDIF
3282       ntdim_2d_xy(1) = ntdim_3d(1)
3283       ntdim_2d_xz(1) = ntdim_3d(1)
3284       ntdim_2d_yz(1) = ntdim_3d(1)
3285
3286    ENDIF
3287
3288!
3289!-- Check, whether a constant diffusion coefficient shall be used
3290    IF ( km_constant /= -1.0_wp )  THEN
3291       IF ( km_constant < 0.0_wp )  THEN
3292          WRITE( message_string, * )  'km_constant = ', km_constant, ' < 0.0'
3293          CALL message( 'check_parameters', 'PA0121', 1, 2, 0, 6, 0 )
3294       ELSE
3295          IF ( prandtl_number < 0.0_wp )  THEN
3296             WRITE( message_string, * )  'prandtl_number = ', prandtl_number,  &
3297                                         ' < 0.0'
3298             CALL message( 'check_parameters', 'PA0122', 1, 2, 0, 6, 0 )
3299          ENDIF
3300          constant_diffusion = .TRUE.
3301
3302          IF ( constant_flux_layer )  THEN
3303             message_string = 'constant_flux_layer is not allowed with fixed ' &
3304                              // 'value of km'
3305             CALL message( 'check_parameters', 'PA0123', 1, 2, 0, 6, 0 )
3306          ENDIF
3307       ENDIF
3308    ENDIF
3309
3310!
3311!-- In case of non-cyclic lateral boundaries and a damping layer for the
3312!-- potential temperature, check the width of the damping layer
3313    IF ( bc_lr /= 'cyclic' ) THEN
3314       IF ( pt_damping_width < 0.0_wp  .OR.                                    &
3315            pt_damping_width > REAL( nx * dx ) )  THEN
3316          message_string = 'pt_damping_width out of range'
3317          CALL message( 'check_parameters', 'PA0124', 1, 2, 0, 6, 0 )
3318       ENDIF
3319    ENDIF
3320
3321    IF ( bc_ns /= 'cyclic' )  THEN
3322       IF ( pt_damping_width < 0.0_wp  .OR.                                    &
3323            pt_damping_width > REAL( ny * dy ) )  THEN
3324          message_string = 'pt_damping_width out of range'
3325          CALL message( 'check_parameters', 'PA0124', 1, 2, 0, 6, 0 )
3326       ENDIF
3327    ENDIF
3328
3329!
3330!-- Check value range for zeta = z/L
3331    IF ( zeta_min >= zeta_max )  THEN
3332       WRITE( message_string, * )  'zeta_min = ', zeta_min, ' must be less ',  &
3333                                   'than zeta_max = ', zeta_max
3334       CALL message( 'check_parameters', 'PA0125', 1, 2, 0, 6, 0 )
3335    ENDIF
3336
3337!
3338!-- Check random generator
3339    IF ( (random_generator /= 'system-specific'      .AND.                     &
3340          random_generator /= 'random-parallel'   )  .AND.                     &
3341          random_generator /= 'numerical-recipes' )  THEN
3342       message_string = 'unknown random generator: random_generator = "' //    &
3343                        TRIM( random_generator ) // '"'
3344       CALL message( 'check_parameters', 'PA0135', 1, 2, 0, 6, 0 )
3345    ENDIF
3346
3347!
3348!-- Determine upper and lower hight level indices for random perturbations
3349    IF ( disturbance_level_b == -9999999.9_wp )  THEN
3350       IF ( ocean )  THEN
3351          disturbance_level_b     = zu((nzt*2)/3)
3352          disturbance_level_ind_b = ( nzt * 2 ) / 3
3353       ELSE
3354          disturbance_level_b     = zu(nzb+3)
3355          disturbance_level_ind_b = nzb + 3
3356       ENDIF
3357    ELSEIF ( disturbance_level_b < zu(3) )  THEN
3358       WRITE( message_string, * )  'disturbance_level_b = ',                   &
3359                           disturbance_level_b, ' must be >= ', zu(3), '(zu(3))'
3360       CALL message( 'check_parameters', 'PA0126', 1, 2, 0, 6, 0 )
3361    ELSEIF ( disturbance_level_b > zu(nzt-2) )  THEN
3362       WRITE( message_string, * )  'disturbance_level_b = ',                   &
3363                   disturbance_level_b, ' must be <= ', zu(nzt-2), '(zu(nzt-2))'
3364       CALL message( 'check_parameters', 'PA0127', 1, 2, 0, 6, 0 )
3365    ELSE
3366       DO  k = 3, nzt-2
3367          IF ( disturbance_level_b <= zu(k) )  THEN
3368             disturbance_level_ind_b = k
3369             EXIT
3370          ENDIF
3371       ENDDO
3372    ENDIF
3373
3374    IF ( disturbance_level_t == -9999999.9_wp )  THEN
3375       IF ( ocean )  THEN
3376          disturbance_level_t     = zu(nzt-3)
3377          disturbance_level_ind_t = nzt - 3
3378       ELSE
3379          disturbance_level_t     = zu(nzt/3)
3380          disturbance_level_ind_t = nzt / 3
3381       ENDIF
3382    ELSEIF ( disturbance_level_t > zu(nzt-2) )  THEN
3383       WRITE( message_string, * )  'disturbance_level_t = ',                   &
3384                   disturbance_level_t, ' must be <= ', zu(nzt-2), '(zu(nzt-2))'
3385       CALL message( 'check_parameters', 'PA0128', 1, 2, 0, 6, 0 )
3386    ELSEIF ( disturbance_level_t < disturbance_level_b )  THEN
3387       WRITE( message_string, * )  'disturbance_level_t = ',                   &
3388                   disturbance_level_t, ' must be >= disturbance_level_b = ',  &
3389                   disturbance_level_b
3390       CALL message( 'check_parameters', 'PA0129', 1, 2, 0, 6, 0 )
3391    ELSE
3392       DO  k = 3, nzt-2
3393          IF ( disturbance_level_t <= zu(k) )  THEN
3394             disturbance_level_ind_t = k
3395             EXIT
3396          ENDIF
3397       ENDDO
3398    ENDIF
3399
3400!
3401!-- Check again whether the levels determined this way are ok.
3402!-- Error may occur at automatic determination and too few grid points in
3403!-- z-direction.
3404    IF ( disturbance_level_ind_t < disturbance_level_ind_b )  THEN
3405       WRITE( message_string, * )  'disturbance_level_ind_t = ',               &
3406                disturbance_level_ind_t, ' must be >= disturbance_level_b = ', &
3407                disturbance_level_b
3408       CALL message( 'check_parameters', 'PA0130', 1, 2, 0, 6, 0 )
3409    ENDIF
3410
3411!
3412!-- Determine the horizontal index range for random perturbations.
3413!-- In case of non-cyclic horizontal boundaries, no perturbations are imposed
3414!-- near the inflow and the perturbation area is further limited to ...(1)
3415!-- after the initial phase of the flow.
3416   
3417    IF ( bc_lr /= 'cyclic' )  THEN
3418       IF ( inflow_disturbance_begin == -1 )  THEN
3419          inflow_disturbance_begin = MIN( 10, nx/2 )
3420       ENDIF
3421       IF ( inflow_disturbance_begin < 0  .OR.  inflow_disturbance_begin > nx )&
3422       THEN
3423          message_string = 'inflow_disturbance_begin out of range'
3424          CALL message( 'check_parameters', 'PA0131', 1, 2, 0, 6, 0 )
3425       ENDIF
3426       IF ( inflow_disturbance_end == -1 )  THEN
3427          inflow_disturbance_end = MIN( 100, 3*nx/4 )
3428       ENDIF
3429       IF ( inflow_disturbance_end < 0  .OR.  inflow_disturbance_end > nx )    &
3430       THEN
3431          message_string = 'inflow_disturbance_end out of range'
3432          CALL message( 'check_parameters', 'PA0132', 1, 2, 0, 6, 0 )
3433       ENDIF
3434    ELSEIF ( bc_ns /= 'cyclic' )  THEN
3435       IF ( inflow_disturbance_begin == -1 )  THEN
3436          inflow_disturbance_begin = MIN( 10, ny/2 )
3437       ENDIF
3438       IF ( inflow_disturbance_begin < 0  .OR.  inflow_disturbance_begin > ny )&
3439       THEN
3440          message_string = 'inflow_disturbance_begin out of range'
3441          CALL message( 'check_parameters', 'PA0131', 1, 2, 0, 6, 0 )
3442       ENDIF
3443       IF ( inflow_disturbance_end == -1 )  THEN
3444          inflow_disturbance_end = MIN( 100, 3*ny/4 )
3445       ENDIF
3446       IF ( inflow_disturbance_end < 0  .OR.  inflow_disturbance_end > ny )    &
3447       THEN
3448          message_string = 'inflow_disturbance_end out of range'
3449          CALL message( 'check_parameters', 'PA0132', 1, 2, 0, 6, 0 )
3450       ENDIF
3451    ENDIF
3452
3453    IF ( random_generator == 'random-parallel' )  THEN
3454       dist_nxl = nxl;  dist_nxr = nxr
3455       dist_nys = nys;  dist_nyn = nyn
3456       IF ( bc_lr == 'radiation/dirichlet' )  THEN
3457          dist_nxr    = MIN( nx - inflow_disturbance_begin, nxr )
3458          dist_nxl(1) = MAX( nx - inflow_disturbance_end, nxl )
3459       ELSEIF ( bc_lr == 'dirichlet/radiation' )  THEN
3460          dist_nxl    = MAX( inflow_disturbance_begin, nxl )
3461          dist_nxr(1) = MIN( inflow_disturbance_end, nxr )
3462       ENDIF
3463       IF ( bc_ns == 'dirichlet/radiation' )  THEN
3464          dist_nyn    = MIN( ny - inflow_disturbance_begin, nyn )
3465          dist_nys(1) = MAX( ny - inflow_disturbance_end, nys )
3466       ELSEIF ( bc_ns == 'radiation/dirichlet' )  THEN
3467          dist_nys    = MAX( inflow_disturbance_begin, nys )
3468          dist_nyn(1) = MIN( inflow_disturbance_end, nyn )
3469       ENDIF
3470    ELSE
3471       dist_nxl = 0;  dist_nxr = nx
3472       dist_nys = 0;  dist_nyn = ny
3473       IF ( bc_lr == 'radiation/dirichlet' )  THEN
3474          dist_nxr    = nx - inflow_disturbance_begin
3475          dist_nxl(1) = nx - inflow_disturbance_end
3476       ELSEIF ( bc_lr == 'dirichlet/radiation' )  THEN
3477          dist_nxl    = inflow_disturbance_begin
3478          dist_nxr(1) = inflow_disturbance_end
3479       ENDIF
3480       IF ( bc_ns == 'dirichlet/radiation' )  THEN
3481          dist_nyn    = ny - inflow_disturbance_begin
3482          dist_nys(1) = ny - inflow_disturbance_end
3483       ELSEIF ( bc_ns == 'radiation/dirichlet' )  THEN
3484          dist_nys    = inflow_disturbance_begin
3485          dist_nyn(1) = inflow_disturbance_end
3486       ENDIF
3487    ENDIF
3488
3489!
3490!-- A turbulent inflow requires Dirichlet conditions at the respective inflow
3491!-- boundary (so far, a turbulent inflow is realized from the left side only)
3492    IF ( turbulent_inflow  .AND.  bc_lr /= 'dirichlet/radiation' )  THEN
3493       message_string = 'turbulent_inflow = .T. requires a Dirichlet ' //      &
3494                        'condition at the inflow boundary'
3495       CALL message( 'check_parameters', 'PA0133', 1, 2, 0, 6, 0 )
3496    ENDIF
3497
3498!
3499!-- Turbulent inflow requires that 3d arrays have been cyclically filled with
3500!-- data from prerun in the first main run
3501    IF ( turbulent_inflow  .AND.  initializing_actions /= 'cyclic_fill'  .AND. &
3502         initializing_actions /= 'read_restart_data' )  THEN
3503       message_string = 'turbulent_inflow = .T. requires ' //                  &
3504                        'initializing_actions = ''cyclic_fill'' '
3505       CALL message( 'check_parameters', 'PA0055', 1, 2, 0, 6, 0 )
3506    ENDIF
3507
3508!
3509!-- In case of turbulent inflow calculate the index of the recycling plane
3510    IF ( turbulent_inflow )  THEN
3511       IF ( recycling_width == 9999999.9_wp )  THEN
3512!
3513!--       Set the default value for the width of the recycling domain
3514          recycling_width = 0.1_wp * nx * dx
3515       ELSE
3516          IF ( recycling_width < dx  .OR.  recycling_width > nx * dx )  THEN
3517             WRITE( message_string, * )  'illegal value for recycling_width:', &
3518                                         ' ', recycling_width
3519             CALL message( 'check_parameters', 'PA0134', 1, 2, 0, 6, 0 )
3520          ENDIF
3521       ENDIF
3522!
3523!--    Calculate the index
3524       recycling_plane = recycling_width / dx
3525!
3526!--    Because the y-shift is done with a distance of INT( npey / 2 ) no shift
3527!--    is possible if there is only one PE in y direction.
3528       IF ( recycling_yshift .AND. pdims(2) < 2 )  THEN
3529          WRITE( message_string, * )  'recycling_yshift = .T. requires more',  &
3530                                      ' than one processor in y direction'
3531          CALL message( 'check_parameters', 'PA0421', 1, 2, 0, 6, 0 )
3532       ENDIF
3533    ENDIF
3534
3535!
3536!-- Determine damping level index for 1D model
3537    IF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
3538       IF ( damp_level_1d == -1.0_wp )  THEN
3539          damp_level_1d     = zu(nzt+1)
3540          damp_level_ind_1d = nzt + 1
3541       ELSEIF ( damp_level_1d < 0.0_wp  .OR.  damp_level_1d > zu(nzt+1) )  THEN
3542          WRITE( message_string, * )  'damp_level_1d = ', damp_level_1d,       &
3543                 ' must be > 0.0 and < ', zu(nzt+1), '(zu(nzt+1))'
3544          CALL message( 'check_parameters', 'PA0136', 1, 2, 0, 6, 0 )
3545       ELSE
3546          DO  k = 1, nzt+1
3547             IF ( damp_level_1d <= zu(k) )  THEN
3548                damp_level_ind_1d = k
3549                EXIT
3550             ENDIF
3551          ENDDO
3552       ENDIF
3553    ENDIF
3554
3555!
3556!-- Check some other 1d-model parameters
3557    IF ( TRIM( mixing_length_1d ) /= 'as_in_3d_model'  .AND.                   &
3558         TRIM( mixing_length_1d ) /= 'blackadar' )  THEN
3559       message_string = 'mixing_length_1d = "' // TRIM( mixing_length_1d ) //  &
3560                        '" is unknown'
3561       CALL message( 'check_parameters', 'PA0137', 1, 2, 0, 6, 0 )
3562    ENDIF
3563    IF ( TRIM( dissipation_1d ) /= 'as_in_3d_model'  .AND.                     &
3564         TRIM( dissipation_1d ) /= 'detering' )  THEN
3565       message_string = 'dissipation_1d = "' // TRIM( dissipation_1d ) //      &
3566                        '" is unknown'
3567       CALL message( 'check_parameters', 'PA0138', 1, 2, 0, 6, 0 )
3568    ENDIF
3569
3570!
3571!-- Set time for the next user defined restart (time_restart is the
3572!-- internal parameter for steering restart events)
3573    IF ( restart_time /= 9999999.9_wp )  THEN
3574       IF ( restart_time > time_since_reference_point )  THEN
3575          time_restart = restart_time
3576       ENDIF
3577    ELSE
3578!
3579!--    In case of a restart run, set internal parameter to default (no restart)
3580!--    if the NAMELIST-parameter restart_time is omitted
3581       time_restart = 9999999.9_wp
3582    ENDIF
3583
3584!
3585!-- Set default value of the time needed to terminate a model run
3586    IF ( termination_time_needed == -1.0_wp )  THEN
3587       IF ( host(1:3) == 'ibm' )  THEN
3588          termination_time_needed = 300.0_wp
3589       ELSE
3590          termination_time_needed = 35.0_wp
3591       ENDIF
3592    ENDIF
3593
3594!
3595!-- Check the time needed to terminate a model run
3596    IF ( host(1:3) == 't3e' )  THEN
3597!
3598!--    Time needed must be at least 30 seconds on all CRAY machines, because
3599!--    MPP_TREMAIN gives the remaining CPU time only in steps of 30 seconds
3600       IF ( termination_time_needed <= 30.0_wp )  THEN
3601          WRITE( message_string, * )  'termination_time_needed = ',            &
3602                 termination_time_needed, ' must be > 30.0 on host "',         &
3603                 TRIM( host ), '"'
3604          CALL message( 'check_parameters', 'PA0139', 1, 2, 0, 6, 0 )
3605       ENDIF
3606    ELSEIF ( host(1:3) == 'ibm' )  THEN
3607!
3608!--    On IBM-regatta machines the time should be at least 300 seconds,
3609!--    because the job time consumed before executing palm (for compiling,
3610!--    copying of files, etc.) has to be regarded
3611       IF ( termination_time_needed < 300.0_wp )  THEN
3612          WRITE( message_string, * )  'termination_time_needed = ',            &
3613                 termination_time_needed, ' should be >= 300.0 on host "',     &
3614                 TRIM( host ), '"'
3615          CALL message( 'check_parameters', 'PA0140', 1, 2, 0, 6, 0 )
3616       ENDIF
3617    ENDIF
3618
3619!
3620!-- Check pressure gradient conditions
3621    IF ( dp_external  .AND.  conserve_volume_flow )  THEN
3622       WRITE( message_string, * )  'Both dp_external and conserve_volume_flo', &
3623            'w are .TRUE. but one of them must be .FALSE.'
3624       CALL message( 'check_parameters', 'PA0150', 1, 2, 0, 6, 0 )
3625    ENDIF
3626    IF ( dp_external )  THEN
3627       IF ( dp_level_b < zu(nzb)  .OR.  dp_level_b > zu(nzt) )  THEN
3628          WRITE( message_string, * )  'dp_level_b = ', dp_level_b, ' is out ', &
3629               ' of range'
3630          CALL message( 'check_parameters', 'PA0151', 1, 2, 0, 6, 0 )
3631       ENDIF
3632       IF ( .NOT. ANY( dpdxy /= 0.0_wp ) )  THEN
3633          WRITE( message_string, * )  'dp_external is .TRUE. but dpdxy is ze', &
3634               'ro, i.e. the external pressure gradient & will not be applied'
3635          CALL message( 'check_parameters', 'PA0152', 0, 1, 0, 6, 0 )
3636       ENDIF
3637    ENDIF
3638    IF ( ANY( dpdxy /= 0.0_wp )  .AND.  .NOT.  dp_external )  THEN
3639       WRITE( message_string, * )  'dpdxy is nonzero but dp_external is ',     &
3640            '.FALSE., i.e. the external pressure gradient & will not be applied'
3641       CALL message( 'check_parameters', 'PA0153', 0, 1, 0, 6, 0 )
3642    ENDIF
3643    IF ( conserve_volume_flow )  THEN
3644       IF ( TRIM( conserve_volume_flow_mode ) == 'default' )  THEN
3645
3646          conserve_volume_flow_mode = 'initial_profiles'
3647
3648       ELSEIF ( TRIM( conserve_volume_flow_mode ) /= 'initial_profiles' .AND.  &
3649            TRIM( conserve_volume_flow_mode ) /= 'inflow_profile' .AND.        &
3650            TRIM( conserve_volume_flow_mode ) /= 'bulk_velocity' )  THEN
3651          WRITE( message_string, * )  'unknown conserve_volume_flow_mode: ',   &
3652               conserve_volume_flow_mode
3653          CALL message( 'check_parameters', 'PA0154', 1, 2, 0, 6, 0 )
3654       ENDIF
3655       IF ( (bc_lr /= 'cyclic'  .OR.  bc_ns /= 'cyclic')  .AND.                &
3656          TRIM( conserve_volume_flow_mode ) == 'bulk_velocity' )  THEN
3657          WRITE( message_string, * )  'non-cyclic boundary conditions ',       &
3658               'require  conserve_volume_flow_mode = ''initial_profiles'''
3659          CALL message( 'check_parameters', 'PA0155', 1, 2, 0, 6, 0 )
3660       ENDIF
3661       IF ( bc_lr == 'cyclic'  .AND.  bc_ns == 'cyclic'  .AND.                 &
3662            TRIM( conserve_volume_flow_mode ) == 'inflow_profile' )  THEN
3663          WRITE( message_string, * )  'cyclic boundary conditions ',           &
3664               'require conserve_volume_flow_mode = ''initial_profiles''',     &
3665               ' or ''bulk_velocity'''
3666          CALL message( 'check_parameters', 'PA0156', 1, 2, 0, 6, 0 )
3667       ENDIF
3668    ENDIF
3669    IF ( ( u_bulk /= 0.0_wp  .OR.  v_bulk /= 0.0_wp )  .AND.                   &
3670         ( .NOT. conserve_volume_flow  .OR.                                    &
3671         TRIM( conserve_volume_flow_mode ) /= 'bulk_velocity' ) )  THEN
3672       WRITE( message_string, * )  'nonzero bulk velocity requires ',          &
3673            'conserve_volume_flow = .T. and ',                                 &
3674            'conserve_volume_flow_mode = ''bulk_velocity'''
3675       CALL message( 'check_parameters', 'PA0157', 1, 2, 0, 6, 0 )
3676    ENDIF
3677
3678!
3679!-- Check particle attributes
3680    IF ( particle_color /= 'none' )  THEN
3681       IF ( particle_color /= 'absuv'  .AND.  particle_color /= 'pt*'  .AND.   &
3682            particle_color /= 'z' )  THEN
3683          message_string = 'illegal value for parameter particle_color: ' //   &
3684                           TRIM( particle_color)
3685          CALL message( 'check_parameters', 'PA0313', 1, 2, 0, 6, 0 )
3686       ELSE
3687          IF ( color_interval(2) <= color_interval(1) )  THEN
3688             message_string = 'color_interval(2) <= color_interval(1)'
3689             CALL message( 'check_parameters', 'PA0315', 1, 2, 0, 6, 0 )
3690          ENDIF
3691       ENDIF
3692    ENDIF
3693
3694    IF ( particle_dvrpsize /= 'none' )  THEN
3695       IF ( particle_dvrpsize /= 'absw' )  THEN
3696          message_string = 'illegal value for parameter particle_dvrpsize:' // &
3697                           ' ' // TRIM( particle_color)
3698          CALL message( 'check_parameters', 'PA0314', 1, 2, 0, 6, 0 )
3699       ELSE
3700          IF ( dvrpsize_interval(2) <= dvrpsize_interval(1) )  THEN
3701             message_string = 'dvrpsize_interval(2) <= dvrpsize_interval(1)'
3702             CALL message( 'check_parameters', 'PA0316', 1, 2, 0, 6, 0 )
3703          ENDIF
3704       ENDIF
3705    ENDIF
3706
3707!
3708!-- Check nudging and large scale forcing from external file
3709    IF ( nudging  .AND.  (  .NOT.  large_scale_forcing ) )  THEN
3710       message_string = 'Nudging requires large_scale_forcing = .T.. &'//      &
3711                        'Surface fluxes and geostrophic wind should be &'//    &
3712                        'prescribed in file LSF_DATA'
3713       CALL message( 'check_parameters', 'PA0374', 1, 2, 0, 6, 0 )
3714    ENDIF
3715
3716    IF ( large_scale_forcing  .AND.  ( bc_lr /= 'cyclic'  .OR.                 &
3717                                    bc_ns /= 'cyclic' ) )  THEN
3718       message_string = 'Non-cyclic lateral boundaries do not allow for &' //  &
3719                        'the usage of large scale forcing from external file.'
3720       CALL message( 'check_parameters', 'PA0375', 1, 2, 0, 6, 0 )
3721    ENDIF
3722
3723    IF ( large_scale_forcing  .AND.  (  .NOT.  humidity ) )  THEN
3724       message_string = 'The usage of large scale forcing from external &'//   & 
3725                        'file LSF_DATA requires humidity = .T..'
3726       CALL message( 'check_parameters', 'PA0376', 1, 2, 0, 6, 0 )
3727    ENDIF
3728
3729    IF ( large_scale_forcing  .AND.  passive_scalar )  THEN
3730       message_string = 'The usage of large scale forcing from external &'//   & 
3731                        'file LSF_DATA is not implemented for passive scalars'
3732       CALL message( 'check_parameters', 'PA0440', 1, 2, 0, 6, 0 )
3733    ENDIF
3734
3735    IF ( large_scale_forcing  .AND.  topography /= 'flat'                      &
3736                              .AND.  .NOT.  urban_surface )  THEN
3737       message_string = 'The usage of large scale forcing from external &'//   & 
3738                        'file LSF_DATA is not implemented for non-flat topography'
3739       CALL message( 'check_parameters', 'PA0377', 1, 2, 0, 6, 0 )
3740    ENDIF
3741
3742    IF ( large_scale_forcing  .AND.  ocean )  THEN
3743       message_string = 'The usage of large scale forcing from external &'//   & 
3744                        'file LSF_DATA is not implemented for ocean runs'
3745       CALL message( 'check_parameters', 'PA0378', 1, 2, 0, 6, 0 )
3746    ENDIF
3747
3748!
3749!-- Prevent empty time records in volume, cross-section and masked data in case
3750!-- of non-parallel netcdf-output in restart runs
3751    IF ( netcdf_data_format < 5 )  THEN
3752       IF ( TRIM( initializing_actions ) == 'read_restart_data' )  THEN
3753          do3d_time_count    = 0
3754          do2d_xy_time_count = 0
3755          do2d_xz_time_count = 0
3756          do2d_yz_time_count = 0
3757          domask_time_count  = 0
3758       ENDIF
3759    ENDIF
3760
3761!
3762!-- Check for valid setting of most_method
3763    IF ( TRIM( most_method ) /= 'circular'  .AND.                              &
3764         TRIM( most_method ) /= 'newton'    .AND.                              &
3765         TRIM( most_method ) /= 'lookup' )  THEN
3766       message_string = 'most_method = "' // TRIM( most_method ) //            &
3767                        '" is unknown'
3768       CALL message( 'check_parameters', 'PA0416', 1, 2, 0, 6, 0 )
3769    ENDIF
3770
3771!
3772!-- Check roughness length, which has to be smaller than dz/2
3773    IF ( ( constant_flux_layer .OR.  &
3774           INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 ) &
3775         .AND. roughness_length > 0.5 * dz )  THEN
3776       message_string = 'roughness_length must be smaller than dz/2'
3777       CALL message( 'check_parameters', 'PA0424', 1, 2, 0, 6, 0 )
3778    ENDIF
3779
3780    CALL location_message( 'finished', .TRUE. )
3781!
3782!-- Check &userpar parameters
3783    CALL user_check_parameters
3784
3785 CONTAINS
3786
3787!------------------------------------------------------------------------------!
3788! Description:
3789! ------------
3790!> Check the length of data output intervals. In case of parallel NetCDF output
3791!> the time levels of the output files need to be fixed. Therefore setting the
3792!> output interval to 0.0s (usually used to output each timestep) is not
3793!> possible as long as a non-fixed timestep is used.
3794!------------------------------------------------------------------------------!
3795
3796    SUBROUTINE check_dt_do( dt_do, dt_do_name )
3797
3798       IMPLICIT NONE
3799
3800       CHARACTER (LEN=*), INTENT (IN) :: dt_do_name !< parin variable name
3801
3802       REAL(wp), INTENT (INOUT)       :: dt_do      !< data output interval
3803
3804       IF ( dt_do == 0.0_wp )  THEN
3805          IF ( dt_fixed )  THEN
3806             WRITE( message_string, '(A,F9.4,A)' )  'Output at every '  //     &
3807                    'timestep is desired (' // dt_do_name // ' = 0.0).&'//     &
3808                    'Setting the output interval to the fixed timestep '//     &
3809                    'dt = ', dt, 's.'
3810             CALL message( 'check_parameters', 'PA0060', 0, 0, 0, 6, 0 )
3811             dt_do = dt
3812          ELSE
3813             message_string = dt_do_name // ' = 0.0 while using a ' //         &
3814                              'variable timestep and parallel netCDF4 ' //     &
3815                              'is not allowed.'
3816             CALL message( 'check_parameters', 'PA0081', 1, 2, 0, 6, 0 )
3817          ENDIF
3818       ENDIF
3819
3820    END SUBROUTINE check_dt_do
3821   
3822   
3823   
3824   
3825!------------------------------------------------------------------------------!
3826! Description:
3827! ------------
3828!> Inititalizes the vertical profiles of scalar quantities.
3829!------------------------------------------------------------------------------!
3830
3831    SUBROUTINE init_vertical_profiles( vertical_gradient_level_ind,            &
3832                                       vertical_gradient_level,                &
3833                                       vertical_gradient,                      &
3834                                       pr_init, surface_value, bc_t_val )
3835
3836
3837       IMPLICIT NONE   
3838   
3839       INTEGER(iwp) ::  i     !< counter
3840       INTEGER(iwp), DIMENSION(1:10) ::  vertical_gradient_level_ind !< vertical grid indices for gradient levels
3841       
3842       REAL(wp)     ::  bc_t_val      !< model top gradient       
3843       REAL(wp)     ::  gradient      !< vertica gradient of the respective quantity
3844       REAL(wp)     ::  surface_value !< surface value of the respecitve quantity
3845
3846       
3847       REAL(wp), DIMENSION(0:nz+1) ::  pr_init                 !< initialisation profile
3848       REAL(wp), DIMENSION(1:10)   ::  vertical_gradient       !< given vertical gradient
3849       REAL(wp), DIMENSION(1:10)   ::  vertical_gradient_level !< given vertical gradient level
3850       
3851       i = 1
3852       gradient = 0.0_wp
3853
3854       IF (  .NOT.  ocean )  THEN
3855
3856          vertical_gradient_level_ind(1) = 0
3857          DO  k = 1, nzt+1
3858             IF ( i < 11 )  THEN
3859                IF ( vertical_gradient_level(i) < zu(k)  .AND.            &
3860                     vertical_gradient_level(i) >= 0.0_wp )  THEN
3861                   gradient = vertical_gradient(i) / 100.0_wp
3862                   vertical_gradient_level_ind(i) = k - 1
3863                   i = i + 1
3864                ENDIF
3865             ENDIF
3866             IF ( gradient /= 0.0_wp )  THEN
3867                IF ( k /= 1 )  THEN
3868                   pr_init(k) = pr_init(k-1) + dzu(k) * gradient
3869                ELSE
3870                   pr_init(k) = pr_init(k-1) + dzu(k) * gradient
3871                ENDIF
3872             ELSE
3873                pr_init(k) = pr_init(k-1)
3874             ENDIF
3875   !
3876   !--       Avoid negative values
3877             IF ( pr_init(k) < 0.0_wp )  THEN
3878                pr_init(k) = 0.0_wp
3879             ENDIF
3880          ENDDO
3881
3882       ELSE
3883
3884          vertical_gradient_level_ind(1) = nzt+1
3885          DO  k = nzt, 0, -1
3886             IF ( i < 11 )  THEN
3887                IF ( vertical_gradient_level(i) > zu(k)  .AND.            &
3888                     vertical_gradient_level(i) <= 0.0_wp )  THEN
3889                   gradient = vertical_gradient(i) / 100.0_wp
3890                   vertical_gradient_level_ind(i) = k + 1
3891                   i = i + 1
3892                ENDIF
3893             ENDIF
3894             IF ( gradient /= 0.0_wp )  THEN
3895                IF ( k /= nzt )  THEN
3896                   pr_init(k) = pr_init(k+1) - dzu(k+1) * gradient
3897                ELSE
3898                   pr_init(k)   = surface_value - 0.5_wp * dzu(k+1) * gradient
3899                   pr_init(k+1) = surface_value + 0.5_wp * dzu(k+1) * gradient
3900                ENDIF
3901             ELSE
3902                pr_init(k) = pr_init(k+1)
3903             ENDIF
3904   !
3905   !--       Avoid negative humidities
3906             IF ( pr_init(k) < 0.0_wp )  THEN
3907                pr_init(k) = 0.0_wp
3908             ENDIF
3909          ENDDO
3910
3911       ENDIF
3912       
3913!
3914!--    In case of no given gradients, choose zero gradient conditions
3915       IF ( vertical_gradient_level(1) == -1.0_wp )  THEN
3916          vertical_gradient_level(1) = 0.0_wp
3917       ENDIF
3918!
3919!--    Store gradient at the top boundary for possible Neumann boundary condition
3920       bc_t_val  = ( pr_init(nzt+1) - pr_init(nzt) ) / dzu(nzt+1)
3921   
3922    END SUBROUTINE init_vertical_profiles
3923   
3924   
3925     
3926!------------------------------------------------------------------------------!
3927! Description:
3928! ------------
3929!> Set the bottom and top boundary conditions for humidity and scalars.
3930!------------------------------------------------------------------------------!
3931
3932    SUBROUTINE set_bc_scalars( sq, bc_b, bc_t, ibc_b, ibc_t, err_nr_b, err_nr_t ) 
3933
3934
3935       IMPLICIT NONE   
3936   
3937       CHARACTER (LEN=1)   ::  sq                       !<
3938       CHARACTER (LEN=*)   ::  bc_b
3939       CHARACTER (LEN=*)   ::  bc_t
3940       CHARACTER (LEN=*)   ::  err_nr_b
3941       CHARACTER (LEN=*)   ::  err_nr_t
3942
3943       INTEGER(iwp)        ::  ibc_b
3944       INTEGER(iwp)        ::  ibc_t
3945
3946!
3947!--    Set Integer flags and check for possilbe errorneous settings for bottom
3948!--    boundary condition
3949       IF ( bc_b == 'dirichlet' )  THEN
3950          ibc_b = 0
3951       ELSEIF ( bc_b == 'neumann' )  THEN
3952          ibc_b = 1
3953       ELSE
3954          message_string = 'unknown boundary condition: bc_' // TRIM( sq ) //  &
3955                           '_b ="' // TRIM( bc_b ) // '"'
3956          CALL message( 'check_parameters', err_nr_b, 1, 2, 0, 6, 0 )
3957       ENDIF
3958!
3959!--    Set Integer flags and check for possilbe errorneous settings for top
3960!--    boundary condition
3961       IF ( bc_t == 'dirichlet' )  THEN
3962          ibc_t = 0
3963       ELSEIF ( bc_t == 'neumann' )  THEN
3964          ibc_t = 1
3965       ELSEIF ( bc_t == 'initial_gradient' )  THEN
3966          ibc_t = 2
3967       ELSEIF ( bc_t == 'nested' )  THEN
3968          ibc_t = 3
3969       ELSE
3970          message_string = 'unknown boundary condition: bc_' // TRIM( sq ) //  &
3971                           '_t ="' // TRIM( bc_t ) // '"'
3972          CALL message( 'check_parameters', err_nr_t, 1, 2, 0, 6, 0 )
3973       ENDIF
3974       
3975   
3976    END SUBROUTINE set_bc_scalars   
3977
3978
3979
3980!------------------------------------------------------------------------------!
3981! Description:
3982! ------------
3983!> Check for consistent settings of bottom boundary conditions for humidity
3984!> and scalars.
3985!------------------------------------------------------------------------------!
3986
3987    SUBROUTINE check_bc_scalars( sq, bc_b, ibc_b,                 &
3988                                 err_nr_1, err_nr_2,       &
3989                                 constant_flux, surface_initial_change )
3990
3991
3992       IMPLICIT NONE   
3993   
3994       CHARACTER (LEN=1)   ::  sq                       !<
3995       CHARACTER (LEN=*)   ::  bc_b
3996       CHARACTER (LEN=*)   ::  err_nr_1
3997       CHARACTER (LEN=*)   ::  err_nr_2
3998       
3999       INTEGER(iwp)        ::  ibc_b
4000       
4001       LOGICAL             ::  constant_flux
4002       
4003       REAL(wp)            ::  surface_initial_change
4004 
4005!
4006!--    A given surface value implies Dirichlet boundary condition for
4007!--    the respective quantity. In this case specification of a constant flux is
4008!--    forbidden. However, an exception is made for large-scale forcing as well
4009!--    as land-surface model.
4010       IF ( .NOT. land_surface  .AND.  .NOT. large_scale_forcing )  THEN
4011          IF ( ibc_b == 0  .AND.  constant_flux )  THEN
4012             message_string = 'boundary condition: bc_' // TRIM( sq ) // '_b ' //  &
4013                              '= "' // TRIM( bc_b ) // '" is not allowed with ' // &
4014                              'prescribed surface flux'
4015             CALL message( 'check_parameters', err_nr_1, 1, 2, 0, 6, 0 )
4016          ENDIF
4017       ENDIF
4018       IF ( constant_flux  .AND.  surface_initial_change /= 0.0_wp )  THEN
4019          WRITE( message_string, * )  'a prescribed surface flux is not allo', &
4020                 'wed with ', sq, '_surface_initial_change (/=0) = ',          &
4021                 surface_initial_change
4022          CALL message( 'check_parameters', err_nr_2, 1, 2, 0, 6, 0 )
4023       ENDIF 
4024       
4025   
4026    END SUBROUTINE check_bc_scalars   
4027   
4028   
4029
4030 END SUBROUTINE check_parameters
Note: See TracBrowser for help on using the repository browser.