Home

Context Navigation

check_parameters.f90 @ 92

Last change on this file since 92 was 90, checked in by raasch, 17 years ago

New:
---
Calculation and output of user-defined profiles. New &userpar parameters data_output_pr_user and max_pr_user.

check_parameters, flow_statistics, modules, parin, read_var_list, user_interface, write_var_list

Changed:

Division through dt_3d replaced by multiplication of the inverse. For performance optimisation, this is done in the loop calculating the divergence instead of using a seperate loop. (pres.f90) var_hom and var_sum renamed pr_palm.

data_output_profiles, flow_statistics, init_3d_model, modules, parin, pres, read_var_list, run_control, time_integration

Errors:

Bugfix: work_fft*_vec removed from some PRIVATE-declarations (poisfft).

Bugfix: field_chr renamed field_char (user_interface).

Bugfix: output of use_upstream_for_tke (header).

header, poisfft, user_interface

Property svn:keywords set to Id

File size: 89.4 KB

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

Note: See TracBrowser for help on using the repository browser.

Context Navigation

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

Download in other formats: