Home

Context Navigation

check_parameters.f90 @ 404

Last change on this file since 404 was 392, checked in by raasch, 15 years ago

New:
---

Adapted for machine lck
(mrun, mbuild, subjob)

bc_lr/bc_ns in most subroutines replaced by LOGICAL variables bc_lr_cyc,
bc_ns_cyc for speed optimization
(check_parameters, diffusion_u, diffusion_v, diffusion_w, modules)

Additional timestep criterion in case of simulations with plant canopy (timestep)

Check for illegal entries in section_xy|xz|yz that exceed nz+1|ny+1|nx+1
(check_parameters)

Clipping of dvrp output implemented. Default colourtable for particles
implemented, particle attributes (color, dvrp_size) can be set with new
parameters particle_color, particle_dvrpsize, color_interval,
dvrpsize_interval (init_dvrp, data_output_dvrp, modules, user_data_output_dvrp).
Slicer attributes (dvrp) are set with new routine set_slicer_attributes_dvrp
and are controlled with existing parameters slicer_range_limits.
(set_slicer_attributes_dvrp)

Ocean atmosphere coupling allows to use independent precursor runs in order
to account for different spin-up times. The time when coupling has to be
started is given by new inipar parameter coupling_start_time. The precursor
ocean run has to be started using new mrun option "-y" in order to add
appendix "_O" to all output files.
(check_for_restart, check_parameters, data_output_2d, data_output_3d,
data_output_profiles, data_output_ptseries, data_output_spectra,
data_output_tseries, header, init_coupling, modules, mrun,
parin, read_var_list, surface_coupler, time_integration, write_var_list)

Polygon reduction for topography and ground plate isosurface. Reduction level
for buildings can be chosen with parameter cluster_size. (init_dvrp)

External pressure gradient (check_parameters, header, init_3d_model, modules,
parin, prognostic_equations, read_var_list, write_var_list)

New topography case 'single_street_canyon' (header, init_grid, modules, parin,
read_var_list, user_check_parameters, user_header, user_init_grid, write_var_list)

Option to predefine a target bulk velocity for conserve_volume_flow
(check_parameters, header, init_3d_model, modules, parin, read_var_list,
write_var_list)

Option for user defined 2D data output in xy cross sections at z=nzb+1
(data_output_2d, user_data_output_2d)

xy cross section output of surface heatfluxes (latent, sensible)
(average_3d_data, check_parameters, data_output_2d, modules, read_3d_binary,
sum_up_3d_data, write_3d_binary)

average_3d_data, check_for_restart, check_parameters, data_output_2d, data_output_3d, data_output_dvrp, data_output_profiles, data_output_ptseries, data_output_spectra, data_output_tseries, init_coupling, init_dvrp, init_grid, init_3d_model, header, mbuild, modules, mrun, package_parin, parin, prognostic_equations, read_3d_binary, read_var_list, subjob, surface_coupler, timestep, time_integration, user_check_parameters, user_data_output_2d, user_data_output_dvrp, user_header, user_init_grid, write_3d_binary, write_var_list

New: set_particle_attributes, set_slicer_attributes_dvrp

Changed:

lcmuk changed to lc to avoid problems with Intel compiler on sgi-ice
(poisfft)

For extended NetCDF files, the updated title attribute includes an update of
time_average_text where appropriate. (netcdf)

In case of restart runs without extension, initial profiles are not written
to NetCDF-file anymore. (data_output_profiles, modules, read_var_list, write_var_list)

Small change in formatting of the message handling routine concering the output in the
job protocoll. (message)

initializing_actions='read_data_for_recycling' renamed to 'cyclic_fill', now
independent of turbulent_inflow (check_parameters, header, init_3d_model)

2 NetCDF error numbers changed. (data_output_3d)

A Link to the website appendix_a.html is printed for further information
about the possible errors. (message)

Temperature gradient criterion for estimating the boundary layer height
replaced by the gradient criterion of Sullivan et al. (1998). (flow_statistics)

NetCDF unit attribute in timeseries output in case of statistic regions added
(netcdf)

Output of NetCDF messages with aid of message handling routine.
(check_open, close_file, data_output_2d, data_output_3d,
data_output_profiles, data_output_ptseries, data_output_spectra,
data_output_tseries, netcdf, output_particles_netcdf)

Output of messages replaced by message handling routine.
(advec_particles, advec_s_bc, buoyancy, calc_spectra, check_for_restart,
check_open, coriolis, cpu_log, data_output_2d, data_output_3d, data_output_dvrp,
data_output_profiles, data_output_spectra, fft_xy, flow_statistics, header,
init_1d_model, init_3d_model, init_dvrp, init_grid, init_particles, init_pegrid,
netcdf, parin, plant_canopy_model, poisfft_hybrid, poismg, read_3d_binary,
read_var_list, surface_coupler, temperton_fft, timestep, user_actions,
user_data_output_dvrp, user_dvrp_coltab, user_init_grid, user_init_plant_canopy,
user_parin, user_read_restart_data, user_spectra )

Maximum number of tails is calculated from maximum number of particles and
skip_particles_for_tail (init_particles)

Value of vertical_particle_advection may differ for each particle group
(advec_particles, header, modules)

First constant in array den also defined as type double. (eqn_state_seawater)

Parameter dvrp_psize moved from particles_par to dvrp_graphics_par. (package_parin)

topography_grid_convention moved from userpar to inipar (check_parameters,
header, parin, read_var_list, user_check_parameters, user_header,
user_init_grid, user_parin, write_var_list)

Default value of grid_matching changed to strict.

Adjustments for runs on lcxt4 (necessary due to an software update on CRAY) and
for coupled runs on ibmy (mrun, subjob)

advec_particles, advec_s_bc, buoyancy, calc_spectra, check_for_restart, check_open, check_parameters, close_file, coriolis, cpu_log, data_output_2d, data_output_3d, data_output_dvrp, data_output_profiles, data_output_ptseries, data_output_spectra, data_output_tseries, eqn_state_seawater, fft_xy, flow_statistics, header, init_1d_model, init_3d_model, init_dvrp, init_grid, init_particles, init_pegrid, message, mrun, netcdf, output_particles_netcdf, package_parin, parin, plant_canopy_model, poisfft, poisfft_hybrid, poismg, read_3d_binary, read_var_list, sort_particles, subjob, user_check_parameters, user_header, user_init_grid, user_parin, surface_coupler, temperton_fft, timestep, user_actions, user_data_output_dvrp, user_dvrp_coltab, user_init_grid, user_init_plant_canopy, user_parin, user_read_restart_data, user_spectra, write_var_list

Errors:

Bugfix: Initial hydrostatic pressure profile in case of ocean runs is now
calculated in 5 iteration steps. (init_ocean)

Bugfix: wrong sign in buoyancy production of ocean part in case of not using
the reference density (only in 3D routine production_e) (production_e)

Bugfix: output of averaged 2d/3d quantities requires that an avaraging
interval has been set, respective error message is included (check_parameters)

Bugfix: Output on unit 14 only if requested by write_binary.
(user_last_actions)

Bugfix to avoid zero division by km_neutral (production_e)

Bugfix for extended NetCDF files: In order to avoid 'data mode' errors if
updated attributes are larger than their original size, NF90_PUT_ATT is called
in 'define mode' enclosed by NF90_REDEF and NF90_ENDDEF calls. This implies a
possible performance loss; an alternative strategy would be to ensure equal
attribute size in a job chain. (netcdf)

Bugfix: correction of initial volume flow for non-flat topography (init_3d_model)
Bugfix: zero initialization of arrays within buildings for 'cyclic_fill' (init_3d_model)

Bugfix: to_be_resorted => s_av for time-averaged scalars (data_output_2d, data_output_3d)

Bugfix: error in formatting the output (message)

Bugfix: avoid that ngp_2dh_s_inner becomes zero (init_3_model)

Typographical error: unit of wpt in dots_unit (modules)

Bugfix: error in check, if particles moved further than one subdomain length.
This check must not be applied for newly released particles. (advec_particles)

Bugfix: several tail counters are initialized, particle_tail_coordinates is
only written to file if its third index is > 0, arrays for tails are allocated
with a minimum size of 10 tails if there is no tail initially (init_particles,
advec_particles)

Bugfix: pressure included for profile output (check_parameters)

Bugfix: Type of count and count_rate changed to default INTEGER on NEC machines
(cpu_log)

Bugfix: output if particle time series only if particle advection is switched
on. (time_integration)

Bugfix: qsws was calculated in case of constant heatflux = .FALSE. (prandtl_fluxes)

Bugfix: averaging along z is not allowed for 2d quantities (e.g. u* and z0) (data_output_2d)

Typographical errors (netcdf)

If the inversion height calculated by the prerun is zero, inflow_damping_height
must be explicitly specified (init_3d_model)

Small bugfix concerning 3d 64bit netcdf output format (header)

Bugfix: dt_fixed removed from the restart file, because otherwise, no change
from a fixed to a variable timestep would be possible in restart runs.
(read_var_list, write_var_list)

Bugfix: initial setting of time_coupling in coupled restart runs (time_integration)

advec_particles, check_parameters, cpu_log, data_output_2d, data_output_3d, header, init_3d_model, init_particles, init_ocean, modules, netcdf, prandtl_fluxes, production_e, read_var_list, time_integration, user_last_actions, write_var_list

Property svn:keywords set to Id

File size: 119.9 KB

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

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

Context Navigation

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

Download in other formats: