source: palm/tags/release-4.0/SOURCE/header.f90 @ 3800

Last change on this file since 3800 was 1497, checked in by maronga, 9 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 85.1 KB
Line 
1 SUBROUTINE header
2
3!--------------------------------------------------------------------------------!
4! This file is part of PALM.
5!
6! PALM is free software: you can redistribute it and/or modify it under the terms
7! of the GNU General Public License as published by the Free Software Foundation,
8! either version 3 of the License, or (at your option) any later version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 1997-2014 Leibniz Universitaet Hannover
18!--------------------------------------------------------------------------------!
19!
20! Current revisions:
21! -----------------
22!
23!
24! Former revisions:
25! -----------------
26! $Id: header.f90 1497 2014-12-02 17:28:07Z forkel $
27!
28! 1496 2014-12-02 17:25:50Z maronga
29! Renamed: "radiation -> "cloud_top_radiation"
30!
31! 1484 2014-10-21 10:53:05Z kanani
32! Changes due to new module structure of the plant canopy model:
33!   module plant_canopy_model_mod and output for new canopy model parameters
34!   (alpha_lad, beta_lad, lai_beta,...) added,
35!   drag_coefficient, leaf_surface_concentration and scalar_exchange_coefficient
36!   renamed to canopy_drag_coeff, leaf_surface_conc and leaf_scalar_exch_coeff,
37!   learde renamed leaf_area_density.
38! Bugfix: DO-WHILE-loop for lad header information additionally restricted
39! by maximum number of gradient levels (currently 10)
40!
41! 1482 2014-10-18 12:34:45Z raasch
42! information about calculated or predefined virtual processor topology adjusted
43!
44! 1468 2014-09-24 14:06:57Z maronga
45! Adapted for use on up to 6-digit processor cores
46!
47! 1429 2014-07-15 12:53:45Z knoop
48! header exended to provide ensemble_member_nr if specified
49!
50! 1376 2014-04-26 11:21:22Z boeske
51! Correction of typos
52!
53! 1365 2014-04-22 15:03:56Z boeske
54! New section 'Large scale forcing and nudging':
55! output of large scale forcing and nudging information,
56! new section for initial profiles created
57!
58! 1359 2014-04-11 17:15:14Z hoffmann
59! dt_sort_particles removed
60!
61! 1353 2014-04-08 15:21:23Z heinze
62! REAL constants provided with KIND-attribute
63!
64! 1327 2014-03-21 11:00:16Z raasch
65! parts concerning iso2d and avs output removed,
66! -netcdf output queries
67!
68! 1324 2014-03-21 09:13:16Z suehring
69! Bugfix: module spectrum added
70!
71! 1322 2014-03-20 16:38:49Z raasch
72! REAL functions provided with KIND-attribute,
73! some REAL constants defined as wp-kind
74!
75! 1320 2014-03-20 08:40:49Z raasch
76! ONLY-attribute added to USE-statements,
77! kind-parameters added to all INTEGER and REAL declaration statements,
78! kinds are defined in new module kinds,
79! revision history before 2012 removed,
80! comment fields (!:) to be used for variable explanations added to
81! all variable declaration statements
82!
83! 1308 2014-03-13 14:58:42Z fricke
84! output of the fixed number of output time levels
85! output_format adjusted for masked data if netcdf_data_format > 5
86!
87! 1299 2014-03-06 13:15:21Z heinze
88! output for using large_scale subsidence in combination
89! with large_scale_forcing
90! reformatting, more detailed explanations
91!
92! 1241 2013-10-30 11:36:58Z heinze
93! output for nudging + large scale forcing from external file
94!
95! 1216 2013-08-26 09:31:42Z raasch
96! output for transpose_compute_overlap
97!
98! 1212 2013-08-15 08:46:27Z raasch
99! output for poisfft_hybrid removed
100!
101! 1179 2013-06-14 05:57:58Z raasch
102! output of reference_state, use_reference renamed use_single_reference_value
103!
104! 1159 2013-05-21 11:58:22Z fricke
105! +use_cmax
106!
107! 1115 2013-03-26 18:16:16Z hoffmann
108! descriptions for Seifert-Beheng-cloud-physics-scheme added
109!
110! 1111 2013-03-08 23:54:10Z raasch
111! output of accelerator board information
112! ibc_p_b = 2 removed
113!
114! 1108 2013-03-05 07:03:32Z raasch
115! bugfix for r1106
116!
117! 1106 2013-03-04 05:31:38Z raasch
118! some format changes for coupled runs
119!
120! 1092 2013-02-02 11:24:22Z raasch
121! unused variables removed
122!
123! 1036 2012-10-22 13:43:42Z raasch
124! code put under GPL (PALM 3.9)
125!
126! 1031 2012-10-19 14:35:30Z raasch
127! output of netCDF data format modified
128!
129! 1015 2012-09-27 09:23:24Z raasch
130! output of Adjustment of mixing length to the Prandtl mixing length at first
131! grid point above ground removed
132!
133! 1003 2012-09-14 14:35:53Z raasch
134! output of information about equal/unequal subdomain size removed
135!
136! 1001 2012-09-13 14:08:46Z raasch
137! all actions concerning leapfrog- and upstream-spline-scheme removed
138!
139! 978 2012-08-09 08:28:32Z fricke
140! -km_damp_max, outflow_damping_width
141! +pt_damping_factor, pt_damping_width
142! +z0h
143!
144! 964 2012-07-26 09:14:24Z raasch
145! output of profil-related quantities removed
146!
147! 940 2012-07-09 14:31:00Z raasch
148! Output in case of simulations for pure neutral stratification (no pt-equation
149! solved)
150!
151! 927 2012-06-06 19:15:04Z raasch
152! output of masking_method for mg-solver
153!
154! 868 2012-03-28 12:21:07Z raasch
155! translation velocity in Galilean transformation changed to 0.6 * ug
156!
157! 833 2012-02-22 08:55:55Z maronga
158! Adjusted format for leaf area density
159!
160! 828 2012-02-21 12:00:36Z raasch
161! output of dissipation_classes + radius_classes
162!
163! 825 2012-02-19 03:03:44Z raasch
164! Output of cloud physics parameters/quantities complemented and restructured
165!
166! Revision 1.1  1997/08/11 06:17:20  raasch
167! Initial revision
168!
169!
170! Description:
171! ------------
172! Writing a header with all important informations about the actual run.
173! This subroutine is called three times, two times at the beginning
174! (writing information on files RUN_CONTROL and HEADER) and one time at the
175! end of the run, then writing additional information about CPU-usage on file
176! header.
177!-----------------------------------------------------------------------------!
178
179    USE arrays_3d,                                                             &
180        ONLY:  pt_init, qsws, q_init, sa_init, shf, ug, vg, w_subs, zu
181       
182    USE control_parameters
183       
184    USE cloud_parameters,                                                      &
185        ONLY:  cp, curvature_solution_effects, c_sedimentation,                &
186               limiter_sedimentation, l_v, nc_const, r_d, ventilation_effect
187       
188    USE cpulog,                                                                &
189        ONLY:  log_point_s
190       
191    USE dvrp_variables,                                                        &
192        ONLY:  use_seperate_pe_for_dvrp_output
193       
194    USE grid_variables,                                                        &
195        ONLY:  dx, dy
196       
197    USE indices,                                                               &
198        ONLY:  mg_loc_ind, nnx, nny, nnz, nx, ny, nxl_mg, nxr_mg, nyn_mg,      &
199               nys_mg, nzt, nzt_mg
200       
201    USE kinds
202   
203    USE model_1d,                                                              &
204        ONLY:  damp_level_ind_1d, dt_pr_1d, dt_run_control_1d, end_time_1d
205       
206    USE particle_attributes,                                                   &
207        ONLY:  bc_par_b, bc_par_lr, bc_par_ns, bc_par_t, collision_kernel,     &
208               density_ratio, dissipation_classes, dt_min_part, dt_prel,       &
209               dt_write_particle_data, end_time_prel,                          &
210               maximum_number_of_tailpoints, maximum_tailpoint_age,            &
211               minimum_tailpoint_distance, number_of_particle_groups,          &
212               particle_advection, particle_advection_start,                   &
213               particles_per_point, pdx, pdy, pdz,  psb, psl, psn, psr, pss,   &
214               pst, radius, radius_classes, random_start_position,             &
215               total_number_of_particles, use_particle_tails,                  &
216               use_sgs_for_particles, total_number_of_tails,                   &
217               vertical_particle_advection, write_particle_statistics
218       
219    USE pegrid
220
221    USE plant_canopy_model_mod,                                                &
222        ONLY:  alpha_lad, beta_lad, calc_beta_lad_profile, canopy_drag_coeff,  &
223               canopy_mode, cthf, lad, lad_surface, lad_vertical_gradient,     &
224               lad_vertical_gradient_level, lad_vertical_gradient_level_ind,   &
225               lai_beta, leaf_scalar_exch_coeff, leaf_surface_conc, pch_index, &
226               plant_canopy
227   
228    USE spectrum,                                                              &
229        ONLY:  comp_spectra_level, data_output_sp, plot_spectra_level,         &
230               spectra_direction
231
232    IMPLICIT NONE
233
234    CHARACTER (LEN=1)  ::  prec                !:
235   
236    CHARACTER (LEN=2)  ::  do2d_mode           !:
237   
238    CHARACTER (LEN=5)  ::  section_chr         !:
239   
240    CHARACTER (LEN=10) ::  coor_chr            !:
241    CHARACTER (LEN=10) ::  host_chr            !:
242   
243    CHARACTER (LEN=16) ::  begin_chr           !:
244   
245    CHARACTER (LEN=26) ::  ver_rev             !:
246   
247    CHARACTER (LEN=40) ::  output_format       !:
248   
249    CHARACTER (LEN=70) ::  char1               !:
250    CHARACTER (LEN=70) ::  char2               !:
251    CHARACTER (LEN=70) ::  dopr_chr            !:
252    CHARACTER (LEN=70) ::  do2d_xy             !:
253    CHARACTER (LEN=70) ::  do2d_xz             !:
254    CHARACTER (LEN=70) ::  do2d_yz             !:
255    CHARACTER (LEN=70) ::  do3d_chr            !:
256    CHARACTER (LEN=70) ::  domask_chr          !:
257    CHARACTER (LEN=70) ::  run_classification  !:
258   
259    CHARACTER (LEN=85) ::  roben               !:
260    CHARACTER (LEN=85) ::  runten              !:
261   
262    CHARACTER (LEN=86) ::  coordinates         !:
263    CHARACTER (LEN=86) ::  gradients           !:
264    CHARACTER (LEN=86) ::  leaf_area_density   !:
265    CHARACTER (LEN=86) ::  slices              !:
266    CHARACTER (LEN=86) ::  temperatures        !:
267    CHARACTER (LEN=86) ::  ugcomponent         !:
268    CHARACTER (LEN=86) ::  vgcomponent         !:
269
270    CHARACTER (LEN=1), DIMENSION(1:3) ::  dir = (/ 'x', 'y', 'z' /)  !:
271
272    INTEGER(iwp) ::  av        !:
273    INTEGER(iwp) ::  bh        !:
274    INTEGER(iwp) ::  blx       !:
275    INTEGER(iwp) ::  bly       !:
276    INTEGER(iwp) ::  bxl       !:
277    INTEGER(iwp) ::  bxr       !:
278    INTEGER(iwp) ::  byn       !:
279    INTEGER(iwp) ::  bys       !:
280    INTEGER(iwp) ::  ch        !:
281    INTEGER(iwp) ::  count     !:
282    INTEGER(iwp) ::  cwx       !:
283    INTEGER(iwp) ::  cwy       !:
284    INTEGER(iwp) ::  cxl       !:
285    INTEGER(iwp) ::  cxr       !:
286    INTEGER(iwp) ::  cyn       !:
287    INTEGER(iwp) ::  cys       !:
288    INTEGER(iwp) ::  dim       !:
289    INTEGER(iwp) ::  i         !:
290    INTEGER(iwp) ::  io        !:
291    INTEGER(iwp) ::  j         !:
292    INTEGER(iwp) ::  k         !:
293    INTEGER(iwp) ::  l         !:
294    INTEGER(iwp) ::  ll        !:
295    INTEGER(iwp) ::  mpi_type  !:
296   
297    REAL(wp) ::  canopy_height                    !: canopy height (in m)
298    REAL(wp) ::  cpuseconds_per_simulated_second  !:
299
300!
301!-- Open the output file. At the end of the simulation, output is directed
302!-- to unit 19.
303    IF ( ( runnr == 0 .OR. force_print_header )  .AND. &
304         .NOT. simulated_time_at_begin /= simulated_time )  THEN
305       io = 15   !  header output on file RUN_CONTROL
306    ELSE
307       io = 19   !  header output on file HEADER
308    ENDIF
309    CALL check_open( io )
310
311!
312!-- At the end of the run, output file (HEADER) will be rewritten with
313!-- new informations
314    IF ( io == 19 .AND. simulated_time_at_begin /= simulated_time ) REWIND( 19 )
315
316!
317!-- Determine kind of model run
318    IF ( TRIM( initializing_actions ) == 'read_restart_data' )  THEN
319       run_classification = '3D - restart run'
320    ELSEIF ( TRIM( initializing_actions ) == 'cyclic_fill' )  THEN
321       run_classification = '3D - run with cyclic fill of 3D - prerun data'
322    ELSEIF ( INDEX( initializing_actions, 'set_constant_profiles' ) /= 0 )  THEN
323       run_classification = '3D - run without 1D - prerun'
324    ELSEIF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
325       run_classification = '3D - run with 1D - prerun'
326    ELSEIF ( INDEX( initializing_actions, 'by_user' ) /=0 )  THEN
327       run_classification = '3D - run initialized by user'
328    ELSE
329       message_string = ' unknown action(s): ' // TRIM( initializing_actions )
330       CALL message( 'header', 'PA0191', 0, 0, 0, 6, 0 )
331    ENDIF
332    IF ( ocean )  THEN
333       run_classification = 'ocean - ' // run_classification
334    ELSE
335       run_classification = 'atmosphere - ' // run_classification
336    ENDIF
337
338!
339!-- Run-identification, date, time, host
340    host_chr = host(1:10)
341    ver_rev = TRIM( version ) // '  ' // TRIM( revision )
342    WRITE ( io, 100 )  ver_rev, TRIM( run_classification )
343    IF ( TRIM( coupling_mode ) /= 'uncoupled' )  THEN
344#if defined( __mpi2 )
345       mpi_type = 2
346#else
347       mpi_type = 1
348#endif
349       WRITE ( io, 101 )  mpi_type, coupling_mode
350    ENDIF
351#if defined( __parallel )
352    IF ( coupling_start_time /= 0.0_wp )  THEN
353       IF ( coupling_start_time > simulated_time_at_begin )  THEN
354          WRITE ( io, 109 )
355       ELSE
356          WRITE ( io, 114 )
357       ENDIF
358    ENDIF
359#endif
360    IF ( ensemble_member_nr /= 0 )  THEN
361       WRITE ( io, 512 )  run_date, run_identifier, run_time, runnr,           &
362                       ADJUSTR( host_chr ), ensemble_member_nr
363    ELSE
364       WRITE ( io, 102 )  run_date, run_identifier, run_time, runnr,           &
365                       ADJUSTR( host_chr )
366    ENDIF
367#if defined( __parallel )
368    IF ( npex == -1  .AND.  npey == -1 )  THEN
369       char1 = 'calculated'
370    ELSE
371       char1 = 'predefined'
372    ENDIF
373    IF ( threads_per_task == 1 )  THEN
374       WRITE ( io, 103 )  numprocs, pdims(1), pdims(2), TRIM( char1 )
375    ELSE
376       WRITE ( io, 104 )  numprocs*threads_per_task, numprocs, &
377                          threads_per_task, pdims(1), pdims(2), TRIM( char1 )
378    ENDIF
379    IF ( num_acc_per_node /= 0 )  WRITE ( io, 117 )  num_acc_per_node   
380    IF ( ( host(1:3) == 'ibm'  .OR.  host(1:3) == 'nec'  .OR.    &
381           host(1:2) == 'lc'   .OR.  host(1:3) == 'dec' )  .AND. &
382         npex == -1  .AND.  pdims(2) == 1 )                      &
383    THEN
384       WRITE ( io, 106 )
385    ELSEIF ( pdims(2) == 1 )  THEN
386       WRITE ( io, 107 )  'x'
387    ELSEIF ( pdims(1) == 1 )  THEN
388       WRITE ( io, 107 )  'y'
389    ENDIF
390    IF ( use_seperate_pe_for_dvrp_output )  WRITE ( io, 105 )
391    IF ( numprocs /= maximum_parallel_io_streams )  THEN
392       WRITE ( io, 108 )  maximum_parallel_io_streams
393    ENDIF
394#else
395    IF ( num_acc_per_node /= 0 )  WRITE ( io, 120 )  num_acc_per_node
396#endif
397    WRITE ( io, 99 )
398
399!
400!-- Numerical schemes
401    WRITE ( io, 110 )
402    IF ( psolver(1:7) == 'poisfft' )  THEN
403       WRITE ( io, 111 )  TRIM( fft_method )
404       IF ( transpose_compute_overlap )  WRITE( io, 115 )
405    ELSEIF ( psolver == 'sor' )  THEN
406       WRITE ( io, 112 )  nsor_ini, nsor, omega_sor
407    ELSEIF ( psolver == 'multigrid' )  THEN
408       WRITE ( io, 135 )  cycle_mg, maximum_grid_level, ngsrb
409       IF ( mg_cycles == -1 )  THEN
410          WRITE ( io, 140 )  residual_limit
411       ELSE
412          WRITE ( io, 141 )  mg_cycles
413       ENDIF
414       IF ( mg_switch_to_pe0_level == 0 )  THEN
415          WRITE ( io, 136 )  nxr_mg(1)-nxl_mg(1)+1, nyn_mg(1)-nys_mg(1)+1, &
416                             nzt_mg(1)
417       ELSEIF (  mg_switch_to_pe0_level /= -1 )  THEN
418          WRITE ( io, 137 )  mg_switch_to_pe0_level,            &
419                             mg_loc_ind(2,0)-mg_loc_ind(1,0)+1, &
420                             mg_loc_ind(4,0)-mg_loc_ind(3,0)+1, &
421                             nzt_mg(mg_switch_to_pe0_level),    &
422                             nxr_mg(1)-nxl_mg(1)+1, nyn_mg(1)-nys_mg(1)+1, &
423                             nzt_mg(1)
424       ENDIF
425       IF ( masking_method )  WRITE ( io, 144 )
426    ENDIF
427    IF ( call_psolver_at_all_substeps  .AND. timestep_scheme(1:5) == 'runge' ) &
428    THEN
429       WRITE ( io, 142 )
430    ENDIF
431
432    IF ( momentum_advec == 'pw-scheme' )  THEN
433       WRITE ( io, 113 )
434    ELSEIF (momentum_advec == 'ws-scheme' )  THEN
435       WRITE ( io, 503 )
436    ENDIF
437    IF ( scalar_advec == 'pw-scheme' )  THEN
438       WRITE ( io, 116 )
439    ELSEIF ( scalar_advec == 'ws-scheme' )  THEN
440       WRITE ( io, 504 )
441    ELSE
442       WRITE ( io, 118 )
443    ENDIF
444
445    WRITE ( io, 139 )  TRIM( loop_optimization )
446
447    IF ( galilei_transformation )  THEN
448       IF ( use_ug_for_galilei_tr )  THEN
449          char1 = '0.6 * geostrophic wind'
450       ELSE
451          char1 = 'mean wind in model domain'
452       ENDIF
453       IF ( simulated_time_at_begin == simulated_time )  THEN
454          char2 = 'at the start of the run'
455       ELSE
456          char2 = 'at the end of the run'
457       ENDIF
458       WRITE ( io, 119 )  TRIM( char1 ), TRIM( char2 ),                        &
459                          advected_distance_x/1000.0_wp,                       &
460                          advected_distance_y/1000.0_wp
461    ENDIF
462    WRITE ( io, 122 )  timestep_scheme
463    IF ( use_upstream_for_tke )  WRITE ( io, 143 )
464    IF ( rayleigh_damping_factor /= 0.0_wp )  THEN
465       IF ( .NOT. ocean )  THEN
466          WRITE ( io, 123 )  'above', rayleigh_damping_height, &
467               rayleigh_damping_factor
468       ELSE
469          WRITE ( io, 123 )  'below', rayleigh_damping_height, &
470               rayleigh_damping_factor
471       ENDIF
472    ENDIF
473    IF ( neutral )  WRITE ( io, 131 )  pt_surface
474    IF ( humidity )  THEN
475       IF ( .NOT. cloud_physics )  THEN
476          WRITE ( io, 129 )
477       ELSE
478          WRITE ( io, 130 )
479       ENDIF
480    ENDIF
481    IF ( passive_scalar )  WRITE ( io, 134 )
482    IF ( conserve_volume_flow )  THEN
483       WRITE ( io, 150 )  conserve_volume_flow_mode
484       IF ( TRIM( conserve_volume_flow_mode ) == 'bulk_velocity' )  THEN
485          WRITE ( io, 151 )  u_bulk, v_bulk
486       ENDIF
487    ELSEIF ( dp_external )  THEN
488       IF ( dp_smooth )  THEN
489          WRITE ( io, 152 )  dpdxy, dp_level_b, ', vertically smoothed.'
490       ELSE
491          WRITE ( io, 152 )  dpdxy, dp_level_b, '.'
492       ENDIF
493    ENDIF
494    WRITE ( io, 99 )
495
496!
497!-- Runtime and timestep informations
498    WRITE ( io, 200 )
499    IF ( .NOT. dt_fixed )  THEN
500       WRITE ( io, 201 )  dt_max, cfl_factor
501    ELSE
502       WRITE ( io, 202 )  dt
503    ENDIF
504    WRITE ( io, 203 )  simulated_time_at_begin, end_time
505
506    IF ( time_restart /= 9999999.9_wp  .AND. &
507         simulated_time_at_begin == simulated_time )  THEN
508       IF ( dt_restart == 9999999.9_wp )  THEN
509          WRITE ( io, 204 )  ' Restart at:       ',time_restart
510       ELSE
511          WRITE ( io, 205 )  ' Restart at:       ',time_restart, dt_restart
512       ENDIF
513    ENDIF
514
515    IF ( simulated_time_at_begin /= simulated_time )  THEN
516       i = MAX ( log_point_s(10)%counts, 1 )
517       IF ( ( simulated_time - simulated_time_at_begin ) == 0.0_wp )  THEN
518          cpuseconds_per_simulated_second = 0.0_wp
519       ELSE
520          cpuseconds_per_simulated_second = log_point_s(10)%sum / &
521                                            ( simulated_time -    &
522                                              simulated_time_at_begin )
523       ENDIF
524       WRITE ( io, 206 )  simulated_time, log_point_s(10)%sum,      &
525                          log_point_s(10)%sum / REAL( i, KIND=wp ), &
526                          cpuseconds_per_simulated_second
527       IF ( time_restart /= 9999999.9_wp  .AND.  time_restart < end_time )  THEN
528          IF ( dt_restart == 9999999.9_wp )  THEN
529             WRITE ( io, 204 )  ' Next restart at:     ',time_restart
530          ELSE
531             WRITE ( io, 205 )  ' Next restart at:     ',time_restart, dt_restart
532          ENDIF
533       ENDIF
534    ENDIF
535
536
537!
538!-- Start time for coupled runs, if independent precursor runs for atmosphere
539!-- and ocean are used or have been used. In this case, coupling_start_time
540!-- defines the time when the coupling is switched on.
541    IF ( coupling_start_time /= 0.0_wp )  THEN
542       WRITE ( io, 207 )  coupling_start_time
543    ENDIF
544
545!
546!-- Computational grid
547    IF ( .NOT. ocean )  THEN
548       WRITE ( io, 250 )  dx, dy, dz, (nx+1)*dx, (ny+1)*dy, zu(nzt+1)
549       IF ( dz_stretch_level_index < nzt+1 )  THEN
550          WRITE ( io, 252 )  dz_stretch_level, dz_stretch_level_index, &
551                             dz_stretch_factor, dz_max
552       ENDIF
553    ELSE
554       WRITE ( io, 250 )  dx, dy, dz, (nx+1)*dx, (ny+1)*dy, zu(0)
555       IF ( dz_stretch_level_index > 0 )  THEN
556          WRITE ( io, 252 )  dz_stretch_level, dz_stretch_level_index, &
557                             dz_stretch_factor, dz_max
558       ENDIF
559    ENDIF
560    WRITE ( io, 254 )  nx, ny, nzt+1, MIN( nnx, nx+1 ), MIN( nny, ny+1 ), &
561                       MIN( nnz+2, nzt+2 )
562    IF ( sloping_surface )  WRITE ( io, 260 )  alpha_surface
563
564!
565!-- Large scale forcing and nudging
566    WRITE ( io, 160 )
567    IF ( large_scale_forcing )  THEN
568       WRITE ( io, 162 )
569       WRITE ( io, 163 )
570
571       IF ( large_scale_subsidence )  THEN
572          IF ( .NOT. use_subsidence_tendencies )  THEN
573             WRITE ( io, 164 )
574          ELSE
575             WRITE ( io, 165 )
576          ENDIF
577       ENDIF
578
579       IF ( bc_pt_b == 'dirichlet' )  THEN
580          WRITE ( io, 180 )
581       ELSEIF ( bc_pt_b == 'neumann' )  THEN
582          WRITE ( io, 181 )
583       ENDIF
584
585       IF ( bc_q_b == 'dirichlet' )  THEN
586          WRITE ( io, 182 )
587       ELSEIF ( bc_q_b == 'neumann' )  THEN
588          WRITE ( io, 183 )
589       ENDIF
590
591       WRITE ( io, 167 )
592       IF ( nudging )  THEN
593          WRITE ( io, 170 )
594       ENDIF
595    ELSE
596       WRITE ( io, 161 )
597       WRITE ( io, 171 )
598    ENDIF
599    IF ( large_scale_subsidence )  THEN
600       WRITE ( io, 168 )
601       WRITE ( io, 169 )
602    ENDIF
603
604!
605!-- Profile for the large scale vertial velocity
606!-- Building output strings, starting with surface value
607    IF ( large_scale_subsidence )  THEN
608       temperatures = '   0.0'
609       gradients = '------'
610       slices = '     0'
611       coordinates = '   0.0'
612       i = 1
613       DO  WHILE ( subs_vertical_gradient_level_i(i) /= -9999 )
614
615          WRITE (coor_chr,'(E10.2,7X)')  &
616                                w_subs(subs_vertical_gradient_level_i(i))
617          temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr )
618
619          WRITE (coor_chr,'(E10.2,7X)')  subs_vertical_gradient(i)
620          gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
621
622          WRITE (coor_chr,'(I10,7X)')  subs_vertical_gradient_level_i(i)
623          slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
624
625          WRITE (coor_chr,'(F10.2,7X)')  subs_vertical_gradient_level(i)
626          coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
627
628          IF ( i == 10 )  THEN
629             EXIT
630          ELSE
631             i = i + 1
632          ENDIF
633
634       ENDDO
635
636 
637       IF ( .NOT. large_scale_forcing )  THEN
638          WRITE ( io, 426 )  TRIM( coordinates ), TRIM( temperatures ), &
639                             TRIM( gradients ), TRIM( slices )
640       ENDIF
641
642
643    ENDIF
644
645!-- Profile of the geostrophic wind (component ug)
646!-- Building output strings
647    WRITE ( ugcomponent, '(F6.2)' )  ug_surface
648    gradients = '------'
649    slices = '     0'
650    coordinates = '   0.0'
651    i = 1
652    DO  WHILE ( ug_vertical_gradient_level_ind(i) /= -9999 )
653     
654       WRITE (coor_chr,'(F6.2,1X)')  ug(ug_vertical_gradient_level_ind(i))
655       ugcomponent = TRIM( ugcomponent ) // '  ' // TRIM( coor_chr )
656
657       WRITE (coor_chr,'(F6.2,1X)')  ug_vertical_gradient(i)
658       gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
659
660       WRITE (coor_chr,'(I6,1X)')  ug_vertical_gradient_level_ind(i)
661       slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
662
663       WRITE (coor_chr,'(F6.1,1X)')  ug_vertical_gradient_level(i)
664       coordinates = TRIM( coordinates ) // '  ' // TRIM( coor_chr )
665
666       IF ( i == 10 )  THEN
667          EXIT
668       ELSE
669          i = i + 1
670       ENDIF
671
672    ENDDO
673
674    IF ( .NOT. large_scale_forcing )  THEN
675       WRITE ( io, 423 )  TRIM( coordinates ), TRIM( ugcomponent ), &
676                          TRIM( gradients ), TRIM( slices )
677    ENDIF
678
679!-- Profile of the geostrophic wind (component vg)
680!-- Building output strings
681    WRITE ( vgcomponent, '(F6.2)' )  vg_surface
682    gradients = '------'
683    slices = '     0'
684    coordinates = '   0.0'
685    i = 1
686    DO  WHILE ( vg_vertical_gradient_level_ind(i) /= -9999 )
687
688       WRITE (coor_chr,'(F6.2,1X)')  vg(vg_vertical_gradient_level_ind(i))
689       vgcomponent = TRIM( vgcomponent ) // '  ' // TRIM( coor_chr )
690
691       WRITE (coor_chr,'(F6.2,1X)')  vg_vertical_gradient(i)
692       gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
693
694       WRITE (coor_chr,'(I6,1X)')  vg_vertical_gradient_level_ind(i)
695       slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
696
697       WRITE (coor_chr,'(F6.1,1X)')  vg_vertical_gradient_level(i)
698       coordinates = TRIM( coordinates ) // '  ' // TRIM( coor_chr )
699
700       IF ( i == 10 )  THEN
701          EXIT
702       ELSE
703          i = i + 1
704       ENDIF
705 
706    ENDDO
707
708    IF ( .NOT. large_scale_forcing )  THEN
709       WRITE ( io, 424 )  TRIM( coordinates ), TRIM( vgcomponent ), &
710                          TRIM( gradients ), TRIM( slices )
711    ENDIF
712
713!
714!-- Topography
715    WRITE ( io, 270 )  topography
716    SELECT CASE ( TRIM( topography ) )
717
718       CASE ( 'flat' )
719          ! no actions necessary
720
721       CASE ( 'single_building' )
722          blx = INT( building_length_x / dx )
723          bly = INT( building_length_y / dy )
724          bh  = INT( building_height / dz )
725
726          IF ( building_wall_left == 9999999.9_wp )  THEN
727             building_wall_left = ( nx + 1 - blx ) / 2 * dx
728          ENDIF
729          bxl = INT ( building_wall_left / dx + 0.5_wp )
730          bxr = bxl + blx
731
732          IF ( building_wall_south == 9999999.9_wp )  THEN
733             building_wall_south = ( ny + 1 - bly ) / 2 * dy
734          ENDIF
735          bys = INT ( building_wall_south / dy + 0.5_wp )
736          byn = bys + bly
737
738          WRITE ( io, 271 )  building_length_x, building_length_y, &
739                             building_height, bxl, bxr, bys, byn
740
741       CASE ( 'single_street_canyon' )
742          ch  = NINT( canyon_height / dz )
743          IF ( canyon_width_x /= 9999999.9_wp )  THEN
744!
745!--          Street canyon in y direction
746             cwx = NINT( canyon_width_x / dx )
747             IF ( canyon_wall_left == 9999999.9_wp )  THEN
748                canyon_wall_left = ( nx + 1 - cwx ) / 2 * dx
749             ENDIF
750             cxl = NINT( canyon_wall_left / dx )
751             cxr = cxl + cwx
752             WRITE ( io, 272 )  'y', canyon_height, ch, 'u', cxl, cxr
753
754          ELSEIF ( canyon_width_y /= 9999999.9_wp )  THEN
755!
756!--          Street canyon in x direction
757             cwy = NINT( canyon_width_y / dy )
758             IF ( canyon_wall_south == 9999999.9_wp )  THEN
759                canyon_wall_south = ( ny + 1 - cwy ) / 2 * dy
760             ENDIF
761             cys = NINT( canyon_wall_south / dy )
762             cyn = cys + cwy
763             WRITE ( io, 272 )  'x', canyon_height, ch, 'v', cys, cyn
764          ENDIF
765
766    END SELECT
767
768    IF ( TRIM( topography ) /= 'flat' )  THEN
769       IF ( TRIM( topography_grid_convention ) == ' ' )  THEN
770          IF ( TRIM( topography ) == 'single_building' .OR.  &
771               TRIM( topography ) == 'single_street_canyon' )  THEN
772             WRITE ( io, 278 )
773          ELSEIF ( TRIM( topography ) == 'read_from_file' )  THEN
774             WRITE ( io, 279 )
775          ENDIF
776       ELSEIF ( TRIM( topography_grid_convention ) == 'cell_edge' )  THEN
777          WRITE ( io, 278 )
778       ELSEIF ( TRIM( topography_grid_convention ) == 'cell_center' )  THEN
779          WRITE ( io, 279 )
780       ENDIF
781    ENDIF
782
783    IF ( plant_canopy )  THEN
784   
785       canopy_height = pch_index * dz
786
787       WRITE ( io, 280 )  canopy_mode, canopy_height, pch_index,               &
788                          canopy_drag_coeff
789       IF ( passive_scalar )  THEN
790          WRITE ( io, 281 )  leaf_scalar_exch_coeff,                           &
791                             leaf_surface_conc
792       ENDIF
793
794!
795!--    Heat flux at the top of vegetation
796       WRITE ( io, 282 )  cthf
797
798!
799!--    Leaf area density profile, calculated either from given vertical
800!--    gradients or from beta probability density function.
801       IF (  .NOT.  calc_beta_lad_profile )  THEN
802
803!--       Building output strings, starting with surface value
804          WRITE ( leaf_area_density, '(F6.4)' )  lad_surface
805          gradients = '------'
806          slices = '     0'
807          coordinates = '   0.0'
808          i = 1
809          DO  WHILE ( i < 11  .AND.  lad_vertical_gradient_level_ind(i) /= -9999 )
810
811             WRITE (coor_chr,'(F7.2)')  lad(lad_vertical_gradient_level_ind(i))
812             leaf_area_density = TRIM( leaf_area_density ) // ' ' // TRIM( coor_chr )
813 
814             WRITE (coor_chr,'(F7.2)')  lad_vertical_gradient(i)
815             gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
816
817             WRITE (coor_chr,'(I7)')  lad_vertical_gradient_level_ind(i)
818             slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
819
820             WRITE (coor_chr,'(F7.1)')  lad_vertical_gradient_level(i)
821             coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
822
823             i = i + 1
824          ENDDO
825
826          WRITE ( io, 283 )  TRIM( coordinates ), TRIM( leaf_area_density ),              &
827                             TRIM( gradients ), TRIM( slices )
828
829       ELSE
830       
831          WRITE ( leaf_area_density, '(F6.4)' )  lad_surface
832          coordinates = '   0.0'
833         
834          DO  k = 1, pch_index
835
836             WRITE (coor_chr,'(F7.2)')  lad(k)
837             leaf_area_density = TRIM( leaf_area_density ) // ' ' // TRIM( coor_chr )
838 
839             WRITE (coor_chr,'(F7.1)')  zu(k)
840             coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
841
842          ENDDO       
843
844          WRITE ( io, 284 ) TRIM( coordinates ), TRIM( leaf_area_density ), alpha_lad,    &
845                            beta_lad, lai_beta
846
847       ENDIF 
848
849    ENDIF
850
851
852!
853!-- Boundary conditions
854    IF ( ibc_p_b == 0 )  THEN
855       runten = 'p(0)     = 0      |'
856    ELSEIF ( ibc_p_b == 1 )  THEN
857       runten = 'p(0)     = p(1)   |'
858    ENDIF
859    IF ( ibc_p_t == 0 )  THEN
860       roben  = 'p(nzt+1) = 0      |'
861    ELSE
862       roben  = 'p(nzt+1) = p(nzt) |'
863    ENDIF
864
865    IF ( ibc_uv_b == 0 )  THEN
866       runten = TRIM( runten ) // ' uv(0)     = -uv(1)                |'
867    ELSE
868       runten = TRIM( runten ) // ' uv(0)     = uv(1)                 |'
869    ENDIF
870    IF ( TRIM( bc_uv_t ) == 'dirichlet_0' )  THEN
871       roben  = TRIM( roben  ) // ' uv(nzt+1) = 0                     |'
872    ELSEIF ( ibc_uv_t == 0 )  THEN
873       roben  = TRIM( roben  ) // ' uv(nzt+1) = ug(nzt+1), vg(nzt+1)  |'
874    ELSE
875       roben  = TRIM( roben  ) // ' uv(nzt+1) = uv(nzt)               |'
876    ENDIF
877
878    IF ( ibc_pt_b == 0 )  THEN
879       runten = TRIM( runten ) // ' pt(0)   = pt_surface'
880    ELSEIF ( ibc_pt_b == 1 )  THEN
881       runten = TRIM( runten ) // ' pt(0)   = pt(1)'
882    ELSEIF ( ibc_pt_b == 2 )  THEN
883       runten = TRIM( runten ) // ' pt(0) = from coupled model'
884    ENDIF
885    IF ( ibc_pt_t == 0 )  THEN
886       roben  = TRIM( roben  ) // ' pt(nzt+1) = pt_top'
887    ELSEIF( ibc_pt_t == 1 )  THEN
888       roben  = TRIM( roben  ) // ' pt(nzt+1) = pt(nzt)'
889    ELSEIF( ibc_pt_t == 2 )  THEN
890       roben  = TRIM( roben  ) // ' pt(nzt+1) = pt(nzt) + dpt/dz_ini'
891
892    ENDIF
893
894    WRITE ( io, 300 )  runten, roben
895
896    IF ( .NOT. constant_diffusion )  THEN
897       IF ( ibc_e_b == 1 )  THEN
898          runten = 'e(0)     = e(1)'
899       ELSE
900          runten = 'e(0)     = e(1) = (u*/0.1)**2'
901       ENDIF
902       roben = 'e(nzt+1) = e(nzt) = e(nzt-1)'
903
904       WRITE ( io, 301 )  'e', runten, roben       
905
906    ENDIF
907
908    IF ( ocean )  THEN
909       runten = 'sa(0)    = sa(1)'
910       IF ( ibc_sa_t == 0 )  THEN
911          roben =  'sa(nzt+1) = sa_surface'
912       ELSE
913          roben =  'sa(nzt+1) = sa(nzt)'
914       ENDIF
915       WRITE ( io, 301 ) 'sa', runten, roben
916    ENDIF
917
918    IF ( humidity )  THEN
919       IF ( ibc_q_b == 0 )  THEN
920          runten = 'q(0)     = q_surface'
921       ELSE
922          runten = 'q(0)     = q(1)'
923       ENDIF
924       IF ( ibc_q_t == 0 )  THEN
925          roben =  'q(nzt)   = q_top'
926       ELSE
927          roben =  'q(nzt)   = q(nzt-1) + dq/dz'
928       ENDIF
929       WRITE ( io, 301 ) 'q', runten, roben
930    ENDIF
931
932    IF ( passive_scalar )  THEN
933       IF ( ibc_q_b == 0 )  THEN
934          runten = 's(0)     = s_surface'
935       ELSE
936          runten = 's(0)     = s(1)'
937       ENDIF
938       IF ( ibc_q_t == 0 )  THEN
939          roben =  's(nzt)   = s_top'
940       ELSE
941          roben =  's(nzt)   = s(nzt-1) + ds/dz'
942       ENDIF
943       WRITE ( io, 301 ) 's', runten, roben
944    ENDIF
945
946    IF ( use_surface_fluxes )  THEN
947       WRITE ( io, 303 )
948       IF ( constant_heatflux )  THEN
949          IF ( large_scale_forcing .AND. lsf_surf )  THEN
950             WRITE ( io, 306 )  shf(0,0)
951          ELSE
952             WRITE ( io, 306 )  surface_heatflux
953          ENDIF
954          IF ( random_heatflux )  WRITE ( io, 307 )
955       ENDIF
956       IF ( humidity  .AND.  constant_waterflux )  THEN
957          IF ( large_scale_forcing .AND. lsf_surf )  THEN
958             WRITE ( io, 311 ) qsws(0,0)
959          ELSE
960             WRITE ( io, 311 ) surface_waterflux
961          ENDIF
962       ENDIF
963       IF ( passive_scalar  .AND.  constant_waterflux )  THEN
964          WRITE ( io, 313 ) surface_waterflux
965       ENDIF
966    ENDIF
967
968    IF ( use_top_fluxes )  THEN
969       WRITE ( io, 304 )
970       IF ( coupling_mode == 'uncoupled' )  THEN
971          WRITE ( io, 320 )  top_momentumflux_u, top_momentumflux_v
972          IF ( constant_top_heatflux )  THEN
973             WRITE ( io, 306 )  top_heatflux
974          ENDIF
975       ELSEIF ( coupling_mode == 'ocean_to_atmosphere' )  THEN
976          WRITE ( io, 316 )
977       ENDIF
978       IF ( ocean  .AND.  constant_top_salinityflux )  THEN
979          WRITE ( io, 309 )  top_salinityflux
980       ENDIF
981       IF ( humidity  .OR.  passive_scalar )  THEN
982          WRITE ( io, 315 )
983       ENDIF
984    ENDIF
985
986    IF ( prandtl_layer )  THEN
987       WRITE ( io, 305 )  (zu(1)-zu(0)), roughness_length, &
988                          z0h_factor*roughness_length, kappa, &
989                          rif_min, rif_max
990       IF ( .NOT. constant_heatflux )  WRITE ( io, 308 )
991       IF ( humidity  .AND.  .NOT. constant_waterflux )  THEN
992          WRITE ( io, 312 )
993       ENDIF
994       IF ( passive_scalar  .AND.  .NOT. constant_waterflux )  THEN
995          WRITE ( io, 314 )
996       ENDIF
997    ELSE
998       IF ( INDEX(initializing_actions, 'set_1d-model_profiles') /= 0 )  THEN
999          WRITE ( io, 310 )  rif_min, rif_max
1000       ENDIF
1001    ENDIF
1002
1003    WRITE ( io, 317 )  bc_lr, bc_ns
1004    IF ( .NOT. bc_lr_cyc  .OR.  .NOT. bc_ns_cyc )  THEN
1005       WRITE ( io, 318 )  use_cmax, pt_damping_width, pt_damping_factor       
1006       IF ( turbulent_inflow )  THEN
1007          WRITE ( io, 319 )  recycling_width, recycling_plane, &
1008                             inflow_damping_height, inflow_damping_width
1009       ENDIF
1010    ENDIF
1011
1012!
1013!-- Initial Profiles
1014    WRITE ( io, 321 )
1015!
1016!-- Initial wind profiles
1017    IF ( u_profile(1) /= 9999999.9_wp )  WRITE ( io, 427 )
1018
1019!
1020!-- Initial temperature profile
1021!-- Building output strings, starting with surface temperature
1022    WRITE ( temperatures, '(F6.2)' )  pt_surface
1023    gradients = '------'
1024    slices = '     0'
1025    coordinates = '   0.0'
1026    i = 1
1027    DO  WHILE ( pt_vertical_gradient_level_ind(i) /= -9999 )
1028
1029       WRITE (coor_chr,'(F7.2)')  pt_init(pt_vertical_gradient_level_ind(i))
1030       temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr )
1031
1032       WRITE (coor_chr,'(F7.2)')  pt_vertical_gradient(i)
1033       gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
1034
1035       WRITE (coor_chr,'(I7)')  pt_vertical_gradient_level_ind(i)
1036       slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
1037
1038       WRITE (coor_chr,'(F7.1)')  pt_vertical_gradient_level(i)
1039       coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
1040
1041       IF ( i == 10 )  THEN
1042          EXIT
1043       ELSE
1044          i = i + 1
1045       ENDIF
1046
1047    ENDDO
1048
1049    IF ( .NOT. nudging )  THEN
1050       WRITE ( io, 420 )  TRIM( coordinates ), TRIM( temperatures ), &
1051                          TRIM( gradients ), TRIM( slices )
1052    ELSE
1053       WRITE ( io, 428 ) 
1054    ENDIF
1055
1056!
1057!-- Initial humidity profile
1058!-- Building output strings, starting with surface humidity
1059    IF ( humidity  .OR.  passive_scalar )  THEN
1060       WRITE ( temperatures, '(E8.1)' )  q_surface
1061       gradients = '--------'
1062       slices = '       0'
1063       coordinates = '     0.0'
1064       i = 1
1065       DO  WHILE ( q_vertical_gradient_level_ind(i) /= -9999 )
1066         
1067          WRITE (coor_chr,'(E8.1,4X)')  q_init(q_vertical_gradient_level_ind(i))
1068          temperatures = TRIM( temperatures ) // '  ' // TRIM( coor_chr )
1069
1070          WRITE (coor_chr,'(E8.1,4X)')  q_vertical_gradient(i)
1071          gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
1072         
1073          WRITE (coor_chr,'(I8,4X)')  q_vertical_gradient_level_ind(i)
1074          slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
1075         
1076          WRITE (coor_chr,'(F8.1,4X)')  q_vertical_gradient_level(i)
1077          coordinates = TRIM( coordinates ) // '  '  // TRIM( coor_chr )
1078
1079          IF ( i == 10 )  THEN
1080             EXIT
1081          ELSE
1082             i = i + 1
1083          ENDIF
1084
1085       ENDDO
1086
1087       IF ( humidity )  THEN
1088          IF ( .NOT. nudging )  THEN
1089             WRITE ( io, 421 )  TRIM( coordinates ), TRIM( temperatures ), &
1090                                TRIM( gradients ), TRIM( slices )
1091          ENDIF
1092       ELSE
1093          WRITE ( io, 422 )  TRIM( coordinates ), TRIM( temperatures ), &
1094                             TRIM( gradients ), TRIM( slices )
1095       ENDIF
1096    ENDIF
1097
1098!
1099!-- Initial salinity profile
1100!-- Building output strings, starting with surface salinity
1101    IF ( ocean )  THEN
1102       WRITE ( temperatures, '(F6.2)' )  sa_surface
1103       gradients = '------'
1104       slices = '     0'
1105       coordinates = '   0.0'
1106       i = 1
1107       DO  WHILE ( sa_vertical_gradient_level_ind(i) /= -9999 )
1108
1109          WRITE (coor_chr,'(F7.2)')  sa_init(sa_vertical_gradient_level_ind(i))
1110          temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr )
1111
1112          WRITE (coor_chr,'(F7.2)')  sa_vertical_gradient(i)
1113          gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
1114
1115          WRITE (coor_chr,'(I7)')  sa_vertical_gradient_level_ind(i)
1116          slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
1117
1118          WRITE (coor_chr,'(F7.1)')  sa_vertical_gradient_level(i)
1119          coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
1120
1121          IF ( i == 10 )  THEN
1122             EXIT
1123          ELSE
1124             i = i + 1
1125          ENDIF
1126
1127       ENDDO
1128
1129       WRITE ( io, 425 )  TRIM( coordinates ), TRIM( temperatures ), &
1130                          TRIM( gradients ), TRIM( slices )
1131    ENDIF
1132
1133
1134!
1135!-- Listing of 1D-profiles
1136    WRITE ( io, 325 )  dt_dopr_listing
1137    IF ( averaging_interval_pr /= 0.0_wp )  THEN
1138       WRITE ( io, 326 )  averaging_interval_pr, dt_averaging_input_pr
1139    ENDIF
1140
1141!
1142!-- DATA output
1143    WRITE ( io, 330 )
1144    IF ( averaging_interval_pr /= 0.0_wp )  THEN
1145       WRITE ( io, 326 )  averaging_interval_pr, dt_averaging_input_pr
1146    ENDIF
1147
1148!
1149!-- 1D-profiles
1150    dopr_chr = 'Profile:'
1151    IF ( dopr_n /= 0 )  THEN
1152       WRITE ( io, 331 )
1153
1154       output_format = ''
1155       output_format = output_format_netcdf
1156       WRITE ( io, 344 )  output_format
1157
1158       DO  i = 1, dopr_n
1159          dopr_chr = TRIM( dopr_chr ) // ' ' // TRIM( data_output_pr(i) ) // ','
1160          IF ( LEN_TRIM( dopr_chr ) >= 60 )  THEN
1161             WRITE ( io, 332 )  dopr_chr
1162             dopr_chr = '       :'
1163          ENDIF
1164       ENDDO
1165
1166       IF ( dopr_chr /= '' )  THEN
1167          WRITE ( io, 332 )  dopr_chr
1168       ENDIF
1169       WRITE ( io, 333 )  dt_dopr, averaging_interval_pr, dt_averaging_input_pr
1170       IF ( skip_time_dopr /= 0.0_wp )  WRITE ( io, 339 )  skip_time_dopr
1171    ENDIF
1172
1173!
1174!-- 2D-arrays
1175    DO  av = 0, 1
1176
1177       i = 1
1178       do2d_xy = ''
1179       do2d_xz = ''
1180       do2d_yz = ''
1181       DO  WHILE ( do2d(av,i) /= ' ' )
1182
1183          l = MAX( 2, LEN_TRIM( do2d(av,i) ) )
1184          do2d_mode = do2d(av,i)(l-1:l)
1185
1186          SELECT CASE ( do2d_mode )
1187             CASE ( 'xy' )
1188                ll = LEN_TRIM( do2d_xy )
1189                do2d_xy = do2d_xy(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
1190             CASE ( 'xz' )
1191                ll = LEN_TRIM( do2d_xz )
1192                do2d_xz = do2d_xz(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
1193             CASE ( 'yz' )
1194                ll = LEN_TRIM( do2d_yz )
1195                do2d_yz = do2d_yz(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
1196          END SELECT
1197
1198          i = i + 1
1199
1200       ENDDO
1201
1202       IF ( ( ( do2d_xy /= ''  .AND.  section(1,1) /= -9999 )  .OR.    &
1203              ( do2d_xz /= ''  .AND.  section(1,2) /= -9999 )  .OR.    &
1204              ( do2d_yz /= ''  .AND.  section(1,3) /= -9999 ) ) )  THEN
1205
1206          IF (  av == 0 )  THEN
1207             WRITE ( io, 334 )  ''
1208          ELSE
1209             WRITE ( io, 334 )  '(time-averaged)'
1210          ENDIF
1211
1212          IF ( do2d_at_begin )  THEN
1213             begin_chr = 'and at the start'
1214          ELSE
1215             begin_chr = ''
1216          ENDIF
1217
1218          output_format = ''
1219          output_format = output_format_netcdf
1220          WRITE ( io, 344 )  output_format
1221
1222          IF ( do2d_xy /= ''  .AND.  section(1,1) /= -9999 )  THEN
1223             i = 1
1224             slices = '/'
1225             coordinates = '/'
1226!
1227!--          Building strings with index and coordinate informations of the
1228!--          slices
1229             DO  WHILE ( section(i,1) /= -9999 )
1230
1231                WRITE (section_chr,'(I5)')  section(i,1)
1232                section_chr = ADJUSTL( section_chr )
1233                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
1234
1235                IF ( section(i,1) == -1 )  THEN
1236                   WRITE (coor_chr,'(F10.1)')  -1.0_wp
1237                ELSE
1238                   WRITE (coor_chr,'(F10.1)')  zu(section(i,1))
1239                ENDIF
1240                coor_chr = ADJUSTL( coor_chr )
1241                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
1242
1243                i = i + 1
1244             ENDDO
1245             IF ( av == 0 )  THEN
1246                WRITE ( io, 335 )  'XY', do2d_xy, dt_do2d_xy, &
1247                                   TRIM( begin_chr ), 'k', TRIM( slices ), &
1248                                   TRIM( coordinates )
1249                IF ( skip_time_do2d_xy /= 0.0_wp )  THEN
1250                   WRITE ( io, 339 )  skip_time_do2d_xy
1251                ENDIF
1252             ELSE
1253                WRITE ( io, 342 )  'XY', do2d_xy, dt_data_output_av, &
1254                                   TRIM( begin_chr ), averaging_interval, &
1255                                   dt_averaging_input, 'k', TRIM( slices ), &
1256                                   TRIM( coordinates )
1257                IF ( skip_time_data_output_av /= 0.0_wp )  THEN
1258                   WRITE ( io, 339 )  skip_time_data_output_av
1259                ENDIF
1260             ENDIF
1261             IF ( netcdf_data_format > 4 )  THEN
1262                WRITE ( io, 352 )  ntdim_2d_xy(av)
1263             ELSE
1264                WRITE ( io, 353 )
1265             ENDIF
1266          ENDIF
1267
1268          IF ( do2d_xz /= ''  .AND.  section(1,2) /= -9999 )  THEN
1269             i = 1
1270             slices = '/'
1271             coordinates = '/'
1272!
1273!--          Building strings with index and coordinate informations of the
1274!--          slices
1275             DO  WHILE ( section(i,2) /= -9999 )
1276
1277                WRITE (section_chr,'(I5)')  section(i,2)
1278                section_chr = ADJUSTL( section_chr )
1279                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
1280
1281                WRITE (coor_chr,'(F10.1)')  section(i,2) * dy
1282                coor_chr = ADJUSTL( coor_chr )
1283                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
1284
1285                i = i + 1
1286             ENDDO
1287             IF ( av == 0 )  THEN
1288                WRITE ( io, 335 )  'XZ', do2d_xz, dt_do2d_xz, &
1289                                   TRIM( begin_chr ), 'j', TRIM( slices ), &
1290                                   TRIM( coordinates )
1291                IF ( skip_time_do2d_xz /= 0.0_wp )  THEN
1292                   WRITE ( io, 339 )  skip_time_do2d_xz
1293                ENDIF
1294             ELSE
1295                WRITE ( io, 342 )  'XZ', do2d_xz, dt_data_output_av, &
1296                                   TRIM( begin_chr ), averaging_interval, &
1297                                   dt_averaging_input, 'j', TRIM( slices ), &
1298                                   TRIM( coordinates )
1299                IF ( skip_time_data_output_av /= 0.0_wp )  THEN
1300                   WRITE ( io, 339 )  skip_time_data_output_av
1301                ENDIF
1302             ENDIF
1303             IF ( netcdf_data_format > 4 )  THEN
1304                WRITE ( io, 352 )  ntdim_2d_xz(av)
1305             ELSE
1306                WRITE ( io, 353 )
1307             ENDIF
1308          ENDIF
1309
1310          IF ( do2d_yz /= ''  .AND.  section(1,3) /= -9999 )  THEN
1311             i = 1
1312             slices = '/'
1313             coordinates = '/'
1314!
1315!--          Building strings with index and coordinate informations of the
1316!--          slices
1317             DO  WHILE ( section(i,3) /= -9999 )
1318
1319                WRITE (section_chr,'(I5)')  section(i,3)
1320                section_chr = ADJUSTL( section_chr )
1321                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
1322
1323                WRITE (coor_chr,'(F10.1)')  section(i,3) * dx
1324                coor_chr = ADJUSTL( coor_chr )
1325                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
1326
1327                i = i + 1
1328             ENDDO
1329             IF ( av == 0 )  THEN
1330                WRITE ( io, 335 )  'YZ', do2d_yz, dt_do2d_yz, &
1331                                   TRIM( begin_chr ), 'i', TRIM( slices ), &
1332                                   TRIM( coordinates )
1333                IF ( skip_time_do2d_yz /= 0.0_wp )  THEN
1334                   WRITE ( io, 339 )  skip_time_do2d_yz
1335                ENDIF
1336             ELSE
1337                WRITE ( io, 342 )  'YZ', do2d_yz, dt_data_output_av, &
1338                                   TRIM( begin_chr ), averaging_interval, &
1339                                   dt_averaging_input, 'i', TRIM( slices ), &
1340                                   TRIM( coordinates )
1341                IF ( skip_time_data_output_av /= 0.0_wp )  THEN
1342                   WRITE ( io, 339 )  skip_time_data_output_av
1343                ENDIF
1344             ENDIF
1345             IF ( netcdf_data_format > 4 )  THEN
1346                WRITE ( io, 352 )  ntdim_2d_yz(av)
1347             ELSE
1348                WRITE ( io, 353 )
1349             ENDIF
1350          ENDIF
1351
1352       ENDIF
1353
1354    ENDDO
1355
1356!
1357!-- 3d-arrays
1358    DO  av = 0, 1
1359
1360       i = 1
1361       do3d_chr = ''
1362       DO  WHILE ( do3d(av,i) /= ' ' )
1363
1364          do3d_chr = TRIM( do3d_chr ) // ' ' // TRIM( do3d(av,i) ) // ','
1365          i = i + 1
1366
1367       ENDDO
1368
1369       IF ( do3d_chr /= '' )  THEN
1370          IF ( av == 0 )  THEN
1371             WRITE ( io, 336 )  ''
1372          ELSE
1373             WRITE ( io, 336 )  '(time-averaged)'
1374          ENDIF
1375
1376          output_format = output_format_netcdf
1377          WRITE ( io, 344 )  output_format
1378
1379          IF ( do3d_at_begin )  THEN
1380             begin_chr = 'and at the start'
1381          ELSE
1382             begin_chr = ''
1383          ENDIF
1384          IF ( av == 0 )  THEN
1385             WRITE ( io, 337 )  do3d_chr, dt_do3d, TRIM( begin_chr ), &
1386                                zu(nz_do3d), nz_do3d
1387          ELSE
1388             WRITE ( io, 343 )  do3d_chr, dt_data_output_av,           &
1389                                TRIM( begin_chr ), averaging_interval, &
1390                                dt_averaging_input, zu(nz_do3d), nz_do3d
1391          ENDIF
1392
1393          IF ( netcdf_data_format > 4 )  THEN
1394             WRITE ( io, 352 )  ntdim_3d(av)
1395          ELSE
1396             WRITE ( io, 353 )
1397          ENDIF
1398
1399          IF ( av == 0 )  THEN
1400             IF ( skip_time_do3d /= 0.0_wp )  THEN
1401                WRITE ( io, 339 )  skip_time_do3d
1402             ENDIF
1403          ELSE
1404             IF ( skip_time_data_output_av /= 0.0_wp )  THEN
1405                WRITE ( io, 339 )  skip_time_data_output_av
1406             ENDIF
1407          ENDIF
1408
1409       ENDIF
1410
1411    ENDDO
1412
1413!
1414!-- masked arrays
1415    IF ( masks > 0 )  WRITE ( io, 345 )  &
1416         mask_scale_x, mask_scale_y, mask_scale_z
1417    DO  mid = 1, masks
1418       DO  av = 0, 1
1419
1420          i = 1
1421          domask_chr = ''
1422          DO  WHILE ( domask(mid,av,i) /= ' ' )
1423             domask_chr = TRIM( domask_chr ) // ' ' //  &
1424                          TRIM( domask(mid,av,i) ) // ','
1425             i = i + 1
1426          ENDDO
1427
1428          IF ( domask_chr /= '' )  THEN
1429             IF ( av == 0 )  THEN
1430                WRITE ( io, 346 )  '', mid
1431             ELSE
1432                WRITE ( io, 346 )  ' (time-averaged)', mid
1433             ENDIF
1434
1435             output_format = output_format_netcdf
1436!--          Parallel output not implemented for mask data, hence
1437!--          output_format must be adjusted.
1438             IF ( netcdf_data_format == 5 ) output_format = 'netCDF4/HDF5'
1439             IF ( netcdf_data_format == 6 ) output_format = 'netCDF4/HDF5 classic'
1440             WRITE ( io, 344 )  output_format
1441
1442             IF ( av == 0 )  THEN
1443                WRITE ( io, 347 )  domask_chr, dt_domask(mid)
1444             ELSE
1445                WRITE ( io, 348 )  domask_chr, dt_data_output_av, &
1446                                   averaging_interval, dt_averaging_input
1447             ENDIF
1448
1449             IF ( av == 0 )  THEN
1450                IF ( skip_time_domask(mid) /= 0.0_wp )  THEN
1451                   WRITE ( io, 339 )  skip_time_domask(mid)
1452                ENDIF
1453             ELSE
1454                IF ( skip_time_data_output_av /= 0.0_wp )  THEN
1455                   WRITE ( io, 339 )  skip_time_data_output_av
1456                ENDIF
1457             ENDIF
1458!
1459!--          output locations
1460             DO  dim = 1, 3
1461                IF ( mask(mid,dim,1) >= 0.0_wp )  THEN
1462                   count = 0
1463                   DO  WHILE ( mask(mid,dim,count+1) >= 0.0_wp )
1464                      count = count + 1
1465                   ENDDO
1466                   WRITE ( io, 349 )  dir(dim), dir(dim), mid, dir(dim), &
1467                                      mask(mid,dim,:count)
1468                ELSEIF ( mask_loop(mid,dim,1) < 0.0_wp .AND.  &
1469                         mask_loop(mid,dim,2) < 0.0_wp .AND.  &
1470                         mask_loop(mid,dim,3) == 0.0_wp )  THEN
1471                   WRITE ( io, 350 )  dir(dim), dir(dim)
1472                ELSEIF ( mask_loop(mid,dim,3) == 0.0_wp )  THEN
1473                   WRITE ( io, 351 )  dir(dim), dir(dim), mid, dir(dim), &
1474                                      mask_loop(mid,dim,1:2)
1475                ELSE
1476                   WRITE ( io, 351 )  dir(dim), dir(dim), mid, dir(dim), &
1477                                      mask_loop(mid,dim,1:3)
1478                ENDIF
1479             ENDDO
1480          ENDIF
1481
1482       ENDDO
1483    ENDDO
1484
1485!
1486!-- Timeseries
1487    IF ( dt_dots /= 9999999.9_wp )  THEN
1488       WRITE ( io, 340 )
1489
1490       output_format = output_format_netcdf
1491       WRITE ( io, 344 )  output_format
1492       WRITE ( io, 341 )  dt_dots
1493    ENDIF
1494
1495#if defined( __dvrp_graphics )
1496!
1497!-- Dvrp-output
1498    IF ( dt_dvrp /= 9999999.9_wp )  THEN
1499       WRITE ( io, 360 )  dt_dvrp, TRIM( dvrp_output ), TRIM( dvrp_host ), &
1500                          TRIM( dvrp_username ), TRIM( dvrp_directory )
1501       i = 1
1502       l = 0
1503       m = 0
1504       DO WHILE ( mode_dvrp(i) /= ' ' )
1505          IF ( mode_dvrp(i)(1:10) == 'isosurface' )  THEN
1506             READ ( mode_dvrp(i), '(10X,I2)' )  j
1507             l = l + 1
1508             IF ( do3d(0,j) /= ' ' )  THEN
1509                WRITE ( io, 361 )  TRIM( do3d(0,j) ), threshold(l), &
1510                                   isosurface_color(:,l)
1511             ENDIF
1512          ELSEIF ( mode_dvrp(i)(1:6) == 'slicer' )  THEN
1513             READ ( mode_dvrp(i), '(6X,I2)' )  j
1514             m = m + 1
1515             IF ( do2d(0,j) /= ' ' )  THEN
1516                WRITE ( io, 362 )  TRIM( do2d(0,j) ), &
1517                                   slicer_range_limits_dvrp(:,m)
1518             ENDIF
1519          ELSEIF ( mode_dvrp(i)(1:9) == 'particles' )  THEN
1520             WRITE ( io, 363 )  dvrp_psize
1521             IF ( particle_dvrpsize /= 'none' )  THEN
1522                WRITE ( io, 364 )  'size', TRIM( particle_dvrpsize ), &
1523                                   dvrpsize_interval
1524             ENDIF
1525             IF ( particle_color /= 'none' )  THEN
1526                WRITE ( io, 364 )  'color', TRIM( particle_color ), &
1527                                   color_interval
1528             ENDIF
1529          ENDIF
1530          i = i + 1
1531       ENDDO
1532
1533       WRITE ( io, 365 )  groundplate_color, superelevation_x, &
1534                          superelevation_y, superelevation, clip_dvrp_l, &
1535                          clip_dvrp_r, clip_dvrp_s, clip_dvrp_n
1536
1537       IF ( TRIM( topography ) /= 'flat' )  THEN
1538          WRITE ( io, 366 )  topography_color
1539          IF ( cluster_size > 1 )  THEN
1540             WRITE ( io, 367 )  cluster_size
1541          ENDIF
1542       ENDIF
1543
1544    ENDIF
1545#endif
1546
1547#if defined( __spectra )
1548!
1549!-- Spectra output
1550    IF ( dt_dosp /= 9999999.9_wp )  THEN
1551       WRITE ( io, 370 )
1552
1553       output_format = output_format_netcdf
1554       WRITE ( io, 344 )  output_format
1555       WRITE ( io, 371 )  dt_dosp
1556       IF ( skip_time_dosp /= 0.0_wp )  WRITE ( io, 339 )  skip_time_dosp
1557       WRITE ( io, 372 )  ( data_output_sp(i), i = 1,10 ),     &
1558                          ( spectra_direction(i), i = 1,10 ),  &
1559                          ( comp_spectra_level(i), i = 1,100 ), &
1560                          ( plot_spectra_level(i), i = 1,100 ), &
1561                          averaging_interval_sp, dt_averaging_input_pr
1562    ENDIF
1563#endif
1564
1565    WRITE ( io, 99 )
1566
1567!
1568!-- Physical quantities
1569    WRITE ( io, 400 )
1570
1571!
1572!-- Geostrophic parameters
1573    WRITE ( io, 410 )  omega, phi, f, fs
1574
1575!
1576!-- Other quantities
1577    WRITE ( io, 411 )  g
1578    WRITE ( io, 412 )  TRIM( reference_state )
1579    IF ( use_single_reference_value )  THEN
1580       IF ( ocean )  THEN
1581          WRITE ( io, 413 )  prho_reference
1582       ELSE
1583          WRITE ( io, 414 )  pt_reference
1584       ENDIF
1585    ENDIF
1586
1587!
1588!-- Cloud physics parameters
1589    IF ( cloud_physics )  THEN
1590       WRITE ( io, 415 )
1591       WRITE ( io, 416 ) surface_pressure, r_d, rho_surface, cp, l_v
1592       IF ( icloud_scheme == 0 )  THEN
1593          WRITE ( io, 510 ) 1.0E-6_wp * nc_const
1594          IF ( precipitation )  WRITE ( io, 511 ) c_sedimentation
1595       ENDIF
1596    ENDIF
1597
1598!
1599!-- Cloud physcis parameters / quantities / numerical methods
1600    WRITE ( io, 430 )
1601    IF ( humidity .AND. .NOT. cloud_physics .AND. .NOT. cloud_droplets)  THEN
1602       WRITE ( io, 431 )
1603    ELSEIF ( humidity  .AND.  cloud_physics )  THEN
1604       WRITE ( io, 432 )
1605       IF ( cloud_top_radiation )  WRITE ( io, 132 )
1606       IF ( icloud_scheme == 1 )  THEN
1607          IF ( precipitation )  WRITE ( io, 133 )
1608       ELSEIF ( icloud_scheme == 0 )  THEN
1609          IF ( drizzle )  WRITE ( io, 506 )
1610          IF ( precipitation )  THEN
1611             WRITE ( io, 505 )
1612             IF ( turbulence )  WRITE ( io, 507 )
1613             IF ( ventilation_effect )  WRITE ( io, 508 )
1614             IF ( limiter_sedimentation )  WRITE ( io, 509 )
1615          ENDIF
1616       ENDIF
1617    ELSEIF ( humidity  .AND.  cloud_droplets )  THEN
1618       WRITE ( io, 433 )
1619       IF ( curvature_solution_effects )  WRITE ( io, 434 )
1620       IF ( collision_kernel /= 'none' )  THEN
1621          WRITE ( io, 435 )  TRIM( collision_kernel )
1622          IF ( collision_kernel(6:9) == 'fast' )  THEN
1623             WRITE ( io, 436 )  radius_classes, dissipation_classes
1624          ENDIF
1625       ELSE
1626          WRITE ( io, 437 )
1627       ENDIF
1628    ENDIF
1629
1630!
1631!-- LES / turbulence parameters
1632    WRITE ( io, 450 )
1633
1634!--
1635! ... LES-constants used must still be added here
1636!--
1637    IF ( constant_diffusion )  THEN
1638       WRITE ( io, 451 )  km_constant, km_constant/prandtl_number, &
1639                          prandtl_number
1640    ENDIF
1641    IF ( .NOT. constant_diffusion)  THEN
1642       IF ( e_init > 0.0_wp )  WRITE ( io, 455 )  e_init
1643       IF ( e_min > 0.0_wp )  WRITE ( io, 454 )  e_min
1644       IF ( wall_adjustment )  WRITE ( io, 453 )  wall_adjustment_factor
1645    ENDIF
1646
1647!
1648!-- Special actions during the run
1649    WRITE ( io, 470 )
1650    IF ( create_disturbances )  THEN
1651       WRITE ( io, 471 )  dt_disturb, disturbance_amplitude,                   &
1652                          zu(disturbance_level_ind_b), disturbance_level_ind_b,&
1653                          zu(disturbance_level_ind_t), disturbance_level_ind_t
1654       IF ( .NOT. bc_lr_cyc  .OR.  .NOT. bc_ns_cyc )  THEN
1655          WRITE ( io, 472 )  inflow_disturbance_begin, inflow_disturbance_end
1656       ELSE
1657          WRITE ( io, 473 )  disturbance_energy_limit
1658       ENDIF
1659       WRITE ( io, 474 )  TRIM( random_generator )
1660    ENDIF
1661    IF ( pt_surface_initial_change /= 0.0_wp )  THEN
1662       WRITE ( io, 475 )  pt_surface_initial_change
1663    ENDIF
1664    IF ( humidity  .AND.  q_surface_initial_change /= 0.0_wp )  THEN
1665       WRITE ( io, 476 )  q_surface_initial_change       
1666    ENDIF
1667    IF ( passive_scalar  .AND.  q_surface_initial_change /= 0.0_wp )  THEN
1668       WRITE ( io, 477 )  q_surface_initial_change       
1669    ENDIF
1670
1671    IF ( particle_advection )  THEN
1672!
1673!--    Particle attributes
1674       WRITE ( io, 480 )  particle_advection_start, dt_prel, bc_par_lr, &
1675                          bc_par_ns, bc_par_b, bc_par_t, particle_maximum_age, &
1676                          end_time_prel
1677       IF ( use_sgs_for_particles )  WRITE ( io, 488 )  dt_min_part
1678       IF ( random_start_position )  WRITE ( io, 481 )
1679       IF ( particles_per_point > 1 )  WRITE ( io, 489 )  particles_per_point
1680       WRITE ( io, 495 )  total_number_of_particles
1681       IF ( use_particle_tails  .AND.  maximum_number_of_tailpoints /= 0 )  THEN
1682          WRITE ( io, 483 )  maximum_number_of_tailpoints
1683          IF ( minimum_tailpoint_distance /= 0 )  THEN
1684             WRITE ( io, 484 )  total_number_of_tails,      &
1685                                minimum_tailpoint_distance, &
1686                                maximum_tailpoint_age
1687          ENDIF
1688       ENDIF
1689       IF ( dt_write_particle_data /= 9999999.9_wp )  THEN
1690          WRITE ( io, 485 )  dt_write_particle_data
1691          IF ( netcdf_data_format > 1 )  THEN
1692             output_format = 'netcdf (64 bit offset) and binary'
1693          ELSE
1694             output_format = 'netcdf and binary'
1695          ENDIF
1696          WRITE ( io, 344 )  output_format
1697       ENDIF
1698       IF ( dt_dopts /= 9999999.9_wp )  WRITE ( io, 494 )  dt_dopts
1699       IF ( write_particle_statistics )  WRITE ( io, 486 )
1700
1701       WRITE ( io, 487 )  number_of_particle_groups
1702
1703       DO  i = 1, number_of_particle_groups
1704          IF ( i == 1  .AND.  density_ratio(i) == 9999999.9_wp )  THEN
1705             WRITE ( io, 490 )  i, 0.0_wp
1706             WRITE ( io, 492 )
1707          ELSE
1708             WRITE ( io, 490 )  i, radius(i)
1709             IF ( density_ratio(i) /= 0.0_wp )  THEN
1710                WRITE ( io, 491 )  density_ratio(i)
1711             ELSE
1712                WRITE ( io, 492 )
1713             ENDIF
1714          ENDIF
1715          WRITE ( io, 493 )  psl(i), psr(i), pss(i), psn(i), psb(i), pst(i), &
1716                             pdx(i), pdy(i), pdz(i)
1717          IF ( .NOT. vertical_particle_advection(i) )  WRITE ( io, 482 )
1718       ENDDO
1719
1720    ENDIF
1721
1722
1723!
1724!-- Parameters of 1D-model
1725    IF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
1726       WRITE ( io, 500 )  end_time_1d, dt_run_control_1d, dt_pr_1d, &
1727                          mixing_length_1d, dissipation_1d
1728       IF ( damp_level_ind_1d /= nzt+1 )  THEN
1729          WRITE ( io, 502 )  zu(damp_level_ind_1d), damp_level_ind_1d
1730       ENDIF
1731    ENDIF
1732
1733!
1734!-- User-defined informations
1735    CALL user_header( io )
1736
1737    WRITE ( io, 99 )
1738
1739!
1740!-- Write buffer contents to disc immediately
1741    CALL local_flush( io )
1742
1743!
1744!-- Here the FORMATs start
1745
1746 99 FORMAT (1X,78('-'))
1747100 FORMAT (/1X,'******************************',4X,44('-')/        &
1748            1X,'* ',A,' *',4X,A/                               &
1749            1X,'******************************',4X,44('-'))
1750101 FORMAT (35X,'coupled run using MPI-',I1,': ',A/ &
1751            35X,42('-'))
1752102 FORMAT (/' Date:                 ',A8,4X,'Run:       ',A20/      &
1753            ' Time:                 ',A8,4X,'Run-No.:   ',I2.2/     &
1754            ' Run on host:        ',A10)
1755#if defined( __parallel )
1756103 FORMAT (' Number of PEs:',10X,I6,4X,'Processor grid (x,y): (',I4,',',I4, &
1757              ')',1X,A)
1758104 FORMAT (' Number of PEs:',10X,I6,4X,'Tasks:',I4,'   threads per task:',I4/ &
1759              35X,'Processor grid (x,y): (',I4,',',I4,')',1X,A)
1760105 FORMAT (35X,'One additional PE is used to handle'/37X,'the dvrp output!')
1761106 FORMAT (35X,'A 1d-decomposition along x is forced'/ &
1762            35X,'because the job is running on an SMP-cluster')
1763107 FORMAT (35X,'A 1d-decomposition along ',A,' is used')
1764108 FORMAT (35X,'Max. # of parallel I/O streams is ',I5)
1765109 FORMAT (35X,'Precursor run for coupled atmos-ocean run'/ &
1766            35X,42('-'))
1767114 FORMAT (35X,'Coupled atmosphere-ocean run following'/ &
1768            35X,'independent precursor runs'/             &
1769            35X,42('-'))
1770117 FORMAT (' Accelerator boards / node:  ',I2)
1771#endif
1772110 FORMAT (/' Numerical Schemes:'/ &
1773             ' -----------------'/)
1774111 FORMAT (' --> Solve perturbation pressure via FFT using ',A,' routines')
1775112 FORMAT (' --> Solve perturbation pressure via SOR-Red/Black-Schema'/ &
1776            '     Iterations (initial/other): ',I3,'/',I3,'  omega = ',F5.3)
1777113 FORMAT (' --> Momentum advection via Piascek-Williams-Scheme (Form C3)', &
1778                  ' or Upstream')
1779115 FORMAT ('     FFT and transpositions are overlapping')
1780116 FORMAT (' --> Scalar advection via Piascek-Williams-Scheme (Form C3)', &
1781                  ' or Upstream')
1782118 FORMAT (' --> Scalar advection via Bott-Chlond-Scheme')
1783119 FORMAT (' --> Galilei-Transform applied to horizontal advection:'/ &
1784            '     translation velocity = ',A/ &
1785            '     distance advected ',A,':  ',F8.3,' km(x)  ',F8.3,' km(y)')
1786120 FORMAT (' Accelerator boards: ',8X,I2)
1787122 FORMAT (' --> Time differencing scheme: ',A)
1788123 FORMAT (' --> Rayleigh-Damping active, starts ',A,' z = ',F8.2,' m'/ &
1789            '     maximum damping coefficient: ',F5.3, ' 1/s')
1790129 FORMAT (' --> Additional prognostic equation for the specific humidity')
1791130 FORMAT (' --> Additional prognostic equation for the total water content')
1792131 FORMAT (' --> No pt-equation solved. Neutral stratification with pt = ', &
1793                  F6.2, ' K assumed')
1794132 FORMAT ('     Parameterization of long-wave radiation processes via'/ &
1795            '     effective emissivity scheme')
1796133 FORMAT ('     Precipitation parameterization via Kessler-Scheme')
1797134 FORMAT (' --> Additional prognostic equation for a passive scalar')
1798135 FORMAT (' --> Solve perturbation pressure via multigrid method (', &
1799                  A,'-cycle)'/ &
1800            '     number of grid levels:                   ',I2/ &
1801            '     Gauss-Seidel red/black iterations:       ',I2)
1802136 FORMAT ('     gridpoints of coarsest subdomain (x,y,z): (',I3,',',I3,',', &
1803                  I3,')')
1804137 FORMAT ('     level data gathered on PE0 at level:     ',I2/ &
1805            '     gridpoints of coarsest subdomain (x,y,z): (',I3,',',I3,',', &
1806                  I3,')'/ &
1807            '     gridpoints of coarsest domain (x,y,z):    (',I3,',',I3,',', &
1808                  I3,')')
1809139 FORMAT (' --> Loop optimization method: ',A)
1810140 FORMAT ('     maximum residual allowed:                ',E10.3)
1811141 FORMAT ('     fixed number of multigrid cycles:        ',I4)
1812142 FORMAT ('     perturbation pressure is calculated at every Runge-Kutta ', &
1813                  'step')
1814143 FORMAT ('     Euler/upstream scheme is used for the SGS turbulent ', &
1815                  'kinetic energy')
1816144 FORMAT ('     masking method is used')
1817150 FORMAT (' --> Volume flow at the right and north boundary will be ', &
1818                  'conserved'/ &
1819            '     using the ',A,' mode')
1820151 FORMAT ('     with u_bulk = ',F7.3,' m/s and v_bulk = ',F7.3,' m/s')
1821152 FORMAT (' --> External pressure gradient directly prescribed by the user:',&
1822           /'     ',2(1X,E12.5),'Pa/m in x/y direction', &
1823           /'     starting from dp_level_b =', F8.3, 'm', A /)
1824160 FORMAT (//' Large scale forcing and nudging:'/ &
1825              ' -------------------------------'/)
1826161 FORMAT (' --> No large scale forcing from external is used (default) ')
1827162 FORMAT (' --> Large scale forcing from external file LSF_DATA is used: ')
1828163 FORMAT ('     - large scale advection tendencies ')
1829164 FORMAT ('     - large scale subsidence velocity w_subs ')
1830165 FORMAT ('     - large scale subsidence tendencies ')
1831167 FORMAT ('     - and geostrophic wind components ug and vg')
1832168 FORMAT (' --> Large-scale vertical motion is used in the ', &
1833                  'prognostic equation(s) for')
1834169 FORMAT ('     the scalar(s) only')
1835170 FORMAT (' --> Nudging is used')
1836171 FORMAT (' --> No nudging is used (default) ')
1837180 FORMAT ('     - prescribed surface values for temperature')
1838181 FORMAT ('     - prescribed surface fluxes for temperature')
1839182 FORMAT ('     - prescribed surface values for humidity')
1840183 FORMAT ('     - prescribed surface fluxes for humidity')
1841200 FORMAT (//' Run time and time step information:'/ &
1842             ' ----------------------------------'/)
1843201 FORMAT ( ' Timestep:             variable     maximum value: ',F6.3,' s', &
1844             '    CFL-factor: ',F4.2)
1845202 FORMAT ( ' Timestep:          dt = ',F6.3,' s'/)
1846203 FORMAT ( ' Start time:          ',F9.3,' s'/ &
1847             ' End time:            ',F9.3,' s')
1848204 FORMAT ( A,F9.3,' s')
1849205 FORMAT ( A,F9.3,' s',5X,'restart every',17X,F9.3,' s')
1850206 FORMAT (/' Time reached:        ',F9.3,' s'/ &
1851             ' CPU-time used:       ',F9.3,' s     per timestep:               ', &
1852               '  ',F9.3,' s'/                                                    &
1853             '                                      per second of simulated tim', &
1854               'e: ',F9.3,' s')
1855207 FORMAT ( ' Coupling start time: ',F9.3,' s')
1856250 FORMAT (//' Computational grid and domain size:'/ &
1857              ' ----------------------------------'// &
1858              ' Grid length:      dx =    ',F7.3,' m    dy =    ',F7.3, &
1859              ' m    dz =    ',F7.3,' m'/ &
1860              ' Domain size:       x = ',F10.3,' m     y = ',F10.3, &
1861              ' m  z(u) = ',F10.3,' m'/)
1862252 FORMAT (' dz constant up to ',F10.3,' m (k=',I4,'), then stretched by', &
1863              ' factor: ',F5.3/ &
1864            ' maximum dz not to be exceeded is dz_max = ',F10.3,' m'/)
1865254 FORMAT (' Number of gridpoints (x,y,z):  (0:',I4,', 0:',I4,', 0:',I4,')'/ &
1866            ' Subdomain size (x,y,z):        (  ',I4,',   ',I4,',   ',I4,')'/)
1867260 FORMAT (/' The model has a slope in x-direction. Inclination angle: ',F6.2,&
1868             ' degrees')
1869270 FORMAT (//' Topography informations:'/ &
1870              ' -----------------------'// &
1871              1X,'Topography: ',A)
1872271 FORMAT (  ' Building size (x/y/z) in m: ',F5.1,' / ',F5.1,' / ',F5.1/ &
1873              ' Horizontal index bounds (l/r/s/n): ',I4,' / ',I4,' / ',I4, &
1874                ' / ',I4)
1875272 FORMAT (  ' Single quasi-2D street canyon of infinite length in ',A, &
1876              ' direction' / &
1877              ' Canyon height: ', F6.2, 'm, ch = ', I4, '.'      / &
1878              ' Canyon position (',A,'-walls): cxl = ', I4,', cxr = ', I4, '.')
1879278 FORMAT (' Topography grid definition convention:'/ &
1880            ' cell edge (staggered grid points'/  &
1881            ' (u in x-direction, v in y-direction))' /)
1882279 FORMAT (' Topography grid definition convention:'/ &
1883            ' cell center (scalar grid points)' /)
1884280 FORMAT (//' Vegetation canopy (drag) model:'/ &
1885              ' ------------------------------'// &
1886              ' Canopy mode: ', A / &
1887              ' Canopy height: ',F6.2,'m (',I4,' grid points)' / &
1888              ' Leaf drag coefficient: ',F6.2 /)
1889281 FORMAT (/ ' Scalar exchange coefficient: ',F6.2 / &
1890              ' Scalar concentration at leaf surfaces in kg/m**3: ',F6.2 /)
1891282 FORMAT (' Predefined constant heatflux at the top of the vegetation: ',F6.2,' K m/s')
1892283 FORMAT (/ ' Characteristic levels of the leaf area density:'// &
1893              ' Height:              ',A,'  m'/ &
1894              ' Leaf area density:   ',A,'  m**2/m**3'/ &
1895              ' Gradient:            ',A,'  m**2/m**4'/ &
1896              ' Gridpoint:           ',A)
1897284 FORMAT (//' Characteristic levels of the leaf area density and coefficients:'// &
1898              ' Height:              ',A,'  m'/ &
1899              ' Leaf area density:   ',A,'  m**2/m**3'/ &
1900              ' Coefficient alpha: ',F6.2 / &
1901              ' Coefficient beta: ',F6.2 / &
1902              ' Leaf area index: ',F6.2,'  m**2/m**2' /)
1903               
1904300 FORMAT (//' Boundary conditions:'/ &
1905             ' -------------------'// &
1906             '                     p                    uv             ', &
1907             '                   pt'// &
1908             ' B. bound.: ',A/ &
1909             ' T. bound.: ',A)
1910301 FORMAT (/'                     ',A// &
1911             ' B. bound.: ',A/ &
1912             ' T. bound.: ',A)
1913303 FORMAT (/' Bottom surface fluxes are used in diffusion terms at k=1')
1914304 FORMAT (/' Top surface fluxes are used in diffusion terms at k=nzt')
1915305 FORMAT (//'    Prandtl-Layer between bottom surface and first ', &
1916               'computational u,v-level:'// &
1917             '       zp = ',F6.2,' m   z0 = ',F6.4,' m   z0h = ',F7.5,&
1918             ' m   kappa = ',F4.2/ &
1919             '       Rif value range:   ',F6.2,' <= rif <=',F6.2)
1920306 FORMAT ('       Predefined constant heatflux:   ',F9.6,' K m/s')
1921307 FORMAT ('       Heatflux has a random normal distribution')
1922308 FORMAT ('       Predefined surface temperature')
1923309 FORMAT ('       Predefined constant salinityflux:   ',F9.6,' psu m/s')
1924310 FORMAT (//'    1D-Model:'// &
1925             '       Rif value range:   ',F6.2,' <= rif <=',F6.2)
1926311 FORMAT ('       Predefined constant humidity flux: ',E10.3,' m/s')
1927312 FORMAT ('       Predefined surface humidity')
1928313 FORMAT ('       Predefined constant scalar flux: ',E10.3,' kg/(m**2 s)')
1929314 FORMAT ('       Predefined scalar value at the surface')
1930315 FORMAT ('       Humidity / scalar flux at top surface is 0.0')
1931316 FORMAT ('       Sensible heatflux and momentum flux from coupled ', &
1932                    'atmosphere model')
1933317 FORMAT (//' Lateral boundaries:'/ &
1934            '       left/right:  ',A/    &
1935            '       north/south: ',A)
1936318 FORMAT (/'       use_cmax: ',L1 / &
1937            '       pt damping layer width = ',F8.2,' m, pt ', &
1938                    'damping factor = ',F6.4)
1939319 FORMAT ('       turbulence recycling at inflow switched on'/ &
1940            '       width of recycling domain: ',F7.1,' m   grid index: ',I4/ &
1941            '       inflow damping height: ',F6.1,' m   width: ',F6.1,' m')
1942320 FORMAT ('       Predefined constant momentumflux:  u: ',F9.6,' m**2/s**2'/ &
1943            '                                          v: ',F9.6,' m**2/s**2')
1944321 FORMAT (//' Initial profiles:'/ &
1945              ' ----------------')
1946325 FORMAT (//' List output:'/ &
1947             ' -----------'//  &
1948            '    1D-Profiles:'/    &
1949            '       Output every             ',F8.2,' s')
1950326 FORMAT ('       Time averaged over       ',F8.2,' s'/ &
1951            '       Averaging input every    ',F8.2,' s')
1952330 FORMAT (//' Data output:'/ &
1953             ' -----------'/)
1954331 FORMAT (/'    1D-Profiles:')
1955332 FORMAT (/'       ',A)
1956333 FORMAT ('       Output every             ',F8.2,' s',/ &
1957            '       Time averaged over       ',F8.2,' s'/ &
1958            '       Averaging input every    ',F8.2,' s')
1959334 FORMAT (/'    2D-Arrays',A,':')
1960335 FORMAT (/'       ',A2,'-cross-section  Arrays: ',A/ &
1961            '       Output every             ',F8.2,' s  ',A/ &
1962            '       Cross sections at ',A1,' = ',A/ &
1963            '       scalar-coordinates:   ',A,' m'/)
1964336 FORMAT (/'    3D-Arrays',A,':')
1965337 FORMAT (/'       Arrays: ',A/ &
1966            '       Output every             ',F8.2,' s  ',A/ &
1967            '       Upper output limit at    ',F8.2,' m  (GP ',I4,')'/)
1968339 FORMAT ('       No output during initial ',F8.2,' s')
1969340 FORMAT (/'    Time series:')
1970341 FORMAT ('       Output every             ',F8.2,' s'/)
1971342 FORMAT (/'       ',A2,'-cross-section  Arrays: ',A/ &
1972            '       Output every             ',F8.2,' s  ',A/ &
1973            '       Time averaged over       ',F8.2,' s'/ &
1974            '       Averaging input every    ',F8.2,' s'/ &
1975            '       Cross sections at ',A1,' = ',A/ &
1976            '       scalar-coordinates:   ',A,' m'/)
1977343 FORMAT (/'       Arrays: ',A/ &
1978            '       Output every             ',F8.2,' s  ',A/ &
1979            '       Time averaged over       ',F8.2,' s'/ &
1980            '       Averaging input every    ',F8.2,' s'/ &
1981            '       Upper output limit at    ',F8.2,' m  (GP ',I4,')'/)
1982344 FORMAT ('       Output format: ',A/)
1983345 FORMAT (/'    Scaling lengths for output locations of all subsequent mask IDs:',/ &
1984            '       mask_scale_x (in x-direction): ',F9.3, ' m',/ &
1985            '       mask_scale_y (in y-direction): ',F9.3, ' m',/ &
1986            '       mask_scale_z (in z-direction): ',F9.3, ' m' )
1987346 FORMAT (/'    Masked data output',A,' for mask ID ',I2, ':')
1988347 FORMAT ('       Variables: ',A/ &
1989            '       Output every             ',F8.2,' s')
1990348 FORMAT ('       Variables: ',A/ &
1991            '       Output every             ',F8.2,' s'/ &
1992            '       Time averaged over       ',F8.2,' s'/ &
1993            '       Averaging input every    ',F8.2,' s')
1994349 FORMAT (/'       Output locations in ',A,'-direction in multiples of ', &
1995            'mask_scale_',A,' predefined by array mask_',I2.2,'_',A,':'/ &
1996            13('       ',8(F8.2,',')/) )
1997350 FORMAT (/'       Output locations in ',A,'-direction: ', &
1998            'all gridpoints along ',A,'-direction (default).' )
1999351 FORMAT (/'       Output locations in ',A,'-direction in multiples of ', &
2000            'mask_scale_',A,' constructed from array mask_',I2.2,'_',A,'_loop:'/ &
2001            '          loop begin:',F8.2,', end:',F8.2,', stride:',F8.2 )
2002352 FORMAT  (/'       Number of output time levels allowed: ',I3 /)
2003353 FORMAT  (/'       Number of output time levels allowed: unlimited' /)
2004#if defined( __dvrp_graphics )
2005360 FORMAT ('    Plot-Sequence with dvrp-software:'/ &
2006            '       Output every      ',F7.1,' s'/ &
2007            '       Output mode:      ',A/ &
2008            '       Host / User:      ',A,' / ',A/ &
2009            '       Directory:        ',A// &
2010            '       The sequence contains:')
2011361 FORMAT (/'       Isosurface of "',A,'"    Threshold value: ', E12.3/ &
2012            '          Isosurface color: (',F4.2,',',F4.2,',',F4.2,') (R,G,B)')
2013362 FORMAT (/'       Slicer plane ',A/ &
2014            '       Slicer limits: [',F6.2,',',F6.2,']')
2015363 FORMAT (/'       Particles'/ &
2016            '          particle size:  ',F7.2,' m')
2017364 FORMAT ('          particle ',A,' controlled by "',A,'" with interval [', &
2018                       F6.2,',',F6.2,']')
2019365 FORMAT (/'       Groundplate color: (',F4.2,',',F4.2,',',F4.2,') (R,G,B)'/ &
2020            '       Superelevation along (x,y,z): (',F4.1,',',F4.1,',',F4.1, &
2021                     ')'/ &
2022            '       Clipping limits: from x = ',F9.1,' m to x = ',F9.1,' m'/ &
2023            '                        from y = ',F9.1,' m to y = ',F9.1,' m')
2024366 FORMAT (/'       Topography color: (',F4.2,',',F4.2,',',F4.2,') (R,G,B)')
2025367 FORMAT ('       Polygon reduction for topography: cluster_size = ', I1)
2026#endif
2027#if defined( __spectra )
2028370 FORMAT ('    Spectra:')
2029371 FORMAT ('       Output every ',F7.1,' s'/)
2030372 FORMAT ('       Arrays:     ', 10(A5,',')/                         &
2031            '       Directions: ', 10(A5,',')/                         &
2032            '       height levels  k = ', 20(I3,',')/                  &
2033            '                          ', 20(I3,',')/                  &
2034            '                          ', 20(I3,',')/                  &
2035            '                          ', 20(I3,',')/                  &
2036            '                          ', 19(I3,','),I3,'.'/           &
2037            '       height levels selected for standard plot:'/        &
2038            '                      k = ', 20(I3,',')/                  &
2039            '                          ', 20(I3,',')/                  &
2040            '                          ', 20(I3,',')/                  &
2041            '                          ', 20(I3,',')/                  &
2042            '                          ', 19(I3,','),I3,'.'/           &
2043            '       Time averaged over ', F7.1, ' s,' /                &
2044            '       Profiles for the time averaging are taken every ', &
2045                    F6.1,' s')
2046#endif
2047400 FORMAT (//' Physical quantities:'/ &
2048              ' -------------------'/)
2049410 FORMAT ('    Angular velocity    :   omega = ',E9.3,' rad/s'/  &
2050            '    Geograph. latitude  :   phi   = ',F4.1,' degr'/   &
2051            '    Coriolis parameter  :   f     = ',F9.6,' 1/s'/    &
2052            '                            f*    = ',F9.6,' 1/s')
2053411 FORMAT (/'    Gravity             :   g     = ',F4.1,' m/s**2')
2054412 FORMAT (/'    Reference state used in buoyancy terms: ',A)
2055413 FORMAT ('       Reference density in buoyancy terms: ',F8.3,' kg/m**3')
2056414 FORMAT ('       Reference temperature in buoyancy terms: ',F8.4,' K')
2057415 FORMAT (/'    Cloud physics parameters:'/ &
2058             '    ------------------------'/)
2059416 FORMAT ('        Surface pressure   :   p_0   = ',F7.2,' hPa'/      &
2060            '        Gas constant       :   R     = ',F5.1,' J/(kg K)'/ &
2061            '        Density of air     :   rho_0 = ',F5.3,' kg/m**3'/  &
2062            '        Specific heat cap. :   c_p   = ',F6.1,' J/(kg K)'/ &
2063            '        Vapourization heat :   L_v   = ',E8.2,' J/kg')
2064420 FORMAT (/'    Characteristic levels of the initial temperature profile:'// &
2065            '       Height:        ',A,'  m'/ &
2066            '       Temperature:   ',A,'  K'/ &
2067            '       Gradient:      ',A,'  K/100m'/ &
2068            '       Gridpoint:     ',A)
2069421 FORMAT (/'    Characteristic levels of the initial humidity profile:'// &
2070            '       Height:      ',A,'  m'/ &
2071            '       Humidity:    ',A,'  kg/kg'/ &
2072            '       Gradient:    ',A,'  (kg/kg)/100m'/ &
2073            '       Gridpoint:   ',A)
2074422 FORMAT (/'    Characteristic levels of the initial scalar profile:'// &
2075            '       Height:                  ',A,'  m'/ &
2076            '       Scalar concentration:    ',A,'  kg/m**3'/ &
2077            '       Gradient:                ',A,'  (kg/m**3)/100m'/ &
2078            '       Gridpoint:               ',A)
2079423 FORMAT (/'    Characteristic levels of the geo. wind component ug:'// &
2080            '       Height:      ',A,'  m'/ &
2081            '       ug:          ',A,'  m/s'/ &
2082            '       Gradient:    ',A,'  1/100s'/ &
2083            '       Gridpoint:   ',A)
2084424 FORMAT (/'    Characteristic levels of the geo. wind component vg:'// &
2085            '       Height:      ',A,'  m'/ &
2086            '       vg:          ',A,'  m/s'/ &
2087            '       Gradient:    ',A,'  1/100s'/ &
2088            '       Gridpoint:   ',A)
2089425 FORMAT (/'    Characteristic levels of the initial salinity profile:'// &
2090            '       Height:     ',A,'  m'/ &
2091            '       Salinity:   ',A,'  psu'/ &
2092            '       Gradient:   ',A,'  psu/100m'/ &
2093            '       Gridpoint:  ',A)
2094426 FORMAT (/'    Characteristic levels of the subsidence/ascent profile:'// &
2095            '       Height:      ',A,'  m'/ &
2096            '       w_subs:      ',A,'  m/s'/ &
2097            '       Gradient:    ',A,'  (m/s)/100m'/ &
2098            '       Gridpoint:   ',A)
2099427 FORMAT (/'    Initial wind profiles (u,v) are interpolated from given'// &
2100                  ' profiles')
2101428 FORMAT (/'    Initial profiles (u, v, pt, q) are taken from file '/ &
2102             '    NUDGING_DATA')
2103430 FORMAT (//' Cloud physics quantities / methods:'/ &
2104              ' ----------------------------------'/)
2105431 FORMAT ('    Humidity is treated as purely passive scalar (no condensati', &
2106                 'on)')
2107432 FORMAT ('    Bulk scheme with liquid water potential temperature and'/ &
2108            '    total water content is used.'/ &
2109            '    Condensation is parameterized via 0% - or 100% scheme.')
2110433 FORMAT ('    Cloud droplets treated explicitly using the Lagrangian part', &
2111                 'icle model')
2112434 FORMAT ('    Curvature and solution effecs are considered for growth of', &
2113                 ' droplets < 1.0E-6 m')
2114435 FORMAT ('    Droplet collision is handled by ',A,'-kernel')
2115436 FORMAT ('       Fast kernel with fixed radius- and dissipation classes ', &
2116                    'are used'/ &
2117            '          number of radius classes:       ',I3,'    interval ', &
2118                       '[1.0E-6,2.0E-4] m'/ &
2119            '          number of dissipation classes:   ',I2,'    interval ', &
2120                       '[0,1000] cm**2/s**3')
2121437 FORMAT ('    Droplet collision is switched off')
2122450 FORMAT (//' LES / Turbulence quantities:'/ &
2123              ' ---------------------------'/)
2124451 FORMAT ('    Diffusion coefficients are constant:'/ &
2125            '    Km = ',F6.2,' m**2/s   Kh = ',F6.2,' m**2/s   Pr = ',F5.2)
2126453 FORMAT ('    Mixing length is limited to ',F4.2,' * z')
2127454 FORMAT ('    TKE is not allowed to fall below ',E9.2,' (m/s)**2')
2128455 FORMAT ('    initial TKE is prescribed as ',E9.2,' (m/s)**2')
2129470 FORMAT (//' Actions during the simulation:'/ &
2130              ' -----------------------------'/)
2131471 FORMAT ('    Disturbance impulse (u,v) every :   ',F6.2,' s'/            &
2132            '    Disturbance amplitude           :     ',F4.2, ' m/s'/       &
2133            '    Lower disturbance level         : ',F8.2,' m (GP ',I4,')'/  &
2134            '    Upper disturbance level         : ',F8.2,' m (GP ',I4,')')
2135472 FORMAT ('    Disturbances continued during the run from i/j =',I4, &
2136                 ' to i/j =',I4)
2137473 FORMAT ('    Disturbances cease as soon as the disturbance energy exceeds',&
2138                 1X,F5.3, ' m**2/s**2')
2139474 FORMAT ('    Random number generator used    : ',A/)
2140475 FORMAT ('    The surface temperature is increased (or decreased, ', &
2141                 'respectively, if'/ &
2142            '    the value is negative) by ',F5.2,' K at the beginning of the',&
2143                 ' 3D-simulation'/)
2144476 FORMAT ('    The surface humidity is increased (or decreased, ',&
2145                 'respectively, if the'/ &
2146            '    value is negative) by ',E8.1,' kg/kg at the beginning of', &
2147                 ' the 3D-simulation'/)
2148477 FORMAT ('    The scalar value is increased at the surface (or decreased, ',&
2149                 'respectively, if the'/ &
2150            '    value is negative) by ',E8.1,' kg/m**3 at the beginning of', &
2151                 ' the 3D-simulation'/)
2152480 FORMAT ('    Particles:'/ &
2153            '    ---------'// &
2154            '       Particle advection is active (switched on at t = ', F7.1, &
2155                    ' s)'/ &
2156            '       Start of new particle generations every  ',F6.1,' s'/ &
2157            '       Boundary conditions: left/right: ', A, ' north/south: ', A/&
2158            '                            bottom:     ', A, ' top:         ', A/&
2159            '       Maximum particle age:                 ',F9.1,' s'/ &
2160            '       Advection stopped at t = ',F9.1,' s'/)
2161481 FORMAT ('       Particles have random start positions'/)
2162482 FORMAT ('          Particles are advected only horizontally'/)
2163483 FORMAT ('       Particles have tails with a maximum of ',I3,' points')
2164484 FORMAT ('            Number of tails of the total domain: ',I10/ &
2165            '            Minimum distance between tailpoints: ',F8.2,' m'/ &
2166            '            Maximum age of the end of the tail:  ',F8.2,' s')
2167485 FORMAT ('       Particle data are written on file every ', F9.1, ' s')
2168486 FORMAT ('       Particle statistics are written on file'/)
2169487 FORMAT ('       Number of particle groups: ',I2/)
2170488 FORMAT ('       SGS velocity components are used for particle advection'/ &
2171            '          minimum timestep for advection: ', F7.5/)
2172489 FORMAT ('       Number of particles simultaneously released at each ', &
2173                    'point: ', I5/)
2174490 FORMAT ('       Particle group ',I2,':'/ &
2175            '          Particle radius: ',E10.3, 'm')
2176491 FORMAT ('          Particle inertia is activated'/ &
2177            '             density_ratio (rho_fluid/rho_particle) = ',F5.3/)
2178492 FORMAT ('          Particles are advected only passively (no inertia)'/)
2179493 FORMAT ('          Boundaries of particle source: x:',F8.1,' - ',F8.1,' m'/&
2180            '                                         y:',F8.1,' - ',F8.1,' m'/&
2181            '                                         z:',F8.1,' - ',F8.1,' m'/&
2182            '          Particle distances:  dx = ',F8.1,' m  dy = ',F8.1, &
2183                       ' m  dz = ',F8.1,' m'/)
2184494 FORMAT ('       Output of particle time series in NetCDF format every ', &
2185                    F8.2,' s'/)
2186495 FORMAT ('       Number of particles in total domain: ',I10/)
2187500 FORMAT (//' 1D-Model parameters:'/                           &
2188              ' -------------------'//                           &
2189            '    Simulation time:                   ',F8.1,' s'/ &
2190            '    Run-controll output every:         ',F8.1,' s'/ &
2191            '    Vertical profile output every:     ',F8.1,' s'/ &
2192            '    Mixing length calculation:         ',A/         &
2193            '    Dissipation calculation:           ',A/)
2194502 FORMAT ('    Damping layer starts from ',F7.1,' m (GP ',I4,')'/)
2195503 FORMAT (' --> Momentum advection via Wicker-Skamarock-Scheme 5th order')
2196504 FORMAT (' --> Scalar advection via Wicker-Skamarock-Scheme 5th order')
2197505 FORMAT ('    Precipitation parameterization via Seifert-Beheng-Scheme')
2198506 FORMAT ('    Drizzle parameterization via Stokes law')
2199507 FORMAT ('    Turbulence effects on precipitation process')
2200508 FORMAT ('    Ventilation effects on evaporation of rain drops')
2201509 FORMAT ('    Slope limiter used for sedimentation process')
2202510 FORMAT ('        Droplet density    :   N_c   = ',F6.1,' 1/cm**3')
2203511 FORMAT ('        Sedimentation Courant number:                  '/&
2204            '                               C_s   = ',F3.1,'        ')
2205512 FORMAT (/' Date:                 ',A8,6X,'Run:       ',A20/      &
2206            ' Time:                 ',A8,6X,'Run-No.:   ',I2.2/     &
2207            ' Run on host:        ',A10,6X,'En-No.:    ',I2.2)
2208
2209 END SUBROUTINE header
Note: See TracBrowser for help on using the repository browser.