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

Last change on this file since 3445 was 3421, checked in by gronemeier, 5 years ago

new surface-data output; renamed output variables (pt to theta, rho_air to rho, rho_ocean to rho_sea_water)

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