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

Last change on this file since 1796 was 1796, checked in by raasch, 8 years ago

last commit documented

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