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

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

last commit documented

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