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

Last change on this file since 39 was 39, checked in by raasch, 16 years ago

comments prepared for 3.1c

  • Property svn:keywords set to Id
File size: 84.5 KB
Line 
1 SUBROUTINE check_parameters
2
3!------------------------------------------------------------------------------!
4! Actual revisions:
5! -----------------
6!
7!
8! Former revisions:
9! -----------------
10! $Id: check_parameters.f90 39 2007-03-01 12:46:59Z raasch $
11!
12! 20 2007-02-26 00:12:32Z raasch
13! Temperature and humidity gradients at top are now calculated for nzt+1,
14! top_heatflux and respective boundary condition bc_pt_t is checked
15!
16! RCS Log replace by Id keyword, revision history cleaned up
17!
18! Revision 1.61  2006/08/04 14:20:25  raasch
19! do2d_unit and do3d_unit now defined as 2d-arrays, check of
20! use_upstream_for_tke, default value for dt_dopts,
21! generation of file header moved from routines palm and header to here
22!
23! Revision 1.1  1997/08/26 06:29:23  raasch
24! Initial revision
25!
26!
27! Description:
28! ------------
29! Check control parameters and deduce further quantities.
30!------------------------------------------------------------------------------!
31
32    USE arrays_3d
33    USE constants
34    USE control_parameters
35    USE grid_variables
36    USE indices
37    USE model_1d
38    USE netcdf_control
39    USE particle_attributes
40    USE pegrid
41    USE profil_parameter
42    USE statistics
43    USE transpose_indices
44
45    IMPLICIT NONE
46
47    CHARACTER (LEN=1)   ::  sq
48    CHARACTER (LEN=6)   ::  var
49    CHARACTER (LEN=7)   ::  unit
50    CHARACTER (LEN=8)   ::  date
51    CHARACTER (LEN=10)  ::  time
52    CHARACTER (LEN=100) ::  action
53
54    INTEGER ::  i, ilen, intervals, iter, j, k, nnxh, nnyh, position, prec
55    LOGICAL ::  found, ldum
56    REAL    ::  gradient, maxn, maxp
57
58
59!
60!-- Warning, if host is not set
61    IF ( host(1:1) == ' ' )  THEN
62       IF ( myid == 0 )  THEN
63          PRINT*, '+++ WARNING: check_parameters:'
64          PRINT*, '    "host" is not set. Please check that environment', &
65                       ' variable "localhost"'
66          PRINT*, '    is set before running PALM'
67       ENDIF
68    ENDIF
69
70!
71!-- Generate the file header which is used as a header for most of PALM's
72!-- output files
73    CALL DATE_AND_TIME( date, time )
74    run_date = date(7:8)//'-'//date(5:6)//'-'//date(3:4)
75    run_time = time(1:2)//':'//time(3:4)//':'//time(5:6)
76
77    WRITE ( run_description_header, '(A,2X,A,A,A,I2.2,2X,A,A,2X,A,1X,A)' )  &
78              TRIM( version ),'run: ', TRIM( run_identifier ), '.', &
79              runnr, 'host: ', TRIM( host ), run_date, run_time
80
81!
82!-- Check topography setting (check for illegal parameter combinations)
83    IF ( topography /= 'flat' )  THEN
84       action = ' '
85       IF ( scalar_advec /= 'pw-scheme' )  THEN
86          WRITE( action, '(A,A)' )  'scalar_advec = ', scalar_advec
87       ENDIF
88       IF ( momentum_advec /= 'pw-scheme' )  THEN
89          WRITE( action, '(A,A)' )  'momentum_advec = ', momentum_advec
90       ENDIF
91       IF ( psolver == 'multigrid'  .OR.  psolver == 'sor' )  THEN
92          WRITE( action, '(A,A)' )  'psolver = ', psolver
93       ENDIF
94       IF ( sloping_surface )  THEN
95          WRITE( action, '(A)' )  'sloping surface = .TRUE.'
96       ENDIF
97       IF ( galilei_transformation )  THEN
98          WRITE( action, '(A)' )  'galilei_transformation = .TRUE.'
99       ENDIF
100       IF ( cloud_physics )  THEN
101          WRITE( action, '(A)' )  'cloud_physics = .TRUE.'
102       ENDIF
103       IF ( cloud_droplets )  THEN
104          WRITE( action, '(A)' )  'cloud_droplets = .TRUE.'
105       ENDIF
106       IF ( moisture )  THEN
107          WRITE( action, '(A)' )  'moisture = .TRUE.'
108       ENDIF
109       IF ( .NOT. prandtl_layer )  THEN
110          WRITE( action, '(A)' )  'prandtl_layer = .FALSE.'
111       ENDIF
112       IF ( action /= ' ' )  THEN
113          IF ( myid == 0 )  THEN
114             PRINT*, '+++ check_parameters:'
115             PRINT*, '    a non-flat topography does not allow ', TRIM( action )
116          ENDIF
117          CALL local_stop
118       ENDIF
119    ENDIF
120!
121!-- Check whether there are any illegal values
122!-- Pressure solver:
123    IF ( psolver /= 'poisfft'  .AND.  psolver /= 'poisfft_hybrid'  .AND. &
124         psolver /= 'sor'  .AND.  psolver /= 'multigrid' )  THEN
125       IF ( myid == 0 )  THEN
126          PRINT*, '+++ check_parameters:'
127          PRINT*, '    unknown solver for perturbation pressure: psolver=', &
128                  psolver
129       ENDIF
130       CALL local_stop
131    ENDIF
132
133#if defined( __parallel )
134    IF ( psolver == 'poisfft_hybrid'  .AND.  pdims(2) /= 1 )  THEN
135       IF ( myid == 0 )  THEN
136          PRINT*, '+++ check_parameters:'
137          PRINT*, '    psolver="', TRIM( psolver ), '" only works for a ', &
138                       '1d domain-decomposition along x'
139          PRINT*, '    please do not set npey/=1 in the parameter file'
140       ENDIF
141       CALL local_stop
142    ENDIF
143    IF ( ( psolver == 'poisfft_hybrid'  .OR.  psolver == 'multigrid' )  .AND.  &
144         ( nxra > nxr  .OR.  nyna > nyn  .OR.  nza > nz ) )  THEN
145       IF ( myid == 0 )  THEN
146          PRINT*, '+++ check_parameters:'
147          PRINT*, '    psolver="', TRIM( psolver ), '" does not work for ', &
148                       'subdomains with unequal size'
149          PRINT*, '    please set grid_matching = ''strict'' in the parameter',&
150                       ' file'
151       ENDIF
152       CALL local_stop
153    ENDIF
154#else
155    IF ( psolver == 'poisfft_hybrid' )  THEN
156       IF ( myid == 0 )  THEN
157          PRINT*, '+++ check_parameters:'
158          PRINT*, '    psolver="', TRIM( psolver ), '" only works for a ', &
159                       'parallel environment'
160       ENDIF
161       CALL local_stop
162    ENDIF
163#endif
164
165    IF ( psolver == 'multigrid' )  THEN
166       IF ( cycle_mg == 'w' )  THEN
167          gamma_mg = 2
168       ELSEIF ( cycle_mg == 'v' )  THEN
169          gamma_mg = 1
170       ELSE
171          IF ( myid == 0 )  THEN
172             PRINT*, '+++ check_parameters:'
173             PRINT*, '    unknown multigrid cycle: cycle_mg=', cycle_mg
174          ENDIF
175          CALL local_stop
176       ENDIF
177    ENDIF
178
179    IF ( fft_method /= 'singleton-algorithm'  .AND.  &
180         fft_method /= 'temperton-algorithm'  .AND.  &
181         fft_method /= 'system-specific' )  THEN
182       IF ( myid == 0 )  THEN
183          PRINT*, '+++ check_parameters:'
184          PRINT*, '    unknown fft-algorithm: fft_method=', fft_method
185       ENDIF
186       CALL local_stop
187    ENDIF
188
189!
190!-- Advection schemes:
191    IF ( momentum_advec /= 'pw-scheme' .AND. momentum_advec /= 'ups-scheme' ) &
192    THEN
193       IF ( myid == 0 )  PRINT*, '+++ check_parameters:  unknown advection ', &
194                                 'scheme: momentum_advec=', momentum_advec
195       CALL local_stop
196    ENDIF
197    IF ( ( momentum_advec == 'ups-scheme'  .OR.  scalar_advec == 'ups-scheme' )&
198                                      .AND.  timestep_scheme /= 'euler' )  THEN
199       IF ( myid == 0 )  PRINT*, '+++ check_parameters:  momentum_advec=', &
200                                 momentum_advec, ' is not allowed with ', &
201                                 'timestep_scheme=', timestep_scheme
202       CALL local_stop
203    ENDIF
204
205    IF ( scalar_advec /= 'pw-scheme'  .AND.  scalar_advec /= 'bc-scheme'  .AND.&
206         scalar_advec /= 'ups-scheme' )  THEN
207       IF ( myid == 0 )  PRINT*, '+++ check_parameters:  unknown advection ', &
208                                 'scheme: scalar_advec=', scalar_advec
209       CALL local_stop
210    ENDIF
211
212    IF ( use_sgs_for_particles  .AND.  .NOT. use_upstream_for_tke )  THEN
213       use_upstream_for_tke = .TRUE.
214       IF ( myid == 0 )  THEN
215          PRINT*, '+++ WARNING: check_parameters:  use_upstream_for_tke set ', &
216                       '.TRUE. because use_sgs_for_particles = .TRUE.'
217       ENDIF
218    ENDIF
219
220    IF ( use_upstream_for_tke  .AND.  timestep_scheme(1:8) == 'leapfrog' )  THEN
221       IF ( myid == 0 )  THEN
222          PRINT*, '+++ check_parameters:  use_upstream_for_tke = .TRUE. ', &
223                       'not allowed with timestep_scheme=', timestep_scheme
224       ENDIF
225       CALL local_stop
226    ENDIF
227
228!
229!-- Timestep schemes:
230    SELECT CASE ( TRIM( timestep_scheme ) )
231
232       CASE ( 'euler' )
233          intermediate_timestep_count_max = 1
234          asselin_filter_factor           = 0.0
235
236       CASE ( 'leapfrog', 'leapfrog+euler' )
237          intermediate_timestep_count_max = 1
238
239       CASE ( 'runge-kutta-2' )
240          intermediate_timestep_count_max = 2
241          asselin_filter_factor           = 0.0
242
243       CASE ( 'runge-kutta-3' )
244          intermediate_timestep_count_max = 3
245          asselin_filter_factor           = 0.0
246
247       CASE DEFAULT
248          IF ( myid == 0 )  PRINT*, '+++ check_parameters:  unknown timestep ',&
249                                    'scheme: timestep_scheme=', timestep_scheme
250          CALL local_stop
251
252    END SELECT
253
254    IF ( scalar_advec /= 'pw-scheme'  .AND.  timestep_scheme(1:5) == 'runge' ) &
255    THEN
256       IF ( myid == 0 )  THEN
257          PRINT*, '+++ check_parameters:  scalar advection scheme "', &
258                                          TRIM( scalar_advec ), '"'
259          PRINT*, '    does not work with timestep_scheme "', &
260                                          TRIM( timestep_scheme ), '"'
261       ENDIF
262       CALL local_stop
263    ENDIF
264
265    IF ( momentum_advec /= 'pw-scheme' .AND. timestep_scheme(1:5) == 'runge' ) &
266    THEN
267       IF ( myid == 0 )  THEN
268          PRINT*, '+++ check_parameters:  momentum advection scheme "', &
269                                          TRIM( momentum_advec ), '"'
270          PRINT*, '    does not work with timestep_scheme "', &
271                                          TRIM( timestep_scheme ), '"'
272       ENDIF
273       CALL local_stop
274    ENDIF
275
276
277    IF ( initializing_actions == ' ' )  THEN
278       IF ( myid == 0 )  THEN
279          PRINT*, '+++ check parameters:'
280          PRINT*, '    no value found for initializing_actions'
281       ENDIF
282       CALL local_stop
283    ENDIF
284
285    IF ( TRIM( initializing_actions ) /= 'read_restart_data' )  THEN
286!
287!--    No model continuation run; several initialising actions are possible
288       action = initializing_actions
289       DO WHILE ( TRIM( action ) /= '' )
290          position = INDEX( action, ' ' )
291          SELECT CASE ( action(1:position-1) )
292
293             CASE ( 'set_constant_profiles', 'set_1d-model_profiles', &
294                    'initialize_vortex',     'initialize_ptanom' )
295                action = action(position+1:)
296
297             CASE DEFAULT
298                IF ( myid == 0 )  PRINT*, '+++ check_parameters: initializi', &
299                               'ng_action unkown or not allowed: action = "', &
300                               TRIM(action), '"'
301                CALL local_stop
302
303          END SELECT
304       ENDDO
305    ENDIF
306    IF ( INDEX( initializing_actions, 'set_constant_profiles' ) /= 0  .AND. &
307         INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
308       IF ( myid == 0 )  PRINT*, '+++ check_parameters: initializing_actions', &
309          '"set_constant_profiles" and "set_1d-model_profiles" are not', &
310          ' allowed simultaneously'
311       CALL local_stop
312    ENDIF
313
314    IF ( cloud_physics  .AND.  .NOT. moisture )  THEN
315       IF ( myid == 0 )  PRINT*, '+++ check_parameters: moisture =', &
316                                 moisture, ' is not allowed with ',  &
317                                 'cloud_physics=', cloud_physics
318       CALL local_stop
319    ENDIF
320
321    IF ( moisture  .AND.  sloping_surface )  THEN
322       IF ( myid == 0 )  PRINT*, '+++ check_parameters: moisture = TRUE', &
323                                 'and hang = TRUE are not',               &
324                                 ' allowed simultaneously' 
325       CALL local_stop       
326    ENDIF
327
328    IF ( moisture  .AND.  scalar_advec == 'ups-scheme' )  THEN
329       IF ( myid == 0 )  PRINT*, '+++ check_parameters: UPS-scheme', &
330                                 'is not implemented for moisture'
331       CALL local_stop       
332    ENDIF
333
334    IF ( passive_scalar  .AND.  moisture )  THEN
335       IF ( myid == 0 )  PRINT*, '+++ check_parameters: moisture = TRUE and', &
336                                 'passive_scalar = TRUE is not allowed ',     &
337                                 'simultaneously'
338       CALL local_stop
339    ENDIF
340
341    IF ( passive_scalar  .AND.  scalar_advec == 'ups-scheme' )  THEN
342       IF ( myid == 0 )  PRINT*, '+++ check_parameters: UPS-scheme', &
343                                 'is not implemented for passive_scalar'
344       CALL local_stop       
345    ENDIF
346
347    IF ( grid_matching /= 'strict'  .AND.  grid_matching /= 'match' )  THEN
348       IF ( myid == 0 )  PRINT*, '+++ check_parameters: illegal value "', &
349                                 TRIM( grid_matching ),                   &
350                                 '" found for parameter grid_matching'
351       CALL local_stop       
352    ENDIF
353
354!
355!-- In case of no model continuation run, check initialising parameters and
356!-- deduce further quantities
357    IF ( TRIM( initializing_actions ) /= 'read_restart_data' )  THEN
358
359!
360!--    Initial profiles for 1D and 3D model, respectively
361       u_init  = ug_surface
362       v_init  = vg_surface
363       pt_init = pt_surface
364       IF ( moisture )        q_init = q_surface
365       IF ( passive_scalar )  q_init = s_surface
366
367!
368!--
369!--    If required, compute initial profile of the geostrophic wind
370!--    (component ug)
371       i = 1
372       gradient = 0.0
373       ug_vertical_gradient_level_ind(1) = 0
374       ug(0) = ug_surface
375       DO  k = 1, nzt+1
376          IF ( ug_vertical_gradient_level(i) < zu(k)  .AND. &
377               ug_vertical_gradient_level(i) >= 0.0 )  THEN
378             gradient = ug_vertical_gradient(i) / 100.0
379             ug_vertical_gradient_level_ind(i) = k - 1
380             i = i + 1
381             IF ( i > 10 )  THEN
382                IF ( myid == 0 )  THEN
383                   PRINT*, '+++ check_parameters: upper bound 10 of array', &
384                           ' "ug_vertical_gradient_level_ind" exceeded'
385                ENDIF
386                CALL local_stop
387             ENDIF
388          ENDIF
389          IF ( gradient /= 0.0 )  THEN
390             IF ( k /= 1 )  THEN
391                ug(k) = ug(k-1) + dzu(k) * gradient
392             ELSE
393                ug(k) = ug_surface + 0.5 * dzu(k) * gradient
394             ENDIF
395          ELSE
396             ug(k) = ug(k-1)
397          ENDIF
398       ENDDO
399
400       u_init = ug
401
402!
403!--    In case of no given gradients for ug, choose a vanishing gradient
404       IF ( ug_vertical_gradient_level(1) == -1.0 )  THEN
405          ug_vertical_gradient_level(1) = 0.0
406       ENDIF 
407
408!
409!--
410!--    If required, compute initial profile of the geostrophic wind
411!--    (component vg)
412       i = 1
413       gradient = 0.0
414       vg_vertical_gradient_level_ind(1) = 0
415       vg(0) = vg_surface
416       DO  k = 1, nzt+1
417          IF ( vg_vertical_gradient_level(i) < zu(k)  .AND. &
418               vg_vertical_gradient_level(i) >= 0.0 )  THEN
419             gradient = vg_vertical_gradient(i) / 100.0
420             vg_vertical_gradient_level_ind(i) = k - 1
421             i = i + 1
422             IF ( i > 10 )  THEN
423                IF ( myid == 0 )  THEN
424                   PRINT*, '+++ check_parameters: upper bound 10 of array', &
425                           ' "vg_vertical_gradient_level_ind" exceeded'
426                ENDIF
427                CALL local_stop
428             ENDIF
429          ENDIF
430          IF ( gradient /= 0.0 )  THEN
431             IF ( k /= 1 )  THEN
432                vg(k) = vg(k-1) + dzu(k) * gradient
433             ELSE
434                vg(k) = vg_surface + 0.5 * dzu(k) * gradient
435             ENDIF
436          ELSE
437             vg(k) = vg(k-1)
438          ENDIF
439       ENDDO
440
441       v_init = vg
442 
443!
444!--    In case of no given gradients for vg, choose a vanishing gradient
445       IF ( vg_vertical_gradient_level(1) == -1.0 )  THEN
446          vg_vertical_gradient_level(1) = 0.0
447       ENDIF
448
449!
450!--    If required, compute initial temperature profile using the given
451!--    temperature gradient
452       i = 1
453       gradient = 0.0
454       pt_vertical_gradient_level_ind(1) = 0
455       DO  k = 1, nzt+1
456          IF ( pt_vertical_gradient_level(i) < zu(k)  .AND. &
457               pt_vertical_gradient_level(i) >= 0.0 )  THEN
458             gradient = pt_vertical_gradient(i) / 100.0
459             pt_vertical_gradient_level_ind(i) = k - 1
460             i = i + 1
461             IF ( i > 10 )  THEN
462                IF ( myid == 0 )  THEN
463                   PRINT*, '+++ check_parameters: upper bound 10 of array', &
464                           ' "pt_vertical_gradient_level_ind" exceeded'
465                ENDIF
466                CALL local_stop
467             ENDIF
468          ENDIF
469          IF ( gradient /= 0.0 )  THEN
470             IF ( k /= 1 )  THEN
471                pt_init(k) = pt_init(k-1) + dzu(k) * gradient
472             ELSE
473                pt_init(k) = pt_init(k-1) + 0.5 * dzu(k) * gradient
474             ENDIF
475          ELSE
476             pt_init(k) = pt_init(k-1)
477          ENDIF
478       ENDDO
479
480!
481!--    In case of no given temperature gradients, choose gradient of neutral
482!--    stratification
483       IF ( pt_vertical_gradient_level(1) == -1.0 )  THEN
484          pt_vertical_gradient_level(1) = 0.0
485       ENDIF
486
487!
488!--    Store temperature gradient at the top boundary for possile Neumann
489!--    boundary condition
490       bc_pt_t_val = ( pt_init(nzt+1) - pt_init(nzt) ) / dzu(nzt+1)
491
492!
493!--    If required, compute initial humidity or scalar profile using the given
494!--    humidity/scalar gradient. In case of scalar transport, initially store
495!--    values of the scalar parameters on humidity parameters
496       IF ( passive_scalar )  THEN
497          bc_q_b                    = bc_s_b
498          bc_q_t                    = bc_s_t
499          q_surface                 = s_surface
500          q_surface_initial_change  = s_surface_initial_change
501          q_vertical_gradient       = s_vertical_gradient
502          q_vertical_gradient_level = s_vertical_gradient_level
503          surface_waterflux         = surface_scalarflux
504       ENDIF
505
506       IF ( moisture  .OR.  passive_scalar )  THEN
507
508          i = 1
509          gradient = 0.0
510          q_vertical_gradient_level_ind(1) = 0
511          DO  k = 1, nzt+1
512             IF ( q_vertical_gradient_level(i) < zu(k)  .AND. &
513                  q_vertical_gradient_level(i) >= 0.0 )  THEN
514                gradient = q_vertical_gradient(i) / 100.0
515                q_vertical_gradient_level_ind(i) = k - 1
516                i = i + 1
517                IF ( i > 10 )  THEN
518                   IF ( myid == 0 )  THEN
519                      PRINT*, '+++ check_parameters: upper bound 10 of arr', &
520                              'ay "q_vertical_gradient_level_ind" exceeded'
521                   ENDIF
522                   CALL local_stop
523                ENDIF
524             ENDIF
525             IF ( gradient /= 0.0 )  THEN
526                IF ( k /= 1 )  THEN
527                   q_init(k) = q_init(k-1) + dzu(k) * gradient
528                ELSE
529                   q_init(k) = q_init(k-1) + 0.5 * dzu(k) * gradient
530                ENDIF
531             ELSE
532                q_init(k) = q_init(k-1)
533             ENDIF
534          ENDDO
535
536!
537!--       In case of no given humidity gradients, choose zero gradient
538!--       conditions
539          IF ( q_vertical_gradient_level(1) == -1.0 )  THEN
540             q_vertical_gradient_level(1) = 0.0
541          ENDIF
542
543!
544!--       Store humidity gradient at the top boundary for possile Neumann
545!--       boundary condition
546          bc_q_t_val = ( q_init(nzt+1) - q_init(nzt) ) / dzu(nzt+1)
547
548       ENDIF
549
550    ENDIF
551
552!
553!-- Compute Coriolis parameter
554    f  = 2.0 * omega * SIN( phi / 180.0 * pi )
555    fs = 2.0 * omega * COS( phi / 180.0 * pi )
556
557!
558!-- In case of a given slope, compute the relevant quantities
559    IF ( alpha_surface /= 0.0 )  THEN
560       IF ( ABS( alpha_surface ) > 90.0 )  THEN
561          IF ( myid == 0 )  PRINT*, '+++ check_parameters: ABS( alpha_surface',&
562                                    '=', alpha_surface, ' ) must be < 90.0'
563          CALL local_stop
564       ENDIF
565       sloping_surface = .TRUE.
566       cos_alpha_surface = COS( alpha_surface / 180.0 * pi )
567       sin_alpha_surface = SIN( alpha_surface / 180.0 * pi )
568    ENDIF
569
570!
571!-- Check time step and cfl_factor
572    IF ( dt /= -1.0 )  THEN
573       IF ( dt <= 0.0  .AND.  dt /= -1.0 )  THEN
574          IF ( myid == 0 )  PRINT*, '+++ check_parameters:  dt=', dt, ' <= 0.0'
575          CALL local_stop
576       ENDIF
577       dt_3d = dt
578       dt_fixed = .TRUE.
579    ENDIF
580
581    IF ( cfl_factor <= 0.0  .OR.  cfl_factor > 1.0 )  THEN
582       IF ( cfl_factor == -1.0 )  THEN
583          IF ( momentum_advec == 'ups-scheme'  .OR.  &
584               scalar_advec == 'ups-scheme' )  THEN
585             cfl_factor = 0.8
586          ELSE
587             IF ( timestep_scheme == 'runge-kutta-2' )  THEN
588                cfl_factor = 0.8
589             ELSEIF ( timestep_scheme == 'runge-kutta-3' )  THEN
590                cfl_factor = 0.9
591             ELSE
592                cfl_factor = 0.1
593             ENDIF
594          ENDIF
595       ELSE
596          IF ( myid == 0 )  THEN
597             PRINT*, '+++ check_parameters: cfl_factor=', cfl_factor, &
598                         ' out of range'
599             PRINT*, '+++                   0.0 < cfl_factor <= 1.0 is required'
600          ENDIF
601          CALL local_stop
602       ENDIF
603    ENDIF
604
605!
606!-- Store simulated time at begin
607    simulated_time_at_begin = simulated_time
608
609!
610!-- Set wind speed in the Galilei-transformed system
611    IF ( galilei_transformation )  THEN
612       IF ( use_ug_for_galilei_tr .AND.                &
613            ug_vertical_gradient_level(1) == 0.0 .AND. & 
614            vg_vertical_gradient_level(1) == 0.0 )  THEN
615          u_gtrans = ug_surface
616          v_gtrans = vg_surface
617       ELSEIF ( use_ug_for_galilei_tr .AND.                &
618                ug_vertical_gradient_level(1) /= 0.0 )  THEN
619          IF ( myid == 0 )  THEN
620             PRINT*, '+++ check_parameters:'
621             PRINT*, '    baroclinicity (ug) not allowed'
622             PRINT*, '    simultaneously with galilei transformation'
623          ENDIF
624          CALL local_stop
625       ELSEIF ( use_ug_for_galilei_tr .AND.                &
626                vg_vertical_gradient_level(1) /= 0.0 )  THEN
627          IF ( myid == 0 )  THEN
628             PRINT*, '+++ check_parameters:'
629             PRINT*, '    baroclinicity (vg) not allowed'
630             PRINT*, '    simultaneously with galilei transformation'
631          ENDIF
632          CALL local_stop
633       ELSE
634          IF ( myid == 0 )  THEN
635             PRINT*, '+++ WARNING: check_parameters:'
636             PRINT*, '    variable translation speed used for galilei-tran' // &
637                          'sformation, which'
638             PRINT*, '    may cause instabilities in stably stratified regions'
639          ENDIF
640       ENDIF
641    ENDIF
642
643!
644!-- In case of using a prandtl-layer, calculated (or prescribed) surface
645!-- fluxes have to be used in the diffusion-terms
646    IF ( prandtl_layer )  use_surface_fluxes = .TRUE.
647
648!
649!-- Check boundary conditions and set internal variables:
650!-- Lateral boundary conditions
651    IF ( bc_lr /= 'cyclic'  .AND.  bc_lr /= 'dirichlet/neumann'  .AND. &
652         bc_lr /= 'neumann/dirichlet' )  THEN
653       IF ( myid == 0 )  THEN
654          PRINT*, '+++ check_parameters:'
655          PRINT*, '    unknown boundary condition: bc_lr = ', bc_lr
656       ENDIF
657       CALL local_stop
658    ENDIF
659    IF ( bc_ns /= 'cyclic'  .AND.  bc_ns /= 'dirichlet/neumann'  .AND. &
660         bc_ns /= 'neumann/dirichlet' )  THEN
661       IF ( myid == 0 )  THEN
662          PRINT*, '+++ check_parameters:'
663          PRINT*, '    unknown boundary condition: bc_ns = ', bc_ns
664       ENDIF
665       CALL local_stop
666    ENDIF
667
668!
669!-- Non-cyclic lateral boundaries require the multigrid method and Piascek-
670!-- Willimas advection scheme. Several schemes and tools do not work with
671!-- non-cyclic boundary conditions.
672    IF ( bc_lr /= 'cyclic'  .OR.  bc_ns /= 'cyclic' )  THEN
673       IF ( psolver /= 'multigrid' )  THEN
674          IF ( myid == 0 )  THEN
675             PRINT*, '+++ check_parameters:'
676             PRINT*, '    non-cyclic lateral boundaries do not allow', &
677                          ' psolver = ', psolver 
678          ENDIF
679          CALL local_stop
680       ENDIF
681       IF ( momentum_advec /= 'pw-scheme' )  THEN
682          IF ( myid == 0 )  THEN
683             PRINT*, '+++ check_parameters:'
684             PRINT*, '    non-cyclic lateral boundaries do not allow', &
685                          ' momentum_advec = ', momentum_advec 
686          ENDIF
687          CALL local_stop
688       ENDIF
689       IF ( scalar_advec /= 'pw-scheme' )  THEN
690          IF ( myid == 0 )  THEN
691             PRINT*, '+++ check_parameters:'
692             PRINT*, '    non-cyclic lateral boundaries do not allow', &
693                          ' scalar_advec = ', scalar_advec 
694          ENDIF
695          CALL local_stop
696       ENDIF
697       IF ( galilei_transformation )  THEN
698          IF ( myid == 0 )  THEN
699             PRINT*, '+++ check_parameters:'
700             PRINT*, '    non-cyclic lateral boundaries do not allow', &
701                          ' galilei_transformation = .T.' 
702          ENDIF
703          CALL local_stop
704       ENDIF
705       IF ( conserve_volume_flow )  THEN
706          IF ( myid == 0 )  THEN
707             PRINT*, '+++ check_parameters:'
708             PRINT*, '    non-cyclic lateral boundaries do not allow', &
709                          ' conserve_volume_flow = .T.' 
710          ENDIF
711          CALL local_stop
712       ENDIF
713    ENDIF
714
715!
716!-- Bottom boundary condition for the turbulent Kinetic energy
717    IF ( bc_e_b == 'neumann' )  THEN
718       ibc_e_b = 1
719       IF ( adjust_mixing_length  .AND.  prandtl_layer )  THEN
720          IF ( myid == 0 )  THEN
721             PRINT*, '+++ WARNING: check_parameters:'
722             PRINT*, '    adjust_mixing_length = TRUE and bc_e_b = ', bc_e_b
723          ENDIF
724       ENDIF
725    ELSEIF ( bc_e_b == '(u*)**2+neumann' )  THEN
726       ibc_e_b = 2
727       IF ( .NOT. adjust_mixing_length  .AND.  prandtl_layer )  THEN
728          IF ( myid == 0 )  THEN
729             PRINT*, '+++ WARNING: check_parameters:'
730             PRINT*, '    adjust_mixing_length = FALSE and bc_e_b = ', bc_e_b
731          ENDIF
732       ENDIF
733       IF ( .NOT. prandtl_layer )  THEN
734          bc_e_b = 'neumann'
735          ibc_e_b = 1
736          IF ( myid == 0 )  THEN
737             PRINT*, '+++ WARNING: check_parameters:'
738             PRINT*, '    boundary condition bc_e_b changed to "', bc_e_b, '"'
739          ENDIF
740       ENDIF
741    ELSE
742       IF ( myid == 0 )  THEN
743          PRINT*, '+++ check_parameters:'
744          PRINT*, '    unknown boundary condition: bc_e_b = ', bc_e_b
745       ENDIF
746       CALL local_stop
747    ENDIF
748
749!
750!-- Boundary conditions for perturbation pressure
751    IF ( bc_p_b == 'dirichlet' )  THEN
752       ibc_p_b = 0
753    ELSEIF ( bc_p_b == 'neumann' )  THEN
754       ibc_p_b = 1
755    ELSEIF ( bc_p_b == 'neumann+inhomo' )  THEN
756       ibc_p_b = 2
757    ELSE
758       IF ( myid == 0 )  THEN
759          PRINT*, '+++ check_parameters:'
760          PRINT*, '    unknown boundary condition: bc_p_b = ', bc_p_b
761       ENDIF
762       CALL local_stop
763    ENDIF
764    IF ( ibc_p_b == 2  .AND.  .NOT. prandtl_layer )  THEN
765       IF ( myid == 0 )  THEN
766          PRINT*, '+++ check_parameters:'
767          PRINT*, '    boundary condition: bc_p_b = ', TRIM( bc_p_b ), &
768                       ' not allowed with'
769          PRINT*, '    prandtl_layer = .FALSE.' 
770       ENDIF
771       CALL local_stop
772    ENDIF
773    IF ( bc_p_t == 'dirichlet' )  THEN
774       ibc_p_t = 0
775    ELSEIF ( bc_p_t == 'neumann' )  THEN
776       ibc_p_t = 1
777    ELSE
778       IF ( myid == 0 )  THEN
779          PRINT*, '+++ check_parameters:'
780          PRINT*, '    unknown boundary condition: bc_p_t = ', bc_p_t
781       ENDIF
782       CALL local_stop
783    ENDIF
784
785!
786!-- Boundary conditions for potential temperature
787    IF ( bc_pt_b == 'dirichlet' )  THEN
788       ibc_pt_b = 0
789    ELSEIF ( bc_pt_b == 'neumann' )  THEN
790       ibc_pt_b = 1
791    ELSE
792       IF ( myid == 0 )  THEN
793          PRINT*, '+++ check_parameters:'
794          PRINT*, '    unknown boundary condition: bc_pt_b = ', bc_pt_b
795       ENDIF
796       CALL local_stop
797    ENDIF
798    IF ( bc_pt_t == 'dirichlet' )  THEN
799       ibc_pt_t = 0
800    ELSEIF ( bc_pt_t == 'neumann' )  THEN
801       ibc_pt_t = 1
802    ELSEIF ( bc_pt_t == 'initial_gradient' )  THEN
803       ibc_pt_t = 2
804    ELSE
805       IF ( myid == 0 )  THEN
806          PRINT*, '+++ check_parameters:'
807          PRINT*, '    unknown boundary condition: bc_pt_t = ', bc_pt_t
808       ENDIF
809       CALL local_stop
810    ENDIF
811
812    IF ( surface_heatflux == 9999999.9 )  constant_heatflux     = .FALSE.
813    IF ( top_heatflux     == 9999999.9 )  constant_top_heatflux = .FALSE.
814
815!
816!-- A given surface temperature implies Dirichlet boundary condition for
817!-- temperature. In this case specification of a constant heat flux is
818!-- forbidden.
819    IF ( ibc_pt_b == 0  .AND.   constant_heatflux  .AND. &
820         surface_heatflux /= 0.0 )  THEN
821       IF ( myid == 0 )  THEN
822          PRINT*, '+++ check_parameters:'
823          PRINT*, '    boundary_condition: bc_pt_b = ', bc_pt_b
824          PRINT*, '    is not allowed with constant_heatflux = .TRUE.'
825       ENDIF
826       CALL local_stop
827    ENDIF
828    IF ( constant_heatflux  .AND.  pt_surface_initial_change /= 0.0 )  THEN
829       IF ( myid == 0 )  THEN
830          PRINT*, '+++ check_parameters: constant_heatflux = .TRUE. is not'
831          PRINT*, '    allowed with pt_surface_initial_change (/=0) = ', &
832                  pt_surface_initial_change
833       ENDIF
834       CALL local_stop
835    ENDIF
836
837!
838!-- A given temperature at the top implies Dirichlet boundary condition for
839!-- temperature. In this case specification of a constant heat flux is
840!-- forbidden.
841    IF ( ibc_pt_t == 0  .AND.   constant_top_heatflux  .AND. &
842         top_heatflux /= 0.0 )  THEN
843       IF ( myid == 0 )  THEN
844          PRINT*, '+++ check_parameters:'
845          PRINT*, '    boundary_condition: bc_pt_t = ', bc_pt_t
846          PRINT*, '    is not allowed with constant_top_heatflux = .TRUE.'
847       ENDIF
848       CALL local_stop
849    ENDIF
850
851!
852!-- In case of moisture or passive scalar, set boundary conditions for total
853!-- water content / scalar
854    IF ( moisture  .OR.  passive_scalar ) THEN
855       IF ( moisture )  THEN
856          sq = 'q'
857       ELSE
858          sq = 's'
859       ENDIF
860       IF ( bc_q_b == 'dirichlet' )  THEN
861          ibc_q_b = 0
862       ELSEIF ( bc_q_b == 'neumann' )  THEN
863          ibc_q_b = 1
864       ELSE
865          IF ( myid == 0 )  THEN
866             PRINT*, '+++ check_parameters:'
867             PRINT*, '    unknown boundary condition: bc_', sq, '_b = ', bc_q_b
868          ENDIF
869          CALL local_stop
870       ENDIF
871       IF ( bc_q_t == 'dirichlet' )  THEN
872          ibc_q_t = 0
873       ELSEIF ( bc_q_t == 'neumann' )  THEN
874          ibc_q_t = 1
875       ELSE
876          IF ( myid == 0 )  THEN
877             PRINT*, '+++ check_parameters:'
878             PRINT*, '    unknown boundary condition: bc_', sq, '_t = ', bc_q_t
879          ENDIF
880          CALL local_stop
881       ENDIF
882
883       IF ( surface_waterflux == 0.0 )  constant_waterflux = .FALSE.
884
885!
886!--    A given surface humidity implies Dirichlet boundary condition for
887!--    moisture. In this case specification of a constant water flux is
888!--    forbidden.
889       IF ( ibc_q_b == 0  .AND.  constant_waterflux )  THEN
890          IF ( myid == 0 )  THEN
891             PRINT*, '+++ check_parameters:'
892             PRINT*, '    boundary_condition: bc_', sq, '_b = ', bc_q_b
893             PRINT*, '    is not allowed with prescribed surface flux'
894          ENDIF
895          CALL local_stop
896       ENDIF
897       IF ( constant_waterflux  .AND.  q_surface_initial_change /= 0.0 )  THEN
898          IF ( myid == 0 )  THEN
899             PRINT*, '+++ check_parameters: a prescribed surface flux is not'
900             PRINT*, '    allowed with ', sq, '_surface_initial_change (/=0)', &
901                     ' = ', q_surface_initial_change
902          ENDIF
903          CALL local_stop
904       ENDIF
905       
906    ENDIF
907
908!
909!-- Boundary conditions for horizontal components of wind speed
910    IF ( bc_uv_b == 'dirichlet' )  THEN
911       ibc_uv_b = 0
912    ELSEIF ( bc_uv_b == 'neumann' )  THEN
913       ibc_uv_b = 1
914       IF ( prandtl_layer )  THEN
915          IF ( myid == 0 )  THEN
916             PRINT*, '+++ check_parameters:'
917             PRINT*, '    boundary condition: bc_uv_b = ', TRIM( bc_uv_b ), &
918                          ' is not allowed with'
919             PRINT*, '    prandtl_layer = .TRUE.' 
920          ENDIF
921          CALL local_stop
922       ENDIF
923    ELSE
924       IF ( myid == 0 )  THEN
925          PRINT*, '+++ check_parameters:'
926          PRINT*, '    unknown boundary condition: bc_uv_b = ', bc_uv_b
927       ENDIF
928       CALL local_stop
929    ENDIF
930    IF ( bc_uv_t == 'dirichlet' )  THEN
931       ibc_uv_t = 0
932    ELSEIF ( bc_uv_t == 'neumann' )  THEN
933       ibc_uv_t = 1
934    ELSE
935       IF ( myid == 0 )  THEN
936          PRINT*, '+++ check_parameters:'
937          PRINT*, '    unknown boundary condition: bc_uv_t = ', bc_uv_t
938       ENDIF
939       CALL local_stop
940    ENDIF
941
942!
943!-- Compute and check, respectively, the Rayleigh Damping parameter
944    IF ( rayleigh_damping_factor == -1.0 )  THEN
945       IF ( momentum_advec == 'ups-scheme' )  THEN
946          rayleigh_damping_factor = 0.01
947       ELSE
948          rayleigh_damping_factor = 0.0
949       ENDIF
950    ELSE
951       IF ( rayleigh_damping_factor < 0.0 .OR. rayleigh_damping_factor > 1.0 ) &
952       THEN
953          IF ( myid == 0 )  THEN
954             PRINT*, '+++ check_parameters:'
955             PRINT*, '    rayleigh_damping_factor = ', rayleigh_damping_factor,&
956                          ' out of range [0.0,1.0]'
957          ENDIF
958          CALL local_stop
959       ENDIF
960    ENDIF
961
962    IF ( rayleigh_damping_height == -1.0 )  THEN
963       rayleigh_damping_height = 0.66666666666 * zu(nzt)
964    ELSE
965       IF ( rayleigh_damping_height < 0.0  .OR. &
966            rayleigh_damping_height > zu(nzt) )  THEN
967          IF ( myid == 0 )  THEN
968             PRINT*, '+++ check_parameters:'
969             PRINT*, '    rayleigh_damping_height = ', rayleigh_damping_height,&
970                          ' out of range [0.0,', zu(nzt), ']'
971          ENDIF
972          CALL local_stop
973       ENDIF
974    ENDIF
975
976!
977!-- Check limiters for Upstream-Spline scheme
978    IF ( overshoot_limit_u < 0.0  .OR.  overshoot_limit_v < 0.0  .OR.  &
979         overshoot_limit_w < 0.0  .OR.  overshoot_limit_pt < 0.0  .OR. &
980         overshoot_limit_e < 0.0 )  THEN
981       IF ( myid == 0 )  THEN
982          PRINT*, '+++ check_parameters:'
983          PRINT*, '    overshoot_limit_... < 0.0 is not allowed'
984       ENDIF
985       CALL local_stop
986    ENDIF
987    IF ( ups_limit_u < 0.0 .OR. ups_limit_v < 0.0 .OR. ups_limit_w < 0.0 .OR. &
988         ups_limit_pt < 0.0 .OR. ups_limit_e < 0.0 )  THEN
989       IF ( myid == 0 )  THEN
990          PRINT*, '+++ check_parameters:'
991          PRINT*, '    ups_limit_... < 0.0 is not allowed'
992       ENDIF
993       CALL local_stop
994    ENDIF
995
996!
997!-- Check number of chosen statistic regions. More than 10 regions are not
998!-- allowed, because so far no more than 10 corresponding output files can
999!-- be opened (cf. check_open)
1000    IF ( statistic_regions > 9  .OR.  statistic_regions < 0 )  THEN
1001       IF ( myid == 0 )  THEN
1002          PRINT*, '+++ check_parameters: Number of statistic_regions = ', &
1003                       statistic_regions+1
1004          PRINT*, '    Only 10 regions are allowed'
1005       ENDIF
1006       CALL local_stop
1007    ENDIF
1008    IF ( normalizing_region > statistic_regions  .OR. &
1009         normalizing_region < 0)  THEN
1010       IF ( myid == 0 )  THEN
1011          PRINT*, '+++ check_parameters: normalizing_region = ', &
1012                       normalizing_region, ' is unknown'
1013          PRINT*, '    Must be <= ', statistic_regions
1014       ENDIF
1015       CALL local_stop
1016    ENDIF
1017
1018!
1019!-- Set the default intervals for data output, if necessary
1020!-- NOTE: dt_dosp has already been set in package_parin
1021    IF ( dt_data_output /= 9999999.9 )  THEN
1022       IF ( dt_dopr           == 9999999.9 )  dt_dopr           = dt_data_output
1023       IF ( dt_dopts          == 9999999.9 )  dt_dopts          = dt_data_output
1024       IF ( dt_do2d_xy        == 9999999.9 )  dt_do2d_xy        = dt_data_output
1025       IF ( dt_do2d_xz        == 9999999.9 )  dt_do2d_xz        = dt_data_output
1026       IF ( dt_do2d_yz        == 9999999.9 )  dt_do2d_yz        = dt_data_output
1027       IF ( dt_do3d           == 9999999.9 )  dt_do3d           = dt_data_output
1028       IF ( dt_data_output_av == 9999999.9 )  dt_data_output_av = dt_data_output
1029    ENDIF
1030
1031!
1032!-- Set the default skip time intervals for data output, if necessary
1033    IF ( skip_time_dopr    == 9999999.9 ) &
1034                                       skip_time_dopr    = skip_time_data_output
1035    IF ( skip_time_dosp    == 9999999.9 ) &
1036                                       skip_time_dosp    = skip_time_data_output
1037    IF ( skip_time_do2d_xy == 9999999.9 ) &
1038                                       skip_time_do2d_xy = skip_time_data_output
1039    IF ( skip_time_do2d_xz == 9999999.9 ) &
1040                                       skip_time_do2d_xz = skip_time_data_output
1041    IF ( skip_time_do2d_yz == 9999999.9 ) &
1042                                       skip_time_do2d_yz = skip_time_data_output
1043    IF ( skip_time_do3d    == 9999999.9 ) &
1044                                       skip_time_do3d    = skip_time_data_output
1045    IF ( skip_time_data_output_av == 9999999.9 ) &
1046                                skip_time_data_output_av = skip_time_data_output
1047
1048!
1049!-- Check the average intervals (first for 3d-data, then for profiles and
1050!-- spectra)
1051    IF ( averaging_interval > dt_data_output_av )  THEN
1052       IF ( myid == 0 )  THEN
1053          PRINT*, '+++ check_parameters: average_interval=',              &
1054                       averaging_interval, ' must be <= dt_data_output=', &
1055                       dt_data_output
1056       ENDIF
1057       CALL local_stop
1058    ENDIF
1059
1060    IF ( averaging_interval_pr == 9999999.9 )  THEN
1061       averaging_interval_pr = averaging_interval
1062    ENDIF
1063
1064    IF ( averaging_interval_pr > dt_dopr )  THEN
1065       IF ( myid == 0 )  THEN
1066          PRINT*, '+++ check_parameters: averaging_interval_pr=', &
1067                       averaging_interval_pr, ' must be <= dt_dopr=', dt_dopr
1068       ENDIF
1069       CALL local_stop
1070    ENDIF
1071
1072    IF ( averaging_interval_sp == 9999999.9 )  THEN
1073       averaging_interval_sp = averaging_interval
1074    ENDIF
1075
1076    IF ( averaging_interval_sp > dt_dosp )  THEN
1077       IF ( myid == 0 )  THEN
1078          PRINT*, '+++ check_parameters: averaging_interval_sp=', &
1079                       averaging_interval_sp, ' must be <= dt_dosp=', dt_dosp
1080       ENDIF
1081       CALL local_stop
1082    ENDIF
1083
1084!
1085!-- Set the default interval for profiles entering the temporal average
1086    IF ( dt_averaging_input_pr == 9999999.9 )  THEN
1087       dt_averaging_input_pr = dt_averaging_input
1088    ENDIF
1089
1090!
1091!-- Set the default interval for the output of timeseries to a reasonable
1092!-- value (tries to minimize the number of calls of flow_statistics)
1093    IF ( dt_dots == 9999999.9 )  THEN
1094       IF ( averaging_interval_pr == 0.0 )  THEN
1095          dt_dots = MIN( dt_run_control, dt_dopr )
1096       ELSE
1097          dt_dots = MIN( dt_run_control, dt_averaging_input_pr )
1098       ENDIF
1099    ENDIF
1100
1101!
1102!-- Check the sample rate for averaging (first for 3d-data, then for profiles)
1103    IF ( dt_averaging_input > averaging_interval )  THEN
1104       IF ( myid == 0 )  THEN
1105          PRINT*, '+++ check_parameters: dt_averaging_input=',                &
1106                       dt_averaging_input, ' must be <= averaging_interval=', &
1107                       averaging_interval
1108       ENDIF
1109       CALL local_stop
1110    ENDIF
1111
1112    IF ( dt_averaging_input_pr > averaging_interval_pr )  THEN
1113       IF ( myid == 0 )  THEN
1114          PRINT*, '+++ check_parameters: dt_averaging_input_pr=', &
1115                       dt_averaging_input_pr,                     &
1116                       ' must be <= averaging_interval_pr=',      &
1117                       averaging_interval_pr
1118       ENDIF
1119       CALL local_stop
1120    ENDIF
1121
1122!
1123!-- Determine the number of output profiles and check whether they are
1124!-- permissible
1125    DO  WHILE ( data_output_pr(dopr_n+1) /= '          ' )
1126
1127       dopr_n = dopr_n + 1
1128       i = dopr_n
1129
1130!
1131!--    Determine internal profile number (for hom, homs)
1132!--    and store height levels
1133       SELECT CASE ( TRIM( data_output_pr(i) ) )
1134
1135          CASE ( 'u', '#u' )
1136             dopr_index(i) = 1
1137             hom(:,2,1,:)  = SPREAD( zu, 2, statistic_regions+1 )
1138             IF ( data_output_pr(i)(1:1) == '#' )  THEN
1139                dopr_initial_index(i) = 5
1140                hom(:,2,5,:)          = SPREAD( zu, 2, statistic_regions+1 )
1141                data_output_pr(i)     = data_output_pr(i)(2:)
1142             ENDIF
1143
1144          CASE ( 'v', '#v' )
1145             dopr_index(i) = 2
1146             hom(:,2,2,:) = SPREAD( zu, 2, statistic_regions+1 )
1147             IF ( data_output_pr(i)(1:1) == '#' )  THEN
1148                dopr_initial_index(i) = 6
1149                hom(:,2,6,:)          = SPREAD( zu, 2, statistic_regions+1 )
1150                data_output_pr(i)     = data_output_pr(i)(2:)
1151             ENDIF
1152
1153          CASE ( 'w' )
1154             dopr_index(i) = 3
1155             hom(:,2,3,:) = SPREAD( zw, 2, statistic_regions+1 )
1156
1157          CASE ( 'pt', '#pt' )
1158             IF ( .NOT. cloud_physics ) THEN
1159                dopr_index(i) = 4
1160                hom(:,2,4,:)  = SPREAD( zu, 2, statistic_regions+1 )
1161                IF ( data_output_pr(i)(1:1) == '#' )  THEN
1162                   dopr_initial_index(i) = 7
1163                   hom(:,2,7,:)          = SPREAD( zu, 2, statistic_regions+1 )
1164                   hom(nzb,2,7,:)        = 0.0    ! weil zu(nzb) negativ ist
1165                   data_output_pr(i)     = data_output_pr(i)(2:)
1166                ENDIF
1167             ELSE
1168                dopr_index(i) = 43
1169                hom(:,2,43,:)  = SPREAD( zu, 2, statistic_regions+1 )
1170                IF ( data_output_pr(i)(1:1) == '#' )  THEN
1171                   dopr_initial_index(i) = 28
1172                   hom(:,2,28,:)         = SPREAD( zu, 2, statistic_regions+1 )
1173                   hom(nzb,2,28,:)       = 0.0    ! weil zu(nzb) negativ ist
1174                   data_output_pr(i)     = data_output_pr(i)(2:)
1175                ENDIF
1176             ENDIF
1177
1178          CASE ( 'e' )
1179             dopr_index(i)  = 8
1180             hom(:,2,8,:)   = SPREAD( zu, 2, statistic_regions+1 )
1181             hom(nzb,2,8,:) = 0.0
1182
1183          CASE ( 'km', '#km' )
1184             dopr_index(i)  = 9
1185             hom(:,2,9,:)   = SPREAD( zu, 2, statistic_regions+1 )
1186             hom(nzb,2,9,:) = 0.0
1187             IF ( data_output_pr(i)(1:1) == '#' )  THEN
1188                dopr_initial_index(i) = 23
1189                hom(:,2,23,:)         = hom(:,2,9,:)
1190                data_output_pr(i)     = data_output_pr(i)(2:)
1191             ENDIF
1192
1193          CASE ( 'kh', '#kh' )
1194             dopr_index(i)   = 10
1195             hom(:,2,10,:)   = SPREAD( zu, 2, statistic_regions+1 )
1196             hom(nzb,2,10,:) = 0.0
1197             IF ( data_output_pr(i)(1:1) == '#' )  THEN
1198                dopr_initial_index(i) = 24
1199                hom(:,2,24,:)         = hom(:,2,10,:)
1200                data_output_pr(i)     = data_output_pr(i)(2:)
1201             ENDIF
1202
1203          CASE ( 'l', '#l' )
1204             dopr_index(i)   = 11
1205             hom(:,2,11,:)   = SPREAD( zu, 2, statistic_regions+1 )
1206             hom(nzb,2,11,:) = 0.0
1207             IF ( data_output_pr(i)(1:1) == '#' )  THEN
1208                dopr_initial_index(i) = 25
1209                hom(:,2,25,:)         = hom(:,2,11,:)
1210                data_output_pr(i)     = data_output_pr(i)(2:)
1211             ENDIF
1212
1213          CASE ( 'w"u"' )
1214             dopr_index(i) = 12
1215             hom(:,2,12,:) = SPREAD( zw, 2, statistic_regions+1 )
1216             IF ( prandtl_layer )  hom(nzb,2,12,:) = zu(1)
1217
1218          CASE ( 'w*u*' )
1219             dopr_index(i) = 13
1220             hom(:,2,13,:) = SPREAD( zw, 2, statistic_regions+1 )
1221
1222          CASE ( 'w"v"' )
1223             dopr_index(i) = 14
1224             hom(:,2,14,:) = SPREAD( zw, 2, statistic_regions+1 )
1225             IF ( prandtl_layer )  hom(nzb,2,14,:) = zu(1)
1226
1227          CASE ( 'w*v*' )
1228             dopr_index(i) = 15
1229             hom(:,2,15,:) = SPREAD( zw, 2, statistic_regions+1 )
1230
1231          CASE ( 'w"pt"' )
1232             dopr_index(i) = 16
1233             hom(:,2,16,:) = SPREAD( zw, 2, statistic_regions+1 )
1234
1235          CASE ( 'w*pt*' )
1236             dopr_index(i) = 17
1237             hom(:,2,17,:) = SPREAD( zw, 2, statistic_regions+1 )
1238
1239          CASE ( 'wpt' )
1240             dopr_index(i) = 18
1241             hom(:,2,18,:) = SPREAD( zw, 2, statistic_regions+1 )
1242
1243          CASE ( 'wu' )
1244             dopr_index(i) = 19
1245             hom(:,2,19,:) = SPREAD( zw, 2, statistic_regions+1 )
1246             IF ( prandtl_layer )  hom(nzb,2,19,:) = zu(1)
1247
1248          CASE ( 'wv' )
1249             dopr_index(i) = 20
1250             hom(:,2,20,:) = SPREAD( zw, 2, statistic_regions+1 )
1251             IF ( prandtl_layer )  hom(nzb,2,20,:) = zu(1)
1252
1253          CASE ( 'w*pt*BC' )
1254             dopr_index(i) = 21
1255             hom(:,2,21,:) = SPREAD( zw, 2, statistic_regions+1 )
1256
1257          CASE ( 'wptBC' )
1258             dopr_index(i) = 22
1259             hom(:,2,22,:) = SPREAD( zw, 2, statistic_regions+1 )
1260
1261          CASE ( 'u*2' )
1262             dopr_index(i) = 30
1263             hom(:,2,30,:) = SPREAD( zu, 2, statistic_regions+1 )
1264
1265          CASE ( 'v*2' )
1266             dopr_index(i) = 31
1267             hom(:,2,31,:) = SPREAD( zu, 2, statistic_regions+1 )
1268
1269          CASE ( 'w*2' )
1270             dopr_index(i) = 32
1271             hom(:,2,32,:) = SPREAD( zw, 2, statistic_regions+1 )
1272
1273          CASE ( 'pt*2' )
1274             dopr_index(i) = 33
1275             hom(:,2,33,:) = SPREAD( zu, 2, statistic_regions+1 )
1276
1277          CASE ( 'e*' )
1278             dopr_index(i) = 34
1279             hom(:,2,34,:) = SPREAD( zu, 2, statistic_regions+1 )
1280
1281          CASE ( 'w*2pt*' )
1282             dopr_index(i) = 35
1283             hom(:,2,35,:) = SPREAD( zw, 2, statistic_regions+1 )
1284
1285          CASE ( 'w*pt*2' )
1286             dopr_index(i) = 36
1287             hom(:,2,36,:) = SPREAD( zw, 2, statistic_regions+1 )
1288
1289          CASE ( 'w*e*' )
1290             dopr_index(i) = 37
1291             hom(:,2,37,:) = SPREAD( zw, 2, statistic_regions+1 )
1292
1293          CASE ( 'w*3' )
1294             dopr_index(i) = 38
1295             hom(:,2,38,:) = SPREAD( zw, 2, statistic_regions+1 )
1296
1297          CASE ( 'Sw' )
1298             dopr_index(i) = 39
1299             hom(:,2,39,:) = SPREAD( zw, 2, statistic_regions+1 )
1300
1301          CASE ( 'q', '#q' )
1302             IF ( .NOT. cloud_physics )  THEN
1303                IF ( myid == 0 )  THEN
1304                   PRINT*, '+++ check_parameters:  data_output_pr = ', &
1305                           data_output_pr(i),                          &
1306                           '    is not implemented for cloud_physics = FALSE'
1307                ENDIF
1308                CALL local_stop
1309             ELSE
1310                dopr_index(i) = 41
1311                hom(:,2,41,:)  = SPREAD( zu, 2, statistic_regions+1 )
1312                IF ( data_output_pr(i)(1:1) == '#' )  THEN
1313                   dopr_initial_index(i) = 26
1314                   hom(:,2,26,:)         = SPREAD( zu, 2, statistic_regions+1 )
1315                   hom(nzb,2,26,:)       = 0.0    ! weil zu(nzb) negativ ist
1316                   data_output_pr(i)     = data_output_pr(i)(2:)
1317                ENDIF
1318             ENDIF
1319
1320          CASE ( 's', '#s' )
1321             IF ( .NOT. passive_scalar )  THEN
1322                IF ( myid == 0 )  THEN
1323                   PRINT*, '+++ check_parameters:  data_output_pr = ', &
1324                           data_output_pr(i),                          &
1325                           '    is not implemented for passive_scalar = FALSE'
1326                ENDIF
1327                CALL local_stop
1328             ELSE
1329                dopr_index(i) = 41
1330                hom(:,2,41,:)  = SPREAD( zu, 2, statistic_regions+1 )
1331                IF ( data_output_pr(i)(1:1) == '#' )  THEN
1332                   dopr_initial_index(i) = 26
1333                   hom(:,2,26,:)         = SPREAD( zu, 2, statistic_regions+1 )
1334                   hom(nzb,2,26,:)       = 0.0    ! weil zu(nzb) negativ ist
1335                   data_output_pr(i)     = data_output_pr(i)(2:)
1336                ENDIF
1337             ENDIF
1338
1339          CASE ( 'qv', '#qv' )
1340             IF ( .NOT. cloud_physics ) THEN
1341                dopr_index(i) = 41
1342                hom(:,2,41,:)  = SPREAD( zu, 2, statistic_regions+1 )
1343                IF ( data_output_pr(i)(1:1) == '#' )  THEN
1344                   dopr_initial_index(i) = 26
1345                   hom(:,2,26,:)         = SPREAD( zu, 2, statistic_regions+1 )
1346                   hom(nzb,2,26,:)       = 0.0    ! weil zu(nzb) negativ ist
1347                   data_output_pr(i)     = data_output_pr(i)(2:)
1348                ENDIF
1349             ELSE
1350                dopr_index(i) = 42
1351                hom(:,2,42,:)  = SPREAD( zu, 2, statistic_regions+1 )
1352                IF ( data_output_pr(i)(1:1) == '#' )  THEN
1353                   dopr_initial_index(i) = 27
1354                   hom(:,2,27,:)         = SPREAD( zu, 2, statistic_regions+1 )
1355                   hom(nzb,2,27,:)       = 0.0    ! weil zu(nzb) negativ ist
1356                   data_output_pr(i)     = data_output_pr(i)(2:)
1357                ENDIF
1358             ENDIF
1359
1360          CASE ( 'lpt', '#lpt' )
1361             IF ( .NOT. cloud_physics ) THEN
1362                IF ( myid == 0 )  THEN
1363                   PRINT*, '+++ check_parameters:  data_output_pr = ', &
1364                           data_output_pr(i),                          &
1365                           '    is not implemented for cloud_physics = FALSE'
1366                ENDIF
1367                CALL local_stop
1368             ELSE
1369                dopr_index(i) = 4
1370                hom(:,2,4,:)  = SPREAD( zu, 2, statistic_regions+1 )
1371                IF ( data_output_pr(i)(1:1) == '#' )  THEN
1372                   dopr_initial_index(i) = 7
1373                   hom(:,2,7,:)          = SPREAD( zu, 2, statistic_regions+1 )
1374                   hom(nzb,2,7,:)        = 0.0    ! weil zu(nzb) negativ ist
1375                   data_output_pr(i)     = data_output_pr(i)(2:)
1376                ENDIF
1377             ENDIF
1378
1379          CASE ( 'vpt', '#vpt' )
1380             dopr_index(i) = 44
1381             hom(:,2,44,:)  = SPREAD( zu, 2, statistic_regions+1 )
1382             IF ( data_output_pr(i)(1:1) == '#' )  THEN
1383                dopr_initial_index(i) = 29
1384                hom(:,2,29,:)         = SPREAD( zu, 2, statistic_regions+1 )
1385                hom(nzb,2,29,:)       = 0.0    ! weil zu(nzb) negativ ist
1386                data_output_pr(i)     = data_output_pr(i)(2:)
1387             ENDIF
1388
1389          CASE ( 'w"vpt"' )
1390             dopr_index(i) = 45
1391             hom(:,2,45,:) = SPREAD( zw, 2, statistic_regions+1 )
1392
1393          CASE ( 'w*vpt*' )
1394             dopr_index(i) = 46
1395             hom(:,2,46,:) = SPREAD( zw, 2, statistic_regions+1 )
1396
1397          CASE ( 'wvpt' )
1398             dopr_index(i) = 47
1399             hom(:,2,47,:) = SPREAD( zw, 2, statistic_regions+1 )
1400
1401          CASE ( 'w"q"' )
1402             IF ( .NOT. cloud_physics ) THEN
1403                IF ( myid == 0 )  THEN
1404                   PRINT*, '+++ check_parameters:  data_output_pr = ', &
1405                           data_output_pr(i),                          &
1406                           '    is not implemented for cloud_physics = FALSE'
1407                ENDIF
1408                CALL local_stop
1409             ELSE
1410                dopr_index(i) = 48
1411                hom(:,2,48,:) = SPREAD( zw, 2, statistic_regions+1 )
1412             ENDIF
1413
1414          CASE ( 'w*q*' )
1415             IF ( .NOT. cloud_physics ) THEN
1416                IF ( myid == 0 )  THEN
1417                   PRINT*, '+++ check_parameters:  data_output_pr = ', &
1418                           data_output_pr(i),                          &
1419                           '    is not implemented for cloud_physics = FALSE'
1420                ENDIF
1421                CALL local_stop
1422             ELSE
1423                dopr_index(i) = 49
1424                hom(:,2,49,:) = SPREAD( zw, 2, statistic_regions+1 )
1425             ENDIF
1426
1427          CASE ( 'wq' )
1428             IF ( .NOT. cloud_physics ) THEN
1429                IF ( myid == 0 )  THEN
1430                   PRINT*, '+++ check_parameters:  data_output_pr = ', &
1431                           data_output_pr(i),                          &
1432                           '    is not implemented for cloud_physics = FALSE'
1433                ENDIF
1434                CALL local_stop
1435             ELSE
1436                dopr_index(i) = 50
1437                hom(:,2,50,:) = SPREAD( zw, 2, statistic_regions+1 )
1438             ENDIF
1439
1440          CASE ( 'w"s"' )
1441             IF ( .NOT. passive_scalar ) THEN
1442                IF ( myid == 0 )  THEN
1443                   PRINT*, '+++ check_parameters:  data_output_pr = ', &
1444                           data_output_pr(i),                          &
1445                           '    is not implemented for passive_scalar = FALSE'
1446                ENDIF
1447                CALL local_stop
1448             ELSE
1449                dopr_index(i) = 48
1450                hom(:,2,48,:) = SPREAD( zw, 2, statistic_regions+1 )
1451             ENDIF
1452
1453          CASE ( 'w*s*' )
1454             IF ( .NOT. passive_scalar ) THEN
1455                IF ( myid == 0 )  THEN
1456                   PRINT*, '+++ check_parameters:  data_output_pr = ', &
1457                           data_output_pr(i),                          &
1458                           '    is not implemented for passive_scalar = FALSE'
1459                ENDIF
1460                CALL local_stop
1461             ELSE
1462                dopr_index(i) = 49
1463                hom(:,2,49,:) = SPREAD( zw, 2, statistic_regions+1 )
1464             ENDIF
1465
1466          CASE ( 'ws' )
1467             IF ( .NOT. passive_scalar ) THEN
1468                IF ( myid == 0 )  THEN
1469                   PRINT*, '+++ check_parameters:  data_output_pr = ', &
1470                           data_output_pr(i),                          &
1471                           '    is not implemented for passive_scalar = FALSE'
1472                ENDIF
1473                CALL local_stop
1474             ELSE
1475                dopr_index(i) = 50
1476                hom(:,2,50,:) = SPREAD( zw, 2, statistic_regions+1 )
1477             ENDIF
1478
1479          CASE ( 'w"qv"' )
1480             IF ( moisture  .AND.  .NOT. cloud_physics ) &
1481             THEN
1482                dopr_index(i) = 48
1483                hom(:,2,48,:) = SPREAD( zw, 2, statistic_regions+1 )
1484             ELSEIF( moisture .AND. cloud_physics ) THEN
1485                dopr_index(i) = 51
1486                hom(:,2,51,:) = SPREAD( zw, 2, statistic_regions+1 )
1487             ELSE
1488                IF ( myid == 0 )  THEN
1489                   PRINT*, '+++ check_parameters:  data_output_pr = ', &
1490                           data_output_pr(i),                          &
1491                           '    is not implemented for cloud_physics = FALSE', &
1492                           '    and                    moisture      = FALSE'
1493                ENDIF
1494                CALL local_stop                   
1495             ENDIF
1496
1497          CASE ( 'w*qv*' )
1498             IF ( moisture  .AND.  .NOT. cloud_physics ) &
1499             THEN
1500                dopr_index(i) = 49
1501                hom(:,2,49,:) = SPREAD( zw, 2, statistic_regions+1 )
1502             ELSEIF( moisture .AND. cloud_physics ) THEN
1503                dopr_index(i) = 52
1504                hom(:,2,52,:) = SPREAD( zw, 2, statistic_regions+1 )
1505             ELSE
1506                IF ( myid == 0 )  THEN
1507                   PRINT*, '+++ check_parameters:  data_output_pr = ',         &
1508                           data_output_pr(i),                                  &
1509                           '    is not implemented for cloud_physics = FALSE', &
1510                           '                       and moisture      = FALSE'
1511                ENDIF
1512                CALL local_stop                   
1513             ENDIF
1514
1515          CASE ( 'wqv' )
1516             IF ( moisture  .AND.  .NOT. cloud_physics ) &
1517             THEN
1518                dopr_index(i) = 50
1519                hom(:,2,50,:) = SPREAD( zw, 2, statistic_regions+1 )
1520             ELSEIF( moisture .AND. cloud_physics ) THEN
1521                dopr_index(i) = 53
1522                hom(:,2,53,:) = SPREAD( zw, 2, statistic_regions+1 )
1523             ELSE
1524                IF ( myid == 0 )  THEN
1525                   PRINT*, '+++ check_parameters:  data_output_pr = ',         &
1526                           data_output_pr(i),                                  &
1527                           '    is not implemented for cloud_physics = FALSE', &
1528                           '                       and moisture      = FALSE'
1529                ENDIF
1530                CALL local_stop                   
1531             ENDIF
1532
1533          CASE ( 'ql' )
1534             IF ( .NOT. cloud_physics  .AND.  .NOT. cloud_droplets )  THEN
1535                IF ( myid == 0 )  THEN
1536                   PRINT*, '+++ check_parameters:  data_output_pr = ', &
1537                           data_output_pr(i),                          &
1538                           '    is not implemented for cloud_physics = FALSE'
1539                ENDIF
1540                CALL local_stop
1541             ELSE
1542                dopr_index(i) = 54
1543                hom(:,2,54,:)  = SPREAD( zu, 2, statistic_regions+1 )
1544             ENDIF
1545
1546          CASE ( 'w*u*u*/dz' )
1547             dopr_index(i) = 55
1548             hom(:,2,55,:) = SPREAD( zu, 2, statistic_regions+1 )
1549
1550          CASE ( 'w*p*/dz' )
1551             dopr_index(i) = 56
1552             hom(:,2,56,:) = SPREAD( zu, 2, statistic_regions+1 )
1553
1554          CASE ( 'w"e/dz' )
1555             dopr_index(i) = 57
1556             hom(:,2,57,:) = SPREAD( zu, 2, statistic_regions+1 )
1557
1558          CASE ( 'u"pt"' )
1559             dopr_index(i) = 58
1560             hom(:,2,58,:) = SPREAD( zu, 2, statistic_regions+1 )
1561
1562          CASE ( 'u*pt*' )
1563             dopr_index(i) = 59
1564             hom(:,2,59,:) = SPREAD( zu, 2, statistic_regions+1 )
1565
1566          CASE ( 'upt_t' )
1567             dopr_index(i) = 60
1568             hom(:,2,60,:) = SPREAD( zu, 2, statistic_regions+1 )
1569
1570          CASE ( 'v"pt"' )
1571             dopr_index(i) = 61
1572             hom(:,2,61,:) = SPREAD( zu, 2, statistic_regions+1 )
1573             
1574          CASE ( 'v*pt*' )
1575             dopr_index(i) = 62
1576             hom(:,2,62,:) = SPREAD( zu, 2, statistic_regions+1 )
1577
1578          CASE ( 'vpt_t' )
1579             dopr_index(i) = 63
1580             hom(:,2,63,:) = SPREAD( zu, 2, statistic_regions+1 )
1581
1582
1583          CASE DEFAULT
1584             IF ( myid == 0 )  THEN
1585                PRINT*, '+++ check_parameters:  unknown output profile:  ', &
1586                        'data_output_pr = ', data_output_pr(i)
1587             ENDIF
1588             CALL local_stop
1589
1590       END SELECT
1591!
1592!--    Check to which of the predefined coordinate systems the profile belongs
1593       DO  k = 1, crmax
1594          IF ( INDEX( cross_profiles(k), ' '//TRIM( data_output_pr(i) )//' ' ) &
1595               /=0 ) &
1596          THEN
1597             dopr_crossindex(i) = k
1598             EXIT
1599          ENDIF
1600       ENDDO
1601!
1602!--    Generate the text for the labels of the PROFIL output file. "-characters
1603!--    must be substituted, otherwise PROFIL would interpret them as TeX
1604!--    control characters
1605       dopr_label(i) = data_output_pr(i)
1606       position = INDEX( dopr_label(i) , '"' )
1607       DO WHILE ( position /= 0 )
1608          dopr_label(i)(position:position) = ''''
1609          position = INDEX( dopr_label(i) , '"' )
1610       ENDDO
1611
1612    ENDDO
1613
1614!
1615!-- y-value range of the coordinate system (PROFIL).
1616!-- x-value range determined in plot_1d.
1617    cross_uymin = 0.0
1618    IF ( z_max_do1d == -1.0 )  THEN
1619       cross_uymax = zu(nzt+1)
1620    ELSEIF ( z_max_do1d < zu(nzb+1)  .OR.  z_max_do1d > zu(nzt+1) )  THEN
1621       IF ( myid == 0 )  PRINT*, '+++ check_parameters:  z_max_do1d=',  &
1622                                 z_max_do1d,' must be >= ', zu(nzb+1),  &
1623                                 ' or <= ', zu(nzt+1)
1624       CALL local_stop
1625    ELSE
1626       cross_uymax = z_max_do1d
1627    ENDIF
1628
1629!
1630!-- Check whether the chosen normalizing factor for the coordinate systems is
1631!-- permissible
1632    DO  i = 1, crmax
1633       SELECT CASE ( TRIM( cross_normalized_x(i) ) )  ! TRIM required on IBM
1634
1635          CASE ( '', 'wpt0', 'ws2', 'tsw2', 'ws3', 'ws2tsw', 'wstsw2' )
1636             j = 0
1637
1638          CASE DEFAULT
1639             IF ( myid == 0 )  THEN
1640                PRINT*, '+++ check_parameters: unknown normalize method'
1641                PRINT*, '    cross_normalized_x="',cross_normalized_x(i),'"'
1642             ENDIF
1643             CALL local_stop
1644
1645       END SELECT
1646       SELECT CASE ( TRIM( cross_normalized_y(i) ) )  ! TRIM required on IBM
1647
1648          CASE ( '', 'z_i' )
1649             j = 0
1650
1651          CASE DEFAULT
1652             IF ( myid == 0 )  THEN
1653                PRINT*, '+++ check_parameters: unknown normalize method'
1654                PRINT*, '    cross_normalized_y="',cross_normalized_y(i),'"'
1655             ENDIF
1656             CALL local_stop
1657
1658       END SELECT
1659    ENDDO
1660!
1661!-- Check normalized y-value range of the coordinate system (PROFIL)
1662    IF ( z_max_do1d_normalized /= -1.0  .AND.  z_max_do1d_normalized <= 0.0 ) &
1663    THEN
1664       IF ( myid == 0 )  PRINT*,'+++ check_parameters:  z_max_do1d_normalize', &
1665                                'd=', z_max_do1d_normalized, ' must be >= 0.0'
1666       CALL local_stop
1667    ENDIF
1668
1669!
1670!-- Determine parameters for time series output and check whether permissible
1671    i = 0
1672    DO  WHILE ( data_output_ts(i+1) /= '          '  .AND.  i+1 <= 100 )
1673
1674       dots_n = dots_n + 1
1675       i = i + 1
1676!
1677!--    Check whether time series is permissible and determine internal number
1678       SELECT CASE ( TRIM( data_output_ts(i) ) )
1679
1680          CASE ( 'E' )
1681             dots_index(i) = 1
1682          CASE ( 'E*' )
1683             dots_index(i) = 2
1684          CASE ( 'dt' )
1685             dots_index(i) = 3
1686          CASE ( 'u*' )
1687             dots_index(i) = 4
1688          CASE ( 'th*' )
1689             dots_index(i) = 5
1690          CASE ( 'umax' )
1691             dots_index(i) = 6
1692          CASE ( 'vmax' )
1693             dots_index(i) = 7
1694          CASE ( 'wmax' )
1695             dots_index(i) = 8
1696          CASE ( 'div_new' )
1697             dots_index(i) = 9
1698          CASE ( 'div_old' )
1699             dots_index(i) = 10
1700          CASE ( 'z_i_wpt' )
1701             dots_index(i) = 11
1702          CASE ( 'z_i_pt' )
1703             dots_index(i) = 12
1704          CASE ( 'w*' )
1705             dots_index(i) = 13
1706          CASE ( 'w"pt"0' )
1707             dots_index(i) = 14
1708          CASE ( 'w"pt"' )
1709             dots_index(i) = 15
1710          CASE ( 'wpt' )
1711             dots_index(i) = 16
1712          CASE ( 'pt(0)' )
1713             dots_index(i) = 17
1714          CASE ( 'pt(zp)' )
1715             dots_index(i) = 18
1716          CASE ( 'splptx'  )
1717             dots_index(i) = 19
1718          CASE ( 'splpty'  )
1719             dots_index(i) = 20
1720          CASE ( 'splptz'  )
1721             dots_index(i) = 21
1722          CASE ( 'L'       )
1723             dots_index(i) = 22
1724
1725          CASE DEFAULT
1726             IF ( myid == 0 )  THEN
1727                PRINT*, '+++ check_parameters:  unknown time series:  ', &
1728                             'data_output_ts = ',&
1729                        data_output_ts(i)
1730             ENDIF
1731             CALL local_stop
1732
1733       END SELECT
1734
1735!
1736!--    Check, to which predefined coordinate system the time series belongs, and
1737!--    store corresponding internal number. Furthermore determine, how many and
1738!--    which graphs are being drawn into the corresponding system
1739       DO  k = 1, crmax
1740          IF ( INDEX( cross_ts_profiles(k), ' ' // TRIM( data_output_ts(i) ) &
1741                      // ' ' ) /=0 )  &
1742          THEN
1743             dots_crossindex(i) = k
1744             cross_ts_number_count(k) = cross_ts_number_count(k) + 1
1745             cross_ts_numbers(cross_ts_number_count(k),k) = dots_index(i)
1746             EXIT
1747          ENDIF
1748       ENDDO
1749
1750    ENDDO
1751
1752!
1753!-- Append user-defined data output variables to the standard data output
1754    IF ( data_output_user(1) /= ' ' )  THEN
1755       i = 1
1756       DO  WHILE ( data_output(i) /= ' '  .AND.  i <= 100 )
1757          i = i + 1
1758       ENDDO
1759       j = 1
1760       DO  WHILE ( data_output_user(j) /= ' '  .AND.  j <= 100 )
1761          IF ( i > 100 )  THEN
1762             IF ( myid == 0 )  THEN
1763                PRINT*, '+++ check_parameters: number of output quantitities', &
1764                             ' given by data_output and data_output_user'
1765                PRINT*, '                      exceeds the limit of 100'
1766             ENDIF
1767             CALL local_stop
1768          ENDIF
1769          data_output(i) = data_output_user(j)
1770          i = i + 1
1771          j = j + 1
1772       ENDDO
1773    ENDIF
1774
1775!
1776!-- Check and set steering parameters for 2d/3d data output and averaging
1777    i   = 1
1778    DO  WHILE ( data_output(i) /= ' '  .AND.  i <= 100 )
1779!
1780!--    Check for data averaging
1781       ilen = LEN_TRIM( data_output(i) )
1782       j = 0                                                 ! no data averaging
1783       IF ( ilen > 3 )  THEN
1784          IF ( data_output(i)(ilen-2:ilen) == '_av' )  THEN
1785             j = 1                                           ! data averaging
1786             data_output(i) = data_output(i)(1:ilen-3)
1787          ENDIF
1788       ENDIF
1789!
1790!--    Check for cross section or volume data
1791       ilen = LEN_TRIM( data_output(i) )
1792       k = 0                                                   ! 3d data
1793       var = data_output(i)(1:ilen)
1794       IF ( ilen > 3 )  THEN
1795          IF ( data_output(i)(ilen-2:ilen) == '_xy'  .OR. &
1796               data_output(i)(ilen-2:ilen) == '_xz'  .OR. &
1797               data_output(i)(ilen-2:ilen) == '_yz' )  THEN
1798             k = 1                                             ! 2d data
1799             var = data_output(i)(1:ilen-3)
1800          ENDIF
1801       ENDIF
1802!
1803!--    Check for allowed value and set units
1804       SELECT CASE ( TRIM( var ) )
1805
1806          CASE ( 'e' )
1807             IF ( constant_diffusion )  THEN
1808                IF ( myid == 0 )  THEN
1809                   PRINT*, '+++ check_parameters: output of "', TRIM( var ), &
1810                                '" requires constant_diffusion = .FALSE.'
1811                ENDIF
1812                CALL local_stop
1813             ENDIF
1814             unit = 'm2/s2'
1815
1816          CASE ( 'pc', 'pr' )
1817             IF ( .NOT. particle_advection )  THEN
1818                IF ( myid == 0 )  THEN
1819                   PRINT*, '+++ check_parameters: output of "', TRIM( var ), &
1820                                '" requires particle package'
1821                   PRINT*, '                      (mrun-option "-p particles")'
1822                ENDIF
1823                CALL local_stop
1824             ENDIF
1825             IF ( TRIM( var ) == 'pc' )  unit = 'number'
1826             IF ( TRIM( var ) == 'pr' )  unit = 'm'
1827
1828          CASE ( 'q', 'vpt' )
1829             IF ( .NOT. moisture )  THEN
1830                IF ( myid == 0 )  THEN
1831                   PRINT*, '+++ check_parameters: output of "', TRIM( var ), &
1832                                '" requires moisture = .TRUE.'
1833                ENDIF
1834                CALL local_stop
1835             ENDIF
1836             IF ( TRIM( var ) == 'q'   )  unit = 'kg/kg'
1837             IF ( TRIM( var ) == 'vpt' )  unit = 'K'
1838
1839          CASE ( 'ql' )
1840             IF ( .NOT. ( cloud_physics  .OR.  cloud_droplets ) )  THEN
1841                IF ( myid == 0 )  THEN
1842                   PRINT*, '+++ check_parameters: output of "', TRIM( var ), &
1843                                '" requires cloud_physics = .TRUE.'
1844                   PRINT*, '                      or cloud_droplets = .TRUE.'
1845                ENDIF
1846                CALL local_stop
1847             ENDIF
1848             unit = 'kg/kg'
1849
1850          CASE ( 'ql_c', 'ql_v', 'ql_vp' )
1851             IF ( .NOT. cloud_droplets )  THEN
1852                IF ( myid == 0 )  THEN
1853                   PRINT*, '+++ check_parameters: output of "', TRIM( var ), &
1854                                '" requires cloud_droplets = .TRUE.'
1855                ENDIF
1856                CALL local_stop
1857             ENDIF
1858             IF ( TRIM( var ) == 'ql_c'  )  unit = 'kg/kg'
1859             IF ( TRIM( var ) == 'ql_v'  )  unit = 'm3'
1860             IF ( TRIM( var ) == 'ql_vp' )  unit = 'none'
1861
1862          CASE ( 'qv' )
1863             IF ( .NOT. cloud_physics )  THEN
1864                IF ( myid == 0 )  THEN
1865                   PRINT*, '+++ check_parameters: output of "', TRIM( var ), &
1866                                '" requires cloud_physics = .TRUE.'
1867                ENDIF
1868                CALL local_stop
1869             ENDIF
1870             unit = 'kg/kg'
1871
1872          CASE ( 's' )
1873             IF ( .NOT. passive_scalar )  THEN
1874                IF ( myid == 0 )  THEN
1875                   PRINT*, '+++ check_parameters: output of "', TRIM( var ), &
1876                                '" requires passive_scalar = .TRUE.'
1877                ENDIF
1878                CALL local_stop
1879             ENDIF
1880             unit = 'conc'
1881
1882          CASE ( 'u*', 't*', 'lwp*' )
1883             IF ( k == 0  .OR.  data_output(i)(ilen-2:ilen) /= '_xy' )  THEN
1884                IF ( myid == 0 )  THEN
1885                   PRINT*, '+++ check_parameters:  illegal value for data_',&
1886                                'output: "', TRIM( var ), '" is only allowed'
1887                   PRINT*, '                       for horizontal cross section'
1888                ENDIF
1889                CALL local_stop
1890             ENDIF
1891             IF ( TRIM( var ) == 'lwp*'  .AND.  .NOT. cloud_physics )  THEN
1892                IF ( myid == 0 )  THEN
1893                   PRINT*, '+++ check_parameters: output of "', TRIM( var ), &
1894                                '" requires cloud_physics = .TRUE.'
1895                ENDIF
1896                CALL local_stop
1897             ENDIF
1898             IF ( TRIM( var ) == 'u*'   )  unit = 'm/s'
1899             IF ( TRIM( var ) == 't*'   )  unit = 'K'
1900             IF ( TRIM( var ) == 'lwp*' )  unit = 'kg/kg*m'
1901
1902          CASE ( 'p', 'pt', 'u', 'v', 'w' )
1903             IF ( TRIM( var ) == 'p'  )  unit = 'Pa'
1904             IF ( TRIM( var ) == 'pt' )  unit = 'K'
1905             IF ( TRIM( var ) == 'u'  )  unit = 'm/s'
1906             IF ( TRIM( var ) == 'v'  )  unit = 'm/s'
1907             IF ( TRIM( var ) == 'w'  )  unit = 'm/s'
1908             CONTINUE
1909
1910          CASE DEFAULT
1911             CALL user_check_data_output( var, unit )
1912
1913             IF ( unit == 'illegal' )  THEN
1914                IF ( myid == 0 )  THEN
1915                   IF ( data_output_user(1) /= ' ' )  THEN
1916                      PRINT*, '+++ check_parameters:  illegal value for data_',&
1917                                   'output or data_output_user: "',            &
1918                                   TRIM( data_output(i) ), '"'
1919                   ELSE
1920                      PRINT*, '+++ check_parameters:  illegal value for data_',&
1921                                   'output: "', TRIM( data_output(i) ), '"'
1922                   ENDIF
1923                ENDIF
1924                CALL local_stop
1925             ENDIF
1926
1927       END SELECT
1928!
1929!--    Set the internal steering parameters appropriately
1930       IF ( k == 0 )  THEN
1931          do3d_no(j)              = do3d_no(j) + 1
1932          do3d(j,do3d_no(j))      = data_output(i)
1933          do3d_unit(j,do3d_no(j)) = unit
1934       ELSE
1935          do2d_no(j)              = do2d_no(j) + 1
1936          do2d(j,do2d_no(j))      = data_output(i)
1937          do2d_unit(j,do2d_no(j)) = unit
1938          IF ( data_output(i)(ilen-2:ilen) == '_xy' )  THEN
1939             data_output_xy(j) = .TRUE.
1940          ENDIF
1941          IF ( data_output(i)(ilen-2:ilen) == '_xz' )  THEN
1942             data_output_xz(j) = .TRUE.
1943          ENDIF
1944          IF ( data_output(i)(ilen-2:ilen) == '_yz' )  THEN
1945             data_output_yz(j) = .TRUE.
1946          ENDIF
1947       ENDIF
1948
1949       IF ( j == 1 )  THEN
1950!
1951!--       Check, if variable is already subject to averaging
1952          found = .FALSE.
1953          DO  k = 1, doav_n
1954             IF ( TRIM( doav(k) ) == TRIM( var ) )  found = .TRUE.
1955          ENDDO
1956
1957          IF ( .NOT. found )  THEN
1958             doav_n = doav_n + 1
1959             doav(doav_n) = var
1960          ENDIF
1961       ENDIF
1962
1963       i = i + 1
1964    ENDDO
1965
1966!
1967!-- Store sectional planes in one shared array
1968    section(:,1) = section_xy
1969    section(:,2) = section_xz
1970    section(:,3) = section_yz
1971
1972!
1973!-- Upper plot limit (grid point value) for 1D profiles
1974    IF ( z_max_do1d == -1.0 )  THEN
1975       nz_do1d = nzt+1
1976    ELSE
1977       DO  k = nzb+1, nzt+1
1978          nz_do1d = k
1979          IF ( zw(k) > z_max_do1d )  EXIT
1980       ENDDO
1981    ENDIF
1982
1983!
1984!-- Upper plot limit for 2D vertical sections
1985    IF ( z_max_do2d == -1.0 )  z_max_do2d = zu(nzt)
1986    IF ( z_max_do2d < zu(nzb+1)  .OR.  z_max_do2d > zu(nzt) )  THEN
1987       IF ( myid == 0 )  PRINT*, '+++ check_parameters:  z_max_do2d=',        &
1988                                 z_max_do2d, ' must be >= ', zu(nzb+1),       &
1989                                 '(zu(nzb+1)) and <= ', zu(nzt), ' (zu(nzt))'
1990       CALL local_stop
1991    ENDIF
1992
1993!
1994!-- Upper plot limit for 3D arrays
1995    IF ( nz_do3d == -9999 )  nz_do3d = nzt + 1
1996
1997!
1998!-- Determine and check accuracy for compressed 3D plot output
1999    IF ( do3d_compress )  THEN
2000!
2001!--    Compression only permissible on T3E machines
2002       IF ( host(1:3) /= 't3e' )  THEN
2003          IF ( myid == 0 )  THEN
2004             PRINT*, '+++ check_parameters: do3d_compress = .TRUE. not allow', &
2005                          'ed on host "', TRIM( host ), '"'
2006          ENDIF
2007          CALL local_stop
2008       ENDIF
2009
2010       i = 1
2011       DO  WHILE ( do3d_comp_prec(i) /= ' ' )
2012
2013          ilen = LEN_TRIM( do3d_comp_prec(i) )
2014          IF ( LLT( do3d_comp_prec(i)(ilen:ilen), '0' ) .OR. &
2015               LGT( do3d_comp_prec(i)(ilen:ilen), '9' ) )  THEN
2016             IF ( myid == 0 )  THEN
2017                PRINT*, '+++ check_parameters: illegal precision: ', &
2018                        'do3d_comp_prec(', i, ')="', TRIM(do3d_comp_prec(i)),'"'
2019             ENDIF
2020             CALL local_stop
2021          ENDIF
2022
2023          prec = IACHAR( do3d_comp_prec(i)(ilen:ilen) ) - IACHAR( '0' )
2024          var = do3d_comp_prec(i)(1:ilen-1)
2025
2026          SELECT CASE ( var )
2027
2028             CASE ( 'u' )
2029                j = 1
2030             CASE ( 'v' )
2031                j = 2
2032             CASE ( 'w' )
2033                j = 3
2034             CASE ( 'p' )
2035                j = 4
2036             CASE ( 'pt' )
2037                j = 5
2038
2039             CASE DEFAULT
2040                IF ( myid == 0 )  THEN
2041                   PRINT*, '+++ check_parameters: unknown variable in ', &
2042                               'assignment'
2043                   PRINT*, '    do3d_comp_prec(', i, ')="', &
2044                           TRIM( do3d_comp_prec(i) ),'"'
2045                ENDIF
2046                CALL local_stop               
2047
2048          END SELECT
2049
2050          plot_3d_precision(j)%precision = prec
2051          i = i + 1
2052
2053       ENDDO
2054    ENDIF
2055
2056!
2057!-- Check the data output format(s)
2058    IF ( data_output_format(1) == ' ' )  THEN
2059!
2060!--    Default value
2061       netcdf_output = .TRUE.
2062    ELSE
2063       i = 1
2064       DO  WHILE ( data_output_format(i) /= ' ' )
2065
2066          SELECT CASE ( data_output_format(i) )
2067
2068             CASE ( 'netcdf' )
2069                netcdf_output = .TRUE.
2070             CASE ( 'iso2d' )
2071                iso2d_output  = .TRUE.
2072             CASE ( 'profil' )
2073                profil_output = .TRUE.
2074             CASE ( 'avs' )
2075                avs_output    = .TRUE.
2076
2077             CASE DEFAULT
2078                IF ( myid == 0 )  THEN
2079                   PRINT*, '+++ check_parameters:'
2080                   PRINT*, '    unknown value for data_output_format "', &
2081                                TRIM( data_output_format(i) ),'"'
2082                ENDIF
2083                CALL local_stop               
2084
2085          END SELECT
2086
2087          i = i + 1
2088          IF ( i > 10 )  EXIT
2089
2090       ENDDO
2091
2092    ENDIF
2093
2094!
2095!-- Check netcdf precison
2096    ldum = .FALSE.
2097    CALL define_netcdf_header( 'ch', ldum, 0 )
2098
2099!
2100!-- Check, whether a constant diffusion coefficient shall be used
2101    IF ( km_constant /= -1.0 )  THEN
2102       IF ( km_constant < 0.0 )  THEN
2103          IF ( myid == 0 )  PRINT*, '+++ check_parameters:  km_constant=', &
2104                                    km_constant, ' < 0.0'
2105          CALL local_stop
2106       ELSE
2107          IF ( prandtl_number < 0.0 )  THEN
2108             IF ( myid == 0 )  PRINT*,'+++ check_parameters:  prandtl_number=',&
2109                                      prandtl_number, ' < 0.0'
2110             CALL local_stop
2111          ENDIF
2112          constant_diffusion = .TRUE.
2113
2114          IF ( prandtl_layer )  THEN
2115             IF ( myid == 0 )  PRINT*, '+++ check_parameters:  prandtl_layer ',&
2116                                       'is not allowed with fixed value of km'
2117             CALL local_stop
2118          ENDIF
2119       ENDIF
2120    ENDIF
2121
2122!
2123!-- In case of non-cyclic lateral boundaries, set the default maximum value
2124!-- for the horizontal diffusivity used within the outflow damping layer,
2125!-- and check/set the width of the damping layer
2126    IF ( bc_lr /= 'cyclic' ) THEN
2127       IF ( km_damp_max == -1.0 )  THEN
2128          km_damp_max = 0.5 * dx
2129       ENDIF
2130       IF ( outflow_damping_width == -1.0 )  THEN
2131          outflow_damping_width = MIN( 20, nx/2 )
2132       ENDIF
2133       IF ( outflow_damping_width <= 0  .OR.  outflow_damping_width > nx )  THEN
2134          IF ( myid == 0 )  PRINT*, '+++ check_parameters: outflow_damping w',&
2135                                    'idth out of range'
2136          CALL local_stop
2137       ENDIF
2138    ENDIF
2139
2140    IF ( bc_ns /= 'cyclic' )  THEN
2141       IF ( km_damp_max == -1.0 )  THEN
2142          km_damp_max = 0.5 * dy
2143       ENDIF
2144       IF ( outflow_damping_width == -1.0 )  THEN
2145          outflow_damping_width = MIN( 20, ny/2 )
2146       ENDIF
2147       IF ( outflow_damping_width <= 0  .OR.  outflow_damping_width > ny )  THEN
2148          IF ( myid == 0 )  PRINT*, '+++ check_parameters: outflow_damping w',&
2149                                    'idth out of range'
2150          CALL local_stop
2151       ENDIF
2152    ENDIF
2153
2154!
2155!-- Check value range for rif
2156    IF ( rif_min >= rif_max )  THEN
2157       IF ( myid == 0 )  PRINT*, '+++ check_parameters:  rif_min=', rif_min, &
2158                                 ' must be less than rif_max=', rif_max
2159       CALL local_stop
2160    ENDIF
2161
2162!
2163!-- Determine upper and lower hight level indices for random perturbations
2164    IF ( disturbance_level_b == -1.0 )  THEN
2165       disturbance_level_b     = zu(nzb+3)
2166       disturbance_level_ind_b = nzb + 3
2167    ELSEIF ( disturbance_level_b < zu(3) )  THEN
2168       IF ( myid == 0 )  PRINT*, '+++ check_parameters:  disturbance_level_b=',&
2169                                 disturbance_level_b, ' must be >= ',zu(3),    &
2170                                 '(zu(3))'
2171       CALL local_stop
2172    ELSEIF ( disturbance_level_b > zu(nzt-2) )  THEN
2173       IF ( myid == 0 )  PRINT*, '+++ check_parameters:  disturbance_level_b=',&
2174                                 disturbance_level_b, ' must be <= ',zu(nzt-2),&
2175                                 '(zu(nzt-2))'
2176       CALL local_stop
2177    ELSE
2178       DO  k = 3, nzt-2
2179          IF ( disturbance_level_b <= zu(k) )  THEN
2180             disturbance_level_ind_b = k
2181             EXIT
2182          ENDIF
2183       ENDDO
2184    ENDIF
2185
2186    IF ( disturbance_level_t == -1.0 )  THEN
2187       disturbance_level_t     = zu(nzt/3)
2188       disturbance_level_ind_t = nzt / 3
2189    ELSEIF ( disturbance_level_t > zu(nzt-2) )  THEN
2190       IF ( myid == 0 )  PRINT*, '+++ check_parameters:  disturbance_level_t=',&
2191                                 disturbance_level_t, ' must be <= ',zu(nzt-2),&
2192                                 '(zu(nzt-2))'
2193       CALL local_stop
2194    ELSEIF ( disturbance_level_t < disturbance_level_b )  THEN
2195       IF ( myid == 0 )  PRINT*, '+++ check_parameters:  disturbance_level_t=',&
2196                                 disturbance_level_t, ' must be >= ',          &
2197                                 'disturbance_level_b=', disturbance_level_b
2198       CALL local_stop
2199    ELSE
2200       DO  k = 3, nzt-2
2201          IF ( disturbance_level_t <= zu(k) )  THEN
2202             disturbance_level_ind_t = k
2203             EXIT
2204          ENDIF
2205       ENDDO
2206    ENDIF
2207
2208!
2209!-- Check again whether the levels determined this way are ok.
2210!-- Error may occur at automatic determination and too few grid points in
2211!-- z-direction.
2212    IF ( disturbance_level_ind_t < disturbance_level_ind_b )  THEN
2213       IF ( myid == 0 )  PRINT*, '+++ check_parameters:  ',                &
2214                                 'disturbance_level_ind_t=',               &
2215                                 disturbance_level_ind_t, ' must be >= ',  &
2216                                 'disturbance_level_ind_b=',               &
2217                                 disturbance_level_ind_b
2218       CALL local_stop
2219    ENDIF
2220
2221!
2222!-- Determine the horizontal index range for random perturbations.
2223!-- In case of non-cyclic horizontal boundaries, no perturbations are imposed
2224!-- near the inflow and the perturbation area is further limited to ...(1)
2225!-- after the initial phase of the flow.
2226    dist_nxl = 0;  dist_nxr = nx
2227    dist_nys = 0;  dist_nyn = ny
2228    IF ( bc_lr /= 'cyclic' )  THEN
2229       IF ( inflow_disturbance_begin == -1 )  THEN
2230          inflow_disturbance_begin = MIN( 10, nx/2 )
2231       ENDIF
2232       IF ( inflow_disturbance_begin < 0  .OR.  inflow_disturbance_begin > nx )&
2233       THEN
2234          IF ( myid == 0 )  PRINT*, '+++ check_parameters: inflow_disturbance',&
2235                                    '_begin out of range'
2236          CALL local_stop
2237       ENDIF
2238       IF ( inflow_disturbance_end == -1 )  THEN
2239          inflow_disturbance_end = MIN( 100, 3*nx/4 )
2240       ENDIF
2241       IF ( inflow_disturbance_end < 0  .OR.  inflow_disturbance_end > nx )    &
2242       THEN
2243          IF ( myid == 0 )  PRINT*, '+++ check_parameters: inflow_disturbance',&
2244                                    '_end out of range'
2245          CALL local_stop
2246       ENDIF
2247    ELSEIF ( bc_ns /= 'cyclic' )  THEN
2248       IF ( inflow_disturbance_begin == -1 )  THEN
2249          inflow_disturbance_begin = MIN( 10, ny/2 )
2250       ENDIF
2251       IF ( inflow_disturbance_begin < 0  .OR.  inflow_disturbance_begin > ny )&
2252       THEN
2253          IF ( myid == 0 )  PRINT*, '+++ check_parameters: inflow_disturbance',&
2254                                    '_begin out of range'
2255          CALL local_stop
2256       ENDIF
2257       IF ( inflow_disturbance_end == -1 )  THEN
2258          inflow_disturbance_end = MIN( 100, 3*ny/4 )
2259       ENDIF
2260       IF ( inflow_disturbance_end < 0  .OR.  inflow_disturbance_end > ny )    &
2261       THEN
2262          IF ( myid == 0 )  PRINT*, '+++ check_parameters: inflow_disturbance',&
2263                                    '_end out of range'
2264          CALL local_stop
2265       ENDIF
2266    ENDIF
2267
2268    IF ( outflow_l )  THEN
2269       dist_nxr    = nx - inflow_disturbance_begin
2270       dist_nxl(1) = nx - inflow_disturbance_end
2271    ELSEIF ( outflow_r )  THEN
2272       dist_nxl    = inflow_disturbance_begin
2273       dist_nxr(1) = inflow_disturbance_end
2274    ELSEIF ( outflow_s )  THEN
2275       dist_nyn    = ny - inflow_disturbance_begin
2276       dist_nys(1) = ny - inflow_disturbance_end
2277    ELSEIF ( outflow_n )  THEN
2278       dist_nys    = inflow_disturbance_begin
2279       dist_nyn(1) = inflow_disturbance_end
2280    ENDIF
2281
2282!
2283!-- Check random generator
2284    IF ( random_generator /= 'system-specific'  .AND. &
2285         random_generator /= 'numerical-recipes' )  THEN
2286       IF ( myid == 0 )  THEN
2287          PRINT*, '+++ check_parameters:'
2288          PRINT*, '    unknown random generator: random_generator=', &
2289                  random_generator
2290       ENDIF
2291       CALL local_stop
2292    ENDIF
2293
2294!
2295!-- Determine damping level index for 1D model
2296    IF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
2297       IF ( damp_level_1d == -1.0 )  THEN
2298          damp_level_1d     = zu(nzt+1)
2299          damp_level_ind_1d = nzt + 1
2300       ELSEIF ( damp_level_1d < 0.0  .OR.  damp_level_1d > zu(nzt+1) )  THEN
2301          IF ( myid == 0 )  PRINT*, '+++ check_parameters:  damp_level_1d=', &
2302                                    damp_level_1d, ' must be > 0.0 and < ',  &
2303                                    'zu(nzt+1)', zu(nzt+1)
2304          CALL local_stop
2305       ELSE
2306          DO  k = 1, nzt+1
2307             IF ( damp_level_1d <= zu(k) )  THEN
2308                damp_level_ind_1d = k
2309                EXIT
2310             ENDIF
2311          ENDDO
2312       ENDIF
2313    ENDIF
2314!
2315!-- Check some other 1d-model parameters
2316    IF ( TRIM( mixing_length_1d ) /= 'as_in_3d_model'  .AND. &
2317         TRIM( mixing_length_1d ) /= 'blackadar' )  THEN
2318       IF ( myid == 0 )  PRINT*, '+++ check_parameters: mixing_length_1d = "', &
2319                                 TRIM( mixing_length_1d ), '" is unknown'
2320       CALL local_stop
2321    ENDIF
2322    IF ( TRIM( dissipation_1d ) /= 'as_in_3d_model'  .AND. &
2323         TRIM( dissipation_1d ) /= 'detering' )  THEN
2324       IF ( myid == 0 )  PRINT*, '+++ check_parameters: dissipation_1d = "', &
2325                                 TRIM( dissipation_1d ), '" is unknown'
2326       CALL local_stop
2327    ENDIF
2328
2329!
2330!-- Set time for the next user defined restart (time_restart is the
2331!-- internal parameter for steering restart events)
2332    IF ( restart_time /= 9999999.9 )  THEN
2333       IF ( restart_time > simulated_time )  time_restart = restart_time
2334    ELSE
2335!
2336!--    In case of a restart run, set internal parameter to default (no restart)
2337!--    if the NAMELIST-parameter restart_time is omitted
2338       time_restart = 9999999.9
2339    ENDIF
2340
2341!
2342!-- Set default value of the time needed to terminate a model run
2343    IF ( termination_time_needed == -1.0 )  THEN
2344       IF ( host(1:3) == 'ibm' )  THEN
2345          termination_time_needed = 300.0
2346       ELSE
2347          termination_time_needed = 35.0
2348       ENDIF
2349    ENDIF
2350
2351!
2352!-- Check the time needed to terminate a model run
2353    IF ( host(1:3) == 't3e' )  THEN
2354!
2355!--    Time needed must be at least 30 seconds on all CRAY machines, because
2356!--    MPP_TREMAIN gives the remaining CPU time only in steps of 30 seconds
2357       IF ( termination_time_needed <= 30.0 )  THEN
2358          IF ( myid == 0 )  THEN
2359             PRINT*, '+++ check_parameters:  termination_time_needed', &
2360                      termination_time_needed
2361             PRINT*, '                       must be > 30.0 on host "', host, &
2362                     '"'
2363          ENDIF
2364          CALL local_stop
2365       ENDIF
2366    ELSEIF ( host(1:3) == 'ibm' )  THEN
2367!
2368!--    On IBM-regatta machines the time should be at least 300 seconds,
2369!--    because the job time consumed before executing palm (for compiling,
2370!--    copying of files, etc.) has to be regarded
2371       IF ( termination_time_needed < 300.0 )  THEN
2372          IF ( myid == 0 )  THEN
2373             PRINT*, '+++ WARNING: check_parameters:  termination_time_',  &
2374                         'needed', termination_time_needed
2375             PRINT*, '                                should be >= 300.0', &
2376                         ' on host "', host, '"'
2377          ENDIF
2378       ENDIF
2379    ENDIF
2380
2381
2382 END SUBROUTINE check_parameters
Note: See TracBrowser for help on using the repository browser.