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

Last change on this file since 1552 was 1552, checked in by maronga, 9 years ago

last commit documented

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