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

Last change on this file since 1805 was 1805, checked in by maronga, 5 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 185.6 KB
Line 
1!> @file check_parameters.f90
2!--------------------------------------------------------------------------------!
3! This file is part of PALM.
4!
5! PALM is free software: you can redistribute it and/or modify it under the 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 1805 2016-04-05 16:32:34Z maronga $
26!
27! 1804 2016-04-05 16:30:18Z maronga
28! Removed code for parameter file check (__check)
29!
30! 1795 2016-03-18 15:00:53Z raasch
31! Bugfix: setting of initial scalar profile in ocean
32!
33! 1788 2016-03-10 11:01:04Z maronga
34! Added check for use of most_method = 'lookup' in combination with water
35! surface presribed in the land surface model. Added output of z0q.
36! Syntax layout improved.
37!
38! 1786 2016-03-08 05:49:27Z raasch
39! cpp-direktives for spectra removed, check of spectra level removed
40!
41! 1783 2016-03-06 18:36:17Z raasch
42! netcdf variables and module name changed,
43! check of netcdf precision removed (is done in the netcdf module)
44!
45! 1764 2016-02-28 12:45:19Z raasch
46! output of nest id in run description header,
47! bugfix: check of validity of lateral boundary conditions moved to parin
48!
49! 1762 2016-02-25 12:31:13Z hellstea
50! Introduction of nested domain feature
51!
52! 1745 2016-02-05 13:06:51Z gronemeier
53! Bugfix: check data output intervals to be /= 0.0 in case of parallel NetCDF4
54!
55! 1701 2015-11-02 07:43:04Z maronga
56! Bugfix: definition of rad_net timeseries was missing
57!
58! 1691 2015-10-26 16:17:44Z maronga
59! Added output of Obukhov length (ol) and radiative heating rates for RRTMG.
60! Added checks for use of radiation / lsm with topography.
61!
62! 1682 2015-10-07 23:56:08Z knoop
63! Code annotations made doxygen readable
64!
65! 1606 2015-06-29 10:43:37Z maronga
66! Added check for use of RRTMG without netCDF.
67!
68! 1587 2015-05-04 14:19:01Z maronga
69! Added check for set of albedos when using RRTMG
70!
71! 1585 2015-04-30 07:05:52Z maronga
72! Added support for RRTMG
73!
74! 1575 2015-03-27 09:56:27Z raasch
75! multigrid_fast added as allowed pressure solver
76!
77! 1573 2015-03-23 17:53:37Z suehring
78! Bugfix: check for advection schemes in case of non-cyclic boundary conditions
79!
80! 1557 2015-03-05 16:43:04Z suehring
81! Added checks for monotonic limiter
82!
83! 1555 2015-03-04 17:44:27Z maronga
84! Added output of r_a and r_s. Renumbering of LSM PA-messages.
85!
86! 1553 2015-03-03 17:33:54Z maronga
87! Removed check for missing soil temperature profile as default values were added.
88!
89! 1551 2015-03-03 14:18:16Z maronga
90! Added various checks for land surface and radiation model. In the course of this
91! action, the length of the variable var has to be increased
92!
93! 1504 2014-12-04 09:23:49Z maronga
94! Bugfix: check for large_scale forcing got lost.
95!
96! 1500 2014-12-03 17:42:41Z maronga
97! Boundary conditions changed to dirichlet for land surface model
98!
99! 1496 2014-12-02 17:25:50Z maronga
100! Added checks for the land surface model
101!
102! 1484 2014-10-21 10:53:05Z kanani
103! Changes due to new module structure of the plant canopy model:
104!   module plant_canopy_model_mod added,
105!   checks regarding usage of new method for leaf area density profile
106!   construction added,
107!   lad-profile construction moved to new subroutine init_plant_canopy within
108!   the module plant_canopy_model_mod,
109!   drag_coefficient renamed to canopy_drag_coeff.
110! Missing KIND-attribute for REAL constant added
111!
112! 1455 2014-08-29 10:47:47Z heinze
113! empty time records in volume, cross-section and masked data output prevented 
114! in case of non-parallel netcdf-output in restart runs
115!
116! 1429 2014-07-15 12:53:45Z knoop
117! run_description_header exended to provide ensemble_member_nr if specified
118!
119! 1425 2014-07-05 10:57:53Z knoop
120! bugfix: perturbation domain modified for parallel random number generator
121!
122! 1402 2014-05-09 14:25:13Z raasch
123! location messages modified
124!
125! 1400 2014-05-09 14:03:54Z knoop
126! Check random generator extended by option random-parallel
127!
128! 1384 2014-05-02 14:31:06Z raasch
129! location messages added
130!
131! 1365 2014-04-22 15:03:56Z boeske
132! Usage of large scale forcing for pt and q enabled:
133! output for profiles of large scale advection (td_lsa_lpt, td_lsa_q),
134! large scale subsidence (td_sub_lpt, td_sub_q)
135! and nudging tendencies (td_nud_lpt, td_nud_q, td_nud_u and td_nud_v) added,
136! check if use_subsidence_tendencies is used correctly
137!
138! 1361 2014-04-16 15:17:48Z hoffmann
139! PA0363 removed
140! PA0362 changed
141!
142! 1359 2014-04-11 17:15:14Z hoffmann
143! Do not allow the execution of PALM with use_particle_tails, since particle
144! tails are currently not supported by our new particle structure.
145!
146! PA0084 not necessary for new particle structure
147!
148! 1353 2014-04-08 15:21:23Z heinze
149! REAL constants provided with KIND-attribute
150!
151! 1330 2014-03-24 17:29:32Z suehring
152! In case of SGS-particle velocity advection of TKE is also allowed with
153! dissipative 5th-order scheme.
154!
155! 1327 2014-03-21 11:00:16Z raasch
156! "baroclinicity" renamed "baroclinity", "ocean version" replaced by "ocean mode"
157! bugfix: duplicate error message 56 removed,
158! check of data_output_format and do3d_compress removed
159!
160! 1322 2014-03-20 16:38:49Z raasch
161! some REAL constants defined as wp-kind
162!
163! 1320 2014-03-20 08:40:49Z raasch
164! Kind-parameters added to all INTEGER and REAL declaration statements,
165! kinds are defined in new module kinds,
166! revision history before 2012 removed,
167! comment fields (!:) to be used for variable explanations added to
168! all variable declaration statements
169!
170! 1308 2014-03-13 14:58:42Z fricke
171! +netcdf_data_format_save
172! Calculate fixed number of output time levels for parallel netcdf output.
173! For masked data, parallel netcdf output is not tested so far, hence
174! netcdf_data_format is switched back to non-paralell output.
175!
176! 1299 2014-03-06 13:15:21Z heinze
177! enable usage of large_scale subsidence in combination with large_scale_forcing
178! output for profile of large scale vertical velocity w_subs added
179!
180! 1276 2014-01-15 13:40:41Z heinze
181! Use LSF_DATA also in case of Dirichlet bottom boundary condition for scalars
182!
183! 1241 2013-10-30 11:36:58Z heinze
184! output for profiles of ug and vg added
185! set surface_heatflux and surface_waterflux also in dependence on
186! large_scale_forcing
187! checks for nudging and large scale forcing from external file
188!
189! 1236 2013-09-27 07:21:13Z raasch
190! check number of spectra levels
191!
192! 1216 2013-08-26 09:31:42Z raasch
193! check for transpose_compute_overlap (temporary)
194!
195! 1214 2013-08-21 12:29:17Z kanani
196! additional check for simultaneous use of vertical grid stretching
197! and particle advection
198!
199! 1212 2013-08-15 08:46:27Z raasch
200! checks for poisfft_hybrid removed
201!
202! 1210 2013-08-14 10:58:20Z raasch
203! check for fftw
204!
205! 1179 2013-06-14 05:57:58Z raasch
206! checks and settings of buoyancy parameters and switches revised,
207! initial profile for rho added to hom (id=77)
208!
209! 1174 2013-05-31 10:28:08Z gryschka
210! Bugfix in computing initial profiles for ug, vg, lad, q in case of Atmosphere
211!
212! 1159 2013-05-21 11:58:22Z fricke
213! bc_lr/ns_dirneu/neudir removed
214!
215! 1115 2013-03-26 18:16:16Z hoffmann
216! unused variables removed
217! drizzle can be used without precipitation
218!
219! 1111 2013-03-08 23:54:10Z raasch
220! ibc_p_b = 2 removed
221!
222! 1103 2013-02-20 02:15:53Z raasch
223! Bugfix: turbulent inflow must not require cyclic fill in restart runs
224!
225! 1092 2013-02-02 11:24:22Z raasch
226! unused variables removed
227!
228! 1069 2012-11-28 16:18:43Z maronga
229! allow usage of topography in combination with cloud physics
230!
231! 1065 2012-11-22 17:42:36Z hoffmann
232! Bugfix: It is not allowed to use cloud_scheme = seifert_beheng without
233!         precipitation in order to save computational resources.
234!
235! 1060 2012-11-21 07:19:51Z raasch
236! additional check for parameter turbulent_inflow
237!
238! 1053 2012-11-13 17:11:03Z hoffmann
239! necessary changes for the new two-moment cloud physics scheme added:
240! - check cloud physics scheme (Kessler or Seifert and Beheng)
241! - plant_canopy is not allowed
242! - currently, only cache loop_optimization is allowed
243! - initial profiles of nr, qr
244! - boundary condition of nr, qr
245! - check output quantities (qr, nr, prr)
246!
247! 1036 2012-10-22 13:43:42Z raasch
248! code put under GPL (PALM 3.9)
249!
250! 1031/1034 2012-10-22 11:32:49Z raasch
251! check of netcdf4 parallel file support
252!
253! 1019 2012-09-28 06:46:45Z raasch
254! non-optimized version of prognostic_equations not allowed any more
255!
256! 1015 2012-09-27 09:23:24Z raasch
257! acc allowed for loop optimization,
258! checks for adjustment of mixing length to the Prandtl mixing length removed
259!
260! 1003 2012-09-14 14:35:53Z raasch
261! checks for cases with unequal subdomain sizes removed
262!
263! 1001 2012-09-13 14:08:46Z raasch
264! all actions concerning leapfrog- and upstream-spline-scheme removed
265!
266! 996 2012-09-07 10:41:47Z raasch
267! little reformatting
268
269! 978 2012-08-09 08:28:32Z fricke
270! setting of bc_lr/ns_dirneu/neudir
271! outflow damping layer removed
272! check for z0h*
273! check for pt_damping_width
274!
275! 964 2012-07-26 09:14:24Z raasch
276! check of old profil-parameters removed
277!
278! 940 2012-07-09 14:31:00Z raasch
279! checks for parameter neutral
280!
281! 924 2012-06-06 07:44:41Z maronga
282! Bugfix: preprocessor directives caused error during compilation
283!
284! 892 2012-05-02 13:51:44Z maronga
285! Bugfix for parameter file check ( excluding __netcdf4 )
286!
287! 866 2012-03-28 06:44:41Z raasch
288! use only 60% of the geostrophic wind as translation speed in case of Galilean
289! transformation and use_ug_for_galilei_tr = .T. in order to mimimize the
290! timestep
291!
292! 861 2012-03-26 14:18:34Z suehring
293! Check for topography and ws-scheme removed.
294! Check for loop_optimization = 'vector' and ws-scheme removed.
295!
296! 845 2012-03-07 10:23:05Z maronga
297! Bugfix: exclude __netcdf4 directive part from namelist file check compilation
298!
299! 828 2012-02-21 12:00:36Z raasch
300! check of collision_kernel extended
301!
302! 825 2012-02-19 03:03:44Z raasch
303! check for collision_kernel and curvature_solution_effects
304!
305! 809 2012-01-30 13:32:58Z maronga
306! Bugfix: replaced .AND. and .NOT. with && and ! in the preprocessor directives
307!
308! 807 2012-01-25 11:53:51Z maronga
309! New cpp directive "__check" implemented which is used by check_namelist_files
310!
311! Revision 1.1  1997/08/26 06:29:23  raasch
312! Initial revision
313!
314!
315! Description:
316! ------------
317!> Check control parameters and deduce further quantities.
318!------------------------------------------------------------------------------!
319 SUBROUTINE check_parameters
320 
321
322    USE arrays_3d
323    USE cloud_parameters
324    USE constants
325    USE control_parameters
326    USE dvrp_variables
327    USE grid_variables
328    USE indices
329    USE land_surface_model_mod
330    USE kinds
331    USE model_1d
332    USE netcdf_interface,                                                      &
333        ONLY:  dopr_unit, do2d_unit, do3d_unit, netcdf_data_format,            &
334               netcdf_data_format_string
335    USE particle_attributes
336    USE pegrid
337    USE plant_canopy_model_mod
338    USE pmc_interface,                                                         &
339        ONLY:  cpl_id, nested_run
340    USE profil_parameter
341    USE radiation_model_mod
342    USE statistics
343    USE subsidence_mod
344    USE statistics
345    USE transpose_indices
346
347
348    IMPLICIT NONE
349
350    CHARACTER (LEN=1)   ::  sq                       !<
351    CHARACTER (LEN=15)  ::  var                      !<
352    CHARACTER (LEN=7)   ::  unit                     !<
353    CHARACTER (LEN=8)   ::  date                     !<
354    CHARACTER (LEN=10)  ::  time                     !<
355    CHARACTER (LEN=10)  ::  ensemble_string          !<
356    CHARACTER (LEN=15)  ::  nest_string              !<
357    CHARACTER (LEN=40)  ::  coupling_string          !<
358    CHARACTER (LEN=100) ::  action                   !<
359
360    INTEGER(iwp) ::  i                               !<
361    INTEGER(iwp) ::  ilen                            !<
362    INTEGER(iwp) ::  iremote = 0                     !<
363    INTEGER(iwp) ::  j                               !<
364    INTEGER(iwp) ::  k                               !<
365    INTEGER(iwp) ::  kk                              !<
366    INTEGER(iwp) ::  netcdf_data_format_save         !<
367    INTEGER(iwp) ::  position                        !<
368    INTEGER(iwp) ::  prec                            !<
369   
370    LOGICAL     ::  found                            !<
371   
372    REAL(wp)    ::  gradient                         !<
373    REAL(wp)    ::  remote = 0.0_wp                  !<
374    REAL(wp)    ::  simulation_time_since_reference  !<
375
376
377    CALL location_message( 'checking parameters', .FALSE. )
378
379!
380!-- Check for overlap combinations, which are not realized yet
381    IF ( transpose_compute_overlap )  THEN
382       IF ( numprocs == 1 )  STOP '+++ transpose-compute-overlap not implemented for single PE runs'
383#if defined( __openacc )
384       STOP '+++ transpose-compute-overlap not implemented for GPU usage'
385#endif
386    ENDIF
387
388!
389!-- Warning, if host is not set
390    IF ( host(1:1) == ' ' )  THEN
391       message_string = '"host" is not set. Please check that environment ' // &
392                        'variable "localhost" & is set before running PALM'
393       CALL message( 'check_parameters', 'PA0001', 0, 0, 0, 6, 0 )
394    ENDIF
395
396!
397!-- Check the coupling mode
398    IF ( coupling_mode /= 'uncoupled'            .AND.  &
399         coupling_mode /= 'atmosphere_to_ocean'  .AND.  &
400         coupling_mode /= 'ocean_to_atmosphere' )  THEN
401       message_string = 'illegal coupling mode: ' // TRIM( coupling_mode )
402       CALL message( 'check_parameters', 'PA0002', 1, 2, 0, 6, 0 )
403    ENDIF
404
405!
406!-- Check dt_coupling, restart_time, dt_restart, end_time, dx, dy, nx and ny
407    IF ( coupling_mode /= 'uncoupled')  THEN
408
409       IF ( dt_coupling == 9999999.9_wp )  THEN
410          message_string = 'dt_coupling is not set but required for coup' //   &
411                           'ling mode "' //  TRIM( coupling_mode ) // '"'
412          CALL message( 'check_parameters', 'PA0003', 1, 2, 0, 6, 0 )
413       ENDIF
414
415#if defined( __parallel )
416
417!
418!--    NOTE: coupled runs have not been implemented in the check_namelist_files
419!--    program.
420!--    check_namelist_files will need the following information of the other
421!--    model (atmosphere/ocean).
422!       dt_coupling = remote
423!       dt_max = remote
424!       restart_time = remote
425!       dt_restart= remote
426!       simulation_time_since_reference = remote
427!       dx = remote
428
429
430       IF ( myid == 0 ) THEN
431          CALL MPI_SEND( dt_coupling, 1, MPI_REAL, target_id, 11, comm_inter,  &
432                         ierr )
433          CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 11, comm_inter,       &
434                         status, ierr )
435       ENDIF
436       CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
437 
438       IF ( dt_coupling /= remote )  THEN
439          WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
440                 '": dt_coupling = ', dt_coupling, '& is not equal to ',       &
441                 'dt_coupling_remote = ', remote
442          CALL message( 'check_parameters', 'PA0004', 1, 2, 0, 6, 0 )
443       ENDIF
444       IF ( dt_coupling <= 0.0_wp )  THEN
445
446          IF ( myid == 0  ) THEN
447             CALL MPI_SEND( dt_max, 1, MPI_REAL, target_id, 19, comm_inter, ierr )
448             CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 19, comm_inter,    &
449                            status, ierr )
450          ENDIF   
451          CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
452     
453          dt_coupling = MAX( dt_max, remote )
454          WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
455                 '": dt_coupling <= 0.0 & is not allowed and is reset to ',    &
456                 'MAX(dt_max(A,O)) = ', dt_coupling
457          CALL message( 'check_parameters', 'PA0005', 0, 1, 0, 6, 0 )
458       ENDIF
459
460       IF ( myid == 0 ) THEN
461          CALL MPI_SEND( restart_time, 1, MPI_REAL, target_id, 12, comm_inter, &
462                         ierr )
463          CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 12, comm_inter,       &
464                         status, ierr )
465       ENDIF
466       CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
467 
468       IF ( restart_time /= remote )  THEN
469          WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
470                 '": restart_time = ', restart_time, '& is not equal to ',     &
471                 'restart_time_remote = ', remote
472          CALL message( 'check_parameters', 'PA0006', 1, 2, 0, 6, 0 )
473       ENDIF
474
475       IF ( myid == 0 ) THEN
476          CALL MPI_SEND( dt_restart, 1, MPI_REAL, target_id, 13, comm_inter,   &
477                         ierr )
478          CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 13, comm_inter,       &
479                         status, ierr )
480       ENDIF   
481       CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
482   
483       IF ( dt_restart /= remote )  THEN
484          WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
485                 '": dt_restart = ', dt_restart, '& is not equal to ',         &
486                 'dt_restart_remote = ', remote
487          CALL message( 'check_parameters', 'PA0007', 1, 2, 0, 6, 0 )
488       ENDIF
489
490       simulation_time_since_reference = end_time - coupling_start_time
491
492       IF  ( myid == 0 ) THEN
493          CALL MPI_SEND( simulation_time_since_reference, 1, MPI_REAL, target_id, &
494                         14, comm_inter, ierr )
495          CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 14, comm_inter,       &
496                         status, ierr )   
497       ENDIF
498       CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
499   
500       IF ( simulation_time_since_reference /= remote )  THEN
501          WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
502                 '": simulation_time_since_reference = ',                      &
503                 simulation_time_since_reference, '& is not equal to ',        &
504                 'simulation_time_since_reference_remote = ', remote
505          CALL message( 'check_parameters', 'PA0008', 1, 2, 0, 6, 0 )
506       ENDIF
507
508       IF ( myid == 0 ) THEN
509          CALL MPI_SEND( dx, 1, MPI_REAL, target_id, 15, comm_inter, ierr )
510          CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 15, comm_inter,       &
511                                                             status, ierr )
512       ENDIF
513       CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
514
515
516       IF ( coupling_mode == 'atmosphere_to_ocean') THEN
517
518          IF ( dx < remote ) THEN
519             WRITE( message_string, * ) 'coupling mode "',                     &
520                   TRIM( coupling_mode ),                                      &
521           '": dx in Atmosphere is not equal to or not larger then dx in ocean'
522             CALL message( 'check_parameters', 'PA0009', 1, 2, 0, 6, 0 )
523          ENDIF
524
525          IF ( (nx_a+1)*dx /= (nx_o+1)*remote )  THEN
526             WRITE( message_string, * ) 'coupling mode "',                     &
527                    TRIM( coupling_mode ),                                     &
528             '": Domain size in x-direction is not equal in ocean and atmosphere'
529             CALL message( 'check_parameters', 'PA0010', 1, 2, 0, 6, 0 )
530          ENDIF
531
532       ENDIF
533
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
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 )
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 ( netcdf_data_format > 2 )  THEN
3859#if defined( __netcdf4 )
3860       CONTINUE
3861#else
3862       message_string = 'netCDF: netCDF4 format requested but no ' //          &
3863                        'cpp-directive __netcdf4 given & switch '  //          &
3864                        'back to 64-bit offset format'
3865       CALL message( 'check_parameters', 'PA0171', 0, 1, 0, 6, 0 )
3866       netcdf_data_format = 2
3867#endif
3868    ENDIF
3869    IF ( netcdf_data_format > 4 )  THEN
3870#if defined( __netcdf4 ) && defined( __netcdf4_parallel )
3871       CONTINUE
3872#else
3873       message_string = 'netCDF: netCDF4 parallel output requested but no ' // &
3874                        'cpp-directive __netcdf4_parallel given & switch '  // &
3875                        'back to netCDF4 non-parallel output'
3876       CALL message( 'check_parameters', 'PA0099', 0, 1, 0, 6, 0 )
3877       netcdf_data_format = netcdf_data_format - 2
3878#endif
3879    ENDIF
3880
3881!
3882!-- Calculate fixed number of output time levels for parallel netcdf output.
3883!-- The time dimension has to be defined as limited for parallel output,
3884!-- because otherwise the I/O performance drops significantly.
3885    IF ( netcdf_data_format > 4 )  THEN
3886
3887!
3888!--    Check if any of the follwoing data output interval is 0.0s, which is
3889!--    not allowed for parallel output.
3890       CALL check_dt_do( dt_do3d,    'dt_do3d'    )
3891       CALL check_dt_do( dt_do2d_xy, 'dt_do2d_xy' )
3892       CALL check_dt_do( dt_do2d_xz, 'dt_do2d_xz' )
3893       CALL check_dt_do( dt_do2d_yz, 'dt_do2d_yz' )
3894
3895       ntdim_3d(0)    = INT( ( end_time - skip_time_do3d ) / dt_do3d )
3896       IF ( do3d_at_begin ) ntdim_3d(0) = ntdim_3d(0) + 1
3897       ntdim_3d(1)    = INT( ( end_time - skip_time_data_output_av )           &
3898                             / dt_data_output_av )
3899       ntdim_2d_xy(0) = INT( ( end_time - skip_time_do2d_xy ) / dt_do2d_xy )
3900       ntdim_2d_xz(0) = INT( ( end_time - skip_time_do2d_xz ) / dt_do2d_xz )
3901       ntdim_2d_yz(0) = INT( ( end_time - skip_time_do2d_yz ) / dt_do2d_yz )
3902       IF ( do2d_at_begin )  THEN
3903          ntdim_2d_xy(0) = ntdim_2d_xy(0) + 1
3904          ntdim_2d_xz(0) = ntdim_2d_xz(0) + 1
3905          ntdim_2d_yz(0) = ntdim_2d_yz(0) + 1
3906       ENDIF
3907       ntdim_2d_xy(1) = ntdim_3d(1)
3908       ntdim_2d_xz(1) = ntdim_3d(1)
3909       ntdim_2d_yz(1) = ntdim_3d(1)
3910
3911    ENDIF
3912
3913!
3914!-- Check, whether a constant diffusion coefficient shall be used
3915    IF ( km_constant /= -1.0_wp )  THEN
3916       IF ( km_constant < 0.0_wp )  THEN
3917          WRITE( message_string, * )  'km_constant = ', km_constant, ' < 0.0'
3918          CALL message( 'check_parameters', 'PA0121', 1, 2, 0, 6, 0 )
3919       ELSE
3920          IF ( prandtl_number < 0.0_wp )  THEN
3921             WRITE( message_string, * )  'prandtl_number = ', prandtl_number,  &
3922                                         ' < 0.0'
3923             CALL message( 'check_parameters', 'PA0122', 1, 2, 0, 6, 0 )
3924          ENDIF
3925          constant_diffusion = .TRUE.
3926
3927          IF ( constant_flux_layer )  THEN
3928             message_string = 'constant_flux_layer is not allowed with fixed ' &
3929                              // 'value of km'
3930             CALL message( 'check_parameters', 'PA0123', 1, 2, 0, 6, 0 )
3931          ENDIF
3932       ENDIF
3933    ENDIF
3934
3935!
3936!-- In case of non-cyclic lateral boundaries and a damping layer for the
3937!-- potential temperature, check the width of the damping layer
3938    IF ( bc_lr /= 'cyclic' ) THEN
3939       IF ( pt_damping_width < 0.0_wp  .OR.                                    &
3940            pt_damping_width > REAL( nx * dx ) )  THEN
3941          message_string = 'pt_damping_width out of range'
3942          CALL message( 'check_parameters', 'PA0124', 1, 2, 0, 6, 0 )
3943       ENDIF
3944    ENDIF
3945
3946    IF ( bc_ns /= 'cyclic' )  THEN
3947       IF ( pt_damping_width < 0.0_wp  .OR.                                    &
3948            pt_damping_width > REAL( ny * dy ) )  THEN
3949          message_string = 'pt_damping_width out of range'
3950          CALL message( 'check_parameters', 'PA0124', 1, 2, 0, 6, 0 )
3951       ENDIF
3952    ENDIF
3953
3954!
3955!-- Check value range for zeta = z/L
3956    IF ( zeta_min >= zeta_max )  THEN
3957       WRITE( message_string, * )  'zeta_min = ', zeta_min, ' must be less ',  &
3958                                   'than zeta_max = ', zeta_max
3959       CALL message( 'check_parameters', 'PA0125', 1, 2, 0, 6, 0 )
3960    ENDIF
3961
3962!
3963!-- Check random generator
3964    IF ( (random_generator /= 'system-specific'      .AND.                     &
3965          random_generator /= 'random-parallel'   )  .AND.                     &
3966          random_generator /= 'numerical-recipes' )  THEN
3967       message_string = 'unknown random generator: random_generator = "' //    &
3968                        TRIM( random_generator ) // '"'
3969       CALL message( 'check_parameters', 'PA0135', 1, 2, 0, 6, 0 )
3970    ENDIF
3971
3972!
3973!-- Determine upper and lower hight level indices for random perturbations
3974    IF ( disturbance_level_b == -9999999.9_wp )  THEN
3975       IF ( ocean )  THEN
3976          disturbance_level_b     = zu((nzt*2)/3)
3977          disturbance_level_ind_b = ( nzt * 2 ) / 3
3978       ELSE
3979          disturbance_level_b     = zu(nzb+3)
3980          disturbance_level_ind_b = nzb + 3
3981       ENDIF
3982    ELSEIF ( disturbance_level_b < zu(3) )  THEN
3983       WRITE( message_string, * )  'disturbance_level_b = ',                   &
3984                           disturbance_level_b, ' must be >= ', zu(3), '(zu(3))'
3985       CALL message( 'check_parameters', 'PA0126', 1, 2, 0, 6, 0 )
3986    ELSEIF ( disturbance_level_b > zu(nzt-2) )  THEN
3987       WRITE( message_string, * )  'disturbance_level_b = ',                   &
3988                   disturbance_level_b, ' must be <= ', zu(nzt-2), '(zu(nzt-2))'
3989       CALL message( 'check_parameters', 'PA0127', 1, 2, 0, 6, 0 )
3990    ELSE
3991       DO  k = 3, nzt-2
3992          IF ( disturbance_level_b <= zu(k) )  THEN
3993             disturbance_level_ind_b = k
3994             EXIT
3995          ENDIF
3996       ENDDO
3997    ENDIF
3998
3999    IF ( disturbance_level_t == -9999999.9_wp )  THEN
4000       IF ( ocean )  THEN
4001          disturbance_level_t     = zu(nzt-3)
4002          disturbance_level_ind_t = nzt - 3
4003       ELSE
4004          disturbance_level_t     = zu(nzt/3)
4005          disturbance_level_ind_t = nzt / 3
4006       ENDIF
4007    ELSEIF ( disturbance_level_t > zu(nzt-2) )  THEN
4008       WRITE( message_string, * )  'disturbance_level_t = ',                   &
4009                   disturbance_level_t, ' must be <= ', zu(nzt-2), '(zu(nzt-2))'
4010       CALL message( 'check_parameters', 'PA0128', 1, 2, 0, 6, 0 )
4011    ELSEIF ( disturbance_level_t < disturbance_level_b )  THEN
4012       WRITE( message_string, * )  'disturbance_level_t = ',                   &
4013                   disturbance_level_t, ' must be >= disturbance_level_b = ',  &
4014                   disturbance_level_b
4015       CALL message( 'check_parameters', 'PA0129', 1, 2, 0, 6, 0 )
4016    ELSE
4017       DO  k = 3, nzt-2
4018          IF ( disturbance_level_t <= zu(k) )  THEN
4019             disturbance_level_ind_t = k
4020             EXIT
4021          ENDIF
4022       ENDDO
4023    ENDIF
4024
4025!
4026!-- Check again whether the levels determined this way are ok.
4027!-- Error may occur at automatic determination and too few grid points in
4028!-- z-direction.
4029    IF ( disturbance_level_ind_t < disturbance_level_ind_b )  THEN
4030       WRITE( message_string, * )  'disturbance_level_ind_t = ',               &
4031                disturbance_level_ind_t, ' must be >= disturbance_level_b = ', &
4032                disturbance_level_b
4033       CALL message( 'check_parameters', 'PA0130', 1, 2, 0, 6, 0 )
4034    ENDIF
4035
4036!
4037!-- Determine the horizontal index range for random perturbations.
4038!-- In case of non-cyclic horizontal boundaries, no perturbations are imposed
4039!-- near the inflow and the perturbation area is further limited to ...(1)
4040!-- after the initial phase of the flow.
4041   
4042    IF ( bc_lr /= 'cyclic' )  THEN
4043       IF ( inflow_disturbance_begin == -1 )  THEN
4044          inflow_disturbance_begin = MIN( 10, nx/2 )
4045       ENDIF
4046       IF ( inflow_disturbance_begin < 0  .OR.  inflow_disturbance_begin > nx )&
4047       THEN
4048          message_string = 'inflow_disturbance_begin out of range'
4049          CALL message( 'check_parameters', 'PA0131', 1, 2, 0, 6, 0 )
4050       ENDIF
4051       IF ( inflow_disturbance_end == -1 )  THEN
4052          inflow_disturbance_end = MIN( 100, 3*nx/4 )
4053       ENDIF
4054       IF ( inflow_disturbance_end < 0  .OR.  inflow_disturbance_end > nx )    &
4055       THEN
4056          message_string = 'inflow_disturbance_end out of range'
4057          CALL message( 'check_parameters', 'PA0132', 1, 2, 0, 6, 0 )
4058       ENDIF
4059    ELSEIF ( bc_ns /= 'cyclic' )  THEN
4060       IF ( inflow_disturbance_begin == -1 )  THEN
4061          inflow_disturbance_begin = MIN( 10, ny/2 )
4062       ENDIF
4063       IF ( inflow_disturbance_begin < 0  .OR.  inflow_disturbance_begin > ny )&
4064       THEN
4065          message_string = 'inflow_disturbance_begin out of range'
4066          CALL message( 'check_parameters', 'PA0131', 1, 2, 0, 6, 0 )
4067       ENDIF
4068       IF ( inflow_disturbance_end == -1 )  THEN
4069          inflow_disturbance_end = MIN( 100, 3*ny/4 )
4070       ENDIF
4071       IF ( inflow_disturbance_end < 0  .OR.  inflow_disturbance_end > ny )    &
4072       THEN
4073          message_string = 'inflow_disturbance_end out of range'
4074          CALL message( 'check_parameters', 'PA0132', 1, 2, 0, 6, 0 )
4075       ENDIF
4076    ENDIF
4077
4078    IF ( random_generator == 'random-parallel' )  THEN
4079       dist_nxl = nxl;  dist_nxr = nxr
4080       dist_nys = nys;  dist_nyn = nyn
4081       IF ( bc_lr == 'radiation/dirichlet' )  THEN
4082          dist_nxr    = MIN( nx - inflow_disturbance_begin, nxr )
4083          dist_nxl(1) = MAX( nx - inflow_disturbance_end, nxl )
4084       ELSEIF ( bc_lr == 'dirichlet/radiation' )  THEN
4085          dist_nxl    = MAX( inflow_disturbance_begin, nxl )
4086          dist_nxr(1) = MIN( inflow_disturbance_end, nxr )
4087       ENDIF
4088       IF ( bc_ns == 'dirichlet/radiation' )  THEN
4089          dist_nyn    = MIN( ny - inflow_disturbance_begin, nyn )
4090          dist_nys(1) = MAX( ny - inflow_disturbance_end, nys )
4091       ELSEIF ( bc_ns == 'radiation/dirichlet' )  THEN
4092          dist_nys    = MAX( inflow_disturbance_begin, nys )
4093          dist_nyn(1) = MIN( inflow_disturbance_end, nyn )
4094       ENDIF
4095    ELSE
4096       dist_nxl = 0;  dist_nxr = nx
4097       dist_nys = 0;  dist_nyn = ny
4098       IF ( bc_lr == 'radiation/dirichlet' )  THEN
4099          dist_nxr    = nx - inflow_disturbance_begin
4100          dist_nxl(1) = nx - inflow_disturbance_end
4101       ELSEIF ( bc_lr == 'dirichlet/radiation' )  THEN
4102          dist_nxl    = inflow_disturbance_begin
4103          dist_nxr(1) = inflow_disturbance_end
4104       ENDIF
4105       IF ( bc_ns == 'dirichlet/radiation' )  THEN
4106          dist_nyn    = ny - inflow_disturbance_begin
4107          dist_nys(1) = ny - inflow_disturbance_end
4108       ELSEIF ( bc_ns == 'radiation/dirichlet' )  THEN
4109          dist_nys    = inflow_disturbance_begin
4110          dist_nyn(1) = inflow_disturbance_end
4111       ENDIF
4112    ENDIF
4113
4114!
4115!-- A turbulent inflow requires Dirichlet conditions at the respective inflow
4116!-- boundary (so far, a turbulent inflow is realized from the left side only)
4117    IF ( turbulent_inflow  .AND.  bc_lr /= 'dirichlet/radiation' )  THEN
4118       message_string = 'turbulent_inflow = .T. requires a Dirichlet ' //      &
4119                        'condition at the inflow boundary'
4120       CALL message( 'check_parameters', 'PA0133', 1, 2, 0, 6, 0 )
4121    ENDIF
4122
4123!
4124!-- Turbulent inflow requires that 3d arrays have been cyclically filled with
4125!-- data from prerun in the first main run
4126    IF ( turbulent_inflow  .AND.  initializing_actions /= 'cyclic_fill'  .AND. &
4127         initializing_actions /= 'read_restart_data' )  THEN
4128       message_string = 'turbulent_inflow = .T. requires ' //                  &
4129                        'initializing_actions = ''cyclic_fill'' '
4130       CALL message( 'check_parameters', 'PA0055', 1, 2, 0, 6, 0 )
4131    ENDIF
4132
4133!
4134!-- In case of turbulent inflow calculate the index of the recycling plane
4135    IF ( turbulent_inflow )  THEN
4136       IF ( recycling_width == 9999999.9_wp )  THEN
4137!
4138!--       Set the default value for the width of the recycling domain
4139          recycling_width = 0.1_wp * nx * dx
4140       ELSE
4141          IF ( recycling_width < dx  .OR.  recycling_width > nx * dx )  THEN
4142             WRITE( message_string, * )  'illegal value for recycling_width:', &
4143                                         ' ', recycling_width
4144             CALL message( 'check_parameters', 'PA0134', 1, 2, 0, 6, 0 )
4145          ENDIF
4146       ENDIF
4147!
4148!--    Calculate the index
4149       recycling_plane = recycling_width / dx
4150    ENDIF
4151
4152!
4153!-- Determine damping level index for 1D model
4154    IF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
4155       IF ( damp_level_1d == -1.0_wp )  THEN
4156          damp_level_1d     = zu(nzt+1)
4157          damp_level_ind_1d = nzt + 1
4158       ELSEIF ( damp_level_1d < 0.0_wp  .OR.  damp_level_1d > zu(nzt+1) )  THEN
4159          WRITE( message_string, * )  'damp_level_1d = ', damp_level_1d,       &
4160                 ' must be > 0.0 and < ', zu(nzt+1), '(zu(nzt+1))'
4161          CALL message( 'check_parameters', 'PA0136', 1, 2, 0, 6, 0 )
4162       ELSE
4163          DO  k = 1, nzt+1
4164             IF ( damp_level_1d <= zu(k) )  THEN
4165                damp_level_ind_1d = k
4166                EXIT
4167             ENDIF
4168          ENDDO
4169       ENDIF
4170    ENDIF
4171
4172!
4173!-- Check some other 1d-model parameters
4174    IF ( TRIM( mixing_length_1d ) /= 'as_in_3d_model'  .AND.                   &
4175         TRIM( mixing_length_1d ) /= 'blackadar' )  THEN
4176       message_string = 'mixing_length_1d = "' // TRIM( mixing_length_1d ) //  &
4177                        '" is unknown'
4178       CALL message( 'check_parameters', 'PA0137', 1, 2, 0, 6, 0 )
4179    ENDIF
4180    IF ( TRIM( dissipation_1d ) /= 'as_in_3d_model'  .AND.                     &
4181         TRIM( dissipation_1d ) /= 'detering' )  THEN
4182       message_string = 'dissipation_1d = "' // TRIM( dissipation_1d ) //      &
4183                        '" is unknown'
4184       CALL message( 'check_parameters', 'PA0138', 1, 2, 0, 6, 0 )
4185    ENDIF
4186
4187!
4188!-- Set time for the next user defined restart (time_restart is the
4189!-- internal parameter for steering restart events)
4190    IF ( restart_time /= 9999999.9_wp )  THEN
4191       IF ( restart_time > time_since_reference_point )  THEN
4192          time_restart = restart_time
4193       ENDIF
4194    ELSE
4195!
4196!--    In case of a restart run, set internal parameter to default (no restart)
4197!--    if the NAMELIST-parameter restart_time is omitted
4198       time_restart = 9999999.9_wp
4199    ENDIF
4200
4201!
4202!-- Set default value of the time needed to terminate a model run
4203    IF ( termination_time_needed == -1.0_wp )  THEN
4204       IF ( host(1:3) == 'ibm' )  THEN
4205          termination_time_needed = 300.0_wp
4206       ELSE
4207          termination_time_needed = 35.0_wp
4208       ENDIF
4209    ENDIF
4210
4211!
4212!-- Check the time needed to terminate a model run
4213    IF ( host(1:3) == 't3e' )  THEN
4214!
4215!--    Time needed must be at least 30 seconds on all CRAY machines, because
4216!--    MPP_TREMAIN gives the remaining CPU time only in steps of 30 seconds
4217       IF ( termination_time_needed <= 30.0_wp )  THEN
4218          WRITE( message_string, * )  'termination_time_needed = ',            &
4219                 termination_time_needed, ' must be > 30.0 on host "',         &
4220                 TRIM( host ), '"'
4221          CALL message( 'check_parameters', 'PA0139', 1, 2, 0, 6, 0 )
4222       ENDIF
4223    ELSEIF ( host(1:3) == 'ibm' )  THEN
4224!
4225!--    On IBM-regatta machines the time should be at least 300 seconds,
4226!--    because the job time consumed before executing palm (for compiling,
4227!--    copying of files, etc.) has to be regarded
4228       IF ( termination_time_needed < 300.0_wp )  THEN
4229          WRITE( message_string, * )  'termination_time_needed = ',            &
4230                 termination_time_needed, ' should be >= 300.0 on host "',     &
4231                 TRIM( host ), '"'
4232          CALL message( 'check_parameters', 'PA0140', 1, 2, 0, 6, 0 )
4233       ENDIF
4234    ENDIF
4235
4236!
4237!-- Check pressure gradient conditions
4238    IF ( dp_external  .AND.  conserve_volume_flow )  THEN
4239       WRITE( message_string, * )  'Both dp_external and conserve_volume_flo', &
4240            'w are .TRUE. but one of them must be .FALSE.'
4241       CALL message( 'check_parameters', 'PA0150', 1, 2, 0, 6, 0 )
4242    ENDIF
4243    IF ( dp_external )  THEN
4244       IF ( dp_level_b < zu(nzb)  .OR.  dp_level_b > zu(nzt) )  THEN
4245          WRITE( message_string, * )  'dp_level_b = ', dp_level_b, ' is out ', &
4246               ' of range'
4247          CALL message( 'check_parameters', 'PA0151', 1, 2, 0, 6, 0 )
4248       ENDIF
4249       IF ( .NOT. ANY( dpdxy /= 0.0_wp ) )  THEN
4250          WRITE( message_string, * )  'dp_external is .TRUE. but dpdxy is ze', &
4251               'ro, i.e. the external pressure gradient & will not be applied'
4252          CALL message( 'check_parameters', 'PA0152', 0, 1, 0, 6, 0 )
4253       ENDIF
4254    ENDIF
4255    IF ( ANY( dpdxy /= 0.0_wp )  .AND.  .NOT.  dp_external )  THEN
4256       WRITE( message_string, * )  'dpdxy is nonzero but dp_external is ',     &
4257            '.FALSE., i.e. the external pressure gradient & will not be applied'
4258       CALL message( 'check_parameters', 'PA0153', 0, 1, 0, 6, 0 )
4259    ENDIF
4260    IF ( conserve_volume_flow )  THEN
4261       IF ( TRIM( conserve_volume_flow_mode ) == 'default' )  THEN
4262
4263          conserve_volume_flow_mode = 'initial_profiles'
4264
4265       ELSEIF ( TRIM( conserve_volume_flow_mode ) /= 'initial_profiles' .AND.  &
4266            TRIM( conserve_volume_flow_mode ) /= 'inflow_profile' .AND.        &
4267            TRIM( conserve_volume_flow_mode ) /= 'bulk_velocity' )  THEN
4268          WRITE( message_string, * )  'unknown conserve_volume_flow_mode: ',   &
4269               conserve_volume_flow_mode
4270          CALL message( 'check_parameters', 'PA0154', 1, 2, 0, 6, 0 )
4271       ENDIF
4272       IF ( (bc_lr /= 'cyclic'  .OR.  bc_ns /= 'cyclic')  .AND.                &
4273          TRIM( conserve_volume_flow_mode ) == 'bulk_velocity' )  THEN
4274          WRITE( message_string, * )  'non-cyclic boundary conditions ',       &
4275               'require  conserve_volume_flow_mode = ''initial_profiles'''
4276          CALL message( 'check_parameters', 'PA0155', 1, 2, 0, 6, 0 )
4277       ENDIF
4278       IF ( bc_lr == 'cyclic'  .AND.  bc_ns == 'cyclic'  .AND.                 &
4279            TRIM( conserve_volume_flow_mode ) == 'inflow_profile' )  THEN
4280          WRITE( message_string, * )  'cyclic boundary conditions ',           &
4281               'require conserve_volume_flow_mode = ''initial_profiles''',     &
4282               ' or ''bulk_velocity'''
4283          CALL message( 'check_parameters', 'PA0156', 1, 2, 0, 6, 0 )
4284       ENDIF
4285    ENDIF
4286    IF ( ( u_bulk /= 0.0_wp  .OR.  v_bulk /= 0.0_wp )  .AND.                   &
4287         ( .NOT. conserve_volume_flow  .OR.                                    &
4288         TRIM( conserve_volume_flow_mode ) /= 'bulk_velocity' ) )  THEN
4289       WRITE( message_string, * )  'nonzero bulk velocity requires ',          &
4290            'conserve_volume_flow = .T. and ',                                 &
4291            'conserve_volume_flow_mode = ''bulk_velocity'''
4292       CALL message( 'check_parameters', 'PA0157', 1, 2, 0, 6, 0 )
4293    ENDIF
4294
4295!
4296!-- Check particle attributes
4297    IF ( particle_color /= 'none' )  THEN
4298       IF ( particle_color /= 'absuv'  .AND.  particle_color /= 'pt*'  .AND.   &
4299            particle_color /= 'z' )  THEN
4300          message_string = 'illegal value for parameter particle_color: ' //   &
4301                           TRIM( particle_color)
4302          CALL message( 'check_parameters', 'PA0313', 1, 2, 0, 6, 0 )
4303       ELSE
4304          IF ( color_interval(2) <= color_interval(1) )  THEN
4305             message_string = 'color_interval(2) <= color_interval(1)'
4306             CALL message( 'check_parameters', 'PA0315', 1, 2, 0, 6, 0 )
4307          ENDIF
4308       ENDIF
4309    ENDIF
4310
4311    IF ( particle_dvrpsize /= 'none' )  THEN
4312       IF ( particle_dvrpsize /= 'absw' )  THEN
4313          message_string = 'illegal value for parameter particle_dvrpsize:' // &
4314                           ' ' // TRIM( particle_color)
4315          CALL message( 'check_parameters', 'PA0314', 1, 2, 0, 6, 0 )
4316       ELSE
4317          IF ( dvrpsize_interval(2) <= dvrpsize_interval(1) )  THEN
4318             message_string = 'dvrpsize_interval(2) <= dvrpsize_interval(1)'
4319             CALL message( 'check_parameters', 'PA0316', 1, 2, 0, 6, 0 )
4320          ENDIF
4321       ENDIF
4322    ENDIF
4323
4324!
4325!-- Check nudging and large scale forcing from external file
4326    IF ( nudging  .AND.  (  .NOT.  large_scale_forcing ) )  THEN
4327       message_string = 'Nudging requires large_scale_forcing = .T.. &'//      &
4328                        'Surface fluxes and geostrophic wind should be &'//    &
4329                        'prescribed in file LSF_DATA'
4330       CALL message( 'check_parameters', 'PA0374', 1, 2, 0, 6, 0 )
4331    ENDIF
4332
4333    IF ( large_scale_forcing  .AND.  ( bc_lr /= 'cyclic'  .OR.                 &
4334                                    bc_ns /= 'cyclic' ) )  THEN
4335       message_string = 'Non-cyclic lateral boundaries do not allow for &' //  &
4336                        'the usage of large scale forcing from external file.'
4337       CALL message( 'check_parameters', 'PA0375', 1, 2, 0, 6, 0 )
4338     ENDIF
4339
4340    IF ( large_scale_forcing  .AND.  (  .NOT.  humidity ) )  THEN
4341       message_string = 'The usage of large scale forcing from external &'//   & 
4342                        'file LSF_DATA requires humidity = .T..'
4343       CALL message( 'check_parameters', 'PA0376', 1, 2, 0, 6, 0 )
4344     ENDIF
4345
4346    IF ( large_scale_forcing  .AND.  topography /= 'flat' )  THEN
4347       message_string = 'The usage of large scale forcing from external &'//   & 
4348                        'file LSF_DATA is not implemented for non-flat topography'
4349       CALL message( 'check_parameters', 'PA0377', 1, 2, 0, 6, 0 )
4350    ENDIF
4351
4352    IF ( large_scale_forcing  .AND.  ocean  )  THEN
4353       message_string = 'The usage of large scale forcing from external &'//   & 
4354                        'file LSF_DATA is not implemented for ocean runs'
4355       CALL message( 'check_parameters', 'PA0378', 1, 2, 0, 6, 0 )
4356    ENDIF
4357
4358    CALL location_message( 'finished', .TRUE. )
4359
4360!
4361!-- Prevent empty time records in volume, cross-section and masked data in case of
4362!-- non-parallel netcdf-output in restart runs
4363    IF ( netcdf_data_format < 5 )  THEN
4364       IF ( TRIM( initializing_actions ) == 'read_restart_data' )  THEN
4365          do3d_time_count    = 0
4366          do2d_xy_time_count = 0
4367          do2d_xz_time_count = 0
4368          do2d_yz_time_count = 0
4369          domask_time_count  = 0
4370       ENDIF
4371    ENDIF
4372
4373!
4374!-- Check for valid setting of most_method
4375    IF ( TRIM( most_method ) /= 'circular'  .AND.                              &
4376         TRIM( most_method ) /= 'newton'    .AND.                              &
4377         TRIM( most_method ) /= 'lookup' )  THEN
4378       message_string = 'most_method = "' // TRIM( most_method ) //            &
4379                        '" is unknown'
4380       CALL message( 'check_parameters', 'PA0416', 1, 2, 0, 6, 0 )
4381    ENDIF
4382
4383!
4384!-- Check &userpar parameters
4385    CALL user_check_parameters
4386
4387 CONTAINS
4388
4389!------------------------------------------------------------------------------!
4390! Description:
4391! ------------
4392!> Check the length of data output intervals. In case of parallel NetCDF output
4393!> the time levels of the output files need to be fixed. Therefore setting the
4394!> output interval to 0.0s (usually used to output each timestep) is not
4395!> possible as long as a non-fixed timestep is used.
4396!------------------------------------------------------------------------------!
4397
4398    SUBROUTINE check_dt_do( dt_do, dt_do_name )
4399
4400       IMPLICIT NONE
4401
4402       CHARACTER (LEN=*), INTENT (IN) :: dt_do_name !< parin variable name
4403
4404       REAL(wp), INTENT (INOUT)       :: dt_do      !< data output interval
4405
4406       IF ( dt_do == 0.0_wp )  THEN
4407          IF ( dt_fixed )  THEN
4408             WRITE( message_string, '(A,F9.4,A)' )  'Output at every '  //     &
4409                    'timestep is desired (' // dt_do_name // ' = 0.0).&'//     &
4410                    'Setting the output interval to the fixed timestep '//     &
4411                    'dt = ', dt, 's.'
4412             CALL message( 'check_parameters', 'PA0060', 0, 0, 0, 6, 0 )
4413             dt_do = dt
4414          ELSE
4415             message_string = dt_do_name // ' = 0.0 while using a ' //         &
4416                              'variable timestep and parallel netCDF4 ' //     &
4417                              'is not allowed.'
4418             CALL message( 'check_parameters', 'PA0081', 1, 2, 0, 6, 0 )
4419          ENDIF
4420       ENDIF
4421
4422    END SUBROUTINE check_dt_do
4423
4424 END SUBROUTINE check_parameters
Note: See TracBrowser for help on using the repository browser.