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

Last change on this file since 3448 was 3448, checked in by kanani, 5 years ago

Implementation of human thermal indices (from branch biomet_p2 at r3444)

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