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

Last change on this file since 1574 was 1574, checked in by suehring, 9 years ago

last commit documented

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