Home

Context Navigation

check_parameters.f90 @ 4239

Last change on this file since 4239 was 4234, checked in by gronemeier, 5 years ago

Check if a cross section is specified if any output cross-section quantity is given (check_parameters)

Property svn:keywords set to Id

Property svn:mergeinfo set to False

/palm/branches/chemistry/SOURCE/check_parameters.f90	2047-3190,3218-3297
/palm/branches/forwind/SOURCE/check_parameters.f90	1564-1913
/palm/branches/mosaik_M2/check_parameters.f90	2360-3471
/palm/branches/palm4u/SOURCE/check_parameters.f90	2540-2692
/palm/branches/rans/SOURCE/check_parameters.f90	2078-3128
/palm/branches/resler/SOURCE/check_parameters.f90	2023-3336
/palm/branches/salsa/SOURCE/check_parameters.f90	2503-3581

File size: 139.1 KB

Line
1	!> @file check_parameters.f90
2	!------------------------------------------------------------------------------!
3	! This file is part of the PALM model system.
4	!
5	! PALM is free software: you can redistribute it and/or modify it under the
6	! terms of the GNU General Public License as published by the Free Software
7	! Foundation, either version 3 of the License, or (at your option) any later
8	! version.
9	!
10	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
13	!
14	! You should have received a copy of the GNU General Public License along with
15	! PALM. If not, see <http://www.gnu.org/licenses/>.
16	!
17	! Copyright 1997-2019 Leibniz Universitaet Hannover
18	!------------------------------------------------------------------------------!
19	!
20	! Current revisions:
21	! -----------------
22	!
23	!
24	! Former revisions:
25	! -----------------
26	! 4172 2019-08-20 11:55:33Z oliver.maas
27	! Check if a cross section is specified if any output cross-section quantity
28	! is given
29	!
30	! 11:55:33Z oliver.maas
31	! Overwrite rotation_angle from namelist by value from static driver
32	!
33	! 11:55:33Z oliver.maas
34	! removed conversion from recycle_absolute_quantities to raq, added check and
35	! error message for correct input of recycling_method_for_thermodynamic_quantities
36	!
37	! 11:55:33Z oliver.maas
38	! Corrected "Former revisions" section
39	!
40	! 11:55:33Z oliver.maas
41	! bugfix error message: replaced PA184 by PA0184
42	!
43	! 11:55:33Z oliver.maas
44	! added conversion from recycle_absolute_quantities to raq for recycling of
45	! absolute quantities and added error message PA184 for not implemented quantities
46	!
47	! 4142 2019-08-05 12:38:31Z suehring
48	! Consider spinup in number of output timesteps for averaged 2D output (merge
49	! from branch resler).
50	!
51	! 4069 2019-07-01 14:05:51Z Giersch
52	! Masked output running index mid has been introduced as a local variable to
53	! avoid runtime error (Loop variable has been modified) in time_integration
54	!
55	! 4048 2019-06-21 21:00:21Z knoop
56	! Moved tcm_check_data_output to module_interface
57	!
58	! 4039 2019-06-18 10:32:41Z suehring
59	! Modularize diagnostic output
60	!
61	! 4017 2019-06-06 12:16:46Z schwenkel
62	! output of turbulence intensity added
63	!
64	! 3933 2019-04-25 12:33:20Z kanani
65	! Alphabetical resorting in CASE, condense settings for theta_2m* into one IF clause
66	!
67	! 3885 2019-04-11 11:29:34Z kanani
68	! Changes related to global restructuring of location messages and introduction
69	! of additional debug messages
70	!
71	! 3766 2019-02-26 16:23:41Z raasch
72	! trim added to avoid truncation compiler warnings
73	!
74	! 3761 2019-02-25 15:31:42Z raasch
75	! unused variables removed
76	!
77	! 3735 2019-02-12 09:52:40Z dom_dwd_user
78	! Passing variable j (averaged output?) to
79	! module_interface.f90:chem_check_data_output.
80	!
81	! 3705 2019-01-29 19:56:39Z suehring
82	! bugfix: renamed thetav_t to vtheta_t
83	!
84	! 3702 2019-01-28 13:19:30Z gronemeier
85	! most_method removed
86	!
87	! 3655 2019-01-07 16:51:22Z knoop
88	! Formatting
89	!
90	! Revision 1.1 1997/08/26 06:29:23 raasch
91	! Initial revision
92	!
93	!
94	! Description:
95	! ------------
96	!> Check control parameters and deduce further quantities.
97	!------------------------------------------------------------------------------!
98	SUBROUTINE check_parameters
99
100
101	USE arrays_3d
102
103	USE basic_constants_and_equations_mod
104
105	USE bulk_cloud_model_mod, &
106	ONLY: bulk_cloud_model
107
108	USE chem_modules
109
110	USE chemistry_model_mod, &
111	ONLY: chem_boundary_conds
112
113	USE control_parameters
114
115	USE grid_variables
116
117	USE kinds
118
119	USE indices
120
121	USE model_1d_mod, &
122	ONLY: damp_level_1d, damp_level_ind_1d
123
124	USE module_interface, &
125	ONLY: module_interface_check_parameters, &
126	module_interface_check_data_output_ts, &
127	module_interface_check_data_output_pr, &
128	module_interface_check_data_output
129
130	USE netcdf_data_input_mod, &
131	ONLY: init_model, input_pids_static, netcdf_data_input_check_dynamic, &
132	netcdf_data_input_check_static
133
134	USE netcdf_interface, &
135	ONLY: dopr_unit, do2d_unit, do3d_unit, netcdf_data_format, &
136	netcdf_data_format_string, dots_unit, heatflux_output_unit, &
137	waterflux_output_unit, momentumflux_output_unit, &
138	dots_max, dots_num, dots_label
139
140	USE particle_attributes, &
141	ONLY: particle_advection, use_sgs_for_particles
142
143	USE pegrid
144
145	USE pmc_interface, &
146	ONLY: cpl_id, nested_run
147
148	USE profil_parameter
149
150	USE statistics
151
152	USE subsidence_mod
153
154	USE transpose_indices
155
156	USE vertical_nesting_mod, &
157	ONLY: vnested, &
158	vnest_check_parameters
159
160
161	IMPLICIT NONE
162
163	CHARACTER (LEN=varnamelength) :: var !< variable name
164	CHARACTER (LEN=7) :: unit !< unit of variable
165	CHARACTER (LEN=8) :: date !< current date string
166	CHARACTER (LEN=10) :: time !< current time string
167	CHARACTER (LEN=20) :: ensemble_string !< string containing number of ensemble member
168	CHARACTER (LEN=15) :: nest_string !< string containing id of nested domain
169	CHARACTER (LEN=40) :: coupling_string !< string containing type of coupling
170	CHARACTER (LEN=100) :: action !< flag string
171
172	INTEGER(iwp) :: i !< loop index
173	INTEGER(iwp) :: ilen !< string length
174	INTEGER(iwp) :: j !< loop index
175	INTEGER(iwp) :: k !< loop index
176	INTEGER(iwp) :: kk !< loop index
177	INTEGER(iwp) :: mid !< masked output running index
178	INTEGER(iwp) :: netcdf_data_format_save !< initial value of netcdf_data_format
179	INTEGER(iwp) :: position !< index position of string
180
181	LOGICAL :: found !< flag, true if output variable is already marked for averaging
182
183	REAL(wp) :: dt_spinup_max !< maximum spinup timestep in nested domains
184	REAL(wp) :: gradient !< local gradient
185	REAL(wp) :: remote = 0.0_wp !< MPI id of remote processor
186	REAL(wp) :: spinup_time_max !< maximum spinup time in nested domains
187	REAL(wp) :: time_to_be_simulated_from_reference_point !< time to be simulated from reference point
188
189
190	CALL location_message( 'checking parameters', 'start' )
191	!
192	!-- At first, check static and dynamic input for consistency
193	CALL netcdf_data_input_check_dynamic
194	CALL netcdf_data_input_check_static
195	!
196	!-- Check for overlap combinations, which are not realized yet
197	IF ( transpose_compute_overlap .AND. numprocs == 1 ) THEN
198	message_string = 'transpose-compute-overlap not implemented for single PE runs'
199	CALL message( 'check_parameters', 'PA0000', 1, 2, 0, 6, 0 )
200	ENDIF
201
202	!
203	!-- Check the coupling mode
204	IF ( coupling_mode /= 'uncoupled' .AND. &
205	coupling_mode /= 'precursor_atmos' .AND. &
206	coupling_mode /= 'precursor_ocean' .AND. &
207	coupling_mode /= 'vnested_crse' .AND. &
208	coupling_mode /= 'vnested_fine' .AND. &
209	coupling_mode /= 'atmosphere_to_ocean' .AND. &
210	coupling_mode /= 'ocean_to_atmosphere' ) THEN
211	message_string = 'illegal coupling mode: ' // TRIM( coupling_mode )
212	CALL message( 'check_parameters', 'PA0002', 1, 2, 0, 6, 0 )
213	ENDIF
214
215	!
216	!-- Check if humidity is set to TRUE in case of the atmospheric run (for coupled runs)
217	IF ( coupling_mode == 'atmosphere_to_ocean' .AND. .NOT. humidity) THEN
218	message_string = ' Humidity has to be set to .T. in the _p3d file ' // &
219	'for coupled runs between ocean and atmosphere.'
220	CALL message( 'check_parameters', 'PA0476', 1, 2, 0, 6, 0 )
221	ENDIF
222
223	!
224	!-- Check dt_coupling, restart_time, dt_restart, end_time, dx, dy, nx and ny
225	IF ( coupling_mode /= 'uncoupled' .AND. &
226	coupling_mode(1:8) /= 'vnested_' .AND. &
227	coupling_mode /= 'precursor_atmos' .AND. &
228	coupling_mode /= 'precursor_ocean' ) THEN
229
230	IF ( dt_coupling == 9999999.9_wp ) THEN
231	message_string = 'dt_coupling is not set but required for coup' // &
232	'ling mode "' // TRIM( coupling_mode ) // '"'
233	CALL message( 'check_parameters', 'PA0003', 1, 2, 0, 6, 0 )
234	ENDIF
235
236	#if defined( __parallel )
237
238
239	IF ( myid == 0 ) THEN
240	CALL MPI_SEND( dt_coupling, 1, MPI_REAL, target_id, 11, comm_inter, &
241	ierr )
242	CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 11, comm_inter, &
243	status, ierr )
244	ENDIF
245	CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
246
247	IF ( dt_coupling /= remote ) THEN
248	WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
249	'": dt_coupling = ', dt_coupling, '& is not equal to ', &
250	'dt_coupling_remote = ', remote
251	CALL message( 'check_parameters', 'PA0004', 1, 2, 0, 6, 0 )
252	ENDIF
253	IF ( dt_coupling <= 0.0_wp ) THEN
254
255	IF ( myid == 0 ) THEN
256	CALL MPI_SEND( dt_max, 1, MPI_REAL, target_id, 19, comm_inter, ierr )
257	CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 19, comm_inter, &
258	status, ierr )
259	ENDIF
260	CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
261
262	dt_coupling = MAX( dt_max, remote )
263	WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
264	'": dt_coupling <= 0.0 & is not allowed and is reset to ', &
265	'MAX(dt_max(A,O)) = ', dt_coupling
266	CALL message( 'check_parameters', 'PA0005', 0, 1, 0, 6, 0 )
267	ENDIF
268
269	IF ( myid == 0 ) THEN
270	CALL MPI_SEND( restart_time, 1, MPI_REAL, target_id, 12, comm_inter, &
271	ierr )
272	CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 12, comm_inter, &
273	status, ierr )
274	ENDIF
275	CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
276
277	IF ( restart_time /= remote ) THEN
278	WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
279	'": restart_time = ', restart_time, '& is not equal to ', &
280	'restart_time_remote = ', remote
281	CALL message( 'check_parameters', 'PA0006', 1, 2, 0, 6, 0 )
282	ENDIF
283
284	IF ( myid == 0 ) THEN
285	CALL MPI_SEND( dt_restart, 1, MPI_REAL, target_id, 13, comm_inter, &
286	ierr )
287	CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 13, comm_inter, &
288	status, ierr )
289	ENDIF
290	CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
291
292	IF ( dt_restart /= remote ) THEN
293	WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
294	'": dt_restart = ', dt_restart, '& is not equal to ', &
295	'dt_restart_remote = ', remote
296	CALL message( 'check_parameters', 'PA0007', 1, 2, 0, 6, 0 )
297	ENDIF
298
299	time_to_be_simulated_from_reference_point = end_time-coupling_start_time
300
301	IF ( myid == 0 ) THEN
302	CALL MPI_SEND( time_to_be_simulated_from_reference_point, 1, &
303	MPI_REAL, target_id, 14, comm_inter, ierr )
304	CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 14, comm_inter, &
305	status, ierr )
306	ENDIF
307	CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
308
309	IF ( time_to_be_simulated_from_reference_point /= remote ) THEN
310	WRITE( message_string, * ) 'coupling mode "', TRIM( coupling_mode ), &
311	'": time_to_be_simulated_from_reference_point = ', &
312	time_to_be_simulated_from_reference_point, '& is not equal ', &
313	'to time_to_be_simulated_from_reference_point_remote = ', &
314	remote
315	CALL message( 'check_parameters', 'PA0008', 1, 2, 0, 6, 0 )
316	ENDIF
317
318	IF ( myid == 0 ) THEN
319	CALL MPI_SEND( dx, 1, MPI_REAL, target_id, 15, comm_inter, ierr )
320	CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 15, comm_inter, &
321	status, ierr )
322	ENDIF
323	CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
324
325
326	IF ( coupling_mode == 'atmosphere_to_ocean') THEN
327
328	IF ( dx < remote ) THEN
329	WRITE( message_string, * ) 'coupling mode "', &
330	TRIM( coupling_mode ), &
331	'": dx in Atmosphere is not equal to or not larger than dx in ocean'
332	CALL message( 'check_parameters', 'PA0009', 1, 2, 0, 6, 0 )
333	ENDIF
334
335	IF ( (nx_a+1)dx /= (nx_o+1)remote ) THEN
336	WRITE( message_string, * ) 'coupling mode "', &
337	TRIM( coupling_mode ), &
338	'": Domain size in x-direction is not equal in ocean and atmosphere'
339	CALL message( 'check_parameters', 'PA0010', 1, 2, 0, 6, 0 )
340	ENDIF
341
342	ENDIF
343
344	IF ( myid == 0) THEN
345	CALL MPI_SEND( dy, 1, MPI_REAL, target_id, 16, comm_inter, ierr )
346	CALL MPI_RECV( remote, 1, MPI_REAL, target_id, 16, comm_inter, &
347	status, ierr )
348	ENDIF
349	CALL MPI_BCAST( remote, 1, MPI_REAL, 0, comm2d, ierr)
350
351	IF ( coupling_mode == 'atmosphere_to_ocean') THEN
352
353	IF ( dy < remote ) THEN
354	WRITE( message_string, * ) 'coupling mode "', &
355	TRIM( coupling_mode ), &
356	'": dy in Atmosphere is not equal to or not larger than dy in ocean'
357	CALL message( 'check_parameters', 'PA0011', 1, 2, 0, 6, 0 )
358	ENDIF
359
360	IF ( (ny_a+1)dy /= (ny_o+1)remote ) THEN
361	WRITE( message_string, * ) 'coupling mode "', &
362	TRIM( coupling_mode ), &
363	'": Domain size in y-direction is not equal in ocean and atmosphere'
364	CALL message( 'check_parameters', 'PA0012', 1, 2, 0, 6, 0 )
365	ENDIF
366
367	IF ( MOD(nx_o+1,nx_a+1) /= 0 ) THEN
368	WRITE( message_string, * ) 'coupling mode "', &
369	TRIM( coupling_mode ), &
370	'": nx+1 in ocean is not divisible by nx+1 in', &
371	' atmosphere without remainder'
372	CALL message( 'check_parameters', 'PA0339', 1, 2, 0, 6, 0 )
373	ENDIF
374
375	IF ( MOD(ny_o+1,ny_a+1) /= 0 ) THEN
376	WRITE( message_string, * ) 'coupling mode "', &
377	TRIM( coupling_mode ), &
378	'": ny+1 in ocean is not divisible by ny+1 in', &
379	' atmosphere without remainder'
380
381	CALL message( 'check_parameters', 'PA0340', 1, 2, 0, 6, 0 )
382	ENDIF
383
384	ENDIF
385	#else
386	WRITE( message_string, * ) 'coupling requires PALM to be compiled with',&
387	' cpp-option "-D__parallel"'
388	CALL message( 'check_parameters', 'PA0141', 1, 2, 0, 6, 0 )
389	#endif
390	ENDIF
391
392	#if defined( __parallel )
393	!
394	!-- Exchange via intercommunicator
395	IF ( coupling_mode == 'atmosphere_to_ocean' .AND. myid == 0 ) THEN
396	CALL MPI_SEND( humidity, 1, MPI_LOGICAL, target_id, 19, comm_inter, &
397	ierr )
398	ELSEIF ( coupling_mode == 'ocean_to_atmosphere' .AND. myid == 0) THEN
399	CALL MPI_RECV( humidity_remote, 1, MPI_LOGICAL, target_id, 19, &
400	comm_inter, status, ierr )
401	ENDIF
402	CALL MPI_BCAST( humidity_remote, 1, MPI_LOGICAL, 0, comm2d, ierr)
403
404	#endif
405
406	!
407	!-- User settings for restart times requires that "restart" has been given as
408	!-- file activation string. Otherwise, binary output would not be saved by
409	!-- palmrun.
410	IF ( ( restart_time /= 9999999.9_wp .OR. dt_restart /= 9999999.9_wp ) &
411	.AND. .NOT. write_binary ) THEN
412	WRITE( message_string, * ) 'manual restart settings requires file ', &
413	'activation string "restart"'
414	CALL message( 'check_parameters', 'PA0001', 1, 2, 0, 6, 0 )
415	ENDIF
416
417
418	!
419	!-- Generate the file header which is used as a header for most of PALM's
420	!-- output files
421	CALL DATE_AND_TIME( date, time, run_zone )
422	run_date = date(1:4)//'-'//date(5:6)//'-'//date(7:8)
423	run_time = time(1:2)//':'//time(3:4)//':'//time(5:6)
424	IF ( coupling_mode == 'uncoupled' ) THEN
425	coupling_string = ''
426	ELSEIF ( coupling_mode == 'vnested_crse' ) THEN
427	coupling_string = ' nested (coarse)'
428	ELSEIF ( coupling_mode == 'vnested_fine' ) THEN
429	coupling_string = ' nested (fine)'
430	ELSEIF ( coupling_mode == 'atmosphere_to_ocean' ) THEN
431	coupling_string = ' coupled (atmosphere)'
432	ELSEIF ( coupling_mode == 'ocean_to_atmosphere' ) THEN
433	coupling_string = ' coupled (ocean)'
434	ENDIF
435	IF ( ensemble_member_nr /= 0 ) THEN
436	WRITE( ensemble_string, '(2X,A,I2.2)' ) 'en-no: ', ensemble_member_nr
437	ELSE
438	ensemble_string = ''
439	ENDIF
440	IF ( nested_run ) THEN
441	WRITE( nest_string, '(2X,A,I2.2)' ) 'nest-id: ', cpl_id
442	ELSE
443	nest_string = ''
444	ENDIF
445
446	WRITE ( run_description_header, &
447	'(A,2X,A,2X,A,A,A,I2.2,A,A,A,2X,A,A,2X,A,1X,A)' ) &
448	TRIM( version ), TRIM( revision ), 'run: ', &
449	TRIM( run_identifier ), '.', runnr, TRIM( coupling_string ), &
450	TRIM( nest_string ), TRIM( ensemble_string), 'host: ', TRIM( host ), &
451	run_date, run_time
452
453	!
454	!-- Check the general loop optimization method
455	SELECT CASE ( TRIM( loop_optimization ) )
456
457	CASE ( 'cache', 'vector' )
458	CONTINUE
459
460	CASE DEFAULT
461	message_string = 'illegal value given for loop_optimization: "' // &
462	TRIM( loop_optimization ) // '"'
463	CALL message( 'check_parameters', 'PA0013', 1, 2, 0, 6, 0 )
464
465	END SELECT
466
467	!
468	!-- Check topography setting (check for illegal parameter combinations)
469	IF ( topography /= 'flat' ) THEN
470	action = ' '
471	IF ( scalar_advec /= 'pw-scheme' .AND. scalar_advec /= 'ws-scheme' &
472	) THEN
473	WRITE( action, '(A,A)' ) 'scalar_advec = ', scalar_advec
474	ENDIF
475	IF ( momentum_advec /= 'pw-scheme' .AND. momentum_advec /= 'ws-scheme' )&
476	THEN
477	WRITE( action, '(A,A)' ) 'momentum_advec = ', momentum_advec
478	ENDIF
479	IF ( psolver == 'sor' ) THEN
480	WRITE( action, '(A,A)' ) 'psolver = ', psolver
481	ENDIF
482	IF ( sloping_surface ) THEN
483	WRITE( action, '(A)' ) 'sloping surface = .TRUE.'
484	ENDIF
485	IF ( galilei_transformation ) THEN
486	WRITE( action, '(A)' ) 'galilei_transformation = .TRUE.'
487	ENDIF
488	IF ( cloud_droplets ) THEN
489	WRITE( action, '(A)' ) 'cloud_droplets = .TRUE.'
490	ENDIF
491	IF ( .NOT. constant_flux_layer ) THEN
492	WRITE( action, '(A)' ) 'constant_flux_layer = .FALSE.'
493	ENDIF
494	IF ( action /= ' ' ) THEN
495	message_string = 'a non-flat topography does not allow ' // &
496	TRIM( action )
497	CALL message( 'check_parameters', 'PA0014', 1, 2, 0, 6, 0 )
498	ENDIF
499
500	ENDIF
501
502	!
503	!-- Check approximation
504	IF ( TRIM( approximation ) /= 'boussinesq' .AND. &
505	TRIM( approximation ) /= 'anelastic' ) THEN
506	message_string = 'unknown approximation: approximation = "' // &
507	TRIM( approximation ) // '"'
508	CALL message( 'check_parameters', 'PA0446', 1, 2, 0, 6, 0 )
509	ENDIF
510
511	!
512	!-- Check approximation requirements
513	IF ( TRIM( approximation ) == 'anelastic' .AND. &
514	TRIM( momentum_advec ) /= 'ws-scheme' ) THEN
515	message_string = 'Anelastic approximation requires ' // &
516	'momentum_advec = "ws-scheme"'
517	CALL message( 'check_parameters', 'PA0447', 1, 2, 0, 6, 0 )
518	ENDIF
519	IF ( TRIM( approximation ) == 'anelastic' .AND. &
520	TRIM( psolver ) == 'multigrid' ) THEN
521	message_string = 'Anelastic approximation currently only supports ' // &
522	'psolver = "poisfft", ' // &
523	'psolver = "sor" and ' // &
524	'psolver = "multigrid_noopt"'
525	CALL message( 'check_parameters', 'PA0448', 1, 2, 0, 6, 0 )
526	ENDIF
527	IF ( TRIM( approximation ) == 'anelastic' .AND. &
528	conserve_volume_flow ) THEN
529	message_string = 'Anelastic approximation is not allowed with ' // &
530	'conserve_volume_flow = .TRUE.'
531	CALL message( 'check_parameters', 'PA0449', 1, 2, 0, 6, 0 )
532	ENDIF
533
534	!
535	!-- Check flux input mode
536	IF ( TRIM( flux_input_mode ) /= 'dynamic' .AND. &
537	TRIM( flux_input_mode ) /= 'kinematic' .AND. &
538	TRIM( flux_input_mode ) /= 'approximation-specific' ) THEN
539	message_string = 'unknown flux input mode: flux_input_mode = "' // &
540	TRIM( flux_input_mode ) // '"'
541	CALL message( 'check_parameters', 'PA0450', 1, 2, 0, 6, 0 )
542	ENDIF
543	!
544	!-- Set flux input mode according to approximation if applicable
545	IF ( TRIM( flux_input_mode ) == 'approximation-specific' ) THEN
546	IF ( TRIM( approximation ) == 'anelastic' ) THEN
547	flux_input_mode = 'dynamic'
548	ELSEIF ( TRIM( approximation ) == 'boussinesq' ) THEN
549	flux_input_mode = 'kinematic'
550	ENDIF
551	ENDIF
552
553	!
554	!-- Check flux output mode
555	IF ( TRIM( flux_output_mode ) /= 'dynamic' .AND. &
556	TRIM( flux_output_mode ) /= 'kinematic' .AND. &
557	TRIM( flux_output_mode ) /= 'approximation-specific' ) THEN
558	message_string = 'unknown flux output mode: flux_output_mode = "' // &
559	TRIM( flux_output_mode ) // '"'
560	CALL message( 'check_parameters', 'PA0451', 1, 2, 0, 6, 0 )
561	ENDIF
562	!
563	!-- Set flux output mode according to approximation if applicable
564	IF ( TRIM( flux_output_mode ) == 'approximation-specific' ) THEN
565	IF ( TRIM( approximation ) == 'anelastic' ) THEN
566	flux_output_mode = 'dynamic'
567	ELSEIF ( TRIM( approximation ) == 'boussinesq' ) THEN
568	flux_output_mode = 'kinematic'
569	ENDIF
570	ENDIF
571
572
573	!
574	!-- When the land- or urban-surface model is used, the flux output must be
575	!-- dynamic.
576	IF ( land_surface .OR. urban_surface ) THEN
577	flux_output_mode = 'dynamic'
578	ENDIF
579
580	!
581	!-- Set the flux output units according to flux_output_mode
582	IF ( TRIM( flux_output_mode ) == 'kinematic' ) THEN
583	heatflux_output_unit = 'K m/s'
584	waterflux_output_unit = 'kg/kg m/s'
585	momentumflux_output_unit = 'm2/s2'
586	ELSEIF ( TRIM( flux_output_mode ) == 'dynamic' ) THEN
587	heatflux_output_unit = 'W/m2'
588	waterflux_output_unit = 'W/m2'
589	momentumflux_output_unit = 'N/m2'
590	ENDIF
591
592	!
593	!-- set time series output units for fluxes
594	dots_unit(14:16) = TRIM( heatflux_output_unit )
595	dots_unit(21) = TRIM( waterflux_output_unit )
596	dots_unit(19:20) = TRIM( momentumflux_output_unit )
597
598	!
599	!-- Add other module specific timeseries
600	CALL module_interface_check_data_output_ts( dots_max, dots_num, dots_label, dots_unit )
601
602	!
603	!-- Check if maximum number of allowed timeseries is exceeded
604	IF ( dots_num > dots_max ) THEN
605	WRITE( message_string, * ) 'number of time series quantities exceeds', &
606	' its maximum of dots_max = ', dots_max, &
607	'&Please increase dots_max in modules.f90.'
608	CALL message( 'init_3d_model', 'PA0194', 1, 2, 0, 6, 0 )
609	ENDIF
610
611	!
612	!-- Check whether there are any illegal values
613	!-- Pressure solver:
614	IF ( psolver /= 'poisfft' .AND. psolver /= 'sor' .AND. &
615	psolver /= 'multigrid' .AND. psolver /= 'multigrid_noopt' ) THEN
616	message_string = 'unknown solver for perturbation pressure: psolver' // &
617	' = "' // TRIM( psolver ) // '"'
618	CALL message( 'check_parameters', 'PA0016', 1, 2, 0, 6, 0 )
619	ENDIF
620
621	IF ( psolver(1:9) == 'multigrid' ) THEN
622	IF ( cycle_mg == 'w' ) THEN
623	gamma_mg = 2
624	ELSEIF ( cycle_mg == 'v' ) THEN
625	gamma_mg = 1
626	ELSE
627	message_string = 'unknown multigrid cycle: cycle_mg = "' // &
628	TRIM( cycle_mg ) // '"'
629	CALL message( 'check_parameters', 'PA0020', 1, 2, 0, 6, 0 )
630	ENDIF
631	ENDIF
632
633	IF ( fft_method /= 'singleton-algorithm' .AND. &
634	fft_method /= 'temperton-algorithm' .AND. &
635	fft_method /= 'fftw' .AND. &
636	fft_method /= 'system-specific' ) THEN
637	message_string = 'unknown fft-algorithm: fft_method = "' // &
638	TRIM( fft_method ) // '"'
639	CALL message( 'check_parameters', 'PA0021', 1, 2, 0, 6, 0 )
640	ENDIF
641
642	IF( momentum_advec == 'ws-scheme' .AND. &
643	.NOT. call_psolver_at_all_substeps ) THEN
644	message_string = 'psolver must be called at each RK3 substep when "'// &
645	TRIM(momentum_advec) // ' "is used for momentum_advec'
646	CALL message( 'check_parameters', 'PA0344', 1, 2, 0, 6, 0 )
647	END IF
648	!
649	!-- Advection schemes:
650	IF ( momentum_advec /= 'pw-scheme' .AND. &
651	momentum_advec /= 'ws-scheme' .AND. &
652	momentum_advec /= 'up-scheme' ) &
653	THEN
654	message_string = 'unknown advection scheme: momentum_advec = "' // &
655	TRIM( momentum_advec ) // '"'
656	CALL message( 'check_parameters', 'PA0022', 1, 2, 0, 6, 0 )
657	ENDIF
658	IF ( ( momentum_advec == 'ws-scheme' .OR. scalar_advec == 'ws-scheme' ) &
659	.AND. ( timestep_scheme == 'euler' .OR. &
660	timestep_scheme == 'runge-kutta-2' ) ) &
661	THEN
662	message_string = 'momentum_advec or scalar_advec = "' &
663	// TRIM( momentum_advec ) // '" is not allowed with ' // &
664	'timestep_scheme = "' // TRIM( timestep_scheme ) // '"'
665	CALL message( 'check_parameters', 'PA0023', 1, 2, 0, 6, 0 )
666	ENDIF
667	IF ( scalar_advec /= 'pw-scheme' .AND. scalar_advec /= 'ws-scheme' .AND. &
668	scalar_advec /= 'bc-scheme' .AND. scalar_advec /= 'up-scheme' ) &
669	THEN
670	message_string = 'unknown advection scheme: scalar_advec = "' // &
671	TRIM( scalar_advec ) // '"'
672	CALL message( 'check_parameters', 'PA0024', 1, 2, 0, 6, 0 )
673	ENDIF
674	IF ( scalar_advec == 'bc-scheme' .AND. loop_optimization == 'cache' ) &
675	THEN
676	message_string = 'advection_scheme scalar_advec = "' &
677	// TRIM( scalar_advec ) // '" not implemented for ' // &
678	'loop_optimization = "' // TRIM( loop_optimization ) // '"'
679	CALL message( 'check_parameters', 'PA0026', 1, 2, 0, 6, 0 )
680	ENDIF
681
682	IF ( use_sgs_for_particles .AND. .NOT. cloud_droplets .AND. &
683	.NOT. use_upstream_for_tke .AND. &
684	scalar_advec /= 'ws-scheme' &
685	) THEN
686	use_upstream_for_tke = .TRUE.
687	message_string = 'use_upstream_for_tke is set to .TRUE. because ' // &
688	'use_sgs_for_particles = .TRUE. ' // &
689	'and scalar_advec /= ws-scheme'
690	CALL message( 'check_parameters', 'PA0025', 0, 1, 0, 6, 0 )
691	ENDIF
692
693	!
694	!-- Set LOGICAL switches to enhance performance
695	IF ( momentum_advec == 'ws-scheme' ) ws_scheme_mom = .TRUE.
696	IF ( scalar_advec == 'ws-scheme' ) ws_scheme_sca = .TRUE.
697
698
699	!
700	!-- Timestep schemes:
701	SELECT CASE ( TRIM( timestep_scheme ) )
702
703	CASE ( 'euler' )
704	intermediate_timestep_count_max = 1
705
706	CASE ( 'runge-kutta-2' )
707	intermediate_timestep_count_max = 2
708
709	CASE ( 'runge-kutta-3' )
710	intermediate_timestep_count_max = 3
711
712	CASE DEFAULT
713	message_string = 'unknown timestep scheme: timestep_scheme = "' // &
714	TRIM( timestep_scheme ) // '"'
715	CALL message( 'check_parameters', 'PA0027', 1, 2, 0, 6, 0 )
716
717	END SELECT
718
719	IF ( (momentum_advec /= 'pw-scheme' .AND. momentum_advec /= 'ws-scheme') &
720	.AND. timestep_scheme(1:5) == 'runge' ) THEN
721	message_string = 'momentum advection scheme "' // &
722	TRIM( momentum_advec ) // '" & does not work with ' // &
723	'timestep_scheme "' // TRIM( timestep_scheme ) // '"'
724	CALL message( 'check_parameters', 'PA0029', 1, 2, 0, 6, 0 )
725	ENDIF
726	!
727	!-- Check for proper settings for microphysics
728	IF ( bulk_cloud_model .AND. cloud_droplets ) THEN
729	message_string = 'bulk_cloud_model = .TRUE. is not allowed with ' // &
730	'cloud_droplets = .TRUE.'
731	CALL message( 'check_parameters', 'PA0442', 1, 2, 0, 6, 0 )
732	ENDIF
733
734	!
735	!-- Initializing actions must have been set by the user
736	IF ( TRIM( initializing_actions ) == '' ) THEN
737	message_string = 'no value specified for initializing_actions'
738	CALL message( 'check_parameters', 'PA0149', 1, 2, 0, 6, 0 )
739	ENDIF
740
741	IF ( TRIM( initializing_actions ) /= 'read_restart_data' .AND. &
742	TRIM( initializing_actions ) /= 'cyclic_fill' ) THEN
743	!
744	!-- No restart run: several initialising actions are possible
745	action = initializing_actions
746	DO WHILE ( TRIM( action ) /= '' )
747	position = INDEX( action, ' ' )
748	SELECT CASE ( action(1:position-1) )
749
750	CASE ( 'set_constant_profiles', 'set_1d-model_profiles', &
751	'by_user', 'initialize_vortex', 'initialize_ptanom', &
752	'initialize_bubble', 'inifor' )
753	action = action(position+1:)
754
755	CASE DEFAULT
756	message_string = 'initializing_action = "' // &
757	TRIM( action ) // '" unknown or not allowed'
758	CALL message( 'check_parameters', 'PA0030', 1, 2, 0, 6, 0 )
759
760	END SELECT
761	ENDDO
762	ENDIF
763
764	IF ( TRIM( initializing_actions ) == 'initialize_vortex' .AND. &
765	conserve_volume_flow ) THEN
766	message_string = 'initializing_actions = "initialize_vortex"' // &
767	' is not allowed with conserve_volume_flow = .T.'
768	CALL message( 'check_parameters', 'PA0343', 1, 2, 0, 6, 0 )
769	ENDIF
770
771
772	IF ( INDEX( initializing_actions, 'set_constant_profiles' ) /= 0 .AND. &
773	INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 ) THEN
774	message_string = 'initializing_actions = "set_constant_profiles"' // &
775	' and "set_1d-model_profiles" are not allowed ' // &
776	'simultaneously'
777	CALL message( 'check_parameters', 'PA0031', 1, 2, 0, 6, 0 )
778	ENDIF
779
780	IF ( INDEX( initializing_actions, 'set_constant_profiles' ) /= 0 .AND. &
781	INDEX( initializing_actions, 'by_user' ) /= 0 ) THEN
782	message_string = 'initializing_actions = "set_constant_profiles"' // &
783	' and "by_user" are not allowed simultaneously'
784	CALL message( 'check_parameters', 'PA0032', 1, 2, 0, 6, 0 )
785	ENDIF
786
787	IF ( INDEX( initializing_actions, 'by_user' ) /= 0 .AND. &
788	INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 ) THEN
789	message_string = 'initializing_actions = "by_user" and ' // &
790	'"set_1d-model_profiles" are not allowed simultaneously'
791	CALL message( 'check_parameters', 'PA0033', 1, 2, 0, 6, 0 )
792	ENDIF
793	!
794	!-- In case of spinup and nested run, spinup end time must be identical
795	!-- in order to have synchronously running simulations.
796	IF ( nested_run ) THEN
797	#if defined( __parallel )
798	CALL MPI_ALLREDUCE( spinup_time, spinup_time_max, 1, MPI_REAL, &
799	MPI_MAX, MPI_COMM_WORLD, ierr )
800	CALL MPI_ALLREDUCE( dt_spinup, dt_spinup_max, 1, MPI_REAL, &
801	MPI_MAX, MPI_COMM_WORLD, ierr )
802
803	IF ( spinup_time /= spinup_time_max .OR. dt_spinup /= dt_spinup_max ) &
804	THEN
805	message_string = 'In case of nesting, spinup_time and ' // &
806	'dt_spinup must be identical in all parent ' // &
807	'and child domains.'
808	CALL message( 'check_parameters', 'PA0489', 3, 2, 0, 6, 0 )
809	ENDIF
810	#endif
811	ENDIF
812
813	IF ( bulk_cloud_model .AND. .NOT. humidity ) THEN
814	WRITE( message_string, * ) 'bulk_cloud_model = ', bulk_cloud_model, &
815	' is not allowed with humidity = ', humidity
816	CALL message( 'check_parameters', 'PA0034', 1, 2, 0, 6, 0 )
817	ENDIF
818
819	IF ( humidity .AND. sloping_surface ) THEN
820	message_string = 'humidity = .TRUE. and sloping_surface = .TRUE. ' // &
821	'are not allowed simultaneously'
822	CALL message( 'check_parameters', 'PA0036', 1, 2, 0, 6, 0 )
823	ENDIF
824
825	!-- Check the module settings
826	CALL module_interface_check_parameters
827
828	!
829	!-- In case of no restart run, check initialising parameters and calculate
830	!-- further quantities
831	IF ( TRIM( initializing_actions ) /= 'read_restart_data' ) THEN
832
833	!
834	!-- Initial profiles for 1D and 3D model, respectively (u,v further below)
835	pt_init = pt_surface
836	IF ( humidity ) q_init = q_surface
837	IF ( passive_scalar ) s_init = s_surface
838
839	!--
840	!-- If required, compute initial profile of the geostrophic wind
841	!-- (component ug)
842	i = 1
843	gradient = 0.0_wp
844
845	IF ( .NOT. ocean_mode ) THEN
846
847	ug_vertical_gradient_level_ind(1) = 0
848	ug(0) = ug_surface
849	DO k = 1, nzt+1
850	IF ( i < 11 ) THEN
851	IF ( ug_vertical_gradient_level(i) < zu(k) .AND. &
852	ug_vertical_gradient_level(i) >= 0.0_wp ) THEN
853	gradient = ug_vertical_gradient(i) / 100.0_wp
854	ug_vertical_gradient_level_ind(i) = k - 1
855	i = i + 1
856	ENDIF
857	ENDIF
858	IF ( gradient /= 0.0_wp ) THEN
859	IF ( k /= 1 ) THEN
860	ug(k) = ug(k-1) + dzu(k) * gradient
861	ELSE
862	ug(k) = ug_surface + dzu(k) * gradient
863	ENDIF
864	ELSE
865	ug(k) = ug(k-1)
866	ENDIF
867	ENDDO
868
869	ELSE
870
871	ug_vertical_gradient_level_ind(1) = nzt+1
872	ug(nzt+1) = ug_surface
873	DO k = nzt, nzb, -1
874	IF ( i < 11 ) THEN
875	IF ( ug_vertical_gradient_level(i) > zu(k) .AND. &
876	ug_vertical_gradient_level(i) <= 0.0_wp ) THEN
877	gradient = ug_vertical_gradient(i) / 100.0_wp
878	ug_vertical_gradient_level_ind(i) = k + 1
879	i = i + 1
880	ENDIF
881	ENDIF
882	IF ( gradient /= 0.0_wp ) THEN
883	IF ( k /= nzt ) THEN
884	ug(k) = ug(k+1) - dzu(k+1) * gradient
885	ELSE
886	ug(k) = ug_surface - 0.5_wp * dzu(k+1) * gradient
887	ug(k+1) = ug_surface + 0.5_wp * dzu(k+1) * gradient
888	ENDIF
889	ELSE
890	ug(k) = ug(k+1)
891	ENDIF
892	ENDDO
893
894	ENDIF
895
896	!
897	!-- In case of no given gradients for ug, choose a zero gradient
898	IF ( ug_vertical_gradient_level(1) == -9999999.9_wp ) THEN
899	ug_vertical_gradient_level(1) = 0.0_wp
900	ENDIF
901
902	!
903	!--
904	!-- If required, compute initial profile of the geostrophic wind
905	!-- (component vg)
906	i = 1
907	gradient = 0.0_wp
908
909	IF ( .NOT. ocean_mode ) THEN
910
911	vg_vertical_gradient_level_ind(1) = 0
912	vg(0) = vg_surface
913	DO k = 1, nzt+1
914	IF ( i < 11 ) THEN
915	IF ( vg_vertical_gradient_level(i) < zu(k) .AND. &
916	vg_vertical_gradient_level(i) >= 0.0_wp ) THEN
917	gradient = vg_vertical_gradient(i) / 100.0_wp
918	vg_vertical_gradient_level_ind(i) = k - 1
919	i = i + 1
920	ENDIF
921	ENDIF
922	IF ( gradient /= 0.0_wp ) THEN
923	IF ( k /= 1 ) THEN
924	vg(k) = vg(k-1) + dzu(k) * gradient
925	ELSE
926	vg(k) = vg_surface + dzu(k) * gradient
927	ENDIF
928	ELSE
929	vg(k) = vg(k-1)
930	ENDIF
931	ENDDO
932
933	ELSE
934
935	vg_vertical_gradient_level_ind(1) = nzt+1
936	vg(nzt+1) = vg_surface
937	DO k = nzt, nzb, -1
938	IF ( i < 11 ) THEN
939	IF ( vg_vertical_gradient_level(i) > zu(k) .AND. &
940	vg_vertical_gradient_level(i) <= 0.0_wp ) THEN
941	gradient = vg_vertical_gradient(i) / 100.0_wp
942	vg_vertical_gradient_level_ind(i) = k + 1
943	i = i + 1
944	ENDIF
945	ENDIF
946	IF ( gradient /= 0.0_wp ) THEN
947	IF ( k /= nzt ) THEN
948	vg(k) = vg(k+1) - dzu(k+1) * gradient
949	ELSE
950	vg(k) = vg_surface - 0.5_wp * dzu(k+1) * gradient
951	vg(k+1) = vg_surface + 0.5_wp * dzu(k+1) * gradient
952	ENDIF
953	ELSE
954	vg(k) = vg(k+1)
955	ENDIF
956	ENDDO
957
958	ENDIF
959
960	!
961	!-- In case of no given gradients for vg, choose a zero gradient
962	IF ( vg_vertical_gradient_level(1) == -9999999.9_wp ) THEN
963	vg_vertical_gradient_level(1) = 0.0_wp
964	ENDIF
965
966	!
967	!-- Let the initial wind profiles be the calculated ug/vg profiles or
968	!-- interpolate them from wind profile data (if given)
969	IF ( u_profile(1) == 9999999.9_wp .AND. v_profile(1) == 9999999.9_wp ) THEN
970
971	u_init = ug
972	v_init = vg
973
974	ELSEIF ( u_profile(1) == 0.0_wp .AND. v_profile(1) == 0.0_wp ) THEN
975
976	IF ( uv_heights(1) /= 0.0_wp ) THEN
977	message_string = 'uv_heights(1) must be 0.0'
978	CALL message( 'check_parameters', 'PA0345', 1, 2, 0, 6, 0 )
979	ENDIF
980
981	IF ( omega /= 0.0_wp ) THEN
982	message_string = 'Coriolis force must be switched off (by setting omega=0.0)' // &
983	' when prescribing the forcing by u_profile and v_profile'
984	CALL message( 'check_parameters', 'PA0347', 1, 2, 0, 6, 0 )
985	ENDIF
986
987	use_prescribed_profile_data = .TRUE.
988
989	kk = 1
990	u_init(0) = 0.0_wp
991	v_init(0) = 0.0_wp
992
993	DO k = 1, nz+1
994
995	IF ( kk < 200 ) THEN
996	DO WHILE ( uv_heights(kk+1) <= zu(k) )
997	kk = kk + 1
998	IF ( kk == 200 ) EXIT
999	ENDDO
1000	ENDIF
1001
1002	IF ( kk < 200 .AND. uv_heights(kk+1) /= 9999999.9_wp ) THEN
1003	u_init(k) = u_profile(kk) + ( zu(k) - uv_heights(kk) ) / &
1004	( uv_heights(kk+1) - uv_heights(kk) ) * &
1005	( u_profile(kk+1) - u_profile(kk) )
1006	v_init(k) = v_profile(kk) + ( zu(k) - uv_heights(kk) ) / &
1007	( uv_heights(kk+1) - uv_heights(kk) ) * &
1008	( v_profile(kk+1) - v_profile(kk) )
1009	ELSE
1010	u_init(k) = u_profile(kk)
1011	v_init(k) = v_profile(kk)
1012	ENDIF
1013
1014	ENDDO
1015
1016	ELSE
1017
1018	message_string = 'u_profile(1) and v_profile(1) must be 0.0'
1019	CALL message( 'check_parameters', 'PA0346', 1, 2, 0, 6, 0 )
1020
1021	ENDIF
1022
1023	!
1024	!-- Compute initial temperature profile using the given temperature gradients
1025	IF ( .NOT. neutral ) THEN
1026	CALL init_vertical_profiles( pt_vertical_gradient_level_ind, &
1027	pt_vertical_gradient_level, &
1028	pt_vertical_gradient, pt_init, &
1029	pt_surface, bc_pt_t_val )
1030	ENDIF
1031	!
1032	!-- Compute initial humidity profile using the given humidity gradients
1033	IF ( humidity ) THEN
1034	CALL init_vertical_profiles( q_vertical_gradient_level_ind, &
1035	q_vertical_gradient_level, &
1036	q_vertical_gradient, q_init, &
1037	q_surface, bc_q_t_val )
1038	ENDIF
1039	!
1040	!-- Compute initial scalar profile using the given scalar gradients
1041	IF ( passive_scalar ) THEN
1042	CALL init_vertical_profiles( s_vertical_gradient_level_ind, &
1043	s_vertical_gradient_level, &
1044	s_vertical_gradient, s_init, &
1045	s_surface, bc_s_t_val )
1046	ENDIF
1047	!
1048	!-- TODO
1049	!-- Compute initial chemistry profile using the given chemical species gradients
1050	!-- Russo: Is done in chem_init --> kanani: Revise
1051
1052	ENDIF
1053
1054	!
1055	!-- Check if the control parameter use_subsidence_tendencies is used correctly
1056	IF ( use_subsidence_tendencies .AND. .NOT. large_scale_subsidence ) THEN
1057	message_string = 'The usage of use_subsidence_tendencies ' // &
1058	'requires large_scale_subsidence = .T..'
1059	CALL message( 'check_parameters', 'PA0396', 1, 2, 0, 6, 0 )
1060	ELSEIF ( use_subsidence_tendencies .AND. .NOT. large_scale_forcing ) THEN
1061	message_string = 'The usage of use_subsidence_tendencies ' // &
1062	'requires large_scale_forcing = .T..'
1063	CALL message( 'check_parameters', 'PA0397', 1, 2, 0, 6, 0 )
1064	ENDIF
1065
1066	!
1067	!-- Initialize large scale subsidence if required
1068	If ( large_scale_subsidence ) THEN
1069	IF ( subs_vertical_gradient_level(1) /= -9999999.9_wp .AND. &
1070	.NOT. large_scale_forcing ) THEN
1071	CALL init_w_subsidence
1072	ENDIF
1073	!
1074	!-- In case large_scale_forcing is used, profiles for subsidence velocity
1075	!-- are read in from file LSF_DATA
1076
1077	IF ( subs_vertical_gradient_level(1) == -9999999.9_wp .AND. &
1078	.NOT. large_scale_forcing ) THEN
1079	message_string = 'There is no default large scale vertical ' // &
1080	'velocity profile set. Specify the subsidence ' // &
1081	'velocity profile via subs_vertical_gradient ' // &
1082	'and subs_vertical_gradient_level.'
1083	CALL message( 'check_parameters', 'PA0380', 1, 2, 0, 6, 0 )
1084	ENDIF
1085	ELSE
1086	IF ( subs_vertical_gradient_level(1) /= -9999999.9_wp ) THEN
1087	message_string = 'Enable usage of large scale subsidence by ' // &
1088	'setting large_scale_subsidence = .T..'
1089	CALL message( 'check_parameters', 'PA0381', 1, 2, 0, 6, 0 )
1090	ENDIF
1091	ENDIF
1092
1093	!
1094	!-- Overwrite parameters from namelist if necessary and compute Coriolis parameter.
1095	!-- @todo - move initialization of f and fs to coriolis_mod.
1096	IF ( input_pids_static ) THEN
1097	latitude = init_model%latitude
1098	longitude = init_model%longitude
1099	rotation_angle = init_model%rotation_angle
1100	ENDIF
1101
1102	f = 2.0_wp * omega * SIN( latitude / 180.0_wp * pi )
1103	fs = 2.0_wp * omega * COS( latitude / 180.0_wp * pi )
1104
1105	!
1106	!-- Check and set buoyancy related parameters and switches
1107	IF ( reference_state == 'horizontal_average' ) THEN
1108	CONTINUE
1109	ELSEIF ( reference_state == 'initial_profile' ) THEN
1110	use_initial_profile_as_reference = .TRUE.
1111	ELSEIF ( reference_state == 'single_value' ) THEN
1112	use_single_reference_value = .TRUE.
1113	IF ( pt_reference == 9999999.9_wp ) pt_reference = pt_surface
1114	vpt_reference = pt_reference * ( 1.0_wp + 0.61_wp * q_surface )
1115	ELSE
1116	message_string = 'illegal value for reference_state: "' // &
1117	TRIM( reference_state ) // '"'
1118	CALL message( 'check_parameters', 'PA0056', 1, 2, 0, 6, 0 )
1119	ENDIF
1120
1121	!
1122	!-- In case of a given slope, compute the relevant quantities
1123	IF ( alpha_surface /= 0.0_wp ) THEN
1124	IF ( ABS( alpha_surface ) > 90.0_wp ) THEN
1125	WRITE( message_string, * ) 'ABS( alpha_surface = ', alpha_surface, &
1126	' ) must be < 90.0'
1127	CALL message( 'check_parameters', 'PA0043', 1, 2, 0, 6, 0 )
1128	ENDIF
1129	sloping_surface = .TRUE.
1130	cos_alpha_surface = COS( alpha_surface / 180.0_wp * pi )
1131	sin_alpha_surface = SIN( alpha_surface / 180.0_wp * pi )
1132	ENDIF
1133
1134	!
1135	!-- Check time step and cfl_factor
1136	IF ( dt /= -1.0_wp ) THEN
1137	IF ( dt <= 0.0_wp ) THEN
1138	WRITE( message_string, * ) 'dt = ', dt , ' <= 0.0'
1139	CALL message( 'check_parameters', 'PA0044', 1, 2, 0, 6, 0 )
1140	ENDIF
1141	dt_3d = dt
1142	dt_fixed = .TRUE.
1143	ENDIF
1144
1145	IF ( cfl_factor <= 0.0_wp .OR. cfl_factor > 1.0_wp ) THEN
1146	IF ( cfl_factor == -1.0_wp ) THEN
1147	IF ( timestep_scheme == 'runge-kutta-2' ) THEN
1148	cfl_factor = 0.8_wp
1149	ELSEIF ( timestep_scheme == 'runge-kutta-3' ) THEN
1150	cfl_factor = 0.9_wp
1151	ELSE
1152	cfl_factor = 0.9_wp
1153	ENDIF
1154	ELSE
1155	WRITE( message_string, * ) 'cfl_factor = ', cfl_factor, &
1156	' out of range &0.0 < cfl_factor <= 1.0 is required'
1157	CALL message( 'check_parameters', 'PA0045', 1, 2, 0, 6, 0 )
1158	ENDIF
1159	ENDIF
1160
1161	!
1162	!-- Store simulated time at begin
1163	simulated_time_at_begin = simulated_time
1164
1165	!
1166	!-- Store reference time for coupled runs and change the coupling flag,
1167	!-- if ...
1168	IF ( simulated_time == 0.0_wp ) THEN
1169	IF ( coupling_start_time == 0.0_wp ) THEN
1170	time_since_reference_point = 0.0_wp
1171	ELSEIF ( time_since_reference_point < 0.0_wp ) THEN
1172	run_coupled = .FALSE.
1173	ENDIF
1174	ENDIF
1175
1176	!
1177	!-- Set wind speed in the Galilei-transformed system
1178	IF ( galilei_transformation ) THEN
1179	IF ( use_ug_for_galilei_tr .AND. &
1180	ug_vertical_gradient_level(1) == 0.0_wp .AND. &
1181	ug_vertical_gradient(1) == 0.0_wp .AND. &
1182	vg_vertical_gradient_level(1) == 0.0_wp .AND. &
1183	vg_vertical_gradient(1) == 0.0_wp ) THEN
1184	u_gtrans = ug_surface * 0.6_wp
1185	v_gtrans = vg_surface * 0.6_wp
1186	ELSEIF ( use_ug_for_galilei_tr .AND. &
1187	( ug_vertical_gradient_level(1) /= 0.0_wp .OR. &
1188	ug_vertical_gradient(1) /= 0.0_wp ) ) THEN
1189	message_string = 'baroclinity (ug) not allowed simultaneously' // &
1190	' with galilei transformation'
1191	CALL message( 'check_parameters', 'PA0046', 1, 2, 0, 6, 0 )
1192	ELSEIF ( use_ug_for_galilei_tr .AND. &
1193	( vg_vertical_gradient_level(1) /= 0.0_wp .OR. &
1194	vg_vertical_gradient(1) /= 0.0_wp ) ) THEN
1195	message_string = 'baroclinity (vg) not allowed simultaneously' // &
1196	' with galilei transformation'
1197	CALL message( 'check_parameters', 'PA0047', 1, 2, 0, 6, 0 )
1198	ELSE
1199	message_string = 'variable translation speed used for Galilei-' // &
1200	'transformation, which may cause & instabilities in stably ' // &
1201	'stratified regions'
1202	CALL message( 'check_parameters', 'PA0048', 0, 1, 0, 6, 0 )
1203	ENDIF
1204	ENDIF
1205
1206	!
1207	!-- In case of using a prandtl-layer, calculated (or prescribed) surface
1208	!-- fluxes have to be used in the diffusion-terms
1209	IF ( constant_flux_layer ) use_surface_fluxes = .TRUE.
1210
1211	!
1212	!-- Check boundary conditions and set internal variables:
1213	!-- Attention: the lateral boundary conditions have been already checked in
1214	!-- parin
1215	!
1216	!-- Non-cyclic lateral boundaries require the multigrid method and Piascek-
1217	!-- Willimas or Wicker - Skamarock advection scheme. Several schemes
1218	!-- and tools do not work with non-cyclic boundary conditions.
1219	IF ( bc_lr /= 'cyclic' .OR. bc_ns /= 'cyclic' ) THEN
1220	IF ( psolver(1:9) /= 'multigrid' ) THEN
1221	message_string = 'non-cyclic lateral boundaries do not allow ' // &
1222	'psolver = "' // TRIM( psolver ) // '"'
1223	CALL message( 'check_parameters', 'PA0051', 1, 2, 0, 6, 0 )
1224	ENDIF
1225	IF ( momentum_advec /= 'pw-scheme' .AND. &
1226	momentum_advec /= 'ws-scheme' ) THEN
1227
1228	message_string = 'non-cyclic lateral boundaries do not allow ' // &
1229	'momentum_advec = "' // TRIM( momentum_advec ) // '"'
1230	CALL message( 'check_parameters', 'PA0052', 1, 2, 0, 6, 0 )
1231	ENDIF
1232	IF ( scalar_advec /= 'pw-scheme' .AND. &
1233	scalar_advec /= 'ws-scheme' ) THEN
1234	message_string = 'non-cyclic lateral boundaries do not allow ' // &
1235	'scalar_advec = "' // TRIM( scalar_advec ) // '"'
1236	CALL message( 'check_parameters', 'PA0053', 1, 2, 0, 6, 0 )
1237	ENDIF
1238	IF ( galilei_transformation ) THEN
1239	message_string = 'non-cyclic lateral boundaries do not allow ' // &
1240	'galilei_transformation = .T.'
1241	CALL message( 'check_parameters', 'PA0054', 1, 2, 0, 6, 0 )
1242	ENDIF
1243	ENDIF
1244
1245	!
1246	!-- Bottom boundary condition for the turbulent Kinetic energy
1247	IF ( bc_e_b == 'neumann' ) THEN
1248	ibc_e_b = 1
1249	ELSEIF ( bc_e_b == '(u)*2+neumann' ) THEN
1250	ibc_e_b = 2
1251	IF ( .NOT. constant_flux_layer ) THEN
1252	bc_e_b = 'neumann'
1253	ibc_e_b = 1
1254	message_string = 'boundary condition bc_e_b changed to "' // &
1255	TRIM( bc_e_b ) // '"'
1256	CALL message( 'check_parameters', 'PA0057', 0, 1, 0, 6, 0 )
1257	ENDIF
1258	ELSE
1259	message_string = 'unknown boundary condition: bc_e_b = "' // &
1260	TRIM( bc_e_b ) // '"'
1261	CALL message( 'check_parameters', 'PA0058', 1, 2, 0, 6, 0 )
1262	ENDIF
1263
1264	!
1265	!-- Boundary conditions for perturbation pressure
1266	IF ( bc_p_b == 'dirichlet' ) THEN
1267	ibc_p_b = 0
1268	ELSEIF ( bc_p_b == 'neumann' ) THEN
1269	ibc_p_b = 1
1270	ELSE
1271	message_string = 'unknown boundary condition: bc_p_b = "' // &
1272	TRIM( bc_p_b ) // '"'
1273	CALL message( 'check_parameters', 'PA0059', 1, 2, 0, 6, 0 )
1274	ENDIF
1275
1276	IF ( bc_p_t == 'dirichlet' ) THEN
1277	ibc_p_t = 0
1278	!-- TO_DO: later set bc_p_t to neumann before, in case of nested domain
1279	ELSEIF ( bc_p_t == 'neumann' .OR. bc_p_t == 'nested' ) THEN
1280	ibc_p_t = 1
1281	ELSE
1282	message_string = 'unknown boundary condition: bc_p_t = "' // &
1283	TRIM( bc_p_t ) // '"'
1284	CALL message( 'check_parameters', 'PA0061', 1, 2, 0, 6, 0 )
1285	ENDIF
1286
1287	!
1288	!-- Boundary conditions for potential temperature
1289	IF ( coupling_mode == 'atmosphere_to_ocean' ) THEN
1290	ibc_pt_b = 2
1291	ELSE
1292	IF ( bc_pt_b == 'dirichlet' ) THEN
1293	ibc_pt_b = 0
1294	ELSEIF ( bc_pt_b == 'neumann' ) THEN
1295	ibc_pt_b = 1
1296	ELSE
1297	message_string = 'unknown boundary condition: bc_pt_b = "' // &
1298	TRIM( bc_pt_b ) // '"'
1299	CALL message( 'check_parameters', 'PA0062', 1, 2, 0, 6, 0 )
1300	ENDIF
1301	ENDIF
1302
1303	IF ( bc_pt_t == 'dirichlet' ) THEN
1304	ibc_pt_t = 0
1305	ELSEIF ( bc_pt_t == 'neumann' ) THEN
1306	ibc_pt_t = 1
1307	ELSEIF ( bc_pt_t == 'initial_gradient' ) THEN
1308	ibc_pt_t = 2
1309	ELSEIF ( bc_pt_t == 'nested' .OR. bc_pt_t == 'nesting_offline' ) THEN
1310	ibc_pt_t = 3
1311	ELSE
1312	message_string = 'unknown boundary condition: bc_pt_t = "' // &
1313	TRIM( bc_pt_t ) // '"'
1314	CALL message( 'check_parameters', 'PA0063', 1, 2, 0, 6, 0 )
1315	ENDIF
1316
1317	IF ( ANY( wall_heatflux /= 0.0_wp ) .AND. &
1318	surface_heatflux == 9999999.9_wp ) THEN
1319	message_string = 'wall_heatflux additionally requires ' // &
1320	'setting of surface_heatflux'
1321	CALL message( 'check_parameters', 'PA0443', 1, 2, 0, 6, 0 )
1322	ENDIF
1323
1324	!
1325	! This IF clause needs revision, got too complex!!
1326	IF ( surface_heatflux == 9999999.9_wp ) THEN
1327	constant_heatflux = .FALSE.
1328	IF ( large_scale_forcing .OR. land_surface .OR. urban_surface ) THEN
1329	IF ( ibc_pt_b == 0 ) THEN
1330	constant_heatflux = .FALSE.
1331	ELSEIF ( ibc_pt_b == 1 ) THEN
1332	constant_heatflux = .TRUE.
1333	surface_heatflux = 0.0_wp
1334	ENDIF
1335	ENDIF
1336	ELSE
1337	constant_heatflux = .TRUE.
1338	ENDIF
1339
1340	IF ( top_heatflux == 9999999.9_wp ) constant_top_heatflux = .FALSE.
1341
1342	IF ( neutral ) THEN
1343
1344	IF ( surface_heatflux /= 0.0_wp .AND. &
1345	surface_heatflux /= 9999999.9_wp ) THEN
1346	message_string = 'heatflux must not be set for pure neutral flow'
1347	CALL message( 'check_parameters', 'PA0351', 1, 2, 0, 6, 0 )
1348	ENDIF
1349
1350	IF ( top_heatflux /= 0.0_wp .AND. top_heatflux /= 9999999.9_wp ) &
1351	THEN
1352	message_string = 'heatflux must not be set for pure neutral flow'
1353	CALL message( 'check_parameters', 'PA0351', 1, 2, 0, 6, 0 )
1354	ENDIF
1355
1356	ENDIF
1357
1358	IF ( top_momentumflux_u /= 9999999.9_wp .AND. &
1359	top_momentumflux_v /= 9999999.9_wp ) THEN
1360	constant_top_momentumflux = .TRUE.
1361	ELSEIF ( .NOT. ( top_momentumflux_u == 9999999.9_wp .AND. &
1362	top_momentumflux_v == 9999999.9_wp ) ) THEN
1363	message_string = 'both, top_momentumflux_u AND top_momentumflux_v ' // &
1364	'must be set'
1365	CALL message( 'check_parameters', 'PA0064', 1, 2, 0, 6, 0 )
1366	ENDIF
1367
1368	!
1369	!-- A given surface temperature implies Dirichlet boundary condition for
1370	!-- temperature. In this case specification of a constant heat flux is
1371	!-- forbidden.
1372	IF ( ibc_pt_b == 0 .AND. constant_heatflux .AND. &
1373	surface_heatflux /= 0.0_wp ) THEN
1374	message_string = 'boundary_condition: bc_pt_b = "' // TRIM( bc_pt_b ) //&
1375	'& is not allowed with constant_heatflux = .TRUE.'
1376	CALL message( 'check_parameters', 'PA0065', 1, 2, 0, 6, 0 )
1377	ENDIF
1378	IF ( constant_heatflux .AND. pt_surface_initial_change /= 0.0_wp ) THEN
1379	WRITE ( message_string, * ) 'constant_heatflux = .TRUE. is not allo', &
1380	'wed with pt_surface_initial_change (/=0) = ', &
1381	pt_surface_initial_change
1382	CALL message( 'check_parameters', 'PA0066', 1, 2, 0, 6, 0 )
1383	ENDIF
1384
1385	!
1386	!-- A given temperature at the top implies Dirichlet boundary condition for
1387	!-- temperature. In this case specification of a constant heat flux is
1388	!-- forbidden.
1389	IF ( ibc_pt_t == 0 .AND. constant_top_heatflux .AND. &
1390	top_heatflux /= 0.0_wp ) THEN
1391	message_string = 'boundary_condition: bc_pt_t = "' // TRIM( bc_pt_t ) //&
1392	'" is not allowed with constant_top_heatflux = .TRUE.'
1393	CALL message( 'check_parameters', 'PA0067', 1, 2, 0, 6, 0 )
1394	ENDIF
1395
1396	!
1397	!-- Set boundary conditions for total water content
1398	IF ( humidity ) THEN
1399
1400	IF ( ANY( wall_humidityflux /= 0.0_wp ) .AND. &
1401	surface_waterflux == 9999999.9_wp ) THEN
1402	message_string = 'wall_humidityflux additionally requires ' // &
1403	'setting of surface_waterflux'
1404	CALL message( 'check_parameters', 'PA0444', 1, 2, 0, 6, 0 )
1405	ENDIF
1406
1407	CALL set_bc_scalars( 'q', bc_q_b, bc_q_t, ibc_q_b, ibc_q_t, &
1408	'PA0071', 'PA0072' )
1409
1410	IF ( surface_waterflux == 9999999.9_wp ) THEN
1411	constant_waterflux = .FALSE.
1412	IF ( large_scale_forcing .OR. land_surface ) THEN
1413	IF ( ibc_q_b == 0 ) THEN
1414	constant_waterflux = .FALSE.
1415	ELSEIF ( ibc_q_b == 1 ) THEN
1416	constant_waterflux = .TRUE.
1417	ENDIF
1418	ENDIF
1419	ELSE
1420	constant_waterflux = .TRUE.
1421	ENDIF
1422
1423	CALL check_bc_scalars( 'q', bc_q_b, ibc_q_b, 'PA0073', 'PA0074', &
1424	constant_waterflux, q_surface_initial_change )
1425
1426	ENDIF
1427
1428	IF ( passive_scalar ) THEN
1429
1430	IF ( ANY( wall_scalarflux /= 0.0_wp ) .AND. &
1431	surface_scalarflux == 9999999.9_wp ) THEN
1432	message_string = 'wall_scalarflux additionally requires ' // &
1433	'setting of surface_scalarflux'
1434	CALL message( 'check_parameters', 'PA0445', 1, 2, 0, 6, 0 )
1435	ENDIF
1436
1437	IF ( surface_scalarflux == 9999999.9_wp ) constant_scalarflux = .FALSE.
1438
1439	CALL set_bc_scalars( 's', bc_s_b, bc_s_t, ibc_s_b, ibc_s_t, &
1440	'PA0071', 'PA0072' )
1441
1442	CALL check_bc_scalars( 's', bc_s_b, ibc_s_b, 'PA0073', 'PA0074', &
1443	constant_scalarflux, s_surface_initial_change )
1444
1445	IF ( top_scalarflux == 9999999.9_wp ) constant_top_scalarflux = .FALSE.
1446	!
1447	!-- A fixed scalar concentration at the top implies Dirichlet boundary
1448	!-- condition for scalar. Hence, in this case specification of a constant
1449	!-- scalar flux is forbidden.
1450	IF ( ( ibc_s_t == 0 .OR. ibc_s_t == 2 ) .AND. constant_top_scalarflux &
1451	.AND. top_scalarflux /= 0.0_wp ) THEN
1452	message_string = 'boundary condition: bc_s_t = "' // &
1453	TRIM( bc_s_t ) // '" is not allowed with ' // &
1454	'top_scalarflux /= 0.0'
1455	CALL message( 'check_parameters', 'PA0441', 1, 2, 0, 6, 0 )
1456	ENDIF
1457	ENDIF
1458
1459	!
1460	!-- Boundary conditions for chemical species
1461	IF ( air_chemistry ) CALL chem_boundary_conds( 'init' )
1462
1463	!
1464	!-- Boundary conditions for horizontal components of wind speed
1465	IF ( bc_uv_b == 'dirichlet' ) THEN
1466	ibc_uv_b = 0
1467	ELSEIF ( bc_uv_b == 'neumann' ) THEN
1468	ibc_uv_b = 1
1469	IF ( constant_flux_layer ) THEN
1470	message_string = 'boundary condition: bc_uv_b = "' // &
1471	TRIM( bc_uv_b ) // '" is not allowed with constant_flux_layer' &
1472	// ' = .TRUE.'
1473	CALL message( 'check_parameters', 'PA0075', 1, 2, 0, 6, 0 )
1474	ENDIF
1475	ELSE
1476	message_string = 'unknown boundary condition: bc_uv_b = "' // &
1477	TRIM( bc_uv_b ) // '"'
1478	CALL message( 'check_parameters', 'PA0076', 1, 2, 0, 6, 0 )
1479	ENDIF
1480	!
1481	!-- In case of coupled simulations u and v at the ground in atmosphere will be
1482	!-- assigned with the u and v values of the ocean surface
1483	IF ( coupling_mode == 'atmosphere_to_ocean' ) THEN
1484	ibc_uv_b = 2
1485	ENDIF
1486
1487	IF ( coupling_mode == 'ocean_to_atmosphere' ) THEN
1488	bc_uv_t = 'neumann'
1489	ibc_uv_t = 1
1490	ELSE
1491	IF ( bc_uv_t == 'dirichlet' .OR. bc_uv_t == 'dirichlet_0' ) THEN
1492	ibc_uv_t = 0
1493	IF ( bc_uv_t == 'dirichlet_0' ) THEN
1494	!
1495	!-- Velocities for the initial u,v-profiles are set zero at the top
1496	!-- in case of dirichlet_0 conditions
1497	u_init(nzt+1) = 0.0_wp
1498	v_init(nzt+1) = 0.0_wp
1499	ENDIF
1500	ELSEIF ( bc_uv_t == 'neumann' ) THEN
1501	ibc_uv_t = 1
1502	ELSEIF ( bc_uv_t == 'nested' .OR. bc_uv_t == 'nesting_offline' ) THEN
1503	ibc_uv_t = 3
1504	ELSE
1505	message_string = 'unknown boundary condition: bc_uv_t = "' // &
1506	TRIM( bc_uv_t ) // '"'
1507	CALL message( 'check_parameters', 'PA0077', 1, 2, 0, 6, 0 )
1508	ENDIF
1509	ENDIF
1510
1511	!
1512	!-- Compute and check, respectively, the Rayleigh Damping parameter
1513	IF ( rayleigh_damping_factor == -1.0_wp ) THEN
1514	rayleigh_damping_factor = 0.0_wp
1515	ELSE
1516	IF ( rayleigh_damping_factor < 0.0_wp .OR. &
1517	rayleigh_damping_factor > 1.0_wp ) THEN
1518	WRITE( message_string, * ) 'rayleigh_damping_factor = ', &
1519	rayleigh_damping_factor, ' out of range [0.0,1.0]'
1520	CALL message( 'check_parameters', 'PA0078', 1, 2, 0, 6, 0 )
1521	ENDIF
1522	ENDIF
1523
1524	IF ( rayleigh_damping_height == -1.0_wp ) THEN
1525	IF ( .NOT. ocean_mode ) THEN
1526	rayleigh_damping_height = 0.66666666666_wp * zu(nzt)
1527	ELSE
1528	rayleigh_damping_height = 0.66666666666_wp * zu(nzb)
1529	ENDIF
1530	ELSE
1531	IF ( .NOT. ocean_mode ) THEN
1532	IF ( rayleigh_damping_height < 0.0_wp .OR. &
1533	rayleigh_damping_height > zu(nzt) ) THEN
1534	WRITE( message_string, * ) 'rayleigh_damping_height = ', &
1535	rayleigh_damping_height, ' out of range [0.0,', zu(nzt), ']'
1536	CALL message( 'check_parameters', 'PA0079', 1, 2, 0, 6, 0 )
1537	ENDIF
1538	ELSE
1539	IF ( rayleigh_damping_height > 0.0_wp .OR. &
1540	rayleigh_damping_height < zu(nzb) ) THEN
1541	WRITE( message_string, * ) 'rayleigh_damping_height = ', &
1542	rayleigh_damping_height, ' out of range [0.0,', zu(nzb), ']'
1543	CALL message( 'check_parameters', 'PA0079', 1, 2, 0, 6, 0 )
1544	ENDIF
1545	ENDIF
1546	ENDIF
1547
1548	!
1549	!-- Check number of chosen statistic regions
1550	IF ( statistic_regions < 0 ) THEN
1551	WRITE ( message_string, * ) 'number of statistic_regions = ', &
1552	statistic_regions+1, ' is not allowed'
1553	CALL message( 'check_parameters', 'PA0082', 1, 2, 0, 6, 0 )
1554	ENDIF
1555	IF ( normalizing_region > statistic_regions .OR. &
1556	normalizing_region < 0) THEN
1557	WRITE ( message_string, * ) 'normalizing_region = ', &
1558	normalizing_region, ' must be >= 0 and <= ',statistic_regions, &
1559	' (value of statistic_regions)'
1560	CALL message( 'check_parameters', 'PA0083', 1, 2, 0, 6, 0 )
1561	ENDIF
1562
1563	!
1564	!-- Set the default intervals for data output, if necessary
1565	!-- NOTE: dt_dosp has already been set in spectra_parin
1566	IF ( dt_data_output /= 9999999.9_wp ) THEN
1567	IF ( dt_dopr == 9999999.9_wp ) dt_dopr = dt_data_output
1568	IF ( dt_dopts == 9999999.9_wp ) dt_dopts = dt_data_output
1569	IF ( dt_do2d_xy == 9999999.9_wp ) dt_do2d_xy = dt_data_output
1570	IF ( dt_do2d_xz == 9999999.9_wp ) dt_do2d_xz = dt_data_output
1571	IF ( dt_do2d_yz == 9999999.9_wp ) dt_do2d_yz = dt_data_output
1572	IF ( dt_do3d == 9999999.9_wp ) dt_do3d = dt_data_output
1573	IF ( dt_data_output_av == 9999999.9_wp ) dt_data_output_av = dt_data_output
1574	DO mid = 1, max_masks
1575	IF ( dt_domask(mid) == 9999999.9_wp ) dt_domask(mid) = dt_data_output
1576	ENDDO
1577	ENDIF
1578
1579	!
1580	!-- Set the default skip time intervals for data output, if necessary
1581	IF ( skip_time_dopr == 9999999.9_wp ) &
1582	skip_time_dopr = skip_time_data_output
1583	IF ( skip_time_do2d_xy == 9999999.9_wp ) &
1584	skip_time_do2d_xy = skip_time_data_output
1585	IF ( skip_time_do2d_xz == 9999999.9_wp ) &
1586	skip_time_do2d_xz = skip_time_data_output
1587	IF ( skip_time_do2d_yz == 9999999.9_wp ) &
1588	skip_time_do2d_yz = skip_time_data_output
1589	IF ( skip_time_do3d == 9999999.9_wp ) &
1590	skip_time_do3d = skip_time_data_output
1591	IF ( skip_time_data_output_av == 9999999.9_wp ) &
1592	skip_time_data_output_av = skip_time_data_output
1593	DO mid = 1, max_masks
1594	IF ( skip_time_domask(mid) == 9999999.9_wp ) &
1595	skip_time_domask(mid) = skip_time_data_output
1596	ENDDO
1597
1598	!
1599	!-- Check the average intervals (first for 3d-data, then for profiles)
1600	IF ( averaging_interval > dt_data_output_av ) THEN
1601	WRITE( message_string, * ) 'averaging_interval = ', &
1602	averaging_interval, ' must be <= dt_data_output_av = ', &
1603	dt_data_output_av
1604	CALL message( 'check_parameters', 'PA0085', 1, 2, 0, 6, 0 )
1605	ENDIF
1606
1607	IF ( averaging_interval_pr == 9999999.9_wp ) THEN
1608	averaging_interval_pr = averaging_interval
1609	ENDIF
1610
1611	IF ( averaging_interval_pr > dt_dopr ) THEN
1612	WRITE( message_string, * ) 'averaging_interval_pr = ', &
1613	averaging_interval_pr, ' must be <= dt_dopr = ', dt_dopr
1614	CALL message( 'check_parameters', 'PA0086', 1, 2, 0, 6, 0 )
1615	ENDIF
1616
1617	!
1618	!-- Set the default interval for profiles entering the temporal average
1619	IF ( dt_averaging_input_pr == 9999999.9_wp ) THEN
1620	dt_averaging_input_pr = dt_averaging_input
1621	ENDIF
1622
1623	!
1624	!-- Set the default interval for the output of timeseries to a reasonable
1625	!-- value (tries to minimize the number of calls of flow_statistics)
1626	IF ( dt_dots == 9999999.9_wp ) THEN
1627	IF ( averaging_interval_pr == 0.0_wp ) THEN
1628	dt_dots = MIN( dt_run_control, dt_dopr )
1629	ELSE
1630	dt_dots = MIN( dt_run_control, dt_averaging_input_pr )
1631	ENDIF
1632	ENDIF
1633
1634	!
1635	!-- Check the sample rate for averaging (first for 3d-data, then for profiles)
1636	IF ( dt_averaging_input > averaging_interval ) THEN
1637	WRITE( message_string, * ) 'dt_averaging_input = ', &
1638	dt_averaging_input, ' must be <= averaging_interval = ', &
1639	averaging_interval
1640	CALL message( 'check_parameters', 'PA0088', 1, 2, 0, 6, 0 )
1641	ENDIF
1642
1643	IF ( dt_averaging_input_pr > averaging_interval_pr ) THEN
1644	WRITE( message_string, * ) 'dt_averaging_input_pr = ', &
1645	dt_averaging_input_pr, ' must be <= averaging_interval_pr = ', &
1646	averaging_interval_pr
1647	CALL message( 'check_parameters', 'PA0089', 1, 2, 0, 6, 0 )
1648	ENDIF
1649
1650	!
1651	!-- Determine the number of output profiles and check whether they are
1652	!-- permissible
1653	DO WHILE ( data_output_pr(dopr_n+1) /= ' ' )
1654
1655	dopr_n = dopr_n + 1
1656	i = dopr_n
1657
1658	!
1659	!-- Determine internal profile number (for hom, homs)
1660	!-- and store height levels
1661	SELECT CASE ( TRIM( data_output_pr(i) ) )
1662
1663	CASE ( 'u', '#u' )
1664	dopr_index(i) = 1
1665	dopr_unit(i) = 'm/s'
1666	hom(:,2,1,:) = SPREAD( zu, 2, statistic_regions+1 )
1667	IF ( data_output_pr(i)(1:1) == '#' ) THEN
1668	dopr_initial_index(i) = 5
1669	hom(:,2,5,:) = SPREAD( zu, 2, statistic_regions+1 )
1670	data_output_pr(i) = data_output_pr(i)(2:)
1671	ENDIF
1672
1673	CASE ( 'v', '#v' )
1674	dopr_index(i) = 2
1675	dopr_unit(i) = 'm/s'
1676	hom(:,2,2,:) = SPREAD( zu, 2, statistic_regions+1 )
1677	IF ( data_output_pr(i)(1:1) == '#' ) THEN
1678	dopr_initial_index(i) = 6
1679	hom(:,2,6,:) = SPREAD( zu, 2, statistic_regions+1 )
1680	data_output_pr(i) = data_output_pr(i)(2:)
1681	ENDIF
1682
1683	CASE ( 'w' )
1684	dopr_index(i) = 3
1685	dopr_unit(i) = 'm/s'
1686	hom(:,2,3,:) = SPREAD( zw, 2, statistic_regions+1 )
1687
1688	CASE ( 'theta', '#theta' )
1689	IF ( .NOT. bulk_cloud_model ) THEN
1690	dopr_index(i) = 4
1691	dopr_unit(i) = 'K'
1692	hom(:,2,4,:) = SPREAD( zu, 2, statistic_regions+1 )
1693	IF ( data_output_pr(i)(1:1) == '#' ) THEN
1694	dopr_initial_index(i) = 7
1695	hom(:,2,7,:) = SPREAD( zu, 2, statistic_regions+1 )
1696	hom(nzb,2,7,:) = 0.0_wp ! because zu(nzb) is negative
1697	data_output_pr(i) = data_output_pr(i)(2:)
1698	ENDIF
1699	ELSE
1700	dopr_index(i) = 43
1701	dopr_unit(i) = 'K'
1702	hom(:,2,43,:) = SPREAD( zu, 2, statistic_regions+1 )
1703	IF ( data_output_pr(i)(1:1) == '#' ) THEN
1704	dopr_initial_index(i) = 28
1705	hom(:,2,28,:) = SPREAD( zu, 2, statistic_regions+1 )
1706	hom(nzb,2,28,:) = 0.0_wp ! because zu(nzb) is negative
1707	data_output_pr(i) = data_output_pr(i)(2:)
1708	ENDIF
1709	ENDIF
1710
1711	CASE ( 'e', '#e' )
1712	dopr_index(i) = 8
1713	dopr_unit(i) = 'm2/s2'
1714	hom(:,2,8,:) = SPREAD( zu, 2, statistic_regions+1 )
1715	hom(nzb,2,8,:) = 0.0_wp
1716	IF ( data_output_pr(i)(1:1) == '#' ) THEN
1717	dopr_initial_index(i) = 8
1718	hom(:,2,8,:) = SPREAD( zu, 2, statistic_regions+1 )
1719	data_output_pr(i) = data_output_pr(i)(2:)
1720	ENDIF
1721
1722	CASE ( 'km', '#km' )
1723	dopr_index(i) = 9
1724	dopr_unit(i) = 'm2/s'
1725	hom(:,2,9,:) = SPREAD( zu, 2, statistic_regions+1 )
1726	hom(nzb,2,9,:) = 0.0_wp
1727	IF ( data_output_pr(i)(1:1) == '#' ) THEN
1728	dopr_initial_index(i) = 23
1729	hom(:,2,23,:) = hom(:,2,9,:)
1730	data_output_pr(i) = data_output_pr(i)(2:)
1731	ENDIF
1732
1733	CASE ( 'kh', '#kh' )
1734	dopr_index(i) = 10
1735	dopr_unit(i) = 'm2/s'
1736	hom(:,2,10,:) = SPREAD( zu, 2, statistic_regions+1 )
1737	hom(nzb,2,10,:) = 0.0_wp
1738	IF ( data_output_pr(i)(1:1) == '#' ) THEN
1739	dopr_initial_index(i) = 24
1740	hom(:,2,24,:) = hom(:,2,10,:)
1741	data_output_pr(i) = data_output_pr(i)(2:)
1742	ENDIF
1743
1744	CASE ( 'l', '#l' )
1745	dopr_index(i) = 11
1746	dopr_unit(i) = 'm'
1747	hom(:,2,11,:) = SPREAD( zu, 2, statistic_regions+1 )
1748	hom(nzb,2,11,:) = 0.0_wp
1749	IF ( data_output_pr(i)(1:1) == '#' ) THEN
1750	dopr_initial_index(i) = 25
1751	hom(:,2,25,:) = hom(:,2,11,:)
1752	data_output_pr(i) = data_output_pr(i)(2:)
1753	ENDIF
1754
1755	CASE ( 'w"u"' )
1756	dopr_index(i) = 12
1757	dopr_unit(i) = TRIM ( momentumflux_output_unit )
1758	hom(:,2,12,:) = SPREAD( zw, 2, statistic_regions+1 )
1759	IF ( constant_flux_layer ) hom(nzb,2,12,:) = zu(1)
1760
1761	CASE ( 'wu' )
1762	dopr_index(i) = 13
1763	dopr_unit(i) = TRIM ( momentumflux_output_unit )
1764	hom(:,2,13,:) = SPREAD( zw, 2, statistic_regions+1 )
1765
1766	CASE ( 'w"v"' )
1767	dopr_index(i) = 14
1768	dopr_unit(i) = TRIM ( momentumflux_output_unit )
1769	hom(:,2,14,:) = SPREAD( zw, 2, statistic_regions+1 )
1770	IF ( constant_flux_layer ) hom(nzb,2,14,:) = zu(1)
1771
1772	CASE ( 'wv' )
1773	dopr_index(i) = 15
1774	dopr_unit(i) = TRIM ( momentumflux_output_unit )
1775	hom(:,2,15,:) = SPREAD( zw, 2, statistic_regions+1 )
1776
1777	CASE ( 'w"theta"' )
1778	dopr_index(i) = 16
1779	dopr_unit(i) = TRIM ( heatflux_output_unit )
1780	hom(:,2,16,:) = SPREAD( zw, 2, statistic_regions+1 )
1781
1782	CASE ( 'wtheta' )
1783	dopr_index(i) = 17
1784	dopr_unit(i) = TRIM ( heatflux_output_unit )
1785	hom(:,2,17,:) = SPREAD( zw, 2, statistic_regions+1 )
1786
1787	CASE ( 'wtheta' )
1788	dopr_index(i) = 18
1789	dopr_unit(i) = TRIM ( heatflux_output_unit )
1790	hom(:,2,18,:) = SPREAD( zw, 2, statistic_regions+1 )
1791
1792	CASE ( 'wu' )
1793	dopr_index(i) = 19
1794	dopr_unit(i) = TRIM ( momentumflux_output_unit )
1795	hom(:,2,19,:) = SPREAD( zw, 2, statistic_regions+1 )
1796	IF ( constant_flux_layer ) hom(nzb,2,19,:) = zu(1)
1797
1798	CASE ( 'wv' )
1799	dopr_index(i) = 20
1800	dopr_unit(i) = TRIM ( momentumflux_output_unit )
1801	hom(:,2,20,:) = SPREAD( zw, 2, statistic_regions+1 )
1802	IF ( constant_flux_layer ) hom(nzb,2,20,:) = zu(1)
1803
1804	CASE ( 'wthetaBC' )
1805	dopr_index(i) = 21
1806	dopr_unit(i) = TRIM ( heatflux_output_unit )
1807	hom(:,2,21,:) = SPREAD( zw, 2, statistic_regions+1 )
1808
1809	CASE ( 'wthetaBC' )
1810	dopr_index(i) = 22
1811	dopr_unit(i) = TRIM ( heatflux_output_unit )
1812	hom(:,2,22,:) = SPREAD( zw, 2, statistic_regions+1 )
1813
1814	CASE ( 'u*2' )
1815	dopr_index(i) = 30
1816	dopr_unit(i) = 'm2/s2'
1817	hom(:,2,30,:) = SPREAD( zu, 2, statistic_regions+1 )
1818
1819	CASE ( 'v*2' )
1820	dopr_index(i) = 31
1821	dopr_unit(i) = 'm2/s2'
1822	hom(:,2,31,:) = SPREAD( zu, 2, statistic_regions+1 )
1823
1824	CASE ( 'w*2' )
1825	dopr_index(i) = 32
1826	dopr_unit(i) = 'm2/s2'
1827	hom(:,2,32,:) = SPREAD( zw, 2, statistic_regions+1 )
1828
1829	CASE ( 'theta*2' )
1830	dopr_index(i) = 33
1831	dopr_unit(i) = 'K2'
1832	hom(:,2,33,:) = SPREAD( zu, 2, statistic_regions+1 )
1833
1834	CASE ( 'e*' )
1835	dopr_index(i) = 34
1836	dopr_unit(i) = 'm2/s2'
1837	hom(:,2,34,:) = SPREAD( zu, 2, statistic_regions+1 )
1838
1839	CASE ( 'w2theta' )
1840	dopr_index(i) = 35
1841	dopr_unit(i) = 'K m2/s2'
1842	hom(:,2,35,:) = SPREAD( zw, 2, statistic_regions+1 )
1843
1844	CASE ( 'wtheta2' )
1845	dopr_index(i) = 36
1846	dopr_unit(i) = 'K2 m/s'
1847	hom(:,2,36,:) = SPREAD( zw, 2, statistic_regions+1 )
1848
1849	CASE ( 'we' )
1850	dopr_index(i) = 37
1851	dopr_unit(i) = 'm3/s3'
1852	hom(:,2,37,:) = SPREAD( zw, 2, statistic_regions+1 )
1853
1854	CASE ( 'w*3' )
1855	dopr_index(i) = 38
1856	dopr_unit(i) = 'm3/s3'
1857	hom(:,2,38,:) = SPREAD( zw, 2, statistic_regions+1 )
1858
1859	CASE ( 'Sw' )
1860	dopr_index(i) = 39
1861	dopr_unit(i) = 'none'
1862	hom(:,2,39,:) = SPREAD( zw, 2, statistic_regions+1 )
1863
1864	CASE ( 'p' )
1865	dopr_index(i) = 40
1866	dopr_unit(i) = 'Pa'
1867	hom(:,2,40,:) = SPREAD( zu, 2, statistic_regions+1 )
1868
1869	CASE ( 'q', '#q' )
1870	IF ( .NOT. humidity ) THEN
1871	message_string = 'data_output_pr = ' // &
1872	TRIM( data_output_pr(i) ) // ' is not imp' // &
1873	'lemented for humidity = .FALSE.'
1874	CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
1875	ELSE
1876	dopr_index(i) = 41
1877	dopr_unit(i) = 'kg/kg'
1878	hom(:,2,41,:) = SPREAD( zu, 2, statistic_regions+1 )
1879	IF ( data_output_pr(i)(1:1) == '#' ) THEN
1880	dopr_initial_index(i) = 26
1881	hom(:,2,26,:) = SPREAD( zu, 2, statistic_regions+1 )
1882	hom(nzb,2,26,:) = 0.0_wp ! because zu(nzb) is negative
1883	data_output_pr(i) = data_output_pr(i)(2:)
1884	ENDIF
1885	ENDIF
1886
1887	CASE ( 's', '#s' )
1888	IF ( .NOT. passive_scalar ) THEN
1889	message_string = 'data_output_pr = ' // &
1890	TRIM( data_output_pr(i) ) // ' is not imp' // &
1891	'lemented for passive_scalar = .FALSE.'
1892	CALL message( 'check_parameters', 'PA0093', 1, 2, 0, 6, 0 )
1893	ELSE
1894	dopr_index(i) = 115
1895	dopr_unit(i) = 'kg/m3'
1896	hom(:,2,115,:) = SPREAD( zu, 2, statistic_regions+1 )
1897	IF ( data_output_pr(i)(1:1) == '#' ) THEN
1898	dopr_initial_index(i) = 121
1899	hom(:,2,121,:) = SPREAD( zu, 2, statistic_regions+1 )
1900	hom(nzb,2,121,:) = 0.0_wp ! because zu(nzb) is negative
1901	data_output_pr(i) = data_output_pr(i)(2:)
1902	ENDIF
1903	ENDIF
1904
1905	CASE ( 'qv', '#qv' )
1906	IF ( .NOT. bulk_cloud_model ) THEN
1907	dopr_index(i) = 41
1908	dopr_unit(i) = 'kg/kg'
1909	hom(:,2,41,:) = SPREAD( zu, 2, statistic_regions+1 )
1910	IF ( data_output_pr(i)(1:1) == '#' ) THEN
1911	dopr_initial_index(i) = 26
1912	hom(:,2,26,:) = SPREAD( zu, 2, statistic_regions+1 )
1913	hom(nzb,2,26,:) = 0.0_wp ! because zu(nzb) is negative
1914	data_output_pr(i) = data_output_pr(i)(2:)
1915	ENDIF
1916	ELSE
1917	dopr_index(i) = 42
1918	dopr_unit(i) = 'kg/kg'
1919	hom(:,2,42,:) = SPREAD( zu, 2, statistic_regions+1 )
1920	IF ( data_output_pr(i)(1:1) == '#' ) THEN
1921	dopr_initial_index(i) = 27
1922	hom(:,2,27,:) = SPREAD( zu, 2, statistic_regions+1 )
1923	hom(nzb,2,27,:) = 0.0_wp ! because zu(nzb) is negative
1924	data_output_pr(i) = data_output_pr(i)(2:)
1925	ENDIF
1926	ENDIF
1927
1928	CASE ( 'thetal', '#thetal' )
1929	IF ( .NOT. bulk_cloud_model ) THEN
1930	message_string = 'data_output_pr = ' // &
1931	TRIM( data_output_pr(i) ) // ' is not imp' // &
1932	'lemented for bulk_cloud_model = .FALSE.'
1933	CALL message( 'check_parameters', 'PA0094', 1, 2, 0, 6, 0 )
1934	ELSE
1935	dopr_index(i) = 4
1936	dopr_unit(i) = 'K'
1937	hom(:,2,4,:) = SPREAD( zu, 2, statistic_regions+1 )
1938	IF ( data_output_pr(i)(1:1) == '#' ) THEN
1939	dopr_initial_index(i) = 7
1940	hom(:,2,7,:) = SPREAD( zu, 2, statistic_regions+1 )
1941	hom(nzb,2,7,:) = 0.0_wp ! because zu(nzb) is negative
1942	data_output_pr(i) = data_output_pr(i)(2:)
1943	ENDIF
1944	ENDIF
1945
1946	CASE ( 'thetav', '#thetav' )
1947	dopr_index(i) = 44
1948	dopr_unit(i) = 'K'
1949	hom(:,2,44,:) = SPREAD( zu, 2, statistic_regions+1 )
1950	IF ( data_output_pr(i)(1:1) == '#' ) THEN
1951	dopr_initial_index(i) = 29
1952	hom(:,2,29,:) = SPREAD( zu, 2, statistic_regions+1 )
1953	hom(nzb,2,29,:) = 0.0_wp ! because zu(nzb) is negative
1954	data_output_pr(i) = data_output_pr(i)(2:)
1955	ENDIF
1956
1957	CASE ( 'w"thetav"' )
1958	dopr_index(i) = 45
1959	dopr_unit(i) = TRIM ( heatflux_output_unit )
1960	hom(:,2,45,:) = SPREAD( zw, 2, statistic_regions+1 )
1961
1962	CASE ( 'wthetav' )
1963	dopr_index(i) = 46
1964	dopr_unit(i) = TRIM ( heatflux_output_unit )
1965	hom(:,2,46,:) = SPREAD( zw, 2, statistic_regions+1 )
1966
1967	CASE ( 'wthetav' )
1968	dopr_index(i) = 47
1969	dopr_unit(i) = TRIM ( heatflux_output_unit )
1970	hom(:,2,47,:) = SPREAD( zw, 2, statistic_regions+1 )
1971
1972	CASE ( 'w"q"' )
1973	IF ( .NOT. humidity ) THEN
1974	message_string = 'data_output_pr = ' // &
1975	TRIM( data_output_pr(i) ) // ' is not imp' // &
1976	'lemented for humidity = .FALSE.'
1977	CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
1978	ELSE
1979	dopr_index(i) = 48
1980	dopr_unit(i) = TRIM ( waterflux_output_unit )
1981	hom(:,2,48,:) = SPREAD( zw, 2, statistic_regions+1 )
1982	ENDIF
1983
1984	CASE ( 'wq' )
1985	IF ( .NOT. humidity ) THEN
1986	message_string = 'data_output_pr = ' // &
1987	TRIM( data_output_pr(i) ) // ' is not imp' // &
1988	'lemented for humidity = .FALSE.'
1989	CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
1990	ELSE
1991	dopr_index(i) = 49
1992	dopr_unit(i) = TRIM ( waterflux_output_unit )
1993	hom(:,2,49,:) = SPREAD( zw, 2, statistic_regions+1 )
1994	ENDIF
1995
1996	CASE ( 'wq' )
1997	IF ( .NOT. humidity ) THEN
1998	message_string = 'data_output_pr = ' // &
1999	TRIM( data_output_pr(i) ) // ' is not imp' // &
2000	'lemented for humidity = .FALSE.'
2001	CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
2002	ELSE
2003	dopr_index(i) = 50
2004	dopr_unit(i) = TRIM ( waterflux_output_unit )
2005	hom(:,2,50,:) = SPREAD( zw, 2, statistic_regions+1 )
2006	ENDIF
2007
2008	CASE ( 'w"s"' )
2009	IF ( .NOT. passive_scalar ) THEN
2010	message_string = 'data_output_pr = ' // &
2011	TRIM( data_output_pr(i) ) // ' is not imp' // &
2012	'lemented for passive_scalar = .FALSE.'
2013	CALL message( 'check_parameters', 'PA0093', 1, 2, 0, 6, 0 )
2014	ELSE
2015	dopr_index(i) = 117
2016	dopr_unit(i) = 'kg/m3 m/s'
2017	hom(:,2,117,:) = SPREAD( zw, 2, statistic_regions+1 )
2018	ENDIF
2019
2020	CASE ( 'ws' )
2021	IF ( .NOT. passive_scalar ) THEN
2022	message_string = 'data_output_pr = ' // &
2023	TRIM( data_output_pr(i) ) // ' is not imp' // &
2024	'lemented for passive_scalar = .FALSE.'
2025	CALL message( 'check_parameters', 'PA0093', 1, 2, 0, 6, 0 )
2026	ELSE
2027	dopr_index(i) = 114
2028	dopr_unit(i) = 'kg/m3 m/s'
2029	hom(:,2,114,:) = SPREAD( zw, 2, statistic_regions+1 )
2030	ENDIF
2031
2032	CASE ( 'ws' )
2033	IF ( .NOT. passive_scalar ) THEN
2034	message_string = 'data_output_pr = ' // &
2035	TRIM( data_output_pr(i) ) // ' is not imp' // &
2036	'lemented for passive_scalar = .FALSE.'
2037	CALL message( 'check_parameters', 'PA0093', 1, 2, 0, 6, 0 )
2038	ELSE
2039	dopr_index(i) = 118
2040	dopr_unit(i) = 'kg/m3 m/s'
2041	hom(:,2,118,:) = SPREAD( zw, 2, statistic_regions+1 )
2042	ENDIF
2043
2044	CASE ( 'w"qv"' )
2045	IF ( humidity .AND. .NOT. bulk_cloud_model ) THEN
2046	dopr_index(i) = 48
2047	dopr_unit(i) = TRIM ( waterflux_output_unit )
2048	hom(:,2,48,:) = SPREAD( zw, 2, statistic_regions+1 )
2049	ELSEIF ( humidity .AND. bulk_cloud_model ) THEN
2050	dopr_index(i) = 51
2051	dopr_unit(i) = TRIM ( waterflux_output_unit )
2052	hom(:,2,51,:) = SPREAD( zw, 2, statistic_regions+1 )
2053	ELSE
2054	message_string = 'data_output_pr = ' // &
2055	TRIM( data_output_pr(i) ) // ' is not imp' // &
2056	'lemented for bulk_cloud_model = .FALSE. ' // &
2057	'and humidity = .FALSE.'
2058	CALL message( 'check_parameters', 'PA0095', 1, 2, 0, 6, 0 )
2059	ENDIF
2060
2061	CASE ( 'wqv' )
2062	IF ( humidity .AND. .NOT. bulk_cloud_model ) THEN
2063	dopr_index(i) = 49
2064	dopr_unit(i) = TRIM ( waterflux_output_unit )
2065	hom(:,2,49,:) = SPREAD( zw, 2, statistic_regions+1 )
2066	ELSEIF( humidity .AND. bulk_cloud_model ) THEN
2067	dopr_index(i) = 52
2068	dopr_unit(i) = TRIM ( waterflux_output_unit )
2069	hom(:,2,52,:) = SPREAD( zw, 2, statistic_regions+1 )
2070	ELSE
2071	message_string = 'data_output_pr = ' // &
2072	TRIM( data_output_pr(i) ) // ' is not imp' // &
2073	'lemented for bulk_cloud_model = .FALSE. ' // &
2074	'and humidity = .FALSE.'
2075	CALL message( 'check_parameters', 'PA0095', 1, 2, 0, 6, 0 )
2076	ENDIF
2077
2078	CASE ( 'wqv' )
2079	IF ( humidity .AND. .NOT. bulk_cloud_model ) THEN
2080	dopr_index(i) = 50
2081	dopr_unit(i) = TRIM ( waterflux_output_unit )
2082	hom(:,2,50,:) = SPREAD( zw, 2, statistic_regions+1 )
2083	ELSEIF ( humidity .AND. bulk_cloud_model ) THEN
2084	dopr_index(i) = 53
2085	dopr_unit(i) = TRIM ( waterflux_output_unit )
2086	hom(:,2,53,:) = SPREAD( zw, 2, statistic_regions+1 )
2087	ELSE
2088	message_string = 'data_output_pr = ' // &
2089	TRIM( data_output_pr(i) ) // ' is not imp' // &
2090	'lemented for bulk_cloud_model = .FALSE. ' // &
2091	'and humidity = .FALSE.'
2092	CALL message( 'check_parameters', 'PA0095', 1, 2, 0, 6, 0 )
2093	ENDIF
2094
2095	CASE ( 'ql' )
2096	IF ( .NOT. bulk_cloud_model .AND. .NOT. cloud_droplets ) THEN
2097	message_string = 'data_output_pr = ' // &
2098	TRIM( data_output_pr(i) ) // ' is not imp' // &
2099	'lemented for bulk_cloud_model = .FALSE. ' // &
2100	'and cloud_droplets = .FALSE.'
2101	CALL message( 'check_parameters', 'PA0096', 1, 2, 0, 6, 0 )
2102	ELSE
2103	dopr_index(i) = 54
2104	dopr_unit(i) = 'kg/kg'
2105	hom(:,2,54,:) = SPREAD( zu, 2, statistic_regions+1 )
2106	ENDIF
2107
2108	CASE ( 'wuu*:dz' )
2109	dopr_index(i) = 55
2110	dopr_unit(i) = 'm2/s3'
2111	hom(:,2,55,:) = SPREAD( zu, 2, statistic_regions+1 )
2112
2113	CASE ( 'wp:dz' )
2114	dopr_index(i) = 56
2115	dopr_unit(i) = 'm2/s3'
2116	hom(:,2,56,:) = SPREAD( zw, 2, statistic_regions+1 )
2117
2118	CASE ( 'w"e:dz' )
2119	dopr_index(i) = 57
2120	dopr_unit(i) = 'm2/s3'
2121	hom(:,2,57,:) = SPREAD( zu, 2, statistic_regions+1 )
2122
2123	CASE ( 'u"theta"' )
2124	dopr_index(i) = 58
2125	dopr_unit(i) = TRIM ( heatflux_output_unit )
2126	hom(:,2,58,:) = SPREAD( zu, 2, statistic_regions+1 )
2127
2128	CASE ( 'utheta' )
2129	dopr_index(i) = 59
2130	dopr_unit(i) = TRIM ( heatflux_output_unit )
2131	hom(:,2,59,:) = SPREAD( zu, 2, statistic_regions+1 )
2132
2133	CASE ( 'utheta_t' )
2134	dopr_index(i) = 60
2135	dopr_unit(i) = TRIM ( heatflux_output_unit )
2136	hom(:,2,60,:) = SPREAD( zu, 2, statistic_regions+1 )
2137
2138	CASE ( 'v"theta"' )
2139	dopr_index(i) = 61
2140	dopr_unit(i) = TRIM ( heatflux_output_unit )
2141	hom(:,2,61,:) = SPREAD( zu, 2, statistic_regions+1 )
2142
2143	CASE ( 'vtheta' )
2144	dopr_index(i) = 62
2145	dopr_unit(i) = TRIM ( heatflux_output_unit )
2146	hom(:,2,62,:) = SPREAD( zu, 2, statistic_regions+1 )
2147
2148	CASE ( 'vtheta_t' )
2149	dopr_index(i) = 63
2150	dopr_unit(i) = TRIM ( heatflux_output_unit )
2151	hom(:,2,63,:) = SPREAD( zu, 2, statistic_regions+1 )
2152
2153	CASE ( 'wp' )
2154	dopr_index(i) = 68
2155	dopr_unit(i) = 'm3/s3'
2156	hom(:,2,68,:) = SPREAD( zu, 2, statistic_regions+1 )
2157
2158	CASE ( 'w"e' )
2159	dopr_index(i) = 69
2160	dopr_unit(i) = 'm3/s3'
2161	hom(:,2,69,:) = SPREAD( zu, 2, statistic_regions+1 )
2162
2163	CASE ( 'q*2' )
2164	IF ( .NOT. humidity ) THEN
2165	message_string = 'data_output_pr = ' // &
2166	TRIM( data_output_pr(i) ) // ' is not imp' // &
2167	'lemented for humidity = .FALSE.'
2168	CALL message( 'check_parameters', 'PA0092', 1, 2, 0, 6, 0 )
2169	ELSE
2170	dopr_index(i) = 70
2171	dopr_unit(i) = 'kg2/kg2'
2172	hom(:,2,70,:) = SPREAD( zu, 2, statistic_regions+1 )
2173	ENDIF
2174
2175	CASE ( 'hyp' )
2176	dopr_index(i) = 72
2177	dopr_unit(i) = 'hPa'
2178	hom(:,2,72,:) = SPREAD( zu, 2, statistic_regions+1 )
2179
2180	CASE ( 'rho' )
2181	dopr_index(i) = 119
2182	dopr_unit(i) = 'kg/m3'
2183	hom(:,2,119,:) = SPREAD( zu, 2, statistic_regions+1 )
2184
2185	CASE ( 'rho_zw' )
2186	dopr_index(i) = 120
2187	dopr_unit(i) = 'kg/m3'
2188	hom(:,2,120,:) = SPREAD( zw, 2, statistic_regions+1 )
2189
2190	CASE ( 'ug' )
2191	dopr_index(i) = 78
2192	dopr_unit(i) = 'm/s'
2193	hom(:,2,78,:) = SPREAD( zu, 2, statistic_regions+1 )
2194
2195	CASE ( 'vg' )
2196	dopr_index(i) = 79
2197	dopr_unit(i) = 'm/s'
2198	hom(:,2,79,:) = SPREAD( zu, 2, statistic_regions+1 )
2199
2200	CASE ( 'w_subs' )
2201	IF ( .NOT. large_scale_subsidence ) THEN
2202	message_string = 'data_output_pr = ' // &
2203	TRIM( data_output_pr(i) ) // ' is not imp' // &
2204	'lemented for large_scale_subsidence = .FALSE.'
2205	CALL message( 'check_parameters', 'PA0382', 1, 2, 0, 6, 0 )
2206	ELSE
2207	dopr_index(i) = 80
2208	dopr_unit(i) = 'm/s'
2209	hom(:,2,80,:) = SPREAD( zu, 2, statistic_regions+1 )
2210	ENDIF
2211
2212	CASE ( 's*2' )
2213	IF ( .NOT. passive_scalar ) THEN
2214	message_string = 'data_output_pr = ' // &
2215	TRIM( data_output_pr(i) ) // ' is not imp' // &
2216	'lemented for passive_scalar = .FALSE.'
2217	CALL message( 'check_parameters', 'PA0185', 1, 2, 0, 6, 0 )
2218	ELSE
2219	dopr_index(i) = 116
2220	dopr_unit(i) = 'kg2/m6'
2221	hom(:,2,116,:) = SPREAD( zu, 2, statistic_regions+1 )
2222	ENDIF
2223
2224	CASE DEFAULT
2225	unit = 'illegal'
2226	!
2227	!-- Check for other modules
2228	CALL module_interface_check_data_output_pr( data_output_pr(i), i, &
2229	unit, dopr_unit(i) )
2230
2231	!
2232	!-- No valid quantity found
2233	IF ( unit == 'illegal' ) THEN
2234	IF ( data_output_pr_user(1) /= ' ' ) THEN
2235	message_string = 'illegal value for data_output_pr or ' // &
2236	'data_output_pr_user = "' // &
2237	TRIM( data_output_pr(i) ) // '"'
2238	CALL message( 'check_parameters', 'PA0097', 1, 2, 0, 6, 0 )
2239	ELSE
2240	message_string = 'illegal value for data_output_pr = "' // &
2241	TRIM( data_output_pr(i) ) // '"'
2242	CALL message( 'check_parameters', 'PA0098', 1, 2, 0, 6, 0 )
2243	ENDIF
2244	ENDIF
2245
2246	END SELECT
2247
2248	ENDDO
2249
2250
2251	!
2252	!-- Append user-defined data output variables to the standard data output
2253	IF ( data_output_user(1) /= ' ' ) THEN
2254	i = 1
2255	DO WHILE ( data_output(i) /= ' ' .AND. i <= 500 )
2256	i = i + 1
2257	ENDDO
2258	j = 1
2259	DO WHILE ( data_output_user(j) /= ' ' .AND. j <= 500 )
2260	IF ( i > 500 ) THEN
2261	message_string = 'number of output quantitities given by data' // &
2262	'_output and data_output_user exceeds the limit of 500'
2263	CALL message( 'check_parameters', 'PA0102', 1, 2, 0, 6, 0 )
2264	ENDIF
2265	data_output(i) = data_output_user(j)
2266	i = i + 1
2267	j = j + 1
2268	ENDDO
2269	ENDIF
2270
2271	!
2272	!-- Check and set steering parameters for 2d/3d data output and averaging
2273	i = 1
2274	DO WHILE ( data_output(i) /= ' ' .AND. i <= 500 )
2275	!
2276	!-- Check for data averaging
2277	ilen = LEN_TRIM( data_output(i) )
2278	j = 0 ! no data averaging
2279	IF ( ilen > 3 ) THEN
2280	IF ( data_output(i)(ilen-2:ilen) == '_av' ) THEN
2281	j = 1 ! data averaging
2282	data_output(i) = data_output(i)(1:ilen-3)
2283	ENDIF
2284	ENDIF
2285	!
2286	!-- Check for cross section or volume data
2287	ilen = LEN_TRIM( data_output(i) )
2288	k = 0 ! 3d data
2289	var = data_output(i)(1:ilen)
2290	IF ( ilen > 3 ) THEN
2291	IF ( data_output(i)(ilen-2:ilen) == '_xy' .OR. &
2292	data_output(i)(ilen-2:ilen) == '_xz' .OR. &
2293	data_output(i)(ilen-2:ilen) == '_yz' ) THEN
2294	k = 1 ! 2d data
2295	var = data_output(i)(1:ilen-3)
2296	ENDIF
2297	ENDIF
2298
2299	!
2300	!-- Check for allowed value and set units
2301	SELECT CASE ( TRIM( var ) )
2302
2303	CASE ( 'e' )
2304	IF ( constant_diffusion ) THEN
2305	message_string = 'output of "' // TRIM( var ) // '" requi' // &
2306	'res constant_diffusion = .FALSE.'
2307	CALL message( 'check_parameters', 'PA0103', 1, 2, 0, 6, 0 )
2308	ENDIF
2309	unit = 'm2/s2'
2310
2311	CASE ( 'thetal' )
2312	IF ( .NOT. bulk_cloud_model ) THEN
2313	message_string = 'output of "' // TRIM( var ) // '" requi' // &
2314	'res bulk_cloud_model = .TRUE.'
2315	CALL message( 'check_parameters', 'PA0108', 1, 2, 0, 6, 0 )
2316	ENDIF
2317	unit = 'K'
2318
2319	CASE ( 'pc', 'pr' )
2320	IF ( .NOT. particle_advection ) THEN
2321	message_string = 'output of "' // TRIM( var ) // '" requir' // &
2322	'es a "particle_parameters"-NAMELIST in the parameter ' // &
2323	'file (PARIN)'
2324	CALL message( 'check_parameters', 'PA0104', 1, 2, 0, 6, 0 )
2325	ENDIF
2326	IF ( TRIM( var ) == 'pc' ) unit = 'number'
2327	IF ( TRIM( var ) == 'pr' ) unit = 'm'
2328
2329	CASE ( 'q', 'thetav' )
2330	IF ( .NOT. humidity ) THEN
2331	message_string = 'output of "' // TRIM( var ) // '" requi' // &
2332	'res humidity = .TRUE.'
2333	CALL message( 'check_parameters', 'PA0105', 1, 2, 0, 6, 0 )
2334	ENDIF
2335	IF ( TRIM( var ) == 'q' ) unit = 'kg/kg'
2336	IF ( TRIM( var ) == 'thetav' ) unit = 'K'
2337
2338	CASE ( 'ql' )
2339	IF ( .NOT. ( bulk_cloud_model .OR. cloud_droplets ) ) THEN
2340	message_string = 'output of "' // TRIM( var ) // '" requi' // &
2341	'res bulk_cloud_model = .TRUE. or cloud_droplets = .TRUE.'
2342	CALL message( 'check_parameters', 'PA0106', 1, 2, 0, 6, 0 )
2343	ENDIF
2344	unit = 'kg/kg'
2345
2346	CASE ( 'ql_c', 'ql_v', 'ql_vp' )
2347	IF ( .NOT. cloud_droplets ) THEN
2348	message_string = 'output of "' // TRIM( var ) // '" requi' // &
2349	'res cloud_droplets = .TRUE.'
2350	CALL message( 'check_parameters', 'PA0107', 1, 2, 0, 6, 0 )
2351	ENDIF
2352	IF ( TRIM( var ) == 'ql_c' ) unit = 'kg/kg'
2353	IF ( TRIM( var ) == 'ql_v' ) unit = 'm3'
2354	IF ( TRIM( var ) == 'ql_vp' ) unit = 'none'
2355
2356	CASE ( 'qv' )
2357	IF ( .NOT. bulk_cloud_model ) THEN
2358	message_string = 'output of "' // TRIM( var ) // '" requi' // &
2359	'res bulk_cloud_model = .TRUE.'
2360	CALL message( 'check_parameters', 'PA0108', 1, 2, 0, 6, 0 )
2361	ENDIF
2362	unit = 'kg/kg'
2363
2364	CASE ( 's' )
2365	IF ( .NOT. passive_scalar ) THEN
2366	message_string = 'output of "' // TRIM( var ) // '" requi' // &
2367	'res passive_scalar = .TRUE.'
2368	CALL message( 'check_parameters', 'PA0110', 1, 2, 0, 6, 0 )
2369	ENDIF
2370	unit = 'kg/m3'
2371
2372	CASE ( 'p', 'theta', 'u', 'v', 'w' )
2373	IF ( TRIM( var ) == 'p' ) unit = 'Pa'
2374	IF ( TRIM( var ) == 'theta' ) unit = 'K'
2375	IF ( TRIM( var ) == 'u' ) unit = 'm/s'
2376	IF ( TRIM( var ) == 'v' ) unit = 'm/s'
2377	IF ( TRIM( var ) == 'w' ) unit = 'm/s'
2378	CONTINUE
2379
2380	CASE ( 'ghf', 'lwp', 'ol', 'qsws', 'r_a*', &
2381	'shf', 'ssws', 't', 'theta_2m', 'tsurf', 'us', &
2382	'z0', 'z0h', 'z0q*' )
2383	IF ( k == 0 .OR. data_output(i)(ilen-2:ilen) /= '_xy' ) THEN
2384	message_string = 'illegal value for data_output: "' // &
2385	TRIM( var ) // '" & only 2d-horizontal ' // &
2386	'cross sections are allowed for this value'
2387	CALL message( 'check_parameters', 'PA0111', 1, 2, 0, 6, 0 )
2388	ENDIF
2389
2390	IF ( TRIM( var ) == 'lwp*' .AND. .NOT. bulk_cloud_model ) THEN
2391	message_string = 'output of "' // TRIM( var ) // '" requi' // &
2392	'res bulk_cloud_model = .TRUE.'
2393	CALL message( 'check_parameters', 'PA0108', 1, 2, 0, 6, 0 )
2394	ENDIF
2395	IF ( TRIM( var ) == 'qsws*' .AND. .NOT. humidity ) THEN
2396	message_string = 'output of "' // TRIM( var ) // '" requi' // &
2397	'res humidity = .TRUE.'
2398	CALL message( 'check_parameters', 'PA0322', 1, 2, 0, 6, 0 )
2399	ENDIF
2400
2401	IF ( TRIM( var ) == 'ghf*' .AND. .NOT. land_surface ) THEN
2402	message_string = 'output of "' // TRIM( var ) // '" requi' // &
2403	'res land_surface = .TRUE.'
2404	CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 )
2405	ENDIF
2406
2407	IF ( ( TRIM( var ) == 'r_a' .OR. TRIM( var ) == 'ghf' ) &
2408	.AND. .NOT. land_surface .AND. .NOT. urban_surface ) &
2409	THEN
2410	message_string = 'output of "' // TRIM( var ) // '" requi' // &
2411	'res land_surface = .TRUE. or ' // &
2412	'urban_surface = .TRUE.'
2413	CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 )
2414	ENDIF
2415
2416	IF ( TRIM( var ) == 'ssws*' .AND. .NOT. passive_scalar ) THEN
2417	message_string = 'output of "' // TRIM( var ) // '" requi' // &
2418	'res passive_scalar = .TRUE.'
2419	CALL message( 'check_parameters', 'PA0361', 1, 2, 0, 6, 0 )
2420	ENDIF
2421	!
2422	!-- Activate calculation of 2m temperature if output is requested
2423	IF ( TRIM( var ) == 'theta_2m*' ) THEN
2424	do_output_at_2m = .TRUE.
2425	unit = 'K'
2426	ENDIF
2427
2428	IF ( TRIM( var ) == 'ghf*' ) unit = 'W/m2'
2429	IF ( TRIM( var ) == 'lwp*' ) unit = 'kg/m2'
2430	IF ( TRIM( var ) == 'ol*' ) unit = 'm'
2431	IF ( TRIM( var ) == 'qsws*' ) unit = 'kgm/kgs'
2432	IF ( TRIM( var ) == 'r_a*' ) unit = 's/m'
2433	IF ( TRIM( var ) == 'shf' ) unit = 'Km/s'
2434	IF ( TRIM( var ) == 'ssws' ) unit = 'kg/m2s'
2435	IF ( TRIM( var ) == 't*' ) unit = 'K'
2436	IF ( TRIM( var ) == 'tsurf*' ) unit = 'K'
2437	IF ( TRIM( var ) == 'us*' ) unit = 'm/s'
2438	IF ( TRIM( var ) == 'z0*' ) unit = 'm'
2439	IF ( TRIM( var ) == 'z0h*' ) unit = 'm'
2440	!
2441	!-- Output of surface latent and sensible heat flux will be in W/m2
2442	!-- in case of natural- and urban-type surfaces, even if
2443	!-- flux_output_mode is set to kinematic units.
2444	IF ( land_surface .OR. urban_surface ) THEN
2445	IF ( TRIM( var ) == 'shf*' ) unit = 'W/m2'
2446	IF ( TRIM( var ) == 'qsws*' ) unit = 'W/m2'
2447	ENDIF
2448
2449	CASE DEFAULT
2450	!
2451	!-- Check for other modules
2452	CALL module_interface_check_data_output( var, unit, i, j, ilen, k )
2453
2454	IF ( unit == 'illegal' ) THEN
2455	IF ( data_output_user(1) /= ' ' ) THEN
2456	message_string = 'illegal value for data_output or ' // &
2457	'data_output_user = "' // TRIM( data_output(i) ) // '"'
2458	CALL message( 'check_parameters', 'PA0114', 1, 2, 0, 6, 0 )
2459	ELSE
2460	message_string = 'illegal value for data_output = "' // &
2461	TRIM( data_output(i) ) // '"'
2462	CALL message( 'check_parameters', 'PA0115', 1, 2, 0, 6, 0 )
2463	ENDIF
2464	ENDIF
2465
2466	END SELECT
2467	!
2468	!-- Set the internal steering parameters appropriately
2469	IF ( k == 0 ) THEN
2470	do3d_no(j) = do3d_no(j) + 1
2471	do3d(j,do3d_no(j)) = data_output(i)
2472	do3d_unit(j,do3d_no(j)) = unit
2473	ELSE
2474	do2d_no(j) = do2d_no(j) + 1
2475	do2d(j,do2d_no(j)) = data_output(i)
2476	do2d_unit(j,do2d_no(j)) = unit
2477	IF ( data_output(i)(ilen-2:ilen) == '_xy' ) THEN
2478	data_output_xy(j) = .TRUE.
2479	ENDIF
2480	IF ( data_output(i)(ilen-2:ilen) == '_xz' ) THEN
2481	data_output_xz(j) = .TRUE.
2482	ENDIF
2483	IF ( data_output(i)(ilen-2:ilen) == '_yz' ) THEN
2484	data_output_yz(j) = .TRUE.
2485	ENDIF
2486	ENDIF
2487
2488	IF ( j == 1 ) THEN
2489	!
2490	!-- Check, if variable is already subject to averaging
2491	found = .FALSE.
2492	DO k = 1, doav_n
2493	IF ( TRIM( doav(k) ) == TRIM( var ) ) found = .TRUE.
2494	ENDDO
2495
2496	IF ( .NOT. found ) THEN
2497	doav_n = doav_n + 1
2498	doav(doav_n) = var
2499	ENDIF
2500	ENDIF
2501
2502	i = i + 1
2503	ENDDO
2504
2505	!
2506	!-- Averaged 2d or 3d output requires that an averaging interval has been set
2507	IF ( doav_n > 0 .AND. averaging_interval == 0.0_wp ) THEN
2508	WRITE( message_string, * ) 'output of averaged quantity "', &
2509	TRIM( doav(1) ), '_av" requires to set a ', &
2510	'non-zero averaging interval'
2511	CALL message( 'check_parameters', 'PA0323', 1, 2, 0, 6, 0 )
2512	ENDIF
2513
2514	!
2515	!-- Check sectional planes and store them in one shared array
2516	IF ( ANY( section_xy > nz + 1 ) ) THEN
2517	WRITE( message_string, * ) 'section_xy must be <= nz + 1 = ', nz + 1
2518	CALL message( 'check_parameters', 'PA0319', 1, 2, 0, 6, 0 )
2519	ENDIF
2520	IF ( ANY( section_xz > ny + 1 ) ) THEN
2521	WRITE( message_string, * ) 'section_xz must be <= ny + 1 = ', ny + 1
2522	CALL message( 'check_parameters', 'PA0320', 1, 2, 0, 6, 0 )
2523	ENDIF
2524	IF ( ANY( section_yz > nx + 1 ) ) THEN
2525	WRITE( message_string, * ) 'section_yz must be <= nx + 1 = ', nx + 1
2526	CALL message( 'check_parameters', 'PA0321', 1, 2, 0, 6, 0 )
2527	ENDIF
2528	section(:,1) = section_xy
2529	section(:,2) = section_xz
2530	section(:,3) = section_yz
2531
2532	IF ( ANY( data_output_xy ) .AND. .NOT. ANY( section(:,1) /= -9999 ) ) THEN
2533	WRITE( message_string, * ) 'section_xy not defined for requested ' // &
2534	'xy-cross section output.&At least one ' // &
2535	'cross section must be given.'
2536	CALL message( 'check_parameters', 'PA0681', 1, 2, 0, 6, 0 )
2537	ENDIF
2538	IF ( ANY( data_output_xz ) .AND. .NOT. ANY( section(:,2) /= -9999 ) ) THEN
2539	WRITE( message_string, * ) 'section_xz not defined for requested ' // &
2540	'xz-cross section output.&At least one ' // &
2541	'cross section must be given.'
2542	CALL message( 'check_parameters', 'PA0681', 1, 2, 0, 6, 0 )
2543	ENDIF
2544	IF ( ANY( data_output_yz ) .AND. .NOT. ANY( section(:,3) /= -9999 ) ) THEN
2545	WRITE( message_string, * ) 'section_yz not defined for requested ' // &
2546	'yz-cross section output.&At least one ' // &
2547	'cross section must be given.'
2548	CALL message( 'check_parameters', 'PA0681', 1, 2, 0, 6, 0 )
2549	ENDIF
2550	!
2551	!-- Upper plot limit for 3D arrays
2552	IF ( nz_do3d == -9999 ) nz_do3d = nzt + 1
2553
2554	!
2555	!-- Set output format string (used in header)
2556	SELECT CASE ( netcdf_data_format )
2557	CASE ( 1 )
2558	netcdf_data_format_string = 'netCDF classic'
2559	CASE ( 2 )
2560	netcdf_data_format_string = 'netCDF 64bit offset'
2561	CASE ( 3 )
2562	netcdf_data_format_string = 'netCDF4/HDF5'
2563	CASE ( 4 )
2564	netcdf_data_format_string = 'netCDF4/HDF5 classic'
2565	CASE ( 5 )
2566	netcdf_data_format_string = 'parallel netCDF4/HDF5'
2567	CASE ( 6 )
2568	netcdf_data_format_string = 'parallel netCDF4/HDF5 classic'
2569
2570	END SELECT
2571
2572	!
2573	!-- Check mask conditions
2574	DO mid = 1, max_masks
2575	IF ( data_output_masks(mid,1) /= ' ' .OR. &
2576	data_output_masks_user(mid,1) /= ' ' ) THEN
2577	masks = masks + 1
2578	ENDIF
2579	ENDDO
2580
2581	IF ( masks < 0 .OR. masks > max_masks ) THEN
2582	WRITE( message_string, * ) 'illegal value: masks must be >= 0 and ', &
2583	'<= ', max_masks, ' (=max_masks)'
2584	CALL message( 'check_parameters', 'PA0325', 1, 2, 0, 6, 0 )
2585	ENDIF
2586	IF ( masks > 0 ) THEN
2587	mask_scale(1) = mask_scale_x
2588	mask_scale(2) = mask_scale_y
2589	mask_scale(3) = mask_scale_z
2590	IF ( ANY( mask_scale <= 0.0_wp ) ) THEN
2591	WRITE( message_string, * ) &
2592	'illegal value: mask_scale_x, mask_scale_y and mask_scale_z', &
2593	'must be > 0.0'
2594	CALL message( 'check_parameters', 'PA0326', 1, 2, 0, 6, 0 )
2595	ENDIF
2596	!
2597	!-- Generate masks for masked data output
2598	!-- Parallel netcdf output is not tested so far for masked data, hence
2599	!-- netcdf_data_format is switched back to non-parallel output.
2600	netcdf_data_format_save = netcdf_data_format
2601	IF ( netcdf_data_format > 4 ) THEN
2602	IF ( netcdf_data_format == 5 ) netcdf_data_format = 3
2603	IF ( netcdf_data_format == 6 ) netcdf_data_format = 4
2604	message_string = 'netCDF file formats '// &
2605	'5 (parallel netCDF 4) and ' // &
2606	'6 (parallel netCDF 4 Classic model) '// &
2607	'& are currently not supported (not yet tested) ' //&
2608	'for masked data. &Using respective non-parallel' //&
2609	' output for masked data.'
2610	CALL message( 'check_parameters', 'PA0383', 0, 0, 0, 6, 0 )
2611	ENDIF
2612	CALL init_masks
2613	netcdf_data_format = netcdf_data_format_save
2614	ENDIF
2615
2616	!
2617	!-- Check the NetCDF data format
2618	IF ( netcdf_data_format > 2 ) THEN
2619	#if defined( __netcdf4 )
2620	CONTINUE
2621	#else
2622	message_string = 'netCDF: netCDF4 format requested but no ' // &
2623	'cpp-directive __netcdf4 given & switch ' // &
2624	'back to 64-bit offset format'
2625	CALL message( 'check_parameters', 'PA0171', 0, 1, 0, 6, 0 )
2626	netcdf_data_format = 2
2627	#endif
2628	ENDIF
2629	IF ( netcdf_data_format > 4 ) THEN
2630	#if defined( __netcdf4 ) && defined( __netcdf4_parallel )
2631	CONTINUE
2632	#else
2633	message_string = 'netCDF: netCDF4 parallel output requested but no ' // &
2634	'cpp-directive __netcdf4_parallel given & switch ' //&
2635	'back to netCDF4 non-parallel output'
2636	CALL message( 'check_parameters', 'PA0099', 0, 1, 0, 6, 0 )
2637	netcdf_data_format = netcdf_data_format - 2
2638	#endif
2639	ENDIF
2640
2641	!
2642	!-- Calculate fixed number of output time levels for parallel netcdf output.
2643	!-- The time dimension has to be defined as limited for parallel output,
2644	!-- because otherwise the I/O performance drops significantly.
2645	IF ( netcdf_data_format > 4 ) THEN
2646
2647	!
2648	!-- Check if any of the follwoing data output interval is 0.0s, which is
2649	!-- not allowed for parallel output.
2650	CALL check_dt_do( dt_do3d, 'dt_do3d' )
2651	CALL check_dt_do( dt_do2d_xy, 'dt_do2d_xy' )
2652	CALL check_dt_do( dt_do2d_xz, 'dt_do2d_xz' )
2653	CALL check_dt_do( dt_do2d_yz, 'dt_do2d_yz' )
2654	CALL check_dt_do( dt_data_output_av, 'dt_data_output_av' )
2655
2656	!-- Set needed time levels (ntdim) to
2657	!-- saved time levels + to be saved time levels.
2658	ntdim_3d(0) = do3d_time_count(0) + CEILING( &
2659	( end_time - MAX( &
2660	MERGE(skip_time_do3d, skip_time_do3d + spinup_time, &
2661	data_output_during_spinup ), &
2662	simulated_time_at_begin ) &
2663	) / dt_do3d )
2664	IF ( do3d_at_begin ) ntdim_3d(0) = ntdim_3d(0) + 1
2665
2666	ntdim_3d(1) = do3d_time_count(1) + CEILING( &
2667	( end_time - MAX( &
2668	MERGE( skip_time_data_output_av, skip_time_data_output_av &
2669	+ spinup_time, data_output_during_spinup ), &
2670	simulated_time_at_begin ) &
2671	) / dt_data_output_av )
2672
2673	ntdim_2d_xy(0) = do2d_xy_time_count(0) + CEILING( &
2674	( end_time - MAX( &
2675	MERGE(skip_time_do2d_xy, skip_time_do2d_xy + spinup_time, &
2676	data_output_during_spinup ), &
2677	simulated_time_at_begin ) &
2678	) / dt_do2d_xy )
2679
2680	ntdim_2d_xz(0) = do2d_xz_time_count(0) + CEILING( &
2681	( end_time - MAX( &
2682	MERGE(skip_time_do2d_xz, skip_time_do2d_xz + spinup_time, &
2683	data_output_during_spinup ), &
2684	simulated_time_at_begin ) &
2685	) / dt_do2d_xz )
2686
2687	ntdim_2d_yz(0) = do2d_yz_time_count(0) + CEILING( &
2688	( end_time - MAX( &
2689	MERGE(skip_time_do2d_yz, skip_time_do2d_yz + spinup_time, &
2690	data_output_during_spinup ), &
2691	simulated_time_at_begin ) &
2692	) / dt_do2d_yz )
2693
2694	IF ( do2d_at_begin ) THEN
2695	ntdim_2d_xy(0) = ntdim_2d_xy(0) + 1
2696	ntdim_2d_xz(0) = ntdim_2d_xz(0) + 1
2697	ntdim_2d_yz(0) = ntdim_2d_yz(0) + 1
2698	ENDIF
2699	!
2700	!-- Please note, for averaged 2D data skip_time_data_output_av is the relavant
2701	!-- output control parameter.
2702	ntdim_2d_xy(1) = do2d_xy_time_count(1) + CEILING( &
2703	( end_time - MAX( MERGE( skip_time_data_output_av, &
2704	skip_time_data_output_av + spinup_time, &
2705	data_output_during_spinup ), &
2706	simulated_time_at_begin ) &
2707	) / dt_data_output_av )
2708
2709	ntdim_2d_xz(1) = do2d_xz_time_count(1) + CEILING( &
2710	( end_time - MAX( MERGE( skip_time_data_output_av, &
2711	skip_time_data_output_av + spinup_time, &
2712	data_output_during_spinup ), &
2713	simulated_time_at_begin ) &
2714	) / dt_data_output_av )
2715
2716	ntdim_2d_yz(1) = do2d_yz_time_count(1) + CEILING( &
2717	( end_time - MAX( MERGE( skip_time_data_output_av, &
2718	skip_time_data_output_av + spinup_time, &
2719	data_output_during_spinup ), &
2720	simulated_time_at_begin ) &
2721	) / dt_data_output_av )
2722
2723	ENDIF
2724
2725	!
2726	!-- Check, whether a constant diffusion coefficient shall be used
2727	IF ( km_constant /= -1.0_wp ) THEN
2728	IF ( km_constant < 0.0_wp ) THEN
2729	WRITE( message_string, * ) 'km_constant = ', km_constant, ' < 0.0'
2730	CALL message( 'check_parameters', 'PA0121', 1, 2, 0, 6, 0 )
2731	ELSE
2732	IF ( prandtl_number < 0.0_wp ) THEN
2733	WRITE( message_string, * ) 'prandtl_number = ', prandtl_number, &
2734	' < 0.0'
2735	CALL message( 'check_parameters', 'PA0122', 1, 2, 0, 6, 0 )
2736	ENDIF
2737	constant_diffusion = .TRUE.
2738
2739	IF ( constant_flux_layer ) THEN
2740	message_string = 'constant_flux_layer is not allowed with fixed ' &
2741	// 'value of km'
2742	CALL message( 'check_parameters', 'PA0123', 1, 2, 0, 6, 0 )
2743	ENDIF
2744	ENDIF
2745	ENDIF
2746
2747	!
2748	!-- In case of non-cyclic lateral boundaries and a damping layer for the
2749	!-- potential temperature, check the width of the damping layer
2750	IF ( bc_lr /= 'cyclic' ) THEN
2751	IF ( pt_damping_width < 0.0_wp .OR. &
2752	pt_damping_width > REAL( (nx+1) * dx ) ) THEN
2753	message_string = 'pt_damping_width out of range'
2754	CALL message( 'check_parameters', 'PA0124', 1, 2, 0, 6, 0 )
2755	ENDIF
2756	ENDIF
2757
2758	IF ( bc_ns /= 'cyclic' ) THEN
2759	IF ( pt_damping_width < 0.0_wp .OR. &
2760	pt_damping_width > REAL( (ny+1) * dy ) ) THEN
2761	message_string = 'pt_damping_width out of range'
2762	CALL message( 'check_parameters', 'PA0124', 1, 2, 0, 6, 0 )
2763	ENDIF
2764	ENDIF
2765
2766	!
2767	!-- Check value range for zeta = z/L
2768	IF ( zeta_min >= zeta_max ) THEN
2769	WRITE( message_string, * ) 'zeta_min = ', zeta_min, ' must be less ', &
2770	'than zeta_max = ', zeta_max
2771	CALL message( 'check_parameters', 'PA0125', 1, 2, 0, 6, 0 )
2772	ENDIF
2773
2774	!
2775	!-- Check random generator
2776	IF ( (random_generator /= 'system-specific' .AND. &
2777	random_generator /= 'random-parallel' ) .AND. &
2778	random_generator /= 'numerical-recipes' ) THEN
2779	message_string = 'unknown random generator: random_generator = "' // &
2780	TRIM( random_generator ) // '"'
2781	CALL message( 'check_parameters', 'PA0135', 1, 2, 0, 6, 0 )
2782	ENDIF
2783
2784	!
2785	!-- Determine upper and lower hight level indices for random perturbations
2786	IF ( disturbance_level_b == -9999999.9_wp ) THEN
2787	IF ( ocean_mode ) THEN
2788	disturbance_level_b = zu((nzt*2)/3)
2789	disturbance_level_ind_b = ( nzt * 2 ) / 3
2790	ELSE
2791	disturbance_level_b = zu(nzb+3)
2792	disturbance_level_ind_b = nzb + 3
2793	ENDIF
2794	ELSEIF ( disturbance_level_b < zu(3) ) THEN
2795	WRITE( message_string, * ) 'disturbance_level_b = ', &
2796	disturbance_level_b, ' must be >= ', zu(3), '(zu(3))'
2797	CALL message( 'check_parameters', 'PA0126', 1, 2, 0, 6, 0 )
2798	ELSEIF ( disturbance_level_b > zu(nzt-2) ) THEN
2799	WRITE( message_string, * ) 'disturbance_level_b = ', &
2800	disturbance_level_b, ' must be <= ', zu(nzt-2), '(zu(nzt-2))'
2801	CALL message( 'check_parameters', 'PA0127', 1, 2, 0, 6, 0 )
2802	ELSE
2803	DO k = 3, nzt-2
2804	IF ( disturbance_level_b <= zu(k) ) THEN
2805	disturbance_level_ind_b = k
2806	EXIT
2807	ENDIF
2808	ENDDO
2809	ENDIF
2810
2811	IF ( disturbance_level_t == -9999999.9_wp ) THEN
2812	IF ( ocean_mode ) THEN
2813	disturbance_level_t = zu(nzt-3)
2814	disturbance_level_ind_t = nzt - 3
2815	ELSE
2816	disturbance_level_t = zu(nzt/3)
2817	disturbance_level_ind_t = nzt / 3
2818	ENDIF
2819	ELSEIF ( disturbance_level_t > zu(nzt-2) ) THEN
2820	WRITE( message_string, * ) 'disturbance_level_t = ', &
2821	disturbance_level_t, ' must be <= ', zu(nzt-2), '(zu(nzt-2))'
2822	CALL message( 'check_parameters', 'PA0128', 1, 2, 0, 6, 0 )
2823	ELSEIF ( disturbance_level_t < disturbance_level_b ) THEN
2824	WRITE( message_string, * ) 'disturbance_level_t = ', &
2825	disturbance_level_t, ' must be >= disturbance_level_b = ', &
2826	disturbance_level_b
2827	CALL message( 'check_parameters', 'PA0129', 1, 2, 0, 6, 0 )
2828	ELSE
2829	DO k = 3, nzt-2
2830	IF ( disturbance_level_t <= zu(k) ) THEN
2831	disturbance_level_ind_t = k
2832	EXIT
2833	ENDIF
2834	ENDDO
2835	ENDIF
2836
2837	!
2838	!-- Check again whether the levels determined this way are ok.
2839	!-- Error may occur at automatic determination and too few grid points in
2840	!-- z-direction.
2841	IF ( disturbance_level_ind_t < disturbance_level_ind_b ) THEN
2842	WRITE( message_string, * ) 'disturbance_level_ind_t = ', &
2843	disturbance_level_ind_t, ' must be >= ', &
2844	'disturbance_level_ind_b = ', disturbance_level_ind_b
2845	CALL message( 'check_parameters', 'PA0130', 1, 2, 0, 6, 0 )
2846	ENDIF
2847
2848	!
2849	!-- Determine the horizontal index range for random perturbations.
2850	!-- In case of non-cyclic horizontal boundaries, no perturbations are imposed
2851	!-- near the inflow and the perturbation area is further limited to ...(1)
2852	!-- after the initial phase of the flow.
2853
2854	IF ( bc_lr /= 'cyclic' ) THEN
2855	IF ( inflow_disturbance_begin == -1 ) THEN
2856	inflow_disturbance_begin = MIN( 10, nx/2 )
2857	ENDIF
2858	IF ( inflow_disturbance_begin < 0 .OR. inflow_disturbance_begin > nx )&
2859	THEN
2860	message_string = 'inflow_disturbance_begin out of range'
2861	CALL message( 'check_parameters', 'PA0131', 1, 2, 0, 6, 0 )
2862	ENDIF
2863	IF ( inflow_disturbance_end == -1 ) THEN
2864	inflow_disturbance_end = MIN( 100, 3*nx/4 )
2865	ENDIF
2866	IF ( inflow_disturbance_end < 0 .OR. inflow_disturbance_end > nx ) &
2867	THEN
2868	message_string = 'inflow_disturbance_end out of range'
2869	CALL message( 'check_parameters', 'PA0132', 1, 2, 0, 6, 0 )
2870	ENDIF
2871	ELSEIF ( bc_ns /= 'cyclic' ) THEN
2872	IF ( inflow_disturbance_begin == -1 ) THEN
2873	inflow_disturbance_begin = MIN( 10, ny/2 )
2874	ENDIF
2875	IF ( inflow_disturbance_begin < 0 .OR. inflow_disturbance_begin > ny )&
2876	THEN
2877	message_string = 'inflow_disturbance_begin out of range'
2878	CALL message( 'check_parameters', 'PA0131', 1, 2, 0, 6, 0 )
2879	ENDIF
2880	IF ( inflow_disturbance_end == -1 ) THEN
2881	inflow_disturbance_end = MIN( 100, 3*ny/4 )
2882	ENDIF
2883	IF ( inflow_disturbance_end < 0 .OR. inflow_disturbance_end > ny ) &
2884	THEN
2885	message_string = 'inflow_disturbance_end out of range'
2886	CALL message( 'check_parameters', 'PA0132', 1, 2, 0, 6, 0 )
2887	ENDIF
2888	ENDIF
2889
2890	IF ( random_generator == 'random-parallel' ) THEN
2891	dist_nxl = nxl; dist_nxr = nxr
2892	dist_nys = nys; dist_nyn = nyn
2893	IF ( bc_lr == 'radiation/dirichlet' ) THEN
2894	dist_nxr = MIN( nx - inflow_disturbance_begin, nxr )
2895	dist_nxl(1) = MAX( nx - inflow_disturbance_end, nxl )
2896	ELSEIF ( bc_lr == 'dirichlet/radiation' ) THEN
2897	dist_nxl = MAX( inflow_disturbance_begin, nxl )
2898	dist_nxr(1) = MIN( inflow_disturbance_end, nxr )
2899	ELSEIF ( bc_lr == 'nested' .OR. bc_lr == 'nesting_offline' ) THEN
2900	dist_nxl = MAX( inflow_disturbance_begin, nxl )
2901	dist_nxr = MIN( nx - inflow_disturbance_begin, nxr )
2902	ENDIF
2903	IF ( bc_ns == 'dirichlet/radiation' ) THEN
2904	dist_nyn = MIN( ny - inflow_disturbance_begin, nyn )
2905	dist_nys(1) = MAX( ny - inflow_disturbance_end, nys )
2906	ELSEIF ( bc_ns == 'radiation/dirichlet' ) THEN
2907	dist_nys = MAX( inflow_disturbance_begin, nys )
2908	dist_nyn(1) = MIN( inflow_disturbance_end, nyn )
2909	ELSEIF ( bc_ns == 'nested' .OR. bc_ns == 'nesting_offline' ) THEN
2910	dist_nys = MAX( inflow_disturbance_begin, nys )
2911	dist_nyn = MIN( ny - inflow_disturbance_begin, nyn )
2912	ENDIF
2913	ELSE
2914	dist_nxl = 0; dist_nxr = nx
2915	dist_nys = 0; dist_nyn = ny
2916	IF ( bc_lr == 'radiation/dirichlet' ) THEN
2917	dist_nxr = nx - inflow_disturbance_begin
2918	dist_nxl(1) = nx - inflow_disturbance_end
2919	ELSEIF ( bc_lr == 'dirichlet/radiation' ) THEN
2920	dist_nxl = inflow_disturbance_begin
2921	dist_nxr(1) = inflow_disturbance_end
2922	ELSEIF ( bc_lr == 'nested' .OR. bc_lr == 'nesting_offline' ) THEN
2923	dist_nxr = nx - inflow_disturbance_begin
2924	dist_nxl = inflow_disturbance_begin
2925	ENDIF
2926	IF ( bc_ns == 'dirichlet/radiation' ) THEN
2927	dist_nyn = ny - inflow_disturbance_begin
2928	dist_nys(1) = ny - inflow_disturbance_end
2929	ELSEIF ( bc_ns == 'radiation/dirichlet' ) THEN
2930	dist_nys = inflow_disturbance_begin
2931	dist_nyn(1) = inflow_disturbance_end
2932	ELSEIF ( bc_ns == 'nested' .OR. bc_ns == 'nesting_offline' ) THEN
2933	dist_nyn = ny - inflow_disturbance_begin
2934	dist_nys = inflow_disturbance_begin
2935	ENDIF
2936	ENDIF
2937
2938	!
2939	!-- A turbulent inflow requires Dirichlet conditions at the respective inflow
2940	!-- boundary (so far, a turbulent inflow is realized from the left side only)
2941	IF ( turbulent_inflow .AND. bc_lr /= 'dirichlet/radiation' ) THEN
2942	message_string = 'turbulent_inflow = .T. requires a Dirichlet ' // &
2943	'condition at the inflow boundary'
2944	CALL message( 'check_parameters', 'PA0133', 1, 2, 0, 6, 0 )
2945	ENDIF
2946
2947	!
2948	!-- Turbulent inflow requires that 3d arrays have been cyclically filled with
2949	!-- data from prerun in the first main run
2950	IF ( turbulent_inflow .AND. initializing_actions /= 'cyclic_fill' .AND. &
2951	initializing_actions /= 'read_restart_data' ) THEN
2952	message_string = 'turbulent_inflow = .T. requires ' // &
2953	'initializing_actions = ''cyclic_fill'' or ' // &
2954	'initializing_actions = ''read_restart_data'' '
2955	CALL message( 'check_parameters', 'PA0055', 1, 2, 0, 6, 0 )
2956	ENDIF
2957
2958	!
2959	!-- In case of turbulent inflow
2960	IF ( turbulent_inflow ) THEN
2961
2962	!
2963	!-- Calculate the index of the recycling plane
2964	IF ( recycling_width <= dx .OR. recycling_width >= nx * dx ) THEN
2965	WRITE( message_string, * ) 'illegal value for recycling_width: ', &
2966	recycling_width
2967	CALL message( 'check_parameters', 'PA0134', 1, 2, 0, 6, 0 )
2968	ENDIF
2969	!
2970	!-- Calculate the index
2971	recycling_plane = recycling_width / dx
2972	!
2973	!-- Because the y-shift is done with a distance of INT( npey / 2 ) no shift
2974	!-- is possible if there is only one PE in y direction.
2975	IF ( recycling_yshift .AND. pdims(2) < 2 ) THEN
2976	WRITE( message_string, * ) 'recycling_yshift = .T. requires more', &
2977	' than one processor in y direction'
2978	CALL message( 'check_parameters', 'PA0421', 1, 2, 0, 6, 0 )
2979	ENDIF
2980
2981	!
2982	!-- Check for correct input of recycling method for thermodynamic quantities
2983	IF ( TRIM( recycling_method_for_thermodynamic_quantities ) /= 'turbulent_fluctuation' .AND. &
2984	TRIM( recycling_method_for_thermodynamic_quantities ) /= 'absolute_value' ) THEN
2985	WRITE( message_string, * ) 'unknown recycling method for thermodynamic quantities: ', &
2986	TRIM( recycling_method_for_thermodynamic_quantities )
2987	CALL message( 'check_parameters', 'PA0184', 1, 2, 0, 6, 0 )
2988	ENDIF
2989
2990	ENDIF
2991
2992
2993	IF ( turbulent_outflow ) THEN
2994	!
2995	!-- Turbulent outflow requires Dirichlet conditions at the respective inflow
2996	!-- boundary (so far, a turbulent outflow is realized at the right side only)
2997	IF ( bc_lr /= 'dirichlet/radiation' ) THEN
2998	message_string = 'turbulent_outflow = .T. requires ' // &
2999	'bc_lr = "dirichlet/radiation"'
3000	CALL message( 'check_parameters', 'PA0038', 1, 2, 0, 6, 0 )
3001	ENDIF
3002	!
3003	!-- The ouflow-source plane must lay inside the model domain
3004	IF ( outflow_source_plane < dx .OR. &
3005	outflow_source_plane > nx * dx ) THEN
3006	WRITE( message_string, * ) 'illegal value for outflow_source'// &
3007	'_plane: ', outflow_source_plane
3008	CALL message( 'check_parameters', 'PA0145', 1, 2, 0, 6, 0 )
3009	ENDIF
3010	ENDIF
3011
3012	!
3013	!-- Determine damping level index for 1D model
3014	IF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 ) THEN
3015	IF ( damp_level_1d == -1.0_wp ) THEN
3016	damp_level_1d = zu(nzt+1)
3017	damp_level_ind_1d = nzt + 1
3018	ELSEIF ( damp_level_1d < 0.0_wp .OR. damp_level_1d > zu(nzt+1) ) THEN
3019	WRITE( message_string, * ) 'damp_level_1d = ', damp_level_1d, &
3020	' must be >= 0.0 and <= ', zu(nzt+1), '(zu(nzt+1))'
3021	CALL message( 'check_parameters', 'PA0136', 1, 2, 0, 6, 0 )
3022	ELSE
3023	DO k = 1, nzt+1
3024	IF ( damp_level_1d <= zu(k) ) THEN
3025	damp_level_ind_1d = k
3026	EXIT
3027	ENDIF
3028	ENDDO
3029	ENDIF
3030	ENDIF
3031
3032	!
3033	!-- Check some other 1d-model parameters
3034	IF ( TRIM( mixing_length_1d ) /= 'as_in_3d_model' .AND. &
3035	TRIM( mixing_length_1d ) /= 'blackadar' ) THEN
3036	message_string = 'mixing_length_1d = "' // TRIM( mixing_length_1d ) // &
3037	'" is unknown'
3038	CALL message( 'check_parameters', 'PA0137', 1, 2, 0, 6, 0 )
3039	ENDIF
3040	IF ( TRIM( dissipation_1d ) /= 'as_in_3d_model' .AND. &
3041	TRIM( dissipation_1d ) /= 'detering' .AND. &
3042	TRIM( dissipation_1d ) /= 'prognostic' ) THEN
3043	message_string = 'dissipation_1d = "' // TRIM( dissipation_1d ) // &
3044	'" is unknown'
3045	CALL message( 'check_parameters', 'PA0138', 1, 2, 0, 6, 0 )
3046	ENDIF
3047	IF ( TRIM( mixing_length_1d ) /= 'as_in_3d_model' .AND. &
3048	TRIM( dissipation_1d ) == 'as_in_3d_model' ) THEN
3049	message_string = 'dissipation_1d = "' // TRIM( dissipation_1d ) // &
3050	'" requires mixing_length_1d = "as_in_3d_model"'
3051	CALL message( 'check_parameters', 'PA0485', 1, 2, 0, 6, 0 )
3052	ENDIF
3053
3054	!
3055	!-- Set time for the next user defined restart (time_restart is the
3056	!-- internal parameter for steering restart events)
3057	IF ( restart_time /= 9999999.9_wp ) THEN
3058	IF ( restart_time > time_since_reference_point ) THEN
3059	time_restart = restart_time
3060	ENDIF
3061	ELSE
3062	!
3063	!-- In case of a restart run, set internal parameter to default (no restart)
3064	!-- if the NAMELIST-parameter restart_time is omitted
3065	time_restart = 9999999.9_wp
3066	ENDIF
3067
3068	!
3069	!-- Check pressure gradient conditions
3070	IF ( dp_external .AND. conserve_volume_flow ) THEN
3071	WRITE( message_string, * ) 'Both dp_external and conserve_volume_flo', &
3072	'w are .TRUE. but one of them must be .FALSE.'
3073	CALL message( 'check_parameters', 'PA0150', 1, 2, 0, 6, 0 )
3074	ENDIF
3075	IF ( dp_external ) THEN
3076	IF ( dp_level_b < zu(nzb) .OR. dp_level_b > zu(nzt) ) THEN
3077	WRITE( message_string, * ) 'dp_level_b = ', dp_level_b, ' is out ', &
3078	' of range [zu(nzb), zu(nzt)]'
3079	CALL message( 'check_parameters', 'PA0151', 1, 2, 0, 6, 0 )
3080	ENDIF
3081	IF ( .NOT. ANY( dpdxy /= 0.0_wp ) ) THEN
3082	WRITE( message_string, * ) 'dp_external is .TRUE. but dpdxy is ze', &
3083	'ro, i.e. the external pressure gradient will not be applied'
3084	CALL message( 'check_parameters', 'PA0152', 0, 1, 0, 6, 0 )
3085	ENDIF
3086	ENDIF
3087	IF ( ANY( dpdxy /= 0.0_wp ) .AND. .NOT. dp_external ) THEN
3088	WRITE( message_string, * ) 'dpdxy is nonzero but dp_external is ', &
3089	'.FALSE., i.e. the external pressure gradient & will not be applied'
3090	CALL message( 'check_parameters', 'PA0153', 0, 1, 0, 6, 0 )
3091	ENDIF
3092	IF ( conserve_volume_flow ) THEN
3093	IF ( TRIM( conserve_volume_flow_mode ) == 'default' ) THEN
3094
3095	conserve_volume_flow_mode = 'initial_profiles'
3096
3097	ELSEIF ( TRIM( conserve_volume_flow_mode ) /= 'initial_profiles' .AND. &
3098	TRIM( conserve_volume_flow_mode ) /= 'bulk_velocity' ) THEN
3099	WRITE( message_string, * ) 'unknown conserve_volume_flow_mode: ', &
3100	conserve_volume_flow_mode
3101	CALL message( 'check_parameters', 'PA0154', 1, 2, 0, 6, 0 )
3102	ENDIF
3103	IF ( (bc_lr /= 'cyclic' .OR. bc_ns /= 'cyclic') .AND. &
3104	TRIM( conserve_volume_flow_mode ) == 'bulk_velocity' ) THEN
3105	WRITE( message_string, * ) 'non-cyclic boundary conditions ', &
3106	'require conserve_volume_flow_mode = ''initial_profiles'''
3107	CALL message( 'check_parameters', 'PA0155', 1, 2, 0, 6, 0 )
3108	ENDIF
3109	ENDIF
3110	IF ( ( u_bulk /= 0.0_wp .OR. v_bulk /= 0.0_wp ) .AND. &
3111	( .NOT. conserve_volume_flow .OR. &
3112	TRIM( conserve_volume_flow_mode ) /= 'bulk_velocity' ) ) THEN
3113	WRITE( message_string, * ) 'nonzero bulk velocity requires ', &
3114	'conserve_volume_flow = .T. and ', &
3115	'conserve_volume_flow_mode = ''bulk_velocity'''
3116	CALL message( 'check_parameters', 'PA0157', 1, 2, 0, 6, 0 )
3117	ENDIF
3118
3119	!
3120	!-- Prevent empty time records in volume, cross-section and masked data in case
3121	!-- of non-parallel netcdf-output in restart runs
3122	IF ( netcdf_data_format < 5 ) THEN
3123	IF ( TRIM( initializing_actions ) == 'read_restart_data' ) THEN
3124	do3d_time_count = 0
3125	do2d_xy_time_count = 0
3126	do2d_xz_time_count = 0
3127	do2d_yz_time_count = 0
3128	domask_time_count = 0
3129	ENDIF
3130	ENDIF
3131
3132
3133	!
3134	!-- Check roughness length, which has to be smaller than dz/2
3135	IF ( ( constant_flux_layer .OR. &
3136	INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 ) &
3137	.AND. roughness_length >= 0.5 * dz(1) ) THEN
3138	message_string = 'roughness_length must be smaller than dz/2'
3139	CALL message( 'check_parameters', 'PA0424', 1, 2, 0, 6, 0 )
3140	ENDIF
3141
3142	!
3143	!-- Vertical nesting: check fine and coarse grid compatibility for data exchange
3144	IF ( vnested ) CALL vnest_check_parameters
3145
3146	!
3147	!-- Check if topography is read from file in case of complex terrain simulations
3148	IF ( complex_terrain .AND. TRIM( topography ) /= 'read_from_file' ) THEN
3149	message_string = 'complex_terrain requires topography' // &
3150	' = ''read_from_file'''
3151	CALL message( 'check_parameters', 'PA0295', 1, 2, 0, 6, 0 )
3152	ENDIF
3153
3154	!
3155	!-- Check if vertical grid stretching is switched off in case of complex
3156	!-- terrain simulations
3157	IF ( complex_terrain .AND. &
3158	ANY( dz_stretch_level_start /= -9999999.9_wp ) ) THEN
3159	message_string = 'Vertical grid stretching is not allowed for ' // &
3160	'complex_terrain = .T.'
3161	CALL message( 'check_parameters', 'PA0473', 1, 2, 0, 6, 0 )
3162	ENDIF
3163
3164	CALL location_message( 'checking parameters', 'finished' )
3165
3166	CONTAINS
3167
3168	!------------------------------------------------------------------------------!
3169	! Description:
3170	! ------------
3171	!> Check the length of data output intervals. In case of parallel NetCDF output
3172	!> the time levels of the output files need to be fixed. Therefore setting the
3173	!> output interval to 0.0s (usually used to output each timestep) is not
3174	!> possible as long as a non-fixed timestep is used.
3175	!------------------------------------------------------------------------------!
3176
3177	SUBROUTINE check_dt_do( dt_do, dt_do_name )
3178
3179	IMPLICIT NONE
3180
3181	CHARACTER (LEN=*), INTENT (IN) :: dt_do_name !< parin variable name
3182
3183	REAL(wp), INTENT (INOUT) :: dt_do !< data output interval
3184
3185	IF ( dt_do == 0.0_wp ) THEN
3186	IF ( dt_fixed ) THEN
3187	WRITE( message_string, '(A,F9.4,A)' ) 'Output at every ' // &
3188	'timestep is wanted (' // dt_do_name // ' = 0.0).&'// &
3189	'The output interval is set to the fixed timestep dt '// &
3190	'= ', dt, 's.'
3191	CALL message( 'check_parameters', 'PA0060', 0, 0, 0, 6, 0 )
3192	dt_do = dt
3193	ELSE
3194	message_string = dt_do_name // ' = 0.0 while using a ' // &
3195	'variable timestep and parallel netCDF4 ' // &
3196	'is not allowed.'
3197	CALL message( 'check_parameters', 'PA0081', 1, 2, 0, 6, 0 )
3198	ENDIF
3199	ENDIF
3200
3201	END SUBROUTINE check_dt_do
3202
3203
3204
3205	!------------------------------------------------------------------------------!
3206	! Description:
3207	! ------------
3208	!> Set the bottom and top boundary conditions for humidity and scalars.
3209	!------------------------------------------------------------------------------!
3210
3211	SUBROUTINE set_bc_scalars( sq, bc_b, bc_t, ibc_b, ibc_t, err_nr_b, err_nr_t )
3212
3213
3214	IMPLICIT NONE
3215
3216	CHARACTER (LEN=1) :: sq !< name of scalar quantity
3217	CHARACTER (LEN=*) :: bc_b !< bottom boundary condition
3218	CHARACTER (LEN=*) :: bc_t !< top boundary condition
3219	CHARACTER (LEN=*) :: err_nr_b !< error number if bottom bc is unknown
3220	CHARACTER (LEN=*) :: err_nr_t !< error number if top bc is unknown
3221
3222	INTEGER(iwp) :: ibc_b !< index for bottom boundary condition
3223	INTEGER(iwp) :: ibc_t !< index for top boundary condition
3224
3225	!
3226	!-- Set Integer flags and check for possilbe errorneous settings for bottom
3227	!-- boundary condition
3228	IF ( bc_b == 'dirichlet' ) THEN
3229	ibc_b = 0
3230	ELSEIF ( bc_b == 'neumann' ) THEN
3231	ibc_b = 1
3232	ELSE
3233	message_string = 'unknown boundary condition: bc_' // TRIM( sq ) // &
3234	'_b ="' // TRIM( bc_b ) // '"'
3235	CALL message( 'check_parameters', err_nr_b, 1, 2, 0, 6, 0 )
3236	ENDIF
3237	!
3238	!-- Set Integer flags and check for possilbe errorneous settings for top
3239	!-- boundary condition
3240	IF ( bc_t == 'dirichlet' ) THEN
3241	ibc_t = 0
3242	ELSEIF ( bc_t == 'neumann' ) THEN
3243	ibc_t = 1
3244	ELSEIF ( bc_t == 'initial_gradient' ) THEN
3245	ibc_t = 2
3246	ELSEIF ( bc_t == 'nested' .OR. bc_t == 'nesting_offline' ) THEN
3247	ibc_t = 3
3248	ELSE
3249	message_string = 'unknown boundary condition: bc_' // TRIM( sq ) // &
3250	'_t ="' // TRIM( bc_t ) // '"'
3251	CALL message( 'check_parameters', err_nr_t, 1, 2, 0, 6, 0 )
3252	ENDIF
3253
3254
3255	END SUBROUTINE set_bc_scalars
3256
3257
3258
3259	!------------------------------------------------------------------------------!
3260	! Description:
3261	! ------------
3262	!> Check for consistent settings of bottom boundary conditions for humidity
3263	!> and scalars.
3264	!------------------------------------------------------------------------------!
3265
3266	SUBROUTINE check_bc_scalars( sq, bc_b, ibc_b, &
3267	err_nr_1, err_nr_2, &
3268	constant_flux, surface_initial_change )
3269
3270
3271	IMPLICIT NONE
3272
3273	CHARACTER (LEN=1) :: sq !< name of scalar quantity
3274	CHARACTER (LEN=*) :: bc_b !< bottom boundary condition
3275	CHARACTER (LEN=*) :: err_nr_1 !< error number of first error
3276	CHARACTER (LEN=*) :: err_nr_2 !< error number of second error
3277
3278	INTEGER(iwp) :: ibc_b !< index of bottom boundary condition
3279
3280	LOGICAL :: constant_flux !< flag for constant-flux layer
3281
3282	REAL(wp) :: surface_initial_change !< value of initial change at the surface
3283
3284	!
3285	!-- A given surface value implies Dirichlet boundary condition for
3286	!-- the respective quantity. In this case specification of a constant flux is
3287	!-- forbidden. However, an exception is made for large-scale forcing as well
3288	!-- as land-surface model.
3289	IF ( .NOT. land_surface .AND. .NOT. large_scale_forcing ) THEN
3290	IF ( ibc_b == 0 .AND. constant_flux ) THEN
3291	message_string = 'boundary condition: bc_' // TRIM( sq ) // &
3292	'_b ' // '= "' // TRIM( bc_b ) // &
3293	'" is not allowed with prescribed surface flux'
3294	CALL message( 'check_parameters', err_nr_1, 1, 2, 0, 6, 0 )
3295	ENDIF
3296	ENDIF
3297	IF ( constant_flux .AND. surface_initial_change /= 0.0_wp ) THEN
3298	WRITE( message_string, * ) 'a prescribed surface flux is not allo', &
3299	'wed with ', sq, '_surface_initial_change (/=0) = ', &
3300	surface_initial_change
3301	CALL message( 'check_parameters', err_nr_2, 1, 2, 0, 6, 0 )
3302	ENDIF
3303
3304
3305	END SUBROUTINE check_bc_scalars
3306
3307
3308
3309	END SUBROUTINE check_parameters

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

Context Navigation

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

Download in other formats: