Home

Context Navigation

check_parameters.f90 @ 1179

Last change on this file since 1179 was 1179, checked in by raasch, 11 years ago

New:
---
Initial profiles can be used as reference state in the buoyancy term. New parameter
reference_state introduced. Calculation and handling of reference state in buoyancy term revised.
binary version for restart files changed from 3.9 to 3.9a (no downward compatibility!),
initial profile for rho added to hom (id=77)

Errors:

small bugfix for background communication (time_integration)

Property svn:keywords set to Id

File size: 139.5 KB

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

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

Context Navigation

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

Download in other formats: