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

Last change on this file since 1808 was 1807, checked in by gronemeier, 8 years ago

last commit documented

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