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

Rev	Line
[1]	1	SUBROUTINE check_parameters
	2
	3	!------------------------------------------------------------------------------!
	4	! Actual revisions:
	5	! -----------------
[392]	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
[388]	13	! Check profiles fpr prho and hyp.
[376]	14	! Bugfix: output of averaged 2d/3d quantities requires that an avaraging
	15	! interval has been set, respective error message is included
[366]	16	! bc_lr_cyc and bc_ns_cyc are set,
[328]	17	! initializing_actions='read_data_for_recycling' renamed to 'cyclic_fill'
[308]	18	! Check for illegal entries in section_xy\|xz\|yz that exceed nz+1\|ny+1\|nx+1
[291]	19	! Coupling with independent precursor runs.
[264]	20	! Check particle_color, particle_dvrpsize, color_interval, dvrpsize_interval
[232]	21	! Bugfix: pressure included for profile output
[240]	22	! Check pressure gradient conditions
[256]	23	! topography_grid_convention moved from user_check_parameters
	24	! 'single_street_canyon'
[367]	25	! Added shf* and qsws* to the list of available output data
[226]	26	!
	27	! 222 2009-01-12 16:04:16Z letzel
[217]	28	! +user_check_parameters
[213]	29	! Output of messages replaced by message handling routine.
[206]	30	! Implementation of an MPI-1 coupling: replaced myid with target_id,
	31	! deleted __mpi2 directives
[222]	32	! Check that PALM is called with mrun -K parallel for coupling
[198]	33	!
	34	! 197 2008-09-16 15:29:03Z raasch
[177]	35	! Bug fix: Construction of vertical profiles when 10 gradients have been
[181]	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.
[153]	39	! Leaf area density (LAD) explicitly set to its surface value at k=0
[151]	40	! Case of reading data for recycling included in initializing_actions,
[197]	41	! check of turbulent_inflow and calculation of recycling_plane.
	42	! q*2 profile added
[139]	43	!
	44	! 138 2007-11-28 10:03:58Z letzel
[138]	45	! Plant canopy added
[132]	46	! Allow new case bc_uv_t = 'dirichlet_0' for channel flow.
[116]	47	! Multigrid solver allows topography, checking of dt_sort_particles
[121]	48	! Bugfix: initializing u_init and v_init in case of ocean runs
[110]	49	!
	50	! 109 2007-08-28 15:26:47Z letzel
[102]	51	! Check coupling_mode and set default (obligatory) values (like boundary
[106]	52	! conditions for temperature and fluxes) in case of coupled runs.
	53	! +profiles for wp and w"e
[104]	54	! Bugfix: Error message concerning output of particle concentration (pc)
	55	! modified
[108]	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.
[109]	60	! Check and, if necessary, set default value for dt_coupling
[98]	61	!
	62	! 97 2007-06-21 08:23:15Z raasch
[97]	63	! Initial salinity profile is calculated, salinity boundary conditions are
	64	! checked,
[94]	65	! z_max_do1d is checked only in case of ocean = .f.,
[97]	66	! +initial temperature and geostrophic velocity profiles for the ocean version,
	67	! use_pt_reference renamed use_reference
[77]	68	!
[90]	69	! 89 2007-05-25 12:08:31Z raasch
	70	! Check for user-defined profiles
	71	!
[77]	72	! 75 2007-03-22 09:54:05Z raasch
[51]	73	! "by_user" allowed as initializing action, -data_output_ts,
[63]	74	! leapfrog with non-flat topography not allowed any more, loop_optimization
[75]	75	! and pt_reference are checked, moisture renamed humidity,
[72]	76	! output of precipitation amount/rate and roughnes length + check
[73]	77	! possible negative humidities are avoided in initial profile,
[75]	78	! dirichlet/neumann changed to dirichlet/radiation, etc.,
	79	! revision added to run_description_header
[1]	80	!
[39]	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	!
[3]	85	! RCS Log replace by Id keyword, revision history cleaned up
	86	!
[1]	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
[264]	104	USE dvrp_variables
[1]	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
[102]	122	CHARACTER (LEN=40) :: coupling_string
[1]	123	CHARACTER (LEN=100) :: action
	124
[108]	125	INTEGER :: i, ilen, intervals, iremote = 0, iter, j, k, nnxh, nnyh, &
	126	position, prec
[1]	127	LOGICAL :: found, ldum
[291]	128	REAL :: gradient, maxn, maxp, remote = 0.0, &
	129	simulation_time_since_reference
[1]	130
	131	!
	132	!-- Warning, if host is not set
	133	IF ( host(1:1) == ' ' ) THEN
[213]	134	message_string = '"host" is not set. Please check that environment ' // &
	135	'variable "localhost" & is set before running PALM'
[226]	136	CALL message( 'check_parameters', 'PA0001', 0, 0, 0, 6, 0 )
[1]	137	ENDIF
	138
	139	!
[102]	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
[213]	144	message_string = 'illegal coupling mode: ' // TRIM( coupling_mode )
[226]	145	CALL message( 'check_parameters', 'PA0002', 1, 2, 0, 6, 0 )
[102]	146	ENDIF
	147
	148	!
[108]	149	!-- Check dt_coupling, restart_time, dt_restart, end_time, dx, dy, nx and ny
	150	IF ( coupling_mode /= 'uncoupled' ) THEN
[213]	151
[108]	152	IF ( dt_coupling == 9999999.9 ) THEN
[213]	153	message_string = 'dt_coupling is not set but required for coup' // &
	154	'ling mode "' // TRIM( coupling_mode ) // '"'
[226]	155	CALL message( 'check_parameters', 'PA0003', 1, 2, 0, 6, 0 )
[108]	156	ENDIF
[213]	157
[206]	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 )
[108]	163	IF ( dt_coupling /= remote ) THEN
[213]	164	WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
	165	'": dt_coupling = ', dt_coupling, '& is not equal to ', &
	166	'dt_coupling_remote = ', remote
[226]	167	CALL message( 'check_parameters', 'PA0004', 1, 2, 0, 6, 0 )
[108]	168	ENDIF
[109]	169	IF ( dt_coupling <= 0.0 ) THEN
[206]	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 )
[109]	173	dt_coupling = MAX( dt_max, remote )
[213]	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
[226]	177	CALL message( 'check_parameters', 'PA0005', 0, 1, 0, 6, 0 )
[109]	178	ENDIF
[213]	179
[206]	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 )
[108]	184	IF ( restart_time /= remote ) THEN
[213]	185	WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
	186	'": restart_time = ', restart_time, '& is not equal to ', &
	187	'restart_time_remote = ', remote
[226]	188	CALL message( 'check_parameters', 'PA0006', 1, 2, 0, 6, 0 )
[108]	189	ENDIF
[213]	190
[206]	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 )
[108]	195	IF ( dt_restart /= remote ) THEN
[213]	196	WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
	197	'": dt_restart = ', dt_restart, '& is not equal to ', &
	198	'dt_restart_remote = ', remote
[226]	199	CALL message( 'check_parameters', 'PA0007', 1, 2, 0, 6, 0 )
[108]	200	ENDIF
[213]	201
[291]	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 )
[206]	205	CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 14, comm_inter, &
	206	status, ierr )
[291]	207	IF ( simulation_time_since_reference /= remote ) THEN
[213]	208	WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
[291]	209	'": simulation_time_since_reference = ', &
	210	simulation_time_since_reference, '& is not equal to ', &
	211	'simulation_time_since_reference_remote = ', remote
[226]	212	CALL message( 'check_parameters', 'PA0008', 1, 2, 0, 6, 0 )
[108]	213	ENDIF
[213]	214
[206]	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 )
[108]	218	IF ( dx /= remote ) THEN
[213]	219	WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
	220	'": dx = ', dx, '& is not equal to dx_remote = ', remote
[226]	221	CALL message( 'check_parameters', 'PA0009', 1, 2, 0, 6, 0 )
[108]	222	ENDIF
[213]	223
[206]	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 )
[108]	227	IF ( dy /= remote ) THEN
[213]	228	WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
	229	'": dy = ', dy, '& is not equal to dy_remote = ', remote
[226]	230	CALL message( 'check_parameters', 'PA0010', 1, 2, 0, 6, 0 )
[108]	231	ENDIF
[213]	232
[206]	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 )
[108]	236	IF ( nx /= iremote ) THEN
[213]	237	WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
	238	'": nx = ', nx, '& is not equal to nx_remote = ', iremote
[226]	239	CALL message( 'check_parameters', 'PA0011', 1, 2, 0, 6, 0 )
[108]	240	ENDIF
[213]	241
[206]	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 )
[108]	245	IF ( ny /= iremote ) THEN
[213]	246	WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
	247	'": ny = ', ny, '& is not equal to ny_remote = ', iremote
[226]	248	CALL message( 'check_parameters', 'PA0012', 1, 2, 0, 6, 0 )
[108]	249	ENDIF
[222]	250	#else
	251	WRITE( message_string, * ) 'coupling requires PALM to be called with', &
	252	' ''mrun -K parallel'''
[226]	253	CALL message( 'check_parameters', 'PA0141', 1, 2, 0, 6, 0 )
[108]	254	#endif
	255	ENDIF
	256
[206]	257	#if defined( __parallel )
[108]	258	!
	259	!-- Exchange via intercommunicator
	260	IF ( coupling_mode == 'atmosphere_to_ocean' ) THEN
[206]	261	CALL MPI_SEND( humidity, 1, MPI_LOGICAL, target_id, 19, comm_inter, &
	262	ierr )
[108]	263	ELSEIF ( coupling_mode == 'ocean_to_atmosphere' ) THEN
[206]	264	CALL MPI_RECV( humidity_remote, 1, MPI_LOGICAL, target_id, 19, &
	265	comm_inter, status, ierr )
[108]	266	ENDIF
	267	#endif
	268
	269
	270	!
[1]	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)
[102]	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
[1]	283
[102]	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
[1]	289
	290	!
[63]	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
[213]	301	message_string = 'illegal value given for loop_optimization: "' // &
	302	TRIM( loop_optimization ) // '"'
[226]	303	CALL message( 'check_parameters', 'PA0013', 1, 2, 0, 6, 0 )
[63]	304	ENDIF
	305
	306	!
[1]	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
[51]	316	IF ( timestep_scheme(1:8) == 'leapfrog' ) THEN
	317	WRITE( action, '(A,A)' ) 'timestep_scheme = ', timestep_scheme
	318	ENDIF
[114]	319	IF ( psolver == 'sor' ) THEN
[1]	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
[75]	334	IF ( humidity ) THEN
	335	WRITE( action, '(A)' ) 'humidity = .TRUE.'
[1]	336	ENDIF
	337	IF ( .NOT. prandtl_layer ) THEN
	338	WRITE( action, '(A)' ) 'prandtl_layer = .FALSE.'
	339	ENDIF
	340	IF ( action /= ' ' ) THEN
[213]	341	message_string = 'a non-flat topography does not allow ' // &
	342	TRIM( action )
[226]	343	CALL message( 'check_parameters', 'PA0014', 1, 2, 0, 6, 0 )
[1]	344	ENDIF
[256]	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
[1]	381	ENDIF
[94]	382
[1]	383	!
[94]	384	!-- Check ocean setting
	385	IF ( ocean ) THEN
[332]	386
[94]	387	action = ' '
	388	IF ( timestep_scheme(1:8) == 'leapfrog' ) THEN
	389	WRITE( action, '(A,A)' ) 'timestep_scheme = ', timestep_scheme
	390	ENDIF
[97]	391	IF ( momentum_advec == 'ups-scheme' ) THEN
	392	WRITE( action, '(A,A)' ) 'momentum_advec = ', momentum_advec
	393	ENDIF
[94]	394	IF ( action /= ' ' ) THEN
[213]	395	message_string = 'ocean = .T. does not allow ' // TRIM( action )
[226]	396	CALL message( 'check_parameters', 'PA0015', 1, 2, 0, 6, 0 )
[94]	397	ENDIF
	398
[332]	399	ELSEIF ( TRIM( coupling_mode ) == 'uncoupled' .AND. &
	400	TRIM( coupling_char ) == '_O' ) THEN
[343]	401
[94]	402	!
[332]	403	!-- Check whether an (uncoupled) atmospheric run has been declared as an
	404	!-- ocean run (this setting is done via mrun-option -y)
[343]	405
[291]	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 )
[332]	409
[291]	410	ENDIF
	411
	412	!
[1]	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
[213]	417	message_string = 'unknown solver for perturbation pressure: psolver' // &
	418	' = "' // TRIM( psolver ) // '"'
[226]	419	CALL message( 'check_parameters', 'PA0016', 1, 2, 0, 6, 0 )
[1]	420	ENDIF
	421
	422	#if defined( __parallel )
	423	IF ( psolver == 'poisfft_hybrid' .AND. pdims(2) /= 1 ) THEN
[213]	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'
[226]	427	CALL message( 'check_parameters', 'PA0017', 1, 2, 0, 6, 0 )
[1]	428	ENDIF
[181]	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
[213]	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'
[226]	436	CALL message( 'check_parameters', 'PA0018', 1, 2, 0, 6, 0 )
[1]	437	ENDIF
	438	#else
	439	IF ( psolver == 'poisfft_hybrid' ) THEN
[213]	440	message_string = 'psolver = "' // TRIM( psolver ) // '" only works' // &
	441	' for a parallel environment'
[226]	442	CALL message( 'check_parameters', 'PA0019', 1, 2, 0, 6, 0 )
[1]	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
[213]	452	message_string = 'unknown multigrid cycle: cycle_mg = "' // &
	453	TRIM( cycle_mg ) // '"'
[226]	454	CALL message( 'check_parameters', 'PA0020', 1, 2, 0, 6, 0 )
[1]	455	ENDIF
	456	ENDIF
	457
	458	IF ( fft_method /= 'singleton-algorithm' .AND. &
	459	fft_method /= 'temperton-algorithm' .AND. &
	460	fft_method /= 'system-specific' ) THEN
[213]	461	message_string = 'unknown fft-algorithm: fft_method = "' // &
	462	TRIM( fft_method ) // '"'
[226]	463	CALL message( 'check_parameters', 'PA0021', 1, 2, 0, 6, 0 )
[1]	464	ENDIF
	465
	466	!
	467	!-- Advection schemes:
	468	IF ( momentum_advec /= 'pw-scheme' .AND. momentum_advec /= 'ups-scheme' ) &
	469	THEN
[214]	470	message_string = 'unknown advection scheme: momentum_advec = "' // &
	471	TRIM( momentum_advec ) // '"'
[226]	472	CALL message( 'check_parameters', 'PA0022', 1, 2, 0, 6, 0 )
[1]	473	ENDIF
	474	IF ( ( momentum_advec == 'ups-scheme' .OR. scalar_advec == 'ups-scheme' )&
	475	.AND. timestep_scheme /= 'euler' ) THEN
[214]	476	message_string = 'momentum_advec = "' // TRIM( momentum_advec ) // &
	477	'" is not allowed with timestep_scheme = "' // &
	478	TRIM( timestep_scheme ) // '"'
[226]	479	CALL message( 'check_parameters', 'PA0023', 1, 2, 0, 6, 0 )
[1]	480	ENDIF
	481
	482	IF ( scalar_advec /= 'pw-scheme' .AND. scalar_advec /= 'bc-scheme' .AND.&
	483	scalar_advec /= 'ups-scheme' ) THEN
[214]	484	message_string = 'unknown advection scheme: scalar_advec = "' // &
	485	TRIM( scalar_advec ) // '"'
[226]	486	CALL message( 'check_parameters', 'PA0024', 1, 2, 0, 6, 0 )
[1]	487	ENDIF
	488
	489	IF ( use_sgs_for_particles .AND. .NOT. use_upstream_for_tke ) THEN
	490	use_upstream_for_tke = .TRUE.
[214]	491	message_string = 'use_upstream_for_tke set .TRUE. because ' // &
	492	'use_sgs_for_particles = .TRUE.'
[226]	493	CALL message( 'check_parameters', 'PA0025', 0, 1, 0, 6, 0 )
[1]	494	ENDIF
	495
	496	IF ( use_upstream_for_tke .AND. timestep_scheme(1:8) == 'leapfrog' ) THEN
[214]	497	message_string = 'use_upstream_for_tke = .TRUE. not allowed with ' // &
	498	'timestep_scheme = "' // TRIM( timestep_scheme ) // '"'
[226]	499	CALL message( 'check_parameters', 'PA0026', 1, 2, 0, 6, 0 )
[1]	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
[214]	522	message_string = 'unknown timestep scheme: timestep_scheme = "' // &
	523	TRIM( timestep_scheme ) // '"'
[226]	524	CALL message( 'check_parameters', 'PA0027', 1, 2, 0, 6, 0 )
[1]	525
	526	END SELECT
	527
[63]	528	IF ( scalar_advec == 'ups-scheme' .AND. timestep_scheme(1:5) == 'runge' )&
[1]	529	THEN
[214]	530	message_string = 'scalar advection scheme "' // TRIM( scalar_advec ) // &
	531	'" & does not work with timestep_scheme "' // &
	532	TRIM( timestep_scheme ) // '"'
[226]	533	CALL message( 'check_parameters', 'PA0028', 1, 2, 0, 6, 0 )
[1]	534	ENDIF
	535
	536	IF ( momentum_advec /= 'pw-scheme' .AND. timestep_scheme(1:5) == 'runge' ) &
	537	THEN
[214]	538	message_string = 'momentum advection scheme "' // &
	539	TRIM( momentum_advec ) // '" & does not work with ' // &
	540	'timestep_scheme "' // TRIM( timestep_scheme ) // '"'
[226]	541	CALL message( 'check_parameters', 'PA0029', 1, 2, 0, 6, 0 )
[1]	542	ENDIF
	543
[147]	544	IF ( TRIM( initializing_actions ) /= 'read_restart_data' .AND. &
[328]	545	TRIM( initializing_actions ) /= 'cyclic_fill' ) THEN
[1]	546	!
[214]	547	!-- No restart run: several initialising actions are possible
[1]	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', &
[46]	554	'by_user', 'initialize_vortex', 'initialize_ptanom' )
[1]	555	action = action(position+1:)
	556
	557	CASE DEFAULT
[214]	558	message_string = 'initializing_action = "' // &
	559	TRIM( action ) // '" unkown or not allowed'
[226]	560	CALL message( 'check_parameters', 'PA0030', 1, 2, 0, 6, 0 )
[1]	561
	562	END SELECT
	563	ENDDO
	564	ENDIF
[214]	565
[1]	566	IF ( INDEX( initializing_actions, 'set_constant_profiles' ) /= 0 .AND. &
	567	INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 ) THEN
[214]	568	message_string = 'initializing_actions = "set_constant_profiles"' // &
	569	' and "set_1d-model_profiles" are not allowed ' // &
	570	'simultaneously'
[226]	571	CALL message( 'check_parameters', 'PA0031', 1, 2, 0, 6, 0 )
[1]	572	ENDIF
[214]	573
[46]	574	IF ( INDEX( initializing_actions, 'set_constant_profiles' ) /= 0 .AND. &
	575	INDEX( initializing_actions, 'by_user' ) /= 0 ) THEN
[214]	576	message_string = 'initializing_actions = "set_constant_profiles"' // &
	577	' and "by_user" are not allowed simultaneously'
[226]	578	CALL message( 'check_parameters', 'PA0032', 1, 2, 0, 6, 0 )
[46]	579	ENDIF
[214]	580
[46]	581	IF ( INDEX( initializing_actions, 'by_user' ) /= 0 .AND. &
	582	INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 ) THEN
[214]	583	message_string = 'initializing_actions = "by_user" and ' // &
	584	'"set_1d-model_profiles" are not allowed simultaneously'
[226]	585	CALL message( 'check_parameters', 'PA0033', 1, 2, 0, 6, 0 )
[46]	586	ENDIF
[1]	587
[75]	588	IF ( cloud_physics .AND. .NOT. humidity ) THEN
[214]	589	WRITE( message_string, * ) 'cloud_physics = ', cloud_physics, ' is ', &
	590	'not allowed with humidity = ', humidity
[226]	591	CALL message( 'check_parameters', 'PA0034', 1, 2, 0, 6, 0 )
[1]	592	ENDIF
	593
[72]	594	IF ( precipitation .AND. .NOT. cloud_physics ) THEN
[214]	595	WRITE( message_string, * ) 'precipitation = ', precipitation, ' is ', &
	596	'not allowed with cloud_physics = ', cloud_physics
[226]	597	CALL message( 'check_parameters', 'PA0035', 1, 2, 0, 6, 0 )
[72]	598	ENDIF
	599
[75]	600	IF ( humidity .AND. sloping_surface ) THEN
[214]	601	message_string = 'humidity = .TRUE. and sloping_surface = .TRUE. ' // &
	602	'are not allowed simultaneously'
[226]	603	CALL message( 'check_parameters', 'PA0036', 1, 2, 0, 6, 0 )
[1]	604	ENDIF
	605
[75]	606	IF ( humidity .AND. scalar_advec == 'ups-scheme' ) THEN
[214]	607	message_string = 'UPS-scheme is not implemented for humidity = .TRUE.'
[226]	608	CALL message( 'check_parameters', 'PA0037', 1, 2, 0, 6, 0 )
[1]	609	ENDIF
	610
[75]	611	IF ( passive_scalar .AND. humidity ) THEN
[214]	612	message_string = 'humidity = .TRUE. and passive_scalar = .TRUE. ' // &
	613	'is not allowed simultaneously'
[226]	614	CALL message( 'check_parameters', 'PA0038', 1, 2, 0, 6, 0 )
[1]	615	ENDIF
	616
	617	IF ( passive_scalar .AND. scalar_advec == 'ups-scheme' ) THEN
[214]	618	message_string = 'UPS-scheme is not implemented for passive_scalar' // &
	619	' = .TRUE.'
[226]	620	CALL message( 'check_parameters', 'PA0039', 1, 2, 0, 6, 0 )
[1]	621	ENDIF
	622
	623	IF ( grid_matching /= 'strict' .AND. grid_matching /= 'match' ) THEN
[214]	624	message_string = 'illegal value "' // TRIM( grid_matching ) // &
	625	'" found for parameter grid_matching'
[226]	626	CALL message( 'check_parameters', 'PA0040', 1, 2, 0, 6, 0 )
[1]	627	ENDIF
	628
[138]	629	IF ( plant_canopy .AND. ( drag_coefficient == 0.0 ) ) THEN
[214]	630	message_string = 'plant_canopy = .TRUE. requires a non-zero drag ' // &
	631	'coefficient & given value is drag_coefficient = 0.0'
[226]	632	CALL message( 'check_parameters', 'PA0041', 1, 2, 0, 6, 0 )
[138]	633	ENDIF
	634
[1]	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
[94]	645	IF ( humidity ) q_init = q_surface
	646	IF ( ocean ) sa_init = sa_surface
	647	IF ( passive_scalar ) q_init = s_surface
[138]	648	IF ( plant_canopy ) lad = 0.0
[1]	649
	650	!
	651	!--
	652	!-- If required, compute initial profile of the geostrophic wind
	653	!-- (component ug)
	654	i = 1
	655	gradient = 0.0
[97]	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
[177]	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
[1]	668	ENDIF
[177]	669	ENDIF
[97]	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
[1]	676	ELSE
[97]	677	ug(k) = ug(k-1)
[1]	678	ENDIF
[97]	679	ENDDO
[1]	680
[97]	681	ELSE
	682
	683	ug_vertical_gradient_level_ind(1) = nzt+1
[121]	684	ug(nzt+1) = ug_surface
[97]	685	DO k = nzt, 0, -1
[177]	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
[97]	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
[1]	708	u_init = ug
	709
	710	!
	711	!-- In case of no given gradients for ug, choose a vanishing gradient
[97]	712	IF ( ug_vertical_gradient_level(1) == -9999999.9 ) THEN
[1]	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
[97]	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
[177]	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
[1]	734	ENDIF
	735	ENDIF
[97]	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
[1]	742	ELSE
[97]	743	vg(k) = vg(k-1)
[1]	744	ENDIF
[97]	745	ENDDO
[1]	746
[97]	747	ELSE
	748
[121]	749	vg_vertical_gradient_level_ind(1) = nzt+1
	750	vg(nzt+1) = vg_surface
[97]	751	DO k = nzt, 0, -1
[177]	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
[97]	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
[1]	774	v_init = vg
	775
	776	!
	777	!-- In case of no given gradients for vg, choose a vanishing gradient
[97]	778	IF ( vg_vertical_gradient_level(1) == -9999999.9 ) THEN
[1]	779	vg_vertical_gradient_level(1) = 0.0
	780	ENDIF
	781
	782	!
[94]	783	!-- Compute initial temperature profile using the given temperature gradients
[1]	784	i = 1
	785	gradient = 0.0
[94]	786
	787	IF ( .NOT. ocean ) THEN
	788
	789	pt_vertical_gradient_level_ind(1) = 0
	790	DO k = 1, nzt+1
[177]	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
[1]	797	ENDIF
	798	ENDIF
[94]	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
[1]	805	ELSE
[94]	806	pt_init(k) = pt_init(k-1)
[1]	807	ENDIF
[94]	808	ENDDO
[1]	809
[94]	810	ELSE
	811
	812	pt_vertical_gradient_level_ind(1) = nzt+1
	813	DO k = nzt, 0, -1
[177]	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
[94]	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
[1]	836	!
	837	!-- In case of no given temperature gradients, choose gradient of neutral
	838	!-- stratification
[94]	839	IF ( pt_vertical_gradient_level(1) == -9999999.9 ) THEN
[1]	840	pt_vertical_gradient_level(1) = 0.0
	841	ENDIF
	842
	843	!
[94]	844	!-- Store temperature gradient at the top boundary for possible Neumann
[1]	845	!-- boundary condition
[19]	846	bc_pt_t_val = ( pt_init(nzt+1) - pt_init(nzt) ) / dzu(nzt+1)
[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
[75]	862	IF ( humidity .OR. passive_scalar ) THEN
[1]	863
	864	i = 1
	865	gradient = 0.0
	866	q_vertical_gradient_level_ind(1) = 0
	867	DO k = 1, nzt+1
[177]	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
[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
[72]	885	!
	886	!-- Avoid negative humidities
	887	IF ( q_init(k) < 0.0 ) THEN
	888	q_init(k) = 0.0
	889	ENDIF
[1]	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
[19]	902	bc_q_t_val = ( q_init(nzt+1) - q_init(nzt) ) / dzu(nzt+1)
[1]	903
	904	ENDIF
	905
[94]	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
[177]	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
[94]	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
[138]	938	!
[388]	939	!-- If required compute the profile of leaf area density used in the plant
	940	!-- canopy model
[138]	941	IF ( plant_canopy ) THEN
	942
	943	i = 1
	944	gradient = 0.0
[1]	945
[138]	946	IF ( .NOT. ocean ) THEN
[153]	947
	948	lad(0) = lad_surface
[138]	949
	950	lad_vertical_gradient_level_ind(1) = 0
	951	DO k = 1, pch_index
[177]	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
[138]	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
[1]	973	!
[388]	974	!-- In case of no given leaf area density gradients, choose a vanishing
	975	!-- gradient
[138]	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	!
[1]	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	!
[97]	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
[57]	993
	994	!
[97]	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	!
[108]	1003	!-- Ocean version must use flux boundary conditions at the top
	1004	IF ( ocean .AND. .NOT. use_top_fluxes ) THEN
[215]	1005	message_string = 'use_top_fluxes must be .TRUE. in ocean version'
[226]	1006	CALL message( 'check_parameters', 'PA0042', 1, 2, 0, 6, 0 )
[108]	1007	ENDIF
[97]	1008
	1009	!
[1]	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
[215]	1013	WRITE( message_string, * ) 'ABS( alpha_surface = ', alpha_surface, &
	1014	' ) must be < 90.0'
[226]	1015	CALL message( 'check_parameters', 'PA0043', 1, 2, 0, 6, 0 )
[1]	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
[215]	1026	WRITE( message_string, * ) 'dt = ', dt , ' <= 0.0'
[226]	1027	CALL message( 'check_parameters', 'PA0044', 1, 2, 0, 6, 0 )
[1]	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
[215]	1048	WRITE( message_string, * ) 'cfl_factor = ', cfl_factor, &
	1049	' out of range & 0.0 < cfl_factor <= 1.0 is required'
[226]	1050	CALL message( 'check_parameters', 'PA0045', 1, 2, 0, 6, 0 )
[1]	1051	ENDIF
	1052	ENDIF
	1053
	1054	!
	1055	!-- Store simulated time at begin
	1056	simulated_time_at_begin = simulated_time
	1057
	1058	!
[291]	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	!
[1]	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
[215]	1079	message_string = 'baroclinicity (ug) not allowed simultaneously' // &
	1080	' with galilei transformation'
[226]	1081	CALL message( 'check_parameters', 'PA0046', 1, 2, 0, 6, 0 )
[1]	1082	ELSEIF ( use_ug_for_galilei_tr .AND. &
	1083	vg_vertical_gradient_level(1) /= 0.0 ) THEN
[215]	1084	message_string = 'baroclinicity (vg) not allowed simultaneously' // &
	1085	' with galilei transformation'
[226]	1086	CALL message( 'check_parameters', 'PA0047', 1, 2, 0, 6, 0 )
[1]	1087	ELSE
[215]	1088	message_string = 'variable translation speed used for galilei-' // &
	1089	'transformation, which may cause & instabilities in stably ' // &
	1090	'stratified regions'
[226]	1091	CALL message( 'check_parameters', 'PA0048', 0, 1, 0, 6, 0 )
[1]	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
[73]	1103	IF ( bc_lr /= 'cyclic' .AND. bc_lr /= 'dirichlet/radiation' .AND. &
	1104	bc_lr /= 'radiation/dirichlet' ) THEN
[215]	1105	message_string = 'unknown boundary condition: bc_lr = "' // &
	1106	TRIM( bc_lr ) // '"'
[226]	1107	CALL message( 'check_parameters', 'PA0049', 1, 2, 0, 6, 0 )
[1]	1108	ENDIF
[73]	1109	IF ( bc_ns /= 'cyclic' .AND. bc_ns /= 'dirichlet/radiation' .AND. &
	1110	bc_ns /= 'radiation/dirichlet' ) THEN
[215]	1111	message_string = 'unknown boundary condition: bc_ns = "' // &
	1112	TRIM( bc_ns ) // '"'
[226]	1113	CALL message( 'check_parameters', 'PA0050', 1, 2, 0, 6, 0 )
[1]	1114	ENDIF
	1115
	1116	!
[366]	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	!
[1]	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
[215]	1127	message_string = 'non-cyclic lateral boundaries do not allow ' // &
	1128	'psolver = "' // TRIM( psolver ) // '"'
[226]	1129	CALL message( 'check_parameters', 'PA0051', 1, 2, 0, 6, 0 )
[1]	1130	ENDIF
	1131	IF ( momentum_advec /= 'pw-scheme' ) THEN
[215]	1132	message_string = 'non-cyclic lateral boundaries do not allow ' // &
	1133	'momentum_advec = "' // TRIM( momentum_advec ) // '"'
[226]	1134	CALL message( 'check_parameters', 'PA0052', 1, 2, 0, 6, 0 )
[1]	1135	ENDIF
	1136	IF ( scalar_advec /= 'pw-scheme' ) THEN
[215]	1137	message_string = 'non-cyclic lateral boundaries do not allow ' // &
	1138	'scalar_advec = "' // TRIM( scalar_advec ) // '"'
[226]	1139	CALL message( 'check_parameters', 'PA0053', 1, 2, 0, 6, 0 )
[1]	1140	ENDIF
	1141	IF ( galilei_transformation ) THEN
[215]	1142	message_string = 'non-cyclic lateral boundaries do not allow ' // &
	1143	'galilei_transformation = .T.'
[226]	1144	CALL message( 'check_parameters', 'PA0054', 1, 2, 0, 6, 0 )
[1]	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
[215]	1153	message_string = 'adjust_mixing_length = TRUE and bc_e_b = "neumann"'
[226]	1154	CALL message( 'check_parameters', 'PA0055', 0, 1, 0, 6, 0 )
[1]	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
[215]	1159	message_string = 'adjust_mixing_length = FALSE and bc_e_b = "' // &
	1160	TRIM( bc_e_b ) // '"'
[226]	1161	CALL message( 'check_parameters', 'PA0056', 0, 1, 0, 6, 0 )
[1]	1162	ENDIF
	1163	IF ( .NOT. prandtl_layer ) THEN
	1164	bc_e_b = 'neumann'
	1165	ibc_e_b = 1
[215]	1166	message_string = 'boundary condition bc_e_b changed to "' // &
	1167	TRIM( bc_e_b ) // '"'
[226]	1168	CALL message( 'check_parameters', 'PA0057', 0, 1, 0, 6, 0 )
[1]	1169	ENDIF
	1170	ELSE
[215]	1171	message_string = 'unknown boundary condition: bc_e_b = "' // &
	1172	TRIM( bc_e_b ) // '"'
[226]	1173	CALL message( 'check_parameters', 'PA0058', 1, 2, 0, 6, 0 )
[1]	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
[215]	1185	message_string = 'unknown boundary condition: bc_p_b = "' // &
	1186	TRIM( bc_p_b ) // '"'
[226]	1187	CALL message( 'check_parameters', 'PA0059', 1, 2, 0, 6, 0 )
[1]	1188	ENDIF
	1189	IF ( ibc_p_b == 2 .AND. .NOT. prandtl_layer ) THEN
[215]	1190	message_string = 'boundary condition: bc_p_b = "' // TRIM( bc_p_b ) // &
	1191	'" not allowed with prandtl_layer = .FALSE.'
[226]	1192	CALL message( 'check_parameters', 'PA0060', 1, 2, 0, 6, 0 )
[1]	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
[215]	1199	message_string = 'unknown boundary condition: bc_p_t = "' // &
	1200	TRIM( bc_p_t ) // '"'
[226]	1201	CALL message( 'check_parameters', 'PA0061', 1, 2, 0, 6, 0 )
[1]	1202	ENDIF
	1203
	1204	!
	1205	!-- Boundary conditions for potential temperature
[102]	1206	IF ( coupling_mode == 'atmosphere_to_ocean' ) THEN
	1207	ibc_pt_b = 2
[1]	1208	ELSE
[102]	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
[215]	1214	message_string = 'unknown boundary condition: bc_pt_b = "' // &
	1215	TRIM( bc_pt_b ) // '"'
[226]	1216	CALL message( 'check_parameters', 'PA0062', 1, 2, 0, 6, 0 )
[1]	1217	ENDIF
	1218	ENDIF
[102]	1219
[1]	1220	IF ( bc_pt_t == 'dirichlet' ) THEN
	1221	ibc_pt_t = 0
	1222	ELSEIF ( bc_pt_t == 'neumann' ) THEN
	1223	ibc_pt_t = 1
[19]	1224	ELSEIF ( bc_pt_t == 'initial_gradient' ) THEN
	1225	ibc_pt_t = 2
[1]	1226	ELSE
[215]	1227	message_string = 'unknown boundary condition: bc_pt_t = "' // &
	1228	TRIM( bc_pt_t ) // '"'
[226]	1229	CALL message( 'check_parameters', 'PA0063', 1, 2, 0, 6, 0 )
[1]	1230	ENDIF
	1231
[20]	1232	IF ( surface_heatflux == 9999999.9 ) constant_heatflux = .FALSE.
	1233	IF ( top_heatflux == 9999999.9 ) constant_top_heatflux = .FALSE.
[103]	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. &
[215]	1238	top_momentumflux_v == 9999999.9 ) ) THEN
	1239	message_string = 'both, top_momentumflux_u AND top_momentumflux_v ' // &
	1240	'must be set'
[226]	1241	CALL message( 'check_parameters', 'PA0064', 1, 2, 0, 6, 0 )
[103]	1242	ENDIF
[1]	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
[215]	1250	message_string = 'boundary_condition: bc_pt_b = "' // TRIM( bc_pt_b ) //&
	1251	'& is not allowed with constant_heatflux = .TRUE.'
[226]	1252	CALL message( 'check_parameters', 'PA0065', 1, 2, 0, 6, 0 )
[1]	1253	ENDIF
	1254	IF ( constant_heatflux .AND. pt_surface_initial_change /= 0.0 ) THEN
[215]	1255	WRITE ( message_string, * ) 'constant_heatflux = .TRUE. is not allo', &
	1256	'wed with pt_surface_initial_change (/=0) = ', &
	1257	pt_surface_initial_change
[226]	1258	CALL message( 'check_parameters', 'PA0066', 1, 2, 0, 6, 0 )
[1]	1259	ENDIF
	1260
	1261	!
[19]	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
[215]	1267	message_string = 'boundary_condition: bc_pt_t = "' // TRIM( bc_pt_t ) //&
	1268	'" is not allowed with constant_top_heatflux = .TRUE.'
[226]	1269	CALL message( 'check_parameters', 'PA0067', 1, 2, 0, 6, 0 )
[19]	1270	ENDIF
	1271
	1272	!
[95]	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
[215]	1280	message_string = 'unknown boundary condition: bc_sa_t = "' // &
	1281	TRIM( bc_sa_t ) // '"'
[226]	1282	CALL message( 'check_parameters', 'PA0068', 1, 2, 0, 6, 0 )
[95]	1283	ENDIF
	1284
	1285	IF ( top_salinityflux == 9999999.9 ) constant_top_salinityflux = .FALSE.
[97]	1286	IF ( ibc_sa_t == 1 .AND. top_salinityflux == 9999999.9 ) THEN
[215]	1287	message_string = 'boundary condition: bc_sa_t = "' // &
	1288	TRIM( bc_sa_t ) // '" requires to set ' // &
	1289	'top_salinityflux'
[226]	1290	CALL message( 'check_parameters', 'PA0069', 1, 2, 0, 6, 0 )
[97]	1291	ENDIF
[95]	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
[215]	1299	message_string = 'boundary condition: bc_sa_t = "' // &
	1300	TRIM( bc_sa_t ) // '" is not allowed with ' // &
	1301	'constant_top_salinityflux = .TRUE.'
[226]	1302	CALL message( 'check_parameters', 'PA0070', 1, 2, 0, 6, 0 )
[95]	1303	ENDIF
	1304
	1305	ENDIF
	1306
	1307	!
[75]	1308	!-- In case of humidity or passive scalar, set boundary conditions for total
[1]	1309	!-- water content / scalar
[75]	1310	IF ( humidity .OR. passive_scalar ) THEN
	1311	IF ( humidity ) THEN
[1]	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
[215]	1321	message_string = 'unknown boundary condition: bc_' // TRIM( sq ) // &
	1322	'_b ="' // TRIM( bc_q_b ) // '"'
[226]	1323	CALL message( 'check_parameters', 'PA0071', 1, 2, 0, 6, 0 )
[1]	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
[215]	1330	message_string = 'unknown boundary condition: bc_' // TRIM( sq ) // &
	1331	'_t ="' // TRIM( bc_q_t ) // '"'
[226]	1332	CALL message( 'check_parameters', 'PA0072', 1, 2, 0, 6, 0 )
[1]	1333	ENDIF
	1334
	1335	IF ( surface_waterflux == 0.0 ) constant_waterflux = .FALSE.
	1336
	1337	!
	1338	!-- A given surface humidity implies Dirichlet boundary condition for
[75]	1339	!-- humidity. In this case specification of a constant water flux is
[1]	1340	!-- forbidden.
	1341	IF ( ibc_q_b == 0 .AND. constant_waterflux ) THEN
[215]	1342	message_string = 'boundary condition: bc_' // TRIM( sq ) // '_b ' // &
	1343	'= "' // TRIM( bc_q_b ) // '" is not allowed wi' // &
	1344	'th prescribed surface flux'
[226]	1345	CALL message( 'check_parameters', 'PA0073', 1, 2, 0, 6, 0 )
[1]	1346	ENDIF
	1347	IF ( constant_waterflux .AND. q_surface_initial_change /= 0.0 ) THEN
[215]	1348	WRITE( message_string, * ) 'a prescribed surface flux is not allo', &
	1349	'wed with ', sq, '_surface_initial_change (/=0) = ', &
	1350	q_surface_initial_change
[226]	1351	CALL message( 'check_parameters', 'PA0074', 1, 2, 0, 6, 0 )
[1]	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
[215]	1363	message_string = 'boundary condition: bc_uv_b = "' // &
	1364	TRIM( bc_uv_b ) // '" is not allowed with prandtl_layer = .TRUE.'
[226]	1365	CALL message( 'check_parameters', 'PA0075', 1, 2, 0, 6, 0 )
[1]	1366	ENDIF
	1367	ELSE
[215]	1368	message_string = 'unknown boundary condition: bc_uv_b = "' // &
	1369	TRIM( bc_uv_b ) // '"'
[226]	1370	CALL message( 'check_parameters', 'PA0076', 1, 2, 0, 6, 0 )
[1]	1371	ENDIF
[215]	1372
[108]	1373	IF ( coupling_mode == 'ocean_to_atmosphere' ) THEN
	1374	bc_uv_t = 'neumann'
[1]	1375	ibc_uv_t = 1
	1376	ELSE
[132]	1377	IF ( bc_uv_t == 'dirichlet' .OR. bc_uv_t == 'dirichlet_0' ) THEN
[108]	1378	ibc_uv_t = 0
	1379	ELSEIF ( bc_uv_t == 'neumann' ) THEN
	1380	ibc_uv_t = 1
	1381	ELSE
[215]	1382	message_string = 'unknown boundary condition: bc_uv_t = "' // &
	1383	TRIM( bc_uv_t ) // '"'
[226]	1384	CALL message( 'check_parameters', 'PA0077', 1, 2, 0, 6, 0 )
[1]	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
[215]	1399	WRITE( message_string, * ) 'rayleigh_damping_factor = ', &
	1400	rayleigh_damping_factor, ' out of range [0.0,1.0]'
[226]	1401	CALL message( 'check_parameters', 'PA0078', 1, 2, 0, 6, 0 )
[1]	1402	ENDIF
	1403	ENDIF
	1404
	1405	IF ( rayleigh_damping_height == -1.0 ) THEN
[108]	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
[1]	1411	ELSE
[108]	1412	IF ( .NOT. ocean ) THEN
	1413	IF ( rayleigh_damping_height < 0.0 .OR. &
	1414	rayleigh_damping_height > zu(nzt) ) THEN
[215]	1415	WRITE( message_string, * ) 'rayleigh_damping_height = ', &
	1416	rayleigh_damping_height, ' out of range [0.0,', zu(nzt), ']'
[226]	1417	CALL message( 'check_parameters', 'PA0079', 1, 2, 0, 6, 0 )
[1]	1418	ENDIF
[108]	1419	ELSE
	1420	IF ( rayleigh_damping_height > 0.0 .OR. &
	1421	rayleigh_damping_height < zu(nzb) ) THEN
[215]	1422	WRITE( message_string, * ) 'rayleigh_damping_height = ', &
	1423	rayleigh_damping_height, ' out of range [0.0,', zu(nzb), ']'
[226]	1424	CALL message( 'check_parameters', 'PA0079', 1, 2, 0, 6, 0 )
[108]	1425	ENDIF
[1]	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
[215]	1434	message_string = 'overshoot_limit_... < 0.0 is not allowed'
[226]	1435	CALL message( 'check_parameters', 'PA0080', 1, 2, 0, 6, 0 )
[1]	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
[215]	1439	message_string = 'ups_limit_... < 0.0 is not allowed'
[226]	1440	CALL message( 'check_parameters', 'PA0081', 1, 2, 0, 6, 0 )
[1]	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
[215]	1448	WRITE ( message_string, * ) 'number of statistic_regions = ', &
	1449	statistic_regions+1, ' but only 10 regions are allowed'
[226]	1450	CALL message( 'check_parameters', 'PA0082', 1, 2, 0, 6, 0 )
[1]	1451	ENDIF
	1452	IF ( normalizing_region > statistic_regions .OR. &
	1453	normalizing_region < 0) THEN
[215]	1454	WRITE ( message_string, * ) 'normalizing_region = ', &
	1455	normalizing_region, ' must be >= 0 and <= ',statistic_regions, &
	1456	' (value of statistic_regions)'
[226]	1457	CALL message( 'check_parameters', 'PA0083', 1, 2, 0, 6, 0 )
[1]	1458	ENDIF
	1459
	1460	!
[116]	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
[215]	1465	message_string = 'dt_sort_particles is reset to 0.0 because of cloud' //&
	1466	'_droplets = .TRUE.'
[226]	1467	CALL message( 'check_parameters', 'PA0084', 0, 1, 0, 6, 0 )
[116]	1468	ENDIF
	1469
	1470	!
[1]	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
[215]	1504	WRITE( message_string, * ) 'averaging_interval = ', &
	1505	averaging_interval, ' must be <= dt_data_output = ', dt_data_output
[226]	1506	CALL message( 'check_parameters', 'PA0085', 1, 2, 0, 6, 0 )
[1]	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
[215]	1514	WRITE( message_string, * ) 'averaging_interval_pr = ', &
	1515	averaging_interval_pr, ' must be <= dt_dopr = ', dt_dopr
[226]	1516	CALL message( 'check_parameters', 'PA0086', 1, 2, 0, 6, 0 )
[1]	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
[215]	1524	WRITE( message_string, * ) 'averaging_interval_sp = ', &
	1525	averaging_interval_sp, ' must be <= dt_dosp = ', dt_dosp
[226]	1526	CALL message( 'check_parameters', 'PA0087', 1, 2, 0, 6, 0 )
[1]	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
[215]	1549	WRITE( message_string, * ) 'dt_averaging_input = ', &
	1550	dt_averaging_input, ' must be <= averaging_interval = ', &
	1551	averaging_interval
[226]	1552	CALL message( 'check_parameters', 'PA0088', 1, 2, 0, 6, 0 )
[1]	1553	ENDIF
	1554
	1555	IF ( dt_averaging_input_pr > averaging_interval_pr ) THEN
[215]	1556	WRITE( message_string, * ) 'dt_averaging_input_pr = ', &
	1557	dt_averaging_input_pr, ' must be <= averaging_interval_pr = ', &
	1558	averaging_interval_pr
[226]	1559	CALL message( 'check_parameters', 'PA0089', 1, 2, 0, 6, 0 )
[1]	1560	ENDIF
	1561
	1562	!
[72]	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
[215]	1569	WRITE( message_string, * ) 'precipitation_amount_interval = ', &
	1570	precipitation_amount_interval, ' must not be larger than ', &
	1571	'dt_do2d_xy = ', dt_do2d_xy
[226]	1572	CALL message( 'check_parameters', 'PA0090', 1, 2, 0, 6, 0 )
[72]	1573	ENDIF
	1574	ENDIF
	1575	ENDIF
	1576
	1577	!
[1]	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
[87]	1592	dopr_unit(i) = 'm/s'
[1]	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
[87]	1602	dopr_unit(i) = 'm/s'
	1603	hom(:,2,2,:) = SPREAD( zu, 2, statistic_regions+1 )
[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
[87]	1612	dopr_unit(i) = 'm/s'
	1613	hom(:,2,3,:) = SPREAD( zw, 2, statistic_regions+1 )
[1]	1614
	1615	CASE ( 'pt', '#pt' )
	1616	IF ( .NOT. cloud_physics ) THEN
	1617	dopr_index(i) = 4
[87]	1618	dopr_unit(i) = 'K'
[1]	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 )
[87]	1623	hom(nzb,2,7,:) = 0.0 ! because zu(nzb) is negative
[1]	1624	data_output_pr(i) = data_output_pr(i)(2:)
	1625	ENDIF
	1626	ELSE
	1627	dopr_index(i) = 43
[87]	1628	dopr_unit(i) = 'K'
[1]	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 )
[87]	1633	hom(nzb,2,28,:) = 0.0 ! because zu(nzb) is negative
[1]	1634	data_output_pr(i) = data_output_pr(i)(2:)
	1635	ENDIF
	1636	ENDIF
	1637
	1638	CASE ( 'e' )
	1639	dopr_index(i) = 8
[87]	1640	dopr_unit(i) = 'm2/s2'
[1]	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
[87]	1646	dopr_unit(i) = 'm2/s'
[1]	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
[87]	1657	dopr_unit(i) = 'm2/s'
[1]	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
[87]	1668	dopr_unit(i) = 'm'
[1]	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
[87]	1679	dopr_unit(i) = 'm2/s2'
[1]	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
[87]	1685	dopr_unit(i) = 'm2/s2'
[1]	1686	hom(:,2,13,:) = SPREAD( zw, 2, statistic_regions+1 )
	1687
	1688	CASE ( 'w"v"' )
	1689	dopr_index(i) = 14
[87]	1690	dopr_unit(i) = 'm2/s2'
[1]	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
[87]	1696	dopr_unit(i) = 'm2/s2'
[1]	1697	hom(:,2,15,:) = SPREAD( zw, 2, statistic_regions+1 )
	1698
	1699	CASE ( 'w"pt"' )
	1700	dopr_index(i) = 16
[87]	1701	dopr_unit(i) = 'K m/s'
[1]	1702	hom(:,2,16,:) = SPREAD( zw, 2, statistic_regions+1 )
	1703
	1704	CASE ( 'wpt' )
	1705	dopr_index(i) = 17
[87]	1706	dopr_unit(i) = 'K m/s'
[1]	1707	hom(:,2,17,:) = SPREAD( zw, 2, statistic_regions+1 )
	1708
	1709	CASE ( 'wpt' )
	1710	dopr_index(i) = 18
[87]	1711	dopr_unit(i) = 'K m/s'
[1]	1712	hom(:,2,18,:) = SPREAD( zw, 2, statistic_regions+1 )
	1713
	1714	CASE ( 'wu' )
	1715	dopr_index(i) = 19
[87]	1716	dopr_unit(i) = 'm2/s2'
[1]	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
[87]	1722	dopr_unit(i) = 'm2/s2'
[1]	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
[87]	1728	dopr_unit(i) = 'K m/s'
[1]	1729	hom(:,2,21,:) = SPREAD( zw, 2, statistic_regions+1 )
	1730
	1731	CASE ( 'wptBC' )
	1732	dopr_index(i) = 22
[87]	1733	dopr_unit(i) = 'K m/s'
[1]	1734	hom(:,2,22,:) = SPREAD( zw, 2, statistic_regions+1 )
	1735
[96]	1736	CASE ( 'sa', '#sa' )
	1737	IF ( .NOT. ocean ) THEN
[215]	1738	message_string = 'data_output_pr = ' // &
	1739	TRIM( data_output_pr(i) ) // ' is not imp' // &
	1740	'lemented for ocean = .FALSE.'
[226]	1741	CALL message( 'check_parameters', 'PA0091', 1, 2, 0, 6, 0 )
[96]	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
[1]	1754	CASE ( 'u*2' )
	1755	dopr_index(i) = 30
[87]	1756	dopr_unit(i) = 'm2/s2'
[1]	1757	hom(:,2,30,:) = SPREAD( zu, 2, statistic_regions+1 )
	1758
	1759	CASE ( 'v*2' )
	1760	dopr_index(i) = 31
[87]	1761	dopr_unit(i) = 'm2/s2'
[1]	1762	hom(:,2,31,:) = SPREAD( zu, 2, statistic_regions+1 )
	1763
	1764	CASE ( 'w*2' )
	1765	dopr_index(i) = 32
[87]	1766	dopr_unit(i) = 'm2/s2'
[1]	1767	hom(:,2,32,:) = SPREAD( zw, 2, statistic_regions+1 )
	1768
	1769	CASE ( 'pt*2' )
	1770	dopr_index(i) = 33
[87]	1771	dopr_unit(i) = 'K2'
[1]	1772	hom(:,2,33,:) = SPREAD( zu, 2, statistic_regions+1 )
	1773
	1774	CASE ( 'e*' )
	1775	dopr_index(i) = 34
[87]	1776	dopr_unit(i) = 'm2/s2'
[1]	1777	hom(:,2,34,:) = SPREAD( zu, 2, statistic_regions+1 )
	1778
	1779	CASE ( 'w2pt' )
	1780	dopr_index(i) = 35
[87]	1781	dopr_unit(i) = 'K m2/s2'
[1]	1782	hom(:,2,35,:) = SPREAD( zw, 2, statistic_regions+1 )
	1783
	1784	CASE ( 'wpt2' )
	1785	dopr_index(i) = 36
[87]	1786	dopr_unit(i) = 'K2 m/s'
[1]	1787	hom(:,2,36,:) = SPREAD( zw, 2, statistic_regions+1 )
	1788
	1789	CASE ( 'we' )
	1790	dopr_index(i) = 37
[87]	1791	dopr_unit(i) = 'm3/s3'
[1]	1792	hom(:,2,37,:) = SPREAD( zw, 2, statistic_regions+1 )
	1793
	1794	CASE ( 'w*3' )
	1795	dopr_index(i) = 38
[87]	1796	dopr_unit(i) = 'm3/s3'
[1]	1797	hom(:,2,38,:) = SPREAD( zw, 2, statistic_regions+1 )
	1798
	1799	CASE ( 'Sw' )
	1800	dopr_index(i) = 39
[89]	1801	dopr_unit(i) = 'none'
[1]	1802	hom(:,2,39,:) = SPREAD( zw, 2, statistic_regions+1 )
	1803
[232]	1804	CASE ( 'p' )
	1805	dopr_index(i) = 40
	1806	dopr_unit(i) = 'Pa'
	1807	hom(:,2,40,:) = SPREAD( zu, 2, statistic_regions+1 )
	1808
[1]	1809	CASE ( 'q', '#q' )
[108]	1810	IF ( .NOT. humidity ) THEN
[215]	1811	message_string = 'data_output_pr = ' // &
	1812	TRIM( data_output_pr(i) ) // ' is not imp' // &
	1813	'lemented for humidity = .FALSE.'
[226]	1814	CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
[1]	1815	ELSE
	1816	dopr_index(i) = 41
[87]	1817	dopr_unit(i) = 'kg/kg'
	1818	hom(:,2,41,:) = SPREAD( zu, 2, statistic_regions+1 )
[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
[215]	1829	message_string = 'data_output_pr = ' // &
	1830	TRIM( data_output_pr(i) ) // ' is not imp' // &
	1831	'lemented for passive_scalar = .FALSE.'
[226]	1832	CALL message( 'check_parameters', 'PA0093', 1, 2, 0, 6, 0 )
[1]	1833	ELSE
	1834	dopr_index(i) = 41
[87]	1835	dopr_unit(i) = 'kg/m3'
	1836	hom(:,2,41,:) = SPREAD( zu, 2, statistic_regions+1 )
[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
[87]	1848	dopr_unit(i) = 'kg/kg'
	1849	hom(:,2,41,:) = SPREAD( zu, 2, statistic_regions+1 )
[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
[87]	1858	dopr_unit(i) = 'kg/kg'
	1859	hom(:,2,42,:) = SPREAD( zu, 2, statistic_regions+1 )
[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
[215]	1870	message_string = 'data_output_pr = ' // &
	1871	TRIM( data_output_pr(i) ) // ' is not imp' // &
	1872	'lemented for cloud_physics = .FALSE.'
[226]	1873	CALL message( 'check_parameters', 'PA0094', 1, 2, 0, 6, 0 )
[1]	1874	ELSE
	1875	dopr_index(i) = 4
[87]	1876	dopr_unit(i) = 'K'
[1]	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
[87]	1888	dopr_unit(i) = 'K'
	1889	hom(:,2,44,:) = SPREAD( zu, 2, statistic_regions+1 )
[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
[87]	1899	dopr_unit(i) = 'K m/s'
[1]	1900	hom(:,2,45,:) = SPREAD( zw, 2, statistic_regions+1 )
	1901
	1902	CASE ( 'wvpt' )
	1903	dopr_index(i) = 46
[87]	1904	dopr_unit(i) = 'K m/s'
[1]	1905	hom(:,2,46,:) = SPREAD( zw, 2, statistic_regions+1 )
	1906
	1907	CASE ( 'wvpt' )
	1908	dopr_index(i) = 47
[87]	1909	dopr_unit(i) = 'K m/s'
[1]	1910	hom(:,2,47,:) = SPREAD( zw, 2, statistic_regions+1 )
	1911
	1912	CASE ( 'w"q"' )
[108]	1913	IF ( .NOT. humidity ) THEN
[215]	1914	message_string = 'data_output_pr = ' // &
	1915	TRIM( data_output_pr(i) ) // ' is not imp' // &
	1916	'lemented for humidity = .FALSE.'
[226]	1917	CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
[1]	1918	ELSE
	1919	dopr_index(i) = 48
[87]	1920	dopr_unit(i) = 'kg/kg m/s'
[1]	1921	hom(:,2,48,:) = SPREAD( zw, 2, statistic_regions+1 )
	1922	ENDIF
	1923
	1924	CASE ( 'wq' )
[108]	1925	IF ( .NOT. humidity ) THEN
[215]	1926	message_string = 'data_output_pr = ' // &
	1927	TRIM( data_output_pr(i) ) // ' is not imp' // &
	1928	'lemented for humidity = .FALSE.'
[226]	1929	CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
[1]	1930	ELSE
	1931	dopr_index(i) = 49
[87]	1932	dopr_unit(i) = 'kg/kg m/s'
[1]	1933	hom(:,2,49,:) = SPREAD( zw, 2, statistic_regions+1 )
	1934	ENDIF
	1935
	1936	CASE ( 'wq' )
[108]	1937	IF ( .NOT. humidity ) THEN
[215]	1938	message_string = 'data_output_pr = ' // &
	1939	TRIM( data_output_pr(i) ) // ' is not imp' // &
	1940	'lemented for humidity = .FALSE.'
[226]	1941	CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
[1]	1942	ELSE
	1943	dopr_index(i) = 50
[87]	1944	dopr_unit(i) = 'kg/kg m/s'
[1]	1945	hom(:,2,50,:) = SPREAD( zw, 2, statistic_regions+1 )
	1946	ENDIF
	1947
	1948	CASE ( 'w"s"' )
	1949	IF ( .NOT. passive_scalar ) THEN
[215]	1950	message_string = 'data_output_pr = ' // &
	1951	TRIM( data_output_pr(i) ) // ' is not imp' // &
	1952	'lemented for passive_scalar = .FALSE.'
[226]	1953	CALL message( 'check_parameters', 'PA0093', 1, 2, 0, 6, 0 )
[1]	1954	ELSE
	1955	dopr_index(i) = 48
[87]	1956	dopr_unit(i) = 'kg/m3 m/s'
[1]	1957	hom(:,2,48,:) = SPREAD( zw, 2, statistic_regions+1 )
	1958	ENDIF
	1959
	1960	CASE ( 'ws' )
	1961	IF ( .NOT. passive_scalar ) THEN
[215]	1962	message_string = 'data_output_pr = ' // &
	1963	TRIM( data_output_pr(i) ) // ' is not imp' // &
	1964	'lemented for passive_scalar = .FALSE.'
[226]	1965	CALL message( 'check_parameters', 'PA0093', 1, 2, 0, 6, 0 )
[1]	1966	ELSE
	1967	dopr_index(i) = 49
[87]	1968	dopr_unit(i) = 'kg/m3 m/s'
[1]	1969	hom(:,2,49,:) = SPREAD( zw, 2, statistic_regions+1 )
	1970	ENDIF
	1971
	1972	CASE ( 'ws' )
	1973	IF ( .NOT. passive_scalar ) THEN
[215]	1974	message_string = 'data_output_pr = ' // &
	1975	TRIM( data_output_pr(i) ) // ' is not imp' // &
	1976	'lemented for passive_scalar = .FALSE.'
[226]	1977	CALL message( 'check_parameters', 'PA0093', 1, 2, 0, 6, 0 )
[1]	1978	ELSE
	1979	dopr_index(i) = 50
[87]	1980	dopr_unit(i) = 'kg/m3 m/s'
[1]	1981	hom(:,2,50,:) = SPREAD( zw, 2, statistic_regions+1 )
	1982	ENDIF
	1983
	1984	CASE ( 'w"qv"' )
[75]	1985	IF ( humidity .AND. .NOT. cloud_physics ) &
[1]	1986	THEN
	1987	dopr_index(i) = 48
[87]	1988	dopr_unit(i) = 'kg/kg m/s'
[1]	1989	hom(:,2,48,:) = SPREAD( zw, 2, statistic_regions+1 )
[75]	1990	ELSEIF( humidity .AND. cloud_physics ) THEN
[1]	1991	dopr_index(i) = 51
[87]	1992	dopr_unit(i) = 'kg/kg m/s'
[1]	1993	hom(:,2,51,:) = SPREAD( zw, 2, statistic_regions+1 )
	1994	ELSE
[215]	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.'
[226]	1999	CALL message( 'check_parameters', 'PA0095', 1, 2, 0, 6, 0 )
[1]	2000	ENDIF
	2001
	2002	CASE ( 'wqv' )
[75]	2003	IF ( humidity .AND. .NOT. cloud_physics ) &
[1]	2004	THEN
	2005	dopr_index(i) = 49
[87]	2006	dopr_unit(i) = 'kg/kg m/s'
[1]	2007	hom(:,2,49,:) = SPREAD( zw, 2, statistic_regions+1 )
[75]	2008	ELSEIF( humidity .AND. cloud_physics ) THEN
[1]	2009	dopr_index(i) = 52
[87]	2010	dopr_unit(i) = 'kg/kg m/s'
[1]	2011	hom(:,2,52,:) = SPREAD( zw, 2, statistic_regions+1 )
	2012	ELSE
[215]	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.'
[226]	2017	CALL message( 'check_parameters', 'PA0095', 1, 2, 0, 6, 0 )
[1]	2018	ENDIF
	2019
	2020	CASE ( 'wqv' )
[75]	2021	IF ( humidity .AND. .NOT. cloud_physics ) &
[1]	2022	THEN
	2023	dopr_index(i) = 50
[87]	2024	dopr_unit(i) = 'kg/kg m/s'
[1]	2025	hom(:,2,50,:) = SPREAD( zw, 2, statistic_regions+1 )
[75]	2026	ELSEIF( humidity .AND. cloud_physics ) THEN
[1]	2027	dopr_index(i) = 53
[87]	2028	dopr_unit(i) = 'kg/kg m/s'
[1]	2029	hom(:,2,53,:) = SPREAD( zw, 2, statistic_regions+1 )
	2030	ELSE
[215]	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.'
[226]	2035	CALL message( 'check_parameters', 'PA0095', 1, 2, 0, 6, 0 )
[1]	2036	ENDIF
	2037
	2038	CASE ( 'ql' )
	2039	IF ( .NOT. cloud_physics .AND. .NOT. cloud_droplets ) THEN
[215]	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.'
[226]	2044	CALL message( 'check_parameters', 'PA0096', 1, 2, 0, 6, 0 )
[1]	2045	ELSE
	2046	dopr_index(i) = 54
[87]	2047	dopr_unit(i) = 'kg/kg'
[1]	2048	hom(:,2,54,:) = SPREAD( zu, 2, statistic_regions+1 )
	2049	ENDIF
	2050
	2051	CASE ( 'wuu*/dz' )
	2052	dopr_index(i) = 55
[87]	2053	dopr_unit(i) = 'm2/s3'
[1]	2054	hom(:,2,55,:) = SPREAD( zu, 2, statistic_regions+1 )
	2055
	2056	CASE ( 'wp/dz' )
	2057	dopr_index(i) = 56
[87]	2058	dopr_unit(i) = 'm2/s3'
[106]	2059	hom(:,2,56,:) = SPREAD( zw, 2, statistic_regions+1 )
[1]	2060
	2061	CASE ( 'w"e/dz' )
	2062	dopr_index(i) = 57
[87]	2063	dopr_unit(i) = 'm2/s3'
[1]	2064	hom(:,2,57,:) = SPREAD( zu, 2, statistic_regions+1 )
	2065
	2066	CASE ( 'u"pt"' )
	2067	dopr_index(i) = 58
[87]	2068	dopr_unit(i) = 'K m/s'
[1]	2069	hom(:,2,58,:) = SPREAD( zu, 2, statistic_regions+1 )
	2070
	2071	CASE ( 'upt' )
	2072	dopr_index(i) = 59
[87]	2073	dopr_unit(i) = 'K m/s'
[1]	2074	hom(:,2,59,:) = SPREAD( zu, 2, statistic_regions+1 )
	2075
	2076	CASE ( 'upt_t' )
	2077	dopr_index(i) = 60
[87]	2078	dopr_unit(i) = 'K m/s'
[1]	2079	hom(:,2,60,:) = SPREAD( zu, 2, statistic_regions+1 )
	2080
	2081	CASE ( 'v"pt"' )
	2082	dopr_index(i) = 61
[87]	2083	dopr_unit(i) = 'K m/s'
[1]	2084	hom(:,2,61,:) = SPREAD( zu, 2, statistic_regions+1 )
	2085
	2086	CASE ( 'vpt' )
	2087	dopr_index(i) = 62
[87]	2088	dopr_unit(i) = 'K m/s'
[1]	2089	hom(:,2,62,:) = SPREAD( zu, 2, statistic_regions+1 )
	2090
	2091	CASE ( 'vpt_t' )
	2092	dopr_index(i) = 63
[87]	2093	dopr_unit(i) = 'K m/s'
[1]	2094	hom(:,2,63,:) = SPREAD( zu, 2, statistic_regions+1 )
	2095
[96]	2096	CASE ( 'rho' )
[388]	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
[1]	2107
[96]	2108	CASE ( 'w"sa"' )
	2109	IF ( .NOT. ocean ) THEN
[215]	2110	message_string = 'data_output_pr = ' // &
	2111	TRIM( data_output_pr(i) ) // ' is not imp' // &
	2112	'lemented for ocean = .FALSE.'
[226]	2113	CALL message( 'check_parameters', 'PA0091', 1, 2, 0, 6, 0 )
[96]	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
[215]	2122	message_string = 'data_output_pr = ' // &
	2123	TRIM( data_output_pr(i) ) // ' is not imp' // &
	2124	'lemented for ocean = .FALSE.'
[226]	2125	CALL message( 'check_parameters', 'PA0091', 1, 2, 0, 6, 0 )
[96]	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
[215]	2134	message_string = 'data_output_pr = ' // &
	2135	TRIM( data_output_pr(i) ) // ' is not imp' // &
	2136	'lemented for ocean = .FALSE.'
[226]	2137	CALL message( 'check_parameters', 'PA0091', 1, 2, 0, 6, 0 )
[96]	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
[106]	2144	CASE ( 'wp' )
	2145	dopr_index(i) = 68
	2146	dopr_unit(i) = 'm3/s3'
	2147	hom(:,2,68,:) = SPREAD( zu, 2, statistic_regions+1 )
[96]	2148
[106]	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
[197]	2154	CASE ( 'q*2' )
	2155	IF ( .NOT. humidity ) THEN
[215]	2156	message_string = 'data_output_pr = ' // &
	2157	TRIM( data_output_pr(i) ) // ' is not imp' // &
	2158	'lemented for humidity = .FALSE.'
[226]	2159	CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
[197]	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
[106]	2165
[388]	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
[1]	2183	CASE DEFAULT
[87]	2184
	2185	CALL user_check_data_output_pr( data_output_pr(i), i, unit )
	2186
	2187	IF ( unit == 'illegal' ) THEN
[215]	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) ) // '"'
[226]	2192	CALL message( 'check_parameters', 'PA0097', 1, 2, 0, 6, 0 )
[215]	2193	ELSE
	2194	message_string = 'illegal value for data_output_pr = "' // &
	2195	TRIM( data_output_pr(i) ) // '"'
[226]	2196	CALL message( 'check_parameters', 'PA0098', 1, 2, 0, 6, 0 )
[87]	2197	ENDIF
[1]	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.
[94]	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
[215]	2232	WRITE( message_string, * ) 'z_max_do1d = ', z_max_do1d, ' must ', &
	2233	'be >= ', zu(nzb+1), ' or <= ', zu(nzt+1)
[226]	2234	CALL message( 'check_parameters', 'PA0099', 1, 2, 0, 6, 0 )
[94]	2235	ELSE
	2236	cross_uymax = z_max_do1d
	2237	ENDIF
[1]	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
[215]	2250	message_string = 'unknown normalization method cross_normali' // &
	2251	'zed_x = "' // TRIM( cross_normalized_x(i) ) // &
	2252	'"'
[226]	2253	CALL message( 'check_parameters', 'PA0100', 1, 2, 0, 6, 0 )
[1]	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
[215]	2262	message_string = 'unknown normalization method cross_normali' // &
	2263	'zed_y = "' // TRIM( cross_normalized_y(i) ) // &
	2264	'"'
[226]	2265	CALL message( 'check_parameters', 'PA0101', 1, 2, 0, 6, 0 )
[1]	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
[215]	2273	WRITE( message_string, * ) 'z_max_do1d_normalized = ', &
	2274	z_max_do1d_normalized, ' must be >= 0.0'
[226]	2275	CALL message( 'check_parameters', 'PA0101', 1, 2, 0, 6, 0 )
[1]	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
[215]	2289	message_string = 'number of output quantitities given by data' // &
	2290	'_output and data_output_user exceeds the limit of 100'
[226]	2291	CALL message( 'check_parameters', 'PA0102', 1, 2, 0, 6, 0 )
[1]	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
[215]	2332	message_string = 'output of "' // TRIM( var ) // '" requi' // &
	2333	'res constant_diffusion = .FALSE.'
[226]	2334	CALL message( 'check_parameters', 'PA0103', 1, 2, 0, 6, 0 )
[1]	2335	ENDIF
	2336	unit = 'm2/s2'
	2337
	2338	CASE ( 'pc', 'pr' )
	2339	IF ( .NOT. particle_advection ) THEN
[215]	2340	message_string = 'output of "' // TRIM( var ) // '" requir' // &
	2341	'es a "particles_par"-NAMELIST in the parameter file (PARIN)'
[226]	2342	CALL message( 'check_parameters', 'PA0104', 1, 2, 0, 6, 0 )
[1]	2343	ENDIF
	2344	IF ( TRIM( var ) == 'pc' ) unit = 'number'
	2345	IF ( TRIM( var ) == 'pr' ) unit = 'm'
	2346
	2347	CASE ( 'q', 'vpt' )
[75]	2348	IF ( .NOT. humidity ) THEN
[215]	2349	message_string = 'output of "' // TRIM( var ) // '" requi' // &
	2350	'res humidity = .TRUE.'
[226]	2351	CALL message( 'check_parameters', 'PA0105', 1, 2, 0, 6, 0 )
[1]	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
[215]	2358	message_string = 'output of "' // TRIM( var ) // '" requi' // &
	2359	'res cloud_physics = .TRUE. or cloud_droplets = .TRUE.'
[226]	2360	CALL message( 'check_parameters', 'PA0106', 1, 2, 0, 6, 0 )
[1]	2361	ENDIF
	2362	unit = 'kg/kg'
	2363
	2364	CASE ( 'ql_c', 'ql_v', 'ql_vp' )
	2365	IF ( .NOT. cloud_droplets ) THEN
[215]	2366	message_string = 'output of "' // TRIM( var ) // '" requi' // &
	2367	'res cloud_droplets = .TRUE.'
[226]	2368	CALL message( 'check_parameters', 'PA0107', 1, 2, 0, 6, 0 )
[1]	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
[215]	2376	message_string = 'output of "' // TRIM( var ) // '" requi' // &
	2377	'res cloud_physics = .TRUE.'
[226]	2378	CALL message( 'check_parameters', 'PA0108', 1, 2, 0, 6, 0 )
[1]	2379	ENDIF
	2380	unit = 'kg/kg'
	2381
[96]	2382	CASE ( 'rho' )
	2383	IF ( .NOT. ocean ) THEN
[215]	2384	message_string = 'output of "' // TRIM( var ) // '" requi' // &
	2385	'res ocean = .TRUE.'
[226]	2386	CALL message( 'check_parameters', 'PA0109', 1, 2, 0, 6, 0 )
[96]	2387	ENDIF
	2388	unit = 'kg/m3'
	2389
[1]	2390	CASE ( 's' )
	2391	IF ( .NOT. passive_scalar ) THEN
[215]	2392	message_string = 'output of "' // TRIM( var ) // '" requi' // &
	2393	'res passive_scalar = .TRUE.'
[226]	2394	CALL message( 'check_parameters', 'PA0110', 1, 2, 0, 6, 0 )
[1]	2395	ENDIF
	2396	unit = 'conc'
	2397
[96]	2398	CASE ( 'sa' )
	2399	IF ( .NOT. ocean ) THEN
[215]	2400	message_string = 'output of "' // TRIM( var ) // '" requi' // &
	2401	'res ocean = .TRUE.'
[226]	2402	CALL message( 'check_parameters', 'PA0109', 1, 2, 0, 6, 0 )
[96]	2403	ENDIF
	2404	unit = 'psu'
	2405
[354]	2406	CASE ( 'u', 't', 'lwp', 'pra', 'prr', 'qsws', 'shf', 'z0' )
[1]	2407	IF ( k == 0 .OR. data_output(i)(ilen-2:ilen) /= '_xy' ) THEN
[215]	2408	message_string = 'illegal value for data_output: "' // &
	2409	TRIM( var ) // '" & only 2d-horizontal ' // &
	2410	'cross sections are allowed for this value'
[226]	2411	CALL message( 'check_parameters', 'PA0111', 1, 2, 0, 6, 0 )
[1]	2412	ENDIF
	2413	IF ( TRIM( var ) == 'lwp*' .AND. .NOT. cloud_physics ) THEN
[215]	2414	message_string = 'output of "' // TRIM( var ) // '" requi' // &
	2415	'res cloud_physics = .TRUE.'
[226]	2416	CALL message( 'check_parameters', 'PA0108', 1, 2, 0, 6, 0 )
[1]	2417	ENDIF
[72]	2418	IF ( TRIM( var ) == 'pra*' .AND. .NOT. precipitation ) THEN
[215]	2419	message_string = 'output of "' // TRIM( var ) // '" requi' // &
	2420	'res precipitation = .TRUE.'
[226]	2421	CALL message( 'check_parameters', 'PA0112', 1, 2, 0, 6, 0 )
[72]	2422	ENDIF
	2423	IF ( TRIM( var ) == 'pra*' .AND. j == 1 ) THEN
[215]	2424	message_string = 'temporal averaging of precipitation ' // &
	2425	'amount "' // TRIM( var ) // '" is not possible'
[226]	2426	CALL message( 'check_parameters', 'PA0113', 1, 2, 0, 6, 0 )
[72]	2427	ENDIF
	2428	IF ( TRIM( var ) == 'prr*' .AND. .NOT. precipitation ) THEN
[215]	2429	message_string = 'output of "' // TRIM( var ) // '" requi' // &
	2430	'res precipitation = .TRUE.'
[226]	2431	CALL message( 'check_parameters', 'PA0112', 1, 2, 0, 6, 0 )
[72]	2432	ENDIF
[354]	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
[72]	2438
[354]	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'
[72]	2447
[1]	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
[215]	2461	IF ( data_output_user(1) /= ' ' ) THEN
	2462	message_string = 'illegal value for data_output or ' // &
	2463	'data_output_user = "' // TRIM( data_output(i) ) // '"'
[226]	2464	CALL message( 'check_parameters', 'PA0114', 1, 2, 0, 6, 0 )
[215]	2465	ELSE
	2466	message_string = 'illegal value for data_output =' // &
	2467	TRIM( data_output(i) ) // '"'
[226]	2468	CALL message( 'check_parameters', 'PA0115', 1, 2, 0, 6, 0 )
[1]	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	!
[376]	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	!
[308]	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
[1]	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
[215]	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))'
[226]	2556	CALL message( 'check_parameters', 'PA0116', 1, 2, 0, 6, 0 )
[1]	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
[215]	2569	message_string = 'do3d_compress = .TRUE. not allowed on host "' // &
	2570	TRIM( host ) // '"'
[226]	2571	CALL message( 'check_parameters', 'PA0117', 1, 2, 0, 6, 0 )
[1]	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
[215]	2580	WRITE( message_string, * ) 'illegal precision: do3d_comp_prec', &
	2581	'(', i, ') = "', TRIM(do3d_comp_prec(i)),'"'
[226]	2582	CALL message( 'check_parameters', 'PA0118', 1, 2, 0, 6, 0 )
[1]	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
[215]	2602	WRITE( message_string, * ) 'unknown variable "', &
	2603	TRIM( do3d_comp_prec(i) ), '" given for do3d_comp_prec(', &
	2604	i, ')'
[226]	2605	CALL message( 'check_parameters', 'PA0119', 1, 2, 0, 6, 0 )
[1]	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
[215]	2637	message_string = 'unknown value for data_output_format "' // &
	2638	TRIM( data_output_format(i) ) // '"'
[226]	2639	CALL message( 'check_parameters', 'PA0120', 1, 2, 0, 6, 0 )
[1]	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
[215]	2659	WRITE( message_string, * ) 'km_constant = ', km_constant, ' < 0.0'
[226]	2660	CALL message( 'check_parameters', 'PA0121', 1, 2, 0, 6, 0 )
[1]	2661	ELSE
	2662	IF ( prandtl_number < 0.0 ) THEN
[215]	2663	WRITE( message_string, * ) 'prandtl_number = ', prandtl_number, &
	2664	' < 0.0'
[226]	2665	CALL message( 'check_parameters', 'PA0122', 1, 2, 0, 6, 0 )
[1]	2666	ENDIF
	2667	constant_diffusion = .TRUE.
	2668
	2669	IF ( prandtl_layer ) THEN
[215]	2670	message_string = 'prandtl_layer is not allowed with fixed ' // &
	2671	'value of km'
[226]	2672	CALL message( 'check_parameters', 'PA0123', 1, 2, 0, 6, 0 )
[1]	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
[215]	2689	message_string = 'outflow_damping width out of range'
[226]	2690	CALL message( 'check_parameters', 'PA0124', 1, 2, 0, 6, 0 )
[1]	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
[215]	2702	message_string = 'outflow_damping width out of range'
[226]	2703	CALL message( 'check_parameters', 'PA0124', 1, 2, 0, 6, 0 )
[1]	2704	ENDIF
	2705	ENDIF
	2706
	2707	!
	2708	!-- Check value range for rif
	2709	IF ( rif_min >= rif_max ) THEN
[215]	2710	WRITE( message_string, * ) 'rif_min = ', rif_min, ' must be less ', &
	2711	'than rif_max = ', rif_max
[226]	2712	CALL message( 'check_parameters', 'PA0125', 1, 2, 0, 6, 0 )
[1]	2713	ENDIF
	2714
	2715	!
	2716	!-- Determine upper and lower hight level indices for random perturbations
[97]	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
[1]	2725	ELSEIF ( disturbance_level_b < zu(3) ) THEN
[215]	2726	WRITE( message_string, * ) 'disturbance_level_b = ', &
	2727	disturbance_level_b, ' must be >= ', zu(3), '(zu(3))'
[226]	2728	CALL message( 'check_parameters', 'PA0126', 1, 2, 0, 6, 0 )
[1]	2729	ELSEIF ( disturbance_level_b > zu(nzt-2) ) THEN
[215]	2730	WRITE( message_string, * ) 'disturbance_level_b = ', &
	2731	disturbance_level_b, ' must be <= ', zu(nzt-2), '(zu(nzt-2))'
[226]	2732	CALL message( 'check_parameters', 'PA0127', 1, 2, 0, 6, 0 )
[1]	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
[97]	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
[1]	2750	ELSEIF ( disturbance_level_t > zu(nzt-2) ) THEN
[215]	2751	WRITE( message_string, * ) 'disturbance_level_t = ', &
	2752	disturbance_level_t, ' must be <= ', zu(nzt-2), '(zu(nzt-2))'
[226]	2753	CALL message( 'check_parameters', 'PA0128', 1, 2, 0, 6, 0 )
[1]	2754	ELSEIF ( disturbance_level_t < disturbance_level_b ) THEN
[215]	2755	WRITE( message_string, * ) 'disturbance_level_t = ', &
	2756	disturbance_level_t, ' must be >= disturbance_level_b = ', &
	2757	disturbance_level_b
[226]	2758	CALL message( 'check_parameters', 'PA0129', 1, 2, 0, 6, 0 )
[1]	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
[215]	2773	WRITE( message_string, * ) 'disturbance_level_ind_t = ', &
	2774	disturbance_level_ind_t, ' must be >= disturbance_level_b = ', &
	2775	disturbance_level_b
[226]	2776	CALL message( 'check_parameters', 'PA0130', 1, 2, 0, 6, 0 )
[1]	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
[215]	2792	message_string = 'inflow_disturbance_begin out of range'
[226]	2793	CALL message( 'check_parameters', 'PA0131', 1, 2, 0, 6, 0 )
[1]	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
[215]	2800	message_string = 'inflow_disturbance_end out of range'
[226]	2801	CALL message( 'check_parameters', 'PA0132', 1, 2, 0, 6, 0 )
[1]	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
[215]	2809	message_string = 'inflow_disturbance_begin out of range'
[226]	2810	CALL message( 'check_parameters', 'PA0131', 1, 2, 0, 6, 0 )
[1]	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
[215]	2817	message_string = 'inflow_disturbance_end out of range'
[226]	2818	CALL message( 'check_parameters', 'PA0132', 1, 2, 0, 6, 0 )
[1]	2819	ENDIF
	2820	ENDIF
	2821
[73]	2822	IF ( bc_lr == 'radiation/dirichlet' ) THEN
[1]	2823	dist_nxr = nx - inflow_disturbance_begin
	2824	dist_nxl(1) = nx - inflow_disturbance_end
[73]	2825	ELSEIF ( bc_lr == 'dirichlet/radiation' ) THEN
[1]	2826	dist_nxl = inflow_disturbance_begin
	2827	dist_nxr(1) = inflow_disturbance_end
[73]	2828	ENDIF
	2829	IF ( bc_ns == 'dirichlet/radiation' ) THEN
[1]	2830	dist_nyn = ny - inflow_disturbance_begin
	2831	dist_nys(1) = ny - inflow_disturbance_end
[73]	2832	ELSEIF ( bc_ns == 'radiation/dirichlet' ) THEN
[1]	2833	dist_nys = inflow_disturbance_begin
	2834	dist_nyn(1) = inflow_disturbance_end
	2835	ENDIF
	2836
	2837	!
[151]	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
[215]	2841	message_string = 'turbulent_inflow = .T. requires a Dirichlet ' // &
	2842	'condition at the inflow boundary'
[226]	2843	CALL message( 'check_parameters', 'PA0133', 1, 2, 0, 6, 0 )
[151]	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
[215]	2855	WRITE( message_string, * ) 'illegal value for recycling_width:', &
	2856	' ', recycling_width
[226]	2857	CALL message( 'check_parameters', 'PA0134', 1, 2, 0, 6, 0 )
[151]	2858	ENDIF
	2859	ENDIF
	2860	!
	2861	!-- Calculate the index
	2862	recycling_plane = recycling_width / dx
	2863	ENDIF
	2864
	2865	!
[1]	2866	!-- Check random generator
	2867	IF ( random_generator /= 'system-specific' .AND. &
	2868	random_generator /= 'numerical-recipes' ) THEN
[215]	2869	message_string = 'unknown random generator: random_generator = "' // &
	2870	TRIM( random_generator ) // '"'
[226]	2871	CALL message( 'check_parameters', 'PA0135', 1, 2, 0, 6, 0 )
[1]	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
[215]	2881	WRITE( message_string, * ) 'damp_level_1d = ', damp_level_1d, &
	2882	' must be > 0.0 and < ', zu(nzt+1), '(zu(nzt+1))'
[226]	2883	CALL message( 'check_parameters', 'PA0136', 1, 2, 0, 6, 0 )
[1]	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
[215]	2893
[1]	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
[215]	2898	message_string = 'mixing_length_1d = "' // TRIM( mixing_length_1d ) // &
	2899	'" is unknown'
[226]	2900	CALL message( 'check_parameters', 'PA0137', 1, 2, 0, 6, 0 )
[1]	2901	ENDIF
	2902	IF ( TRIM( dissipation_1d ) /= 'as_in_3d_model' .AND. &
	2903	TRIM( dissipation_1d ) /= 'detering' ) THEN
[215]	2904	message_string = 'dissipation_1d = "' // TRIM( dissipation_1d ) // &
	2905	'" is unknown'
[226]	2906	CALL message( 'check_parameters', 'PA0138', 1, 2, 0, 6, 0 )
[1]	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
[291]	2913	IF ( restart_time > time_since_reference_point ) THEN
	2914	time_restart = restart_time
	2915	ENDIF
[1]	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
[215]	2940	WRITE( message_string, * ) 'termination_time_needed = ', &
	2941	termination_time_needed, ' must be > 30.0 on host "', &
	2942	TRIM( host ), '"'
[226]	2943	CALL message( 'check_parameters', 'PA0139', 1, 2, 0, 6, 0 )
[1]	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
[215]	2951	WRITE( message_string, * ) 'termination_time_needed = ', &
	2952	termination_time_needed, ' should be >= 300.0 on host "', &
	2953	TRIM( host ), '"'
[226]	2954	CALL message( 'check_parameters', 'PA0140', 1, 2, 0, 6, 0 )
[1]	2955	ENDIF
	2956	ENDIF
	2957
[217]	2958	!
[240]	2959	!-- Check pressure gradient conditions
	2960	IF ( dp_external .AND. conserve_volume_flow ) THEN
[388]	2961	WRITE( message_string, * ) 'Both dp_external and conserve_volume_flo', &
	2962	'w are .TRUE. but one of them must be .FALSE.'
[240]	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
[388]	2972	WRITE( message_string, * ) 'dp_external is .TRUE. but dpdxy is ze', &
	2973	'ro, i.e. the external pressure gradient & will not be applied'
[240]	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
[241]	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
[240]	3018
	3019	!
[264]	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	!
[217]	3049	!-- Check &userpar parameters
	3050	CALL user_check_parameters
[1]	3051
[217]	3052
[1]	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: