source: palm/trunk/SOURCE/header.f90 @ 1856

Last change on this file since 1856 was 1852, checked in by hoffmann, 8 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 86.0 KB
Line 
1!> @file header.f90
2!--------------------------------------------------------------------------------!
3! This file is part of PALM.
4!
5! PALM is free software: you can redistribute it and/or modify it under the terms
6! of the GNU General Public License as published by the Free Software Foundation,
7! either version 3 of the License, or (at your option) any later version.
8!
9! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
10! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
11! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
12!
13! You should have received a copy of the GNU General Public License along with
14! PALM. If not, see <http://www.gnu.org/licenses/>.
15!
16! Copyright 1997-2016 Leibniz Universitaet Hannover
17!--------------------------------------------------------------------------------!
18!
19! Current revisions:
20! -----------------
21!
22!
23! Former revisions:
24! -----------------
25! $Id: header.f90 1852 2016-04-08 14:07:36Z maronga $
26!
27! 1849 2016-04-08 11:33:18Z hoffmann
28! Adapted for modularization of microphysics
29!
30! 1833 2016-04-07 14:23:03Z raasch
31! spectrum renamed spectra_mod, output of spectra related quantities moved to
32! spectra_mod
33!
34! 1831 2016-04-07 13:15:51Z hoffmann
35! turbulence renamed collision_turbulence,
36! drizzle renamed cloud_water_sedimentation
37!
38! 1826 2016-04-07 12:01:39Z maronga
39! Moved radiation model header output to the respective module.
40! Moved canopy model header output to the respective module.
41!
42! 1822 2016-04-07 07:49:42Z hoffmann
43! Tails removed. icloud_scheme replaced by microphysics_*
44!
45! 1817 2016-04-06 15:44:20Z maronga
46! Moved land_surface_model header output to the respective module.
47!
48! 1808 2016-04-05 19:44:00Z raasch
49! routine local_flush replaced by FORTRAN statement
50!
51! 1797 2016-03-21 16:50:28Z raasch
52! output of nesting datatransfer mode
53!
54! 1791 2016-03-11 10:41:25Z raasch
55! output of nesting informations of all domains
56!
57! 1788 2016-03-10 11:01:04Z maronga
58! Parameter dewfall removed
59!
60! 1786 2016-03-08 05:49:27Z raasch
61! cpp-direktives for spectra removed
62!
63! 1783 2016-03-06 18:36:17Z raasch
64! netcdf module and variable names changed, output of netcdf_deflate
65!
66! 1764 2016-02-28 12:45:19Z raasch
67! output of nesting informations
68!
69! 1697 2015-10-28 17:14:10Z raasch
70! small E- and F-FORMAT changes to avoid informative compiler messages about
71! insufficient field width
72!
73! 1691 2015-10-26 16:17:44Z maronga
74! Renamed prandtl_layer to constant_flux_layer, renames rif_min/rif_max to
75! zeta_min/zeta_max.
76!
77! 1682 2015-10-07 23:56:08Z knoop
78! Code annotations made doxygen readable
79!
80! 1675 2015-10-02 08:28:59Z gronemeier
81! Bugfix: Definition of topography grid levels
82!
83! 1660 2015-09-21 08:15:16Z gronemeier
84! Bugfix: Definition of building/street canyon height if vertical grid stretching
85!         starts below the maximum topography height.
86!
87! 1590 2015-05-08 13:56:27Z maronga
88! Bugfix: Added TRIM statements for character strings for LSM and radiation code
89!
90! 1585 2015-04-30 07:05:52Z maronga
91! Further output for radiation model(s).
92!
93! 1575 2015-03-27 09:56:27Z raasch
94! adjustments for psolver-queries, output of seed_follows_topography
95!
96! 1560 2015-03-06 10:48:54Z keck
97! output for recycling y shift
98!
99! 1557 2015-03-05 16:43:04Z suehring
100! output for monotonic limiter
101!
102! 1551 2015-03-03 14:18:16Z maronga
103! Added informal output for land surface model and radiation model. Removed typo.
104!
105! 1496 2014-12-02 17:25:50Z maronga
106! Renamed: "radiation -> "cloud_top_radiation"
107!
108! 1484 2014-10-21 10:53:05Z kanani
109! Changes due to new module structure of the plant canopy model:
110!   module plant_canopy_model_mod and output for new canopy model parameters
111!   (alpha_lad, beta_lad, lai_beta,...) added,
112!   drag_coefficient, leaf_surface_concentration and scalar_exchange_coefficient
113!   renamed to canopy_drag_coeff, leaf_surface_conc and leaf_scalar_exch_coeff,
114!   learde renamed leaf_area_density.
115! Bugfix: DO-WHILE-loop for lad header information additionally restricted
116! by maximum number of gradient levels (currently 10)
117!
118! 1482 2014-10-18 12:34:45Z raasch
119! information about calculated or predefined virtual processor topology adjusted
120!
121! 1468 2014-09-24 14:06:57Z maronga
122! Adapted for use on up to 6-digit processor cores
123!
124! 1429 2014-07-15 12:53:45Z knoop
125! header exended to provide ensemble_member_nr if specified
126!
127! 1376 2014-04-26 11:21:22Z boeske
128! Correction of typos
129!
130! 1365 2014-04-22 15:03:56Z boeske
131! New section 'Large scale forcing and nudging':
132! output of large scale forcing and nudging information,
133! new section for initial profiles created
134!
135! 1359 2014-04-11 17:15:14Z hoffmann
136! dt_sort_particles removed
137!
138! 1353 2014-04-08 15:21:23Z heinze
139! REAL constants provided with KIND-attribute
140!
141! 1327 2014-03-21 11:00:16Z raasch
142! parts concerning iso2d and avs output removed,
143! -netcdf output queries
144!
145! 1324 2014-03-21 09:13:16Z suehring
146! Bugfix: module spectrum added
147!
148! 1322 2014-03-20 16:38:49Z raasch
149! REAL functions provided with KIND-attribute,
150! some REAL constants defined as wp-kind
151!
152! 1320 2014-03-20 08:40:49Z raasch
153! ONLY-attribute added to USE-statements,
154! kind-parameters added to all INTEGER and REAL declaration statements,
155! kinds are defined in new module kinds,
156! revision history before 2012 removed,
157! comment fields (!:) to be used for variable explanations added to
158! all variable declaration statements
159!
160! 1308 2014-03-13 14:58:42Z fricke
161! output of the fixed number of output time levels
162! output_format adjusted for masked data if netcdf_data_format > 5
163!
164! 1299 2014-03-06 13:15:21Z heinze
165! output for using large_scale subsidence in combination
166! with large_scale_forcing
167! reformatting, more detailed explanations
168!
169! 1241 2013-10-30 11:36:58Z heinze
170! output for nudging + large scale forcing from external file
171!
172! 1216 2013-08-26 09:31:42Z raasch
173! output for transpose_compute_overlap
174!
175! 1212 2013-08-15 08:46:27Z raasch
176! output for poisfft_hybrid removed
177!
178! 1179 2013-06-14 05:57:58Z raasch
179! output of reference_state, use_reference renamed use_single_reference_value
180!
181! 1159 2013-05-21 11:58:22Z fricke
182! +use_cmax
183!
184! 1115 2013-03-26 18:16:16Z hoffmann
185! descriptions for Seifert-Beheng-cloud-physics-scheme added
186!
187! 1111 2013-03-08 23:54:10Z raasch
188! output of accelerator board information
189! ibc_p_b = 2 removed
190!
191! 1108 2013-03-05 07:03:32Z raasch
192! bugfix for r1106
193!
194! 1106 2013-03-04 05:31:38Z raasch
195! some format changes for coupled runs
196!
197! 1092 2013-02-02 11:24:22Z raasch
198! unused variables removed
199!
200! 1036 2012-10-22 13:43:42Z raasch
201! code put under GPL (PALM 3.9)
202!
203! 1031 2012-10-19 14:35:30Z raasch
204! output of netCDF data format modified
205!
206! 1015 2012-09-27 09:23:24Z raasch
207! output of Adjustment of mixing length to the Prandtl mixing length at first
208! grid point above ground removed
209!
210! 1003 2012-09-14 14:35:53Z raasch
211! output of information about equal/unequal subdomain size removed
212!
213! 1001 2012-09-13 14:08:46Z raasch
214! all actions concerning leapfrog- and upstream-spline-scheme removed
215!
216! 978 2012-08-09 08:28:32Z fricke
217! -km_damp_max, outflow_damping_width
218! +pt_damping_factor, pt_damping_width
219! +z0h
220!
221! 964 2012-07-26 09:14:24Z raasch
222! output of profil-related quantities removed
223!
224! 940 2012-07-09 14:31:00Z raasch
225! Output in case of simulations for pure neutral stratification (no pt-equation
226! solved)
227!
228! 927 2012-06-06 19:15:04Z raasch
229! output of masking_method for mg-solver
230!
231! 868 2012-03-28 12:21:07Z raasch
232! translation velocity in Galilean transformation changed to 0.6 * ug
233!
234! 833 2012-02-22 08:55:55Z maronga
235! Adjusted format for leaf area density
236!
237! 828 2012-02-21 12:00:36Z raasch
238! output of dissipation_classes + radius_classes
239!
240! 825 2012-02-19 03:03:44Z raasch
241! Output of cloud physics parameters/quantities complemented and restructured
242!
243! Revision 1.1  1997/08/11 06:17:20  raasch
244! Initial revision
245!
246!
247! Description:
248! ------------
249!> Writing a header with all important information about the current run.
250!> This subroutine is called three times, two times at the beginning
251!> (writing information on files RUN_CONTROL and HEADER) and one time at the
252!> end of the run, then writing additional information about CPU-usage on file
253!> header.
254!-----------------------------------------------------------------------------!
255 SUBROUTINE header
256 
257
258    USE arrays_3d,                                                             &
259        ONLY:  pt_init, qsws, q_init, sa_init, shf, ug, vg, w_subs, zu, zw
260       
261    USE control_parameters
262       
263    USE cloud_parameters,                                                      &
264        ONLY:  cp, l_v, r_d
265
266    USE cpulog,                                                                &
267        ONLY:  log_point_s
268       
269    USE dvrp_variables,                                                        &
270        ONLY:  use_seperate_pe_for_dvrp_output
271       
272    USE grid_variables,                                                        &
273        ONLY:  dx, dy
274       
275    USE indices,                                                               &
276        ONLY:  mg_loc_ind, nnx, nny, nnz, nx, ny, nxl_mg, nxr_mg, nyn_mg,      &
277               nys_mg, nzt, nzt_mg
278       
279    USE kinds
280 
281    USE land_surface_model_mod,                                                &
282        ONLY: land_surface, lsm_header
283
284    USE microphysics_mod,                                                      &
285        ONLY:  cloud_water_sedimentation, collision_turbulence,                &
286               c_sedimentation, limiter_sedimentation, nc_const,               &
287               ventilation_effect
288
289    USE model_1d,                                                              &
290        ONLY:  damp_level_ind_1d, dt_pr_1d, dt_run_control_1d, end_time_1d
291       
292    USE netcdf_interface,                                                      &
293        ONLY:  netcdf_data_format, netcdf_data_format_string, netcdf_deflate
294
295    USE particle_attributes,                                                   &
296        ONLY:  bc_par_b, bc_par_lr, bc_par_ns, bc_par_t, collision_kernel,     &
297               curvature_solution_effects,                                     &
298               density_ratio, dissipation_classes, dt_min_part, dt_prel,       &
299               dt_write_particle_data, end_time_prel,                          &
300               number_of_particle_groups, particle_advection,                  &
301               particle_advection_start,                                       &
302               particles_per_point, pdx, pdy, pdz,  psb, psl, psn, psr, pss,   &
303               pst, radius, radius_classes, random_start_position,             &
304               seed_follows_topography,                                        &
305               total_number_of_particles, use_sgs_for_particles,               &
306               vertical_particle_advection, write_particle_statistics
307       
308    USE pegrid
309
310    USE plant_canopy_model_mod,                                                &
311        ONLY:  pcm_header, plant_canopy
312
313    USE pmc_handle_communicator,                                               &
314        ONLY:  pmc_get_model_info
315
316    USE pmc_interface,                                                         &
317        ONLY:  nested_run, nesting_datatransfer_mode, nesting_mode
318
319    USE radiation_model_mod,                                                   &
320        ONLY:  radiation, radiation_header
321   
322    USE spectra_mod,                                                           &
323        ONLY:  calculate_spectra, spectra_header
324
325    IMPLICIT NONE
326
327    CHARACTER (LEN=1)  ::  prec                !<
328   
329    CHARACTER (LEN=2)  ::  do2d_mode           !<
330   
331    CHARACTER (LEN=5)  ::  section_chr         !<
332   
333    CHARACTER (LEN=10) ::  coor_chr            !<
334    CHARACTER (LEN=10) ::  host_chr            !<
335   
336    CHARACTER (LEN=16) ::  begin_chr           !<
337   
338    CHARACTER (LEN=26) ::  ver_rev             !<
339
340    CHARACTER (LEN=32) ::  cpl_name            !<
341   
342    CHARACTER (LEN=40) ::  output_format       !<
343   
344    CHARACTER (LEN=70) ::  char1               !<
345    CHARACTER (LEN=70) ::  char2               !<
346    CHARACTER (LEN=70) ::  dopr_chr            !<
347    CHARACTER (LEN=70) ::  do2d_xy             !<
348    CHARACTER (LEN=70) ::  do2d_xz             !<
349    CHARACTER (LEN=70) ::  do2d_yz             !<
350    CHARACTER (LEN=70) ::  do3d_chr            !<
351    CHARACTER (LEN=70) ::  domask_chr          !<
352    CHARACTER (LEN=70) ::  run_classification  !<
353   
354    CHARACTER (LEN=85) ::  r_upper             !<
355    CHARACTER (LEN=85) ::  r_lower             !<
356   
357    CHARACTER (LEN=86) ::  coordinates         !<
358    CHARACTER (LEN=86) ::  gradients           !<
359    CHARACTER (LEN=86) ::  slices              !<
360    CHARACTER (LEN=86) ::  temperatures        !<
361    CHARACTER (LEN=86) ::  ugcomponent         !<
362    CHARACTER (LEN=86) ::  vgcomponent         !<
363
364    CHARACTER (LEN=1), DIMENSION(1:3) ::  dir = (/ 'x', 'y', 'z' /)  !<
365
366    INTEGER(iwp) ::  av             !<
367    INTEGER(iwp) ::  bh             !<
368    INTEGER(iwp) ::  blx            !<
369    INTEGER(iwp) ::  bly            !<
370    INTEGER(iwp) ::  bxl            !<
371    INTEGER(iwp) ::  bxr            !<
372    INTEGER(iwp) ::  byn            !<
373    INTEGER(iwp) ::  bys            !<
374    INTEGER(iwp) ::  ch             !<
375    INTEGER(iwp) ::  count          !<
376    INTEGER(iwp) ::  cpl_parent_id  !<
377    INTEGER(iwp) ::  cwx            !<
378    INTEGER(iwp) ::  cwy            !<
379    INTEGER(iwp) ::  cxl            !<
380    INTEGER(iwp) ::  cxr            !<
381    INTEGER(iwp) ::  cyn            !<
382    INTEGER(iwp) ::  cys            !<
383    INTEGER(iwp) ::  dim            !<
384    INTEGER(iwp) ::  i              !<
385    INTEGER(iwp) ::  io             !<
386    INTEGER(iwp) ::  j              !<
387    INTEGER(iwp) ::  k              !<
388    INTEGER(iwp) ::  l              !<
389    INTEGER(iwp) ::  ll             !<
390    INTEGER(iwp) ::  mpi_type       !<
391    INTEGER(iwp) ::  my_cpl_id      !<
392    INTEGER(iwp) ::  n              !<
393    INTEGER(iwp) ::  ncpl           !<
394    INTEGER(iwp) ::  npe_total      !<
395   
396
397    REAL(wp) ::  cpuseconds_per_simulated_second  !<
398    REAL(wp) ::  lower_left_coord_x               !< x-coordinate of nest domain
399    REAL(wp) ::  lower_left_coord_y               !< y-coordinate of nest domain
400
401!
402!-- Open the output file. At the end of the simulation, output is directed
403!-- to unit 19.
404    IF ( ( runnr == 0 .OR. force_print_header )  .AND. &
405         .NOT. simulated_time_at_begin /= simulated_time )  THEN
406       io = 15   !  header output on file RUN_CONTROL
407    ELSE
408       io = 19   !  header output on file HEADER
409    ENDIF
410    CALL check_open( io )
411
412!
413!-- At the end of the run, output file (HEADER) will be rewritten with
414!-- new information
415    IF ( io == 19 .AND. simulated_time_at_begin /= simulated_time ) REWIND( 19 )
416
417!
418!-- Determine kind of model run
419    IF ( TRIM( initializing_actions ) == 'read_restart_data' )  THEN
420       run_classification = 'restart run'
421    ELSEIF ( TRIM( initializing_actions ) == 'cyclic_fill' )  THEN
422       run_classification = 'run with cyclic fill of 3D - prerun data'
423    ELSEIF ( INDEX( initializing_actions, 'set_constant_profiles' ) /= 0 )  THEN
424       run_classification = 'run without 1D - prerun'
425    ELSEIF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
426       run_classification = 'run with 1D - prerun'
427    ELSEIF ( INDEX( initializing_actions, 'by_user' ) /=0 )  THEN
428       run_classification = 'run initialized by user'
429    ELSE
430       message_string = ' unknown action(s): ' // TRIM( initializing_actions )
431       CALL message( 'header', 'PA0191', 0, 0, 0, 6, 0 )
432    ENDIF
433    IF ( nested_run )  run_classification = 'nested ' // run_classification
434    IF ( ocean )  THEN
435       run_classification = 'ocean - ' // run_classification
436    ELSE
437       run_classification = 'atmosphere - ' // run_classification
438    ENDIF
439
440!
441!-- Run-identification, date, time, host
442    host_chr = host(1:10)
443    ver_rev = TRIM( version ) // '  ' // TRIM( revision )
444    WRITE ( io, 100 )  ver_rev, TRIM( run_classification )
445    IF ( TRIM( coupling_mode ) /= 'uncoupled' )  THEN
446#if defined( __mpi2 )
447       mpi_type = 2
448#else
449       mpi_type = 1
450#endif
451       WRITE ( io, 101 )  mpi_type, coupling_mode
452    ENDIF
453#if defined( __parallel )
454    IF ( coupling_start_time /= 0.0_wp )  THEN
455       IF ( coupling_start_time > simulated_time_at_begin )  THEN
456          WRITE ( io, 109 )
457       ELSE
458          WRITE ( io, 114 )
459       ENDIF
460    ENDIF
461#endif
462    IF ( ensemble_member_nr /= 0 )  THEN
463       WRITE ( io, 512 )  run_date, run_identifier, run_time, runnr,           &
464                       ADJUSTR( host_chr ), ensemble_member_nr
465    ELSE
466       WRITE ( io, 102 )  run_date, run_identifier, run_time, runnr,           &
467                       ADJUSTR( host_chr )
468    ENDIF
469#if defined( __parallel )
470    IF ( npex == -1  .AND.  npey == -1 )  THEN
471       char1 = 'calculated'
472    ELSE
473       char1 = 'predefined'
474    ENDIF
475    IF ( threads_per_task == 1 )  THEN
476       WRITE ( io, 103 )  numprocs, pdims(1), pdims(2), TRIM( char1 )
477    ELSE
478       WRITE ( io, 104 )  numprocs*threads_per_task, numprocs, &
479                          threads_per_task, pdims(1), pdims(2), TRIM( char1 )
480    ENDIF
481    IF ( num_acc_per_node /= 0 )  WRITE ( io, 117 )  num_acc_per_node   
482    IF ( ( host(1:3) == 'ibm'  .OR.  host(1:3) == 'nec'  .OR.    &
483           host(1:2) == 'lc'   .OR.  host(1:3) == 'dec' )  .AND. &
484         npex == -1  .AND.  pdims(2) == 1 )                      &
485    THEN
486       WRITE ( io, 106 )
487    ELSEIF ( pdims(2) == 1 )  THEN
488       WRITE ( io, 107 )  'x'
489    ELSEIF ( pdims(1) == 1 )  THEN
490       WRITE ( io, 107 )  'y'
491    ENDIF
492    IF ( use_seperate_pe_for_dvrp_output )  WRITE ( io, 105 )
493    IF ( numprocs /= maximum_parallel_io_streams )  THEN
494       WRITE ( io, 108 )  maximum_parallel_io_streams
495    ENDIF
496#else
497    IF ( num_acc_per_node /= 0 )  WRITE ( io, 120 )  num_acc_per_node
498#endif
499
500!
501!-- Nesting informations
502    IF ( nested_run )  THEN
503
504       WRITE ( io, 600 )  TRIM( nesting_mode ),                                &
505                          TRIM( nesting_datatransfer_mode )
506       CALL pmc_get_model_info( ncpl = ncpl, cpl_id = my_cpl_id )
507
508       DO  n = 1, ncpl
509          CALL pmc_get_model_info( request_for_cpl_id = n, cpl_name = cpl_name,&
510                                   cpl_parent_id = cpl_parent_id,              &
511                                   lower_left_x = lower_left_coord_x,          &
512                                   lower_left_y = lower_left_coord_y,          &
513                                   npe_total = npe_total )
514          IF ( n == my_cpl_id )  THEN
515             char1 = '*'
516          ELSE
517             char1 = ' '
518          ENDIF
519          WRITE ( io, 601 )  TRIM( char1 ), n, cpl_parent_id, npe_total,       &
520                             lower_left_coord_x, lower_left_coord_y,           &
521                             TRIM( cpl_name )
522       ENDDO
523    ENDIF
524    WRITE ( io, 99 )
525
526!
527!-- Numerical schemes
528    WRITE ( io, 110 )
529    IF ( psolver(1:7) == 'poisfft' )  THEN
530       WRITE ( io, 111 )  TRIM( fft_method )
531       IF ( transpose_compute_overlap )  WRITE( io, 115 )
532    ELSEIF ( psolver == 'sor' )  THEN
533       WRITE ( io, 112 )  nsor_ini, nsor, omega_sor
534    ELSEIF ( psolver(1:9) == 'multigrid' )  THEN
535       WRITE ( io, 135 )  TRIM(psolver), cycle_mg, maximum_grid_level, ngsrb
536       IF ( mg_cycles == -1 )  THEN
537          WRITE ( io, 140 )  residual_limit
538       ELSE
539          WRITE ( io, 141 )  mg_cycles
540       ENDIF
541       IF ( mg_switch_to_pe0_level == 0 )  THEN
542          WRITE ( io, 136 )  nxr_mg(1)-nxl_mg(1)+1, nyn_mg(1)-nys_mg(1)+1, &
543                             nzt_mg(1)
544       ELSEIF (  mg_switch_to_pe0_level /= -1 )  THEN
545          WRITE ( io, 137 )  mg_switch_to_pe0_level,            &
546                             mg_loc_ind(2,0)-mg_loc_ind(1,0)+1, &
547                             mg_loc_ind(4,0)-mg_loc_ind(3,0)+1, &
548                             nzt_mg(mg_switch_to_pe0_level),    &
549                             nxr_mg(1)-nxl_mg(1)+1, nyn_mg(1)-nys_mg(1)+1, &
550                             nzt_mg(1)
551       ENDIF
552       IF ( masking_method )  WRITE ( io, 144 )
553    ENDIF
554    IF ( call_psolver_at_all_substeps  .AND. timestep_scheme(1:5) == 'runge' ) &
555    THEN
556       WRITE ( io, 142 )
557    ENDIF
558
559    IF ( momentum_advec == 'pw-scheme' )  THEN
560       WRITE ( io, 113 )
561    ELSEIF (momentum_advec == 'ws-scheme' )  THEN
562       WRITE ( io, 503 )
563    ENDIF
564    IF ( scalar_advec == 'pw-scheme' )  THEN
565       WRITE ( io, 116 )
566    ELSEIF ( scalar_advec == 'ws-scheme' )  THEN
567       WRITE ( io, 504 )
568    ELSEIF ( scalar_advec == 'ws-scheme-mono' )  THEN
569       WRITE ( io, 513 )
570    ELSE
571       WRITE ( io, 118 )
572    ENDIF
573
574    WRITE ( io, 139 )  TRIM( loop_optimization )
575
576    IF ( galilei_transformation )  THEN
577       IF ( use_ug_for_galilei_tr )  THEN
578          char1 = '0.6 * geostrophic wind'
579       ELSE
580          char1 = 'mean wind in model domain'
581       ENDIF
582       IF ( simulated_time_at_begin == simulated_time )  THEN
583          char2 = 'at the start of the run'
584       ELSE
585          char2 = 'at the end of the run'
586       ENDIF
587       WRITE ( io, 119 )  TRIM( char1 ), TRIM( char2 ),                        &
588                          advected_distance_x/1000.0_wp,                       &
589                          advected_distance_y/1000.0_wp
590    ENDIF
591    WRITE ( io, 122 )  timestep_scheme
592    IF ( use_upstream_for_tke )  WRITE ( io, 143 )
593    IF ( rayleigh_damping_factor /= 0.0_wp )  THEN
594       IF ( .NOT. ocean )  THEN
595          WRITE ( io, 123 )  'above', rayleigh_damping_height, &
596               rayleigh_damping_factor
597       ELSE
598          WRITE ( io, 123 )  'below', rayleigh_damping_height, &
599               rayleigh_damping_factor
600       ENDIF
601    ENDIF
602    IF ( neutral )  WRITE ( io, 131 )  pt_surface
603    IF ( humidity )  THEN
604       IF ( .NOT. cloud_physics )  THEN
605          WRITE ( io, 129 )
606       ELSE
607          WRITE ( io, 130 )
608       ENDIF
609    ENDIF
610    IF ( passive_scalar )  WRITE ( io, 134 )
611    IF ( conserve_volume_flow )  THEN
612       WRITE ( io, 150 )  conserve_volume_flow_mode
613       IF ( TRIM( conserve_volume_flow_mode ) == 'bulk_velocity' )  THEN
614          WRITE ( io, 151 )  u_bulk, v_bulk
615       ENDIF
616    ELSEIF ( dp_external )  THEN
617       IF ( dp_smooth )  THEN
618          WRITE ( io, 152 )  dpdxy, dp_level_b, ', vertically smoothed.'
619       ELSE
620          WRITE ( io, 152 )  dpdxy, dp_level_b, '.'
621       ENDIF
622    ENDIF
623    WRITE ( io, 99 )
624
625!
626!-- Runtime and timestep information
627    WRITE ( io, 200 )
628    IF ( .NOT. dt_fixed )  THEN
629       WRITE ( io, 201 )  dt_max, cfl_factor
630    ELSE
631       WRITE ( io, 202 )  dt
632    ENDIF
633    WRITE ( io, 203 )  simulated_time_at_begin, end_time
634
635    IF ( time_restart /= 9999999.9_wp  .AND. &
636         simulated_time_at_begin == simulated_time )  THEN
637       IF ( dt_restart == 9999999.9_wp )  THEN
638          WRITE ( io, 204 )  ' Restart at:       ',time_restart
639       ELSE
640          WRITE ( io, 205 )  ' Restart at:       ',time_restart, dt_restart
641       ENDIF
642    ENDIF
643
644    IF ( simulated_time_at_begin /= simulated_time )  THEN
645       i = MAX ( log_point_s(10)%counts, 1 )
646       IF ( ( simulated_time - simulated_time_at_begin ) == 0.0_wp )  THEN
647          cpuseconds_per_simulated_second = 0.0_wp
648       ELSE
649          cpuseconds_per_simulated_second = log_point_s(10)%sum / &
650                                            ( simulated_time -    &
651                                              simulated_time_at_begin )
652       ENDIF
653       WRITE ( io, 206 )  simulated_time, log_point_s(10)%sum,      &
654                          log_point_s(10)%sum / REAL( i, KIND=wp ), &
655                          cpuseconds_per_simulated_second
656       IF ( time_restart /= 9999999.9_wp  .AND.  time_restart < end_time )  THEN
657          IF ( dt_restart == 9999999.9_wp )  THEN
658             WRITE ( io, 204 )  ' Next restart at:     ',time_restart
659          ELSE
660             WRITE ( io, 205 )  ' Next restart at:     ',time_restart, dt_restart
661          ENDIF
662       ENDIF
663    ENDIF
664
665
666!
667!-- Start time for coupled runs, if independent precursor runs for atmosphere
668!-- and ocean are used or have been used. In this case, coupling_start_time
669!-- defines the time when the coupling is switched on.
670    IF ( coupling_start_time /= 0.0_wp )  THEN
671       WRITE ( io, 207 )  coupling_start_time
672    ENDIF
673
674!
675!-- Computational grid
676    IF ( .NOT. ocean )  THEN
677       WRITE ( io, 250 )  dx, dy, dz, (nx+1)*dx, (ny+1)*dy, zu(nzt+1)
678       IF ( dz_stretch_level_index < nzt+1 )  THEN
679          WRITE ( io, 252 )  dz_stretch_level, dz_stretch_level_index, &
680                             dz_stretch_factor, dz_max
681       ENDIF
682    ELSE
683       WRITE ( io, 250 )  dx, dy, dz, (nx+1)*dx, (ny+1)*dy, zu(0)
684       IF ( dz_stretch_level_index > 0 )  THEN
685          WRITE ( io, 252 )  dz_stretch_level, dz_stretch_level_index, &
686                             dz_stretch_factor, dz_max
687       ENDIF
688    ENDIF
689    WRITE ( io, 254 )  nx, ny, nzt+1, MIN( nnx, nx+1 ), MIN( nny, ny+1 ), &
690                       MIN( nnz+2, nzt+2 )
691    IF ( sloping_surface )  WRITE ( io, 260 )  alpha_surface
692
693!
694!-- Large scale forcing and nudging
695    WRITE ( io, 160 )
696    IF ( large_scale_forcing )  THEN
697       WRITE ( io, 162 )
698       WRITE ( io, 163 )
699
700       IF ( large_scale_subsidence )  THEN
701          IF ( .NOT. use_subsidence_tendencies )  THEN
702             WRITE ( io, 164 )
703          ELSE
704             WRITE ( io, 165 )
705          ENDIF
706       ENDIF
707
708       IF ( bc_pt_b == 'dirichlet' )  THEN
709          WRITE ( io, 180 )
710       ELSEIF ( bc_pt_b == 'neumann' )  THEN
711          WRITE ( io, 181 )
712       ENDIF
713
714       IF ( bc_q_b == 'dirichlet' )  THEN
715          WRITE ( io, 182 )
716       ELSEIF ( bc_q_b == 'neumann' )  THEN
717          WRITE ( io, 183 )
718       ENDIF
719
720       WRITE ( io, 167 )
721       IF ( nudging )  THEN
722          WRITE ( io, 170 )
723       ENDIF
724    ELSE
725       WRITE ( io, 161 )
726       WRITE ( io, 171 )
727    ENDIF
728    IF ( large_scale_subsidence )  THEN
729       WRITE ( io, 168 )
730       WRITE ( io, 169 )
731    ENDIF
732
733!
734!-- Profile for the large scale vertial velocity
735!-- Building output strings, starting with surface value
736    IF ( large_scale_subsidence )  THEN
737       temperatures = '   0.0'
738       gradients = '------'
739       slices = '     0'
740       coordinates = '   0.0'
741       i = 1
742       DO  WHILE ( subs_vertical_gradient_level_i(i) /= -9999 )
743
744          WRITE (coor_chr,'(E10.2,7X)')  &
745                                w_subs(subs_vertical_gradient_level_i(i))
746          temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr )
747
748          WRITE (coor_chr,'(E10.2,7X)')  subs_vertical_gradient(i)
749          gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
750
751          WRITE (coor_chr,'(I10,7X)')  subs_vertical_gradient_level_i(i)
752          slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
753
754          WRITE (coor_chr,'(F10.2,7X)')  subs_vertical_gradient_level(i)
755          coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
756
757          IF ( i == 10 )  THEN
758             EXIT
759          ELSE
760             i = i + 1
761          ENDIF
762
763       ENDDO
764
765 
766       IF ( .NOT. large_scale_forcing )  THEN
767          WRITE ( io, 426 )  TRIM( coordinates ), TRIM( temperatures ), &
768                             TRIM( gradients ), TRIM( slices )
769       ENDIF
770
771
772    ENDIF
773
774!-- Profile of the geostrophic wind (component ug)
775!-- Building output strings
776    WRITE ( ugcomponent, '(F6.2)' )  ug_surface
777    gradients = '------'
778    slices = '     0'
779    coordinates = '   0.0'
780    i = 1
781    DO  WHILE ( ug_vertical_gradient_level_ind(i) /= -9999 )
782     
783       WRITE (coor_chr,'(F6.2,1X)')  ug(ug_vertical_gradient_level_ind(i))
784       ugcomponent = TRIM( ugcomponent ) // '  ' // TRIM( coor_chr )
785
786       WRITE (coor_chr,'(F6.2,1X)')  ug_vertical_gradient(i)
787       gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
788
789       WRITE (coor_chr,'(I6,1X)')  ug_vertical_gradient_level_ind(i)
790       slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
791
792       WRITE (coor_chr,'(F6.1,1X)')  ug_vertical_gradient_level(i)
793       coordinates = TRIM( coordinates ) // '  ' // TRIM( coor_chr )
794
795       IF ( i == 10 )  THEN
796          EXIT
797       ELSE
798          i = i + 1
799       ENDIF
800
801    ENDDO
802
803    IF ( .NOT. large_scale_forcing )  THEN
804       WRITE ( io, 423 )  TRIM( coordinates ), TRIM( ugcomponent ), &
805                          TRIM( gradients ), TRIM( slices )
806    ENDIF
807
808!-- Profile of the geostrophic wind (component vg)
809!-- Building output strings
810    WRITE ( vgcomponent, '(F6.2)' )  vg_surface
811    gradients = '------'
812    slices = '     0'
813    coordinates = '   0.0'
814    i = 1
815    DO  WHILE ( vg_vertical_gradient_level_ind(i) /= -9999 )
816
817       WRITE (coor_chr,'(F6.2,1X)')  vg(vg_vertical_gradient_level_ind(i))
818       vgcomponent = TRIM( vgcomponent ) // '  ' // TRIM( coor_chr )
819
820       WRITE (coor_chr,'(F6.2,1X)')  vg_vertical_gradient(i)
821       gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
822
823       WRITE (coor_chr,'(I6,1X)')  vg_vertical_gradient_level_ind(i)
824       slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
825
826       WRITE (coor_chr,'(F6.1,1X)')  vg_vertical_gradient_level(i)
827       coordinates = TRIM( coordinates ) // '  ' // TRIM( coor_chr )
828
829       IF ( i == 10 )  THEN
830          EXIT
831       ELSE
832          i = i + 1
833       ENDIF
834 
835    ENDDO
836
837    IF ( .NOT. large_scale_forcing )  THEN
838       WRITE ( io, 424 )  TRIM( coordinates ), TRIM( vgcomponent ), &
839                          TRIM( gradients ), TRIM( slices )
840    ENDIF
841
842!
843!-- Topography
844    WRITE ( io, 270 )  topography
845    SELECT CASE ( TRIM( topography ) )
846
847       CASE ( 'flat' )
848          ! no actions necessary
849
850       CASE ( 'single_building' )
851          blx = INT( building_length_x / dx )
852          bly = INT( building_length_y / dy )
853          bh  = MINLOC( ABS( zw - building_height ), 1 ) - 1
854          IF ( ABS( zw(bh  ) - building_height ) == &
855               ABS( zw(bh+1) - building_height )    )  bh = bh + 1
856
857          IF ( building_wall_left == 9999999.9_wp )  THEN
858             building_wall_left = ( nx + 1 - blx ) / 2 * dx
859          ENDIF
860          bxl = INT ( building_wall_left / dx + 0.5_wp )
861          bxr = bxl + blx
862
863          IF ( building_wall_south == 9999999.9_wp )  THEN
864             building_wall_south = ( ny + 1 - bly ) / 2 * dy
865          ENDIF
866          bys = INT ( building_wall_south / dy + 0.5_wp )
867          byn = bys + bly
868
869          WRITE ( io, 271 )  building_length_x, building_length_y, &
870                             building_height, bxl, bxr, bys, byn
871
872       CASE ( 'single_street_canyon' )
873          ch  = MINLOC( ABS( zw - canyon_height ), 1 ) - 1
874          IF ( ABS( zw(ch  ) - canyon_height ) == &
875               ABS( zw(ch+1) - canyon_height )    )  ch = ch + 1
876          IF ( canyon_width_x /= 9999999.9_wp )  THEN
877!
878!--          Street canyon in y direction
879             cwx = NINT( canyon_width_x / dx )
880             IF ( canyon_wall_left == 9999999.9_wp )  THEN
881                canyon_wall_left = ( nx + 1 - cwx ) / 2 * dx
882             ENDIF
883             cxl = NINT( canyon_wall_left / dx )
884             cxr = cxl + cwx
885             WRITE ( io, 272 )  'y', canyon_height, ch, 'u', cxl, cxr
886
887          ELSEIF ( canyon_width_y /= 9999999.9_wp )  THEN
888!
889!--          Street canyon in x direction
890             cwy = NINT( canyon_width_y / dy )
891             IF ( canyon_wall_south == 9999999.9_wp )  THEN
892                canyon_wall_south = ( ny + 1 - cwy ) / 2 * dy
893             ENDIF
894             cys = NINT( canyon_wall_south / dy )
895             cyn = cys + cwy
896             WRITE ( io, 272 )  'x', canyon_height, ch, 'v', cys, cyn
897          ENDIF
898
899    END SELECT
900
901    IF ( TRIM( topography ) /= 'flat' )  THEN
902       IF ( TRIM( topography_grid_convention ) == ' ' )  THEN
903          IF ( TRIM( topography ) == 'single_building' .OR.  &
904               TRIM( topography ) == 'single_street_canyon' )  THEN
905             WRITE ( io, 278 )
906          ELSEIF ( TRIM( topography ) == 'read_from_file' )  THEN
907             WRITE ( io, 279 )
908          ENDIF
909       ELSEIF ( TRIM( topography_grid_convention ) == 'cell_edge' )  THEN
910          WRITE ( io, 278 )
911       ELSEIF ( TRIM( topography_grid_convention ) == 'cell_center' )  THEN
912          WRITE ( io, 279 )
913       ENDIF
914    ENDIF
915
916    IF ( plant_canopy )  CALL pcm_header ( io )
917
918    IF ( land_surface )  CALL lsm_header ( io )
919
920    IF ( radiation )  CALL radiation_header ( io )
921
922!
923!-- Boundary conditions
924    IF ( ibc_p_b == 0 )  THEN
925       r_lower = 'p(0)     = 0      |'
926    ELSEIF ( ibc_p_b == 1 )  THEN
927       r_lower = 'p(0)     = p(1)   |'
928    ENDIF
929    IF ( ibc_p_t == 0 )  THEN
930       r_upper  = 'p(nzt+1) = 0      |'
931    ELSE
932       r_upper  = 'p(nzt+1) = p(nzt) |'
933    ENDIF
934
935    IF ( ibc_uv_b == 0 )  THEN
936       r_lower = TRIM( r_lower ) // ' uv(0)     = -uv(1)                |'
937    ELSE
938       r_lower = TRIM( r_lower ) // ' uv(0)     = uv(1)                 |'
939    ENDIF
940    IF ( TRIM( bc_uv_t ) == 'dirichlet_0' )  THEN
941       r_upper  = TRIM( r_upper  ) // ' uv(nzt+1) = 0                     |'
942    ELSEIF ( ibc_uv_t == 0 )  THEN
943       r_upper  = TRIM( r_upper  ) // ' uv(nzt+1) = ug(nzt+1), vg(nzt+1)  |'
944    ELSE
945       r_upper  = TRIM( r_upper  ) // ' uv(nzt+1) = uv(nzt)               |'
946    ENDIF
947
948    IF ( ibc_pt_b == 0 )  THEN
949       IF ( land_surface )  THEN
950          r_lower = TRIM( r_lower ) // ' pt(0)     = from soil model'
951       ELSE
952          r_lower = TRIM( r_lower ) // ' pt(0)     = pt_surface'
953       ENDIF
954    ELSEIF ( ibc_pt_b == 1 )  THEN
955       r_lower = TRIM( r_lower ) // ' pt(0)     = pt(1)'
956    ELSEIF ( ibc_pt_b == 2 )  THEN
957       r_lower = TRIM( r_lower ) // ' pt(0)     = from coupled model'
958    ENDIF
959    IF ( ibc_pt_t == 0 )  THEN
960       r_upper  = TRIM( r_upper  ) // ' pt(nzt+1) = pt_top'
961    ELSEIF( ibc_pt_t == 1 )  THEN
962       r_upper  = TRIM( r_upper  ) // ' pt(nzt+1) = pt(nzt)'
963    ELSEIF( ibc_pt_t == 2 )  THEN
964       r_upper  = TRIM( r_upper  ) // ' pt(nzt+1) = pt(nzt) + dpt/dz_ini'
965
966    ENDIF
967
968    WRITE ( io, 300 )  r_lower, r_upper
969
970    IF ( .NOT. constant_diffusion )  THEN
971       IF ( ibc_e_b == 1 )  THEN
972          r_lower = 'e(0)     = e(1)'
973       ELSE
974          r_lower = 'e(0)     = e(1) = (u*/0.1)**2'
975       ENDIF
976       r_upper = 'e(nzt+1) = e(nzt) = e(nzt-1)'
977
978       WRITE ( io, 301 )  'e', r_lower, r_upper       
979
980    ENDIF
981
982    IF ( ocean )  THEN
983       r_lower = 'sa(0)    = sa(1)'
984       IF ( ibc_sa_t == 0 )  THEN
985          r_upper =  'sa(nzt+1) = sa_surface'
986       ELSE
987          r_upper =  'sa(nzt+1) = sa(nzt)'
988       ENDIF
989       WRITE ( io, 301 ) 'sa', r_lower, r_upper
990    ENDIF
991
992    IF ( humidity )  THEN
993       IF ( ibc_q_b == 0 )  THEN
994          IF ( land_surface )  THEN
995             r_lower = 'q(0)     = from soil model'
996          ELSE
997             r_lower = 'q(0)     = q_surface'
998          ENDIF
999
1000       ELSE
1001          r_lower = 'q(0)     = q(1)'
1002       ENDIF
1003       IF ( ibc_q_t == 0 )  THEN
1004          r_upper =  'q(nzt)   = q_top'
1005       ELSE
1006          r_upper =  'q(nzt)   = q(nzt-1) + dq/dz'
1007       ENDIF
1008       WRITE ( io, 301 ) 'q', r_lower, r_upper
1009    ENDIF
1010
1011    IF ( passive_scalar )  THEN
1012       IF ( ibc_q_b == 0 )  THEN
1013          r_lower = 's(0)     = s_surface'
1014       ELSE
1015          r_lower = 's(0)     = s(1)'
1016       ENDIF
1017       IF ( ibc_q_t == 0 )  THEN
1018          r_upper =  's(nzt)   = s_top'
1019       ELSE
1020          r_upper =  's(nzt)   = s(nzt-1) + ds/dz'
1021       ENDIF
1022       WRITE ( io, 301 ) 's', r_lower, r_upper
1023    ENDIF
1024
1025    IF ( use_surface_fluxes )  THEN
1026       WRITE ( io, 303 )
1027       IF ( constant_heatflux )  THEN
1028          IF ( large_scale_forcing .AND. lsf_surf )  THEN
1029             WRITE ( io, 306 )  shf(0,0)
1030          ELSE
1031             WRITE ( io, 306 )  surface_heatflux
1032          ENDIF
1033          IF ( random_heatflux )  WRITE ( io, 307 )
1034       ENDIF
1035       IF ( humidity  .AND.  constant_waterflux )  THEN
1036          IF ( large_scale_forcing .AND. lsf_surf )  THEN
1037             WRITE ( io, 311 ) qsws(0,0)
1038          ELSE
1039             WRITE ( io, 311 ) surface_waterflux
1040          ENDIF
1041       ENDIF
1042       IF ( passive_scalar  .AND.  constant_waterflux )  THEN
1043          WRITE ( io, 313 ) surface_waterflux
1044       ENDIF
1045    ENDIF
1046
1047    IF ( use_top_fluxes )  THEN
1048       WRITE ( io, 304 )
1049       IF ( coupling_mode == 'uncoupled' )  THEN
1050          WRITE ( io, 320 )  top_momentumflux_u, top_momentumflux_v
1051          IF ( constant_top_heatflux )  THEN
1052             WRITE ( io, 306 )  top_heatflux
1053          ENDIF
1054       ELSEIF ( coupling_mode == 'ocean_to_atmosphere' )  THEN
1055          WRITE ( io, 316 )
1056       ENDIF
1057       IF ( ocean  .AND.  constant_top_salinityflux )  THEN
1058          WRITE ( io, 309 )  top_salinityflux
1059       ENDIF
1060       IF ( humidity  .OR.  passive_scalar )  THEN
1061          WRITE ( io, 315 )
1062       ENDIF
1063    ENDIF
1064
1065    IF ( constant_flux_layer )  THEN
1066       WRITE ( io, 305 )  (zu(1)-zu(0)), roughness_length,                     &
1067                          z0h_factor*roughness_length, kappa,                  &
1068                          zeta_min, zeta_max
1069       IF ( .NOT. constant_heatflux )  WRITE ( io, 308 )
1070       IF ( humidity  .AND.  .NOT. constant_waterflux )  THEN
1071          WRITE ( io, 312 )
1072       ENDIF
1073       IF ( passive_scalar  .AND.  .NOT. constant_waterflux )  THEN
1074          WRITE ( io, 314 )
1075       ENDIF
1076    ELSE
1077       IF ( INDEX(initializing_actions, 'set_1d-model_profiles') /= 0 )  THEN
1078          WRITE ( io, 310 )  zeta_min, zeta_max
1079       ENDIF
1080    ENDIF
1081
1082    WRITE ( io, 317 )  bc_lr, bc_ns
1083    IF ( .NOT. bc_lr_cyc  .OR.  .NOT. bc_ns_cyc )  THEN
1084       WRITE ( io, 318 )  use_cmax, pt_damping_width, pt_damping_factor       
1085       IF ( turbulent_inflow )  THEN
1086          IF ( .NOT. recycling_yshift ) THEN
1087             WRITE ( io, 319 )  recycling_width, recycling_plane, &
1088                                inflow_damping_height, inflow_damping_width
1089          ELSE
1090             WRITE ( io, 322 )  recycling_width, recycling_plane, &
1091                                inflow_damping_height, inflow_damping_width
1092          END IF
1093       ENDIF
1094    ENDIF
1095
1096!
1097!-- Initial Profiles
1098    WRITE ( io, 321 )
1099!
1100!-- Initial wind profiles
1101    IF ( u_profile(1) /= 9999999.9_wp )  WRITE ( io, 427 )
1102
1103!
1104!-- Initial temperature profile
1105!-- Building output strings, starting with surface temperature
1106    WRITE ( temperatures, '(F6.2)' )  pt_surface
1107    gradients = '------'
1108    slices = '     0'
1109    coordinates = '   0.0'
1110    i = 1
1111    DO  WHILE ( pt_vertical_gradient_level_ind(i) /= -9999 )
1112
1113       WRITE (coor_chr,'(F7.2)')  pt_init(pt_vertical_gradient_level_ind(i))
1114       temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr )
1115
1116       WRITE (coor_chr,'(F7.2)')  pt_vertical_gradient(i)
1117       gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
1118
1119       WRITE (coor_chr,'(I7)')  pt_vertical_gradient_level_ind(i)
1120       slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
1121
1122       WRITE (coor_chr,'(F7.1)')  pt_vertical_gradient_level(i)
1123       coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
1124
1125       IF ( i == 10 )  THEN
1126          EXIT
1127       ELSE
1128          i = i + 1
1129       ENDIF
1130
1131    ENDDO
1132
1133    IF ( .NOT. nudging )  THEN
1134       WRITE ( io, 420 )  TRIM( coordinates ), TRIM( temperatures ), &
1135                          TRIM( gradients ), TRIM( slices )
1136    ELSE
1137       WRITE ( io, 428 ) 
1138    ENDIF
1139
1140!
1141!-- Initial humidity profile
1142!-- Building output strings, starting with surface humidity
1143    IF ( humidity  .OR.  passive_scalar )  THEN
1144       WRITE ( temperatures, '(E8.1)' )  q_surface
1145       gradients = '--------'
1146       slices = '       0'
1147       coordinates = '     0.0'
1148       i = 1
1149       DO  WHILE ( q_vertical_gradient_level_ind(i) /= -9999 )
1150         
1151          WRITE (coor_chr,'(E8.1,4X)')  q_init(q_vertical_gradient_level_ind(i))
1152          temperatures = TRIM( temperatures ) // '  ' // TRIM( coor_chr )
1153
1154          WRITE (coor_chr,'(E8.1,4X)')  q_vertical_gradient(i)
1155          gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
1156         
1157          WRITE (coor_chr,'(I8,4X)')  q_vertical_gradient_level_ind(i)
1158          slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
1159         
1160          WRITE (coor_chr,'(F8.1,4X)')  q_vertical_gradient_level(i)
1161          coordinates = TRIM( coordinates ) // '  '  // TRIM( coor_chr )
1162
1163          IF ( i == 10 )  THEN
1164             EXIT
1165          ELSE
1166             i = i + 1
1167          ENDIF
1168
1169       ENDDO
1170
1171       IF ( humidity )  THEN
1172          IF ( .NOT. nudging )  THEN
1173             WRITE ( io, 421 )  TRIM( coordinates ), TRIM( temperatures ), &
1174                                TRIM( gradients ), TRIM( slices )
1175          ENDIF
1176       ELSE
1177          WRITE ( io, 422 )  TRIM( coordinates ), TRIM( temperatures ), &
1178                             TRIM( gradients ), TRIM( slices )
1179       ENDIF
1180    ENDIF
1181
1182!
1183!-- Initial salinity profile
1184!-- Building output strings, starting with surface salinity
1185    IF ( ocean )  THEN
1186       WRITE ( temperatures, '(F6.2)' )  sa_surface
1187       gradients = '------'
1188       slices = '     0'
1189       coordinates = '   0.0'
1190       i = 1
1191       DO  WHILE ( sa_vertical_gradient_level_ind(i) /= -9999 )
1192
1193          WRITE (coor_chr,'(F7.2)')  sa_init(sa_vertical_gradient_level_ind(i))
1194          temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr )
1195
1196          WRITE (coor_chr,'(F7.2)')  sa_vertical_gradient(i)
1197          gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
1198
1199          WRITE (coor_chr,'(I7)')  sa_vertical_gradient_level_ind(i)
1200          slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
1201
1202          WRITE (coor_chr,'(F7.1)')  sa_vertical_gradient_level(i)
1203          coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
1204
1205          IF ( i == 10 )  THEN
1206             EXIT
1207          ELSE
1208             i = i + 1
1209          ENDIF
1210
1211       ENDDO
1212
1213       WRITE ( io, 425 )  TRIM( coordinates ), TRIM( temperatures ), &
1214                          TRIM( gradients ), TRIM( slices )
1215    ENDIF
1216
1217
1218!
1219!-- Listing of 1D-profiles
1220    WRITE ( io, 325 )  dt_dopr_listing
1221    IF ( averaging_interval_pr /= 0.0_wp )  THEN
1222       WRITE ( io, 326 )  averaging_interval_pr, dt_averaging_input_pr
1223    ENDIF
1224
1225!
1226!-- DATA output
1227    WRITE ( io, 330 )
1228    IF ( averaging_interval_pr /= 0.0_wp )  THEN
1229       WRITE ( io, 326 )  averaging_interval_pr, dt_averaging_input_pr
1230    ENDIF
1231
1232!
1233!-- 1D-profiles
1234    dopr_chr = 'Profile:'
1235    IF ( dopr_n /= 0 )  THEN
1236       WRITE ( io, 331 )
1237
1238       output_format = ''
1239       output_format = netcdf_data_format_string
1240       IF ( netcdf_deflate == 0 )  THEN
1241          WRITE ( io, 344 )  output_format
1242       ELSE
1243          WRITE ( io, 354 )  TRIM( output_format ), netcdf_deflate
1244       ENDIF
1245
1246       DO  i = 1, dopr_n
1247          dopr_chr = TRIM( dopr_chr ) // ' ' // TRIM( data_output_pr(i) ) // ','
1248          IF ( LEN_TRIM( dopr_chr ) >= 60 )  THEN
1249             WRITE ( io, 332 )  dopr_chr
1250             dopr_chr = '       :'
1251          ENDIF
1252       ENDDO
1253
1254       IF ( dopr_chr /= '' )  THEN
1255          WRITE ( io, 332 )  dopr_chr
1256       ENDIF
1257       WRITE ( io, 333 )  dt_dopr, averaging_interval_pr, dt_averaging_input_pr
1258       IF ( skip_time_dopr /= 0.0_wp )  WRITE ( io, 339 )  skip_time_dopr
1259    ENDIF
1260
1261!
1262!-- 2D-arrays
1263    DO  av = 0, 1
1264
1265       i = 1
1266       do2d_xy = ''
1267       do2d_xz = ''
1268       do2d_yz = ''
1269       DO  WHILE ( do2d(av,i) /= ' ' )
1270
1271          l = MAX( 2, LEN_TRIM( do2d(av,i) ) )
1272          do2d_mode = do2d(av,i)(l-1:l)
1273
1274          SELECT CASE ( do2d_mode )
1275             CASE ( 'xy' )
1276                ll = LEN_TRIM( do2d_xy )
1277                do2d_xy = do2d_xy(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
1278             CASE ( 'xz' )
1279                ll = LEN_TRIM( do2d_xz )
1280                do2d_xz = do2d_xz(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
1281             CASE ( 'yz' )
1282                ll = LEN_TRIM( do2d_yz )
1283                do2d_yz = do2d_yz(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
1284          END SELECT
1285
1286          i = i + 1
1287
1288       ENDDO
1289
1290       IF ( ( ( do2d_xy /= ''  .AND.  section(1,1) /= -9999 )  .OR.    &
1291              ( do2d_xz /= ''  .AND.  section(1,2) /= -9999 )  .OR.    &
1292              ( do2d_yz /= ''  .AND.  section(1,3) /= -9999 ) ) )  THEN
1293
1294          IF (  av == 0 )  THEN
1295             WRITE ( io, 334 )  ''
1296          ELSE
1297             WRITE ( io, 334 )  '(time-averaged)'
1298          ENDIF
1299
1300          IF ( do2d_at_begin )  THEN
1301             begin_chr = 'and at the start'
1302          ELSE
1303             begin_chr = ''
1304          ENDIF
1305
1306          output_format = ''
1307          output_format = netcdf_data_format_string
1308          IF ( netcdf_deflate == 0 )  THEN
1309             WRITE ( io, 344 )  output_format
1310          ELSE
1311             WRITE ( io, 354 )  TRIM( output_format ), netcdf_deflate
1312          ENDIF
1313
1314          IF ( do2d_xy /= ''  .AND.  section(1,1) /= -9999 )  THEN
1315             i = 1
1316             slices = '/'
1317             coordinates = '/'
1318!
1319!--          Building strings with index and coordinate information of the
1320!--          slices
1321             DO  WHILE ( section(i,1) /= -9999 )
1322
1323                WRITE (section_chr,'(I5)')  section(i,1)
1324                section_chr = ADJUSTL( section_chr )
1325                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
1326
1327                IF ( section(i,1) == -1 )  THEN
1328                   WRITE (coor_chr,'(F10.1)')  -1.0_wp
1329                ELSE
1330                   WRITE (coor_chr,'(F10.1)')  zu(section(i,1))
1331                ENDIF
1332                coor_chr = ADJUSTL( coor_chr )
1333                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
1334
1335                i = i + 1
1336             ENDDO
1337             IF ( av == 0 )  THEN
1338                WRITE ( io, 335 )  'XY', do2d_xy, dt_do2d_xy, &
1339                                   TRIM( begin_chr ), 'k', TRIM( slices ), &
1340                                   TRIM( coordinates )
1341                IF ( skip_time_do2d_xy /= 0.0_wp )  THEN
1342                   WRITE ( io, 339 )  skip_time_do2d_xy
1343                ENDIF
1344             ELSE
1345                WRITE ( io, 342 )  'XY', do2d_xy, dt_data_output_av, &
1346                                   TRIM( begin_chr ), averaging_interval, &
1347                                   dt_averaging_input, 'k', TRIM( slices ), &
1348                                   TRIM( coordinates )
1349                IF ( skip_time_data_output_av /= 0.0_wp )  THEN
1350                   WRITE ( io, 339 )  skip_time_data_output_av
1351                ENDIF
1352             ENDIF
1353             IF ( netcdf_data_format > 4 )  THEN
1354                WRITE ( io, 352 )  ntdim_2d_xy(av)
1355             ELSE
1356                WRITE ( io, 353 )
1357             ENDIF
1358          ENDIF
1359
1360          IF ( do2d_xz /= ''  .AND.  section(1,2) /= -9999 )  THEN
1361             i = 1
1362             slices = '/'
1363             coordinates = '/'
1364!
1365!--          Building strings with index and coordinate information of the
1366!--          slices
1367             DO  WHILE ( section(i,2) /= -9999 )
1368
1369                WRITE (section_chr,'(I5)')  section(i,2)
1370                section_chr = ADJUSTL( section_chr )
1371                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
1372
1373                WRITE (coor_chr,'(F10.1)')  section(i,2) * dy
1374                coor_chr = ADJUSTL( coor_chr )
1375                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
1376
1377                i = i + 1
1378             ENDDO
1379             IF ( av == 0 )  THEN
1380                WRITE ( io, 335 )  'XZ', do2d_xz, dt_do2d_xz, &
1381                                   TRIM( begin_chr ), 'j', TRIM( slices ), &
1382                                   TRIM( coordinates )
1383                IF ( skip_time_do2d_xz /= 0.0_wp )  THEN
1384                   WRITE ( io, 339 )  skip_time_do2d_xz
1385                ENDIF
1386             ELSE
1387                WRITE ( io, 342 )  'XZ', do2d_xz, dt_data_output_av, &
1388                                   TRIM( begin_chr ), averaging_interval, &
1389                                   dt_averaging_input, 'j', TRIM( slices ), &
1390                                   TRIM( coordinates )
1391                IF ( skip_time_data_output_av /= 0.0_wp )  THEN
1392                   WRITE ( io, 339 )  skip_time_data_output_av
1393                ENDIF
1394             ENDIF
1395             IF ( netcdf_data_format > 4 )  THEN
1396                WRITE ( io, 352 )  ntdim_2d_xz(av)
1397             ELSE
1398                WRITE ( io, 353 )
1399             ENDIF
1400          ENDIF
1401
1402          IF ( do2d_yz /= ''  .AND.  section(1,3) /= -9999 )  THEN
1403             i = 1
1404             slices = '/'
1405             coordinates = '/'
1406!
1407!--          Building strings with index and coordinate information of the
1408!--          slices
1409             DO  WHILE ( section(i,3) /= -9999 )
1410
1411                WRITE (section_chr,'(I5)')  section(i,3)
1412                section_chr = ADJUSTL( section_chr )
1413                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
1414
1415                WRITE (coor_chr,'(F10.1)')  section(i,3) * dx
1416                coor_chr = ADJUSTL( coor_chr )
1417                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
1418
1419                i = i + 1
1420             ENDDO
1421             IF ( av == 0 )  THEN
1422                WRITE ( io, 335 )  'YZ', do2d_yz, dt_do2d_yz, &
1423                                   TRIM( begin_chr ), 'i', TRIM( slices ), &
1424                                   TRIM( coordinates )
1425                IF ( skip_time_do2d_yz /= 0.0_wp )  THEN
1426                   WRITE ( io, 339 )  skip_time_do2d_yz
1427                ENDIF
1428             ELSE
1429                WRITE ( io, 342 )  'YZ', do2d_yz, dt_data_output_av, &
1430                                   TRIM( begin_chr ), averaging_interval, &
1431                                   dt_averaging_input, 'i', TRIM( slices ), &
1432                                   TRIM( coordinates )
1433                IF ( skip_time_data_output_av /= 0.0_wp )  THEN
1434                   WRITE ( io, 339 )  skip_time_data_output_av
1435                ENDIF
1436             ENDIF
1437             IF ( netcdf_data_format > 4 )  THEN
1438                WRITE ( io, 352 )  ntdim_2d_yz(av)
1439             ELSE
1440                WRITE ( io, 353 )
1441             ENDIF
1442          ENDIF
1443
1444       ENDIF
1445
1446    ENDDO
1447
1448!
1449!-- 3d-arrays
1450    DO  av = 0, 1
1451
1452       i = 1
1453       do3d_chr = ''
1454       DO  WHILE ( do3d(av,i) /= ' ' )
1455
1456          do3d_chr = TRIM( do3d_chr ) // ' ' // TRIM( do3d(av,i) ) // ','
1457          i = i + 1
1458
1459       ENDDO
1460
1461       IF ( do3d_chr /= '' )  THEN
1462          IF ( av == 0 )  THEN
1463             WRITE ( io, 336 )  ''
1464          ELSE
1465             WRITE ( io, 336 )  '(time-averaged)'
1466          ENDIF
1467
1468          output_format = netcdf_data_format_string
1469          IF ( netcdf_deflate == 0 )  THEN
1470             WRITE ( io, 344 )  output_format
1471          ELSE
1472             WRITE ( io, 354 )  TRIM( output_format ), netcdf_deflate
1473          ENDIF
1474
1475          IF ( do3d_at_begin )  THEN
1476             begin_chr = 'and at the start'
1477          ELSE
1478             begin_chr = ''
1479          ENDIF
1480          IF ( av == 0 )  THEN
1481             WRITE ( io, 337 )  do3d_chr, dt_do3d, TRIM( begin_chr ), &
1482                                zu(nz_do3d), nz_do3d
1483          ELSE
1484             WRITE ( io, 343 )  do3d_chr, dt_data_output_av,           &
1485                                TRIM( begin_chr ), averaging_interval, &
1486                                dt_averaging_input, zu(nz_do3d), nz_do3d
1487          ENDIF
1488
1489          IF ( netcdf_data_format > 4 )  THEN
1490             WRITE ( io, 352 )  ntdim_3d(av)
1491          ELSE
1492             WRITE ( io, 353 )
1493          ENDIF
1494
1495          IF ( av == 0 )  THEN
1496             IF ( skip_time_do3d /= 0.0_wp )  THEN
1497                WRITE ( io, 339 )  skip_time_do3d
1498             ENDIF
1499          ELSE
1500             IF ( skip_time_data_output_av /= 0.0_wp )  THEN
1501                WRITE ( io, 339 )  skip_time_data_output_av
1502             ENDIF
1503          ENDIF
1504
1505       ENDIF
1506
1507    ENDDO
1508
1509!
1510!-- masked arrays
1511    IF ( masks > 0 )  WRITE ( io, 345 )  &
1512         mask_scale_x, mask_scale_y, mask_scale_z
1513    DO  mid = 1, masks
1514       DO  av = 0, 1
1515
1516          i = 1
1517          domask_chr = ''
1518          DO  WHILE ( domask(mid,av,i) /= ' ' )
1519             domask_chr = TRIM( domask_chr ) // ' ' //  &
1520                          TRIM( domask(mid,av,i) ) // ','
1521             i = i + 1
1522          ENDDO
1523
1524          IF ( domask_chr /= '' )  THEN
1525             IF ( av == 0 )  THEN
1526                WRITE ( io, 346 )  '', mid
1527             ELSE
1528                WRITE ( io, 346 )  ' (time-averaged)', mid
1529             ENDIF
1530
1531             output_format = netcdf_data_format_string
1532!--          Parallel output not implemented for mask data, hence
1533!--          output_format must be adjusted.
1534             IF ( netcdf_data_format == 5 ) output_format = 'netCDF4/HDF5'
1535             IF ( netcdf_data_format == 6 ) output_format = 'netCDF4/HDF5 classic'
1536             IF ( netcdf_deflate == 0 )  THEN
1537                WRITE ( io, 344 )  output_format
1538             ELSE
1539                WRITE ( io, 354 )  TRIM( output_format ), netcdf_deflate
1540             ENDIF
1541
1542             IF ( av == 0 )  THEN
1543                WRITE ( io, 347 )  domask_chr, dt_domask(mid)
1544             ELSE
1545                WRITE ( io, 348 )  domask_chr, dt_data_output_av, &
1546                                   averaging_interval, dt_averaging_input
1547             ENDIF
1548
1549             IF ( av == 0 )  THEN
1550                IF ( skip_time_domask(mid) /= 0.0_wp )  THEN
1551                   WRITE ( io, 339 )  skip_time_domask(mid)
1552                ENDIF
1553             ELSE
1554                IF ( skip_time_data_output_av /= 0.0_wp )  THEN
1555                   WRITE ( io, 339 )  skip_time_data_output_av
1556                ENDIF
1557             ENDIF
1558!
1559!--          output locations
1560             DO  dim = 1, 3
1561                IF ( mask(mid,dim,1) >= 0.0_wp )  THEN
1562                   count = 0
1563                   DO  WHILE ( mask(mid,dim,count+1) >= 0.0_wp )
1564                      count = count + 1
1565                   ENDDO
1566                   WRITE ( io, 349 )  dir(dim), dir(dim), mid, dir(dim), &
1567                                      mask(mid,dim,:count)
1568                ELSEIF ( mask_loop(mid,dim,1) < 0.0_wp .AND.  &
1569                         mask_loop(mid,dim,2) < 0.0_wp .AND.  &
1570                         mask_loop(mid,dim,3) == 0.0_wp )  THEN
1571                   WRITE ( io, 350 )  dir(dim), dir(dim)
1572                ELSEIF ( mask_loop(mid,dim,3) == 0.0_wp )  THEN
1573                   WRITE ( io, 351 )  dir(dim), dir(dim), mid, dir(dim), &
1574                                      mask_loop(mid,dim,1:2)
1575                ELSE
1576                   WRITE ( io, 351 )  dir(dim), dir(dim), mid, dir(dim), &
1577                                      mask_loop(mid,dim,1:3)
1578                ENDIF
1579             ENDDO
1580          ENDIF
1581
1582       ENDDO
1583    ENDDO
1584
1585!
1586!-- Timeseries
1587    IF ( dt_dots /= 9999999.9_wp )  THEN
1588       WRITE ( io, 340 )
1589
1590       output_format = netcdf_data_format_string
1591       IF ( netcdf_deflate == 0 )  THEN
1592          WRITE ( io, 344 )  output_format
1593       ELSE
1594          WRITE ( io, 354 )  TRIM( output_format ), netcdf_deflate
1595       ENDIF
1596       WRITE ( io, 341 )  dt_dots
1597    ENDIF
1598
1599#if defined( __dvrp_graphics )
1600!
1601!-- Dvrp-output
1602    IF ( dt_dvrp /= 9999999.9_wp )  THEN
1603       WRITE ( io, 360 )  dt_dvrp, TRIM( dvrp_output ), TRIM( dvrp_host ), &
1604                          TRIM( dvrp_username ), TRIM( dvrp_directory )
1605       i = 1
1606       l = 0
1607       m = 0
1608       DO WHILE ( mode_dvrp(i) /= ' ' )
1609          IF ( mode_dvrp(i)(1:10) == 'isosurface' )  THEN
1610             READ ( mode_dvrp(i), '(10X,I2)' )  j
1611             l = l + 1
1612             IF ( do3d(0,j) /= ' ' )  THEN
1613                WRITE ( io, 361 )  TRIM( do3d(0,j) ), threshold(l), &
1614                                   isosurface_color(:,l)
1615             ENDIF
1616          ELSEIF ( mode_dvrp(i)(1:6) == 'slicer' )  THEN
1617             READ ( mode_dvrp(i), '(6X,I2)' )  j
1618             m = m + 1
1619             IF ( do2d(0,j) /= ' ' )  THEN
1620                WRITE ( io, 362 )  TRIM( do2d(0,j) ), &
1621                                   slicer_range_limits_dvrp(:,m)
1622             ENDIF
1623          ENDIF
1624          i = i + 1
1625       ENDDO
1626
1627       WRITE ( io, 365 )  groundplate_color, superelevation_x, &
1628                          superelevation_y, superelevation, clip_dvrp_l, &
1629                          clip_dvrp_r, clip_dvrp_s, clip_dvrp_n
1630
1631       IF ( TRIM( topography ) /= 'flat' )  THEN
1632          WRITE ( io, 366 )  topography_color
1633          IF ( cluster_size > 1 )  THEN
1634             WRITE ( io, 367 )  cluster_size
1635          ENDIF
1636       ENDIF
1637
1638    ENDIF
1639#endif
1640
1641!
1642!-- Output of spectra related quantities
1643    IF ( calculate_spectra )  CALL spectra_header( io )
1644
1645    WRITE ( io, 99 )
1646
1647!
1648!-- Physical quantities
1649    WRITE ( io, 400 )
1650
1651!
1652!-- Geostrophic parameters
1653    WRITE ( io, 410 )  phi, omega, f, fs
1654
1655!
1656!-- Other quantities
1657    WRITE ( io, 411 )  g
1658
1659    WRITE ( io, 412 )  TRIM( reference_state )
1660    IF ( use_single_reference_value )  THEN
1661       IF ( ocean )  THEN
1662          WRITE ( io, 413 )  prho_reference
1663       ELSE
1664          WRITE ( io, 414 )  pt_reference
1665       ENDIF
1666    ENDIF
1667
1668!
1669!-- Cloud physics parameters
1670    IF ( cloud_physics )  THEN
1671       WRITE ( io, 415 )
1672       WRITE ( io, 416 ) surface_pressure, r_d, rho_surface, cp, l_v
1673       IF ( microphysics_seifert )  THEN
1674          WRITE ( io, 510 ) 1.0E-6_wp * nc_const
1675          WRITE ( io, 511 ) c_sedimentation
1676       ENDIF
1677    ENDIF
1678
1679!
1680!-- Cloud physcis parameters / quantities / numerical methods
1681    WRITE ( io, 430 )
1682    IF ( humidity .AND. .NOT. cloud_physics .AND. .NOT. cloud_droplets)  THEN
1683       WRITE ( io, 431 )
1684    ELSEIF ( humidity  .AND.  cloud_physics )  THEN
1685       WRITE ( io, 432 )
1686       IF ( cloud_top_radiation )  WRITE ( io, 132 )
1687       IF ( microphysics_kessler )  THEN
1688          WRITE ( io, 133 )
1689       ELSEIF ( microphysics_seifert )  THEN
1690          IF ( cloud_water_sedimentation )  WRITE ( io, 506 )
1691          WRITE ( io, 505 )
1692          IF ( collision_turbulence )  WRITE ( io, 507 )
1693          IF ( ventilation_effect )  WRITE ( io, 508 )
1694          IF ( limiter_sedimentation )  WRITE ( io, 509 )
1695       ENDIF
1696    ELSEIF ( humidity  .AND.  cloud_droplets )  THEN
1697       WRITE ( io, 433 )
1698       IF ( curvature_solution_effects )  WRITE ( io, 434 )
1699       IF ( collision_kernel /= 'none' )  THEN
1700          WRITE ( io, 435 )  TRIM( collision_kernel )
1701          IF ( collision_kernel(6:9) == 'fast' )  THEN
1702             WRITE ( io, 436 )  radius_classes, dissipation_classes
1703          ENDIF
1704       ELSE
1705          WRITE ( io, 437 )
1706       ENDIF
1707    ENDIF
1708
1709!
1710!-- LES / turbulence parameters
1711    WRITE ( io, 450 )
1712
1713!--
1714! ... LES-constants used must still be added here
1715!--
1716    IF ( constant_diffusion )  THEN
1717       WRITE ( io, 451 )  km_constant, km_constant/prandtl_number, &
1718                          prandtl_number
1719    ENDIF
1720    IF ( .NOT. constant_diffusion)  THEN
1721       IF ( e_init > 0.0_wp )  WRITE ( io, 455 )  e_init
1722       IF ( e_min > 0.0_wp )  WRITE ( io, 454 )  e_min
1723       IF ( wall_adjustment )  WRITE ( io, 453 )  wall_adjustment_factor
1724    ENDIF
1725
1726!
1727!-- Special actions during the run
1728    WRITE ( io, 470 )
1729    IF ( create_disturbances )  THEN
1730       WRITE ( io, 471 )  dt_disturb, disturbance_amplitude,                   &
1731                          zu(disturbance_level_ind_b), disturbance_level_ind_b,&
1732                          zu(disturbance_level_ind_t), disturbance_level_ind_t
1733       IF ( .NOT. bc_lr_cyc  .OR.  .NOT. bc_ns_cyc )  THEN
1734          WRITE ( io, 472 )  inflow_disturbance_begin, inflow_disturbance_end
1735       ELSE
1736          WRITE ( io, 473 )  disturbance_energy_limit
1737       ENDIF
1738       WRITE ( io, 474 )  TRIM( random_generator )
1739    ENDIF
1740    IF ( pt_surface_initial_change /= 0.0_wp )  THEN
1741       WRITE ( io, 475 )  pt_surface_initial_change
1742    ENDIF
1743    IF ( humidity  .AND.  q_surface_initial_change /= 0.0_wp )  THEN
1744       WRITE ( io, 476 )  q_surface_initial_change       
1745    ENDIF
1746    IF ( passive_scalar  .AND.  q_surface_initial_change /= 0.0_wp )  THEN
1747       WRITE ( io, 477 )  q_surface_initial_change       
1748    ENDIF
1749
1750    IF ( particle_advection )  THEN
1751!
1752!--    Particle attributes
1753       WRITE ( io, 480 )  particle_advection_start, dt_prel, bc_par_lr, &
1754                          bc_par_ns, bc_par_b, bc_par_t, particle_maximum_age, &
1755                          end_time_prel
1756       IF ( use_sgs_for_particles )  WRITE ( io, 488 )  dt_min_part
1757       IF ( random_start_position )  WRITE ( io, 481 )
1758       IF ( seed_follows_topography )  WRITE ( io, 496 )
1759       IF ( particles_per_point > 1 )  WRITE ( io, 489 )  particles_per_point
1760       WRITE ( io, 495 )  total_number_of_particles
1761       IF ( dt_write_particle_data /= 9999999.9_wp )  THEN
1762          WRITE ( io, 485 )  dt_write_particle_data
1763          IF ( netcdf_data_format > 1 )  THEN
1764             output_format = 'netcdf (64 bit offset) and binary'
1765          ELSE
1766             output_format = 'netcdf and binary'
1767          ENDIF
1768          IF ( netcdf_deflate == 0 )  THEN
1769             WRITE ( io, 344 )  output_format
1770          ELSE
1771             WRITE ( io, 354 )  TRIM( output_format ), netcdf_deflate
1772          ENDIF
1773       ENDIF
1774       IF ( dt_dopts /= 9999999.9_wp )  WRITE ( io, 494 )  dt_dopts
1775       IF ( write_particle_statistics )  WRITE ( io, 486 )
1776
1777       WRITE ( io, 487 )  number_of_particle_groups
1778
1779       DO  i = 1, number_of_particle_groups
1780          IF ( i == 1  .AND.  density_ratio(i) == 9999999.9_wp )  THEN
1781             WRITE ( io, 490 )  i, 0.0_wp
1782             WRITE ( io, 492 )
1783          ELSE
1784             WRITE ( io, 490 )  i, radius(i)
1785             IF ( density_ratio(i) /= 0.0_wp )  THEN
1786                WRITE ( io, 491 )  density_ratio(i)
1787             ELSE
1788                WRITE ( io, 492 )
1789             ENDIF
1790          ENDIF
1791          WRITE ( io, 493 )  psl(i), psr(i), pss(i), psn(i), psb(i), pst(i), &
1792                             pdx(i), pdy(i), pdz(i)
1793          IF ( .NOT. vertical_particle_advection(i) )  WRITE ( io, 482 )
1794       ENDDO
1795
1796    ENDIF
1797
1798
1799!
1800!-- Parameters of 1D-model
1801    IF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
1802       WRITE ( io, 500 )  end_time_1d, dt_run_control_1d, dt_pr_1d, &
1803                          mixing_length_1d, dissipation_1d
1804       IF ( damp_level_ind_1d /= nzt+1 )  THEN
1805          WRITE ( io, 502 )  zu(damp_level_ind_1d), damp_level_ind_1d
1806       ENDIF
1807    ENDIF
1808
1809!
1810!-- User-defined information
1811    CALL user_header( io )
1812
1813    WRITE ( io, 99 )
1814
1815!
1816!-- Write buffer contents to disc immediately
1817    FLUSH( io )
1818
1819!
1820!-- Here the FORMATs start
1821
1822 99 FORMAT (1X,78('-'))
1823100 FORMAT (/1X,'******************************',4X,44('-')/        &
1824            1X,'* ',A,' *',4X,A/                               &
1825            1X,'******************************',4X,44('-'))
1826101 FORMAT (35X,'coupled run using MPI-',I1,': ',A/ &
1827            35X,42('-'))
1828102 FORMAT (/' Date:                 ',A8,4X,'Run:       ',A20/      &
1829            ' Time:                 ',A8,4X,'Run-No.:   ',I2.2/     &
1830            ' Run on host:        ',A10)
1831#if defined( __parallel )
1832103 FORMAT (' Number of PEs:',10X,I6,4X,'Processor grid (x,y): (',I4,',',I4, &
1833              ')',1X,A)
1834104 FORMAT (' Number of PEs:',10X,I6,4X,'Tasks:',I4,'   threads per task:',I4/ &
1835              35X,'Processor grid (x,y): (',I4,',',I4,')',1X,A)
1836105 FORMAT (35X,'One additional PE is used to handle'/37X,'the dvrp output!')
1837106 FORMAT (35X,'A 1d-decomposition along x is forced'/ &
1838            35X,'because the job is running on an SMP-cluster')
1839107 FORMAT (35X,'A 1d-decomposition along ',A,' is used')
1840108 FORMAT (35X,'Max. # of parallel I/O streams is ',I5)
1841109 FORMAT (35X,'Precursor run for coupled atmos-ocean run'/ &
1842            35X,42('-'))
1843114 FORMAT (35X,'Coupled atmosphere-ocean run following'/ &
1844            35X,'independent precursor runs'/             &
1845            35X,42('-'))
1846117 FORMAT (' Accelerator boards / node:  ',I2)
1847#endif
1848110 FORMAT (/' Numerical Schemes:'/ &
1849             ' -----------------'/)
1850111 FORMAT (' --> Solve perturbation pressure via FFT using ',A,' routines')
1851112 FORMAT (' --> Solve perturbation pressure via SOR-Red/Black-Schema'/ &
1852            '     Iterations (initial/other): ',I3,'/',I3,'  omega =',F6.3)
1853113 FORMAT (' --> Momentum advection via Piascek-Williams-Scheme (Form C3)', &
1854                  ' or Upstream')
1855115 FORMAT ('     FFT and transpositions are overlapping')
1856116 FORMAT (' --> Scalar advection via Piascek-Williams-Scheme (Form C3)', &
1857                  ' or Upstream')
1858118 FORMAT (' --> Scalar advection via Bott-Chlond-Scheme')
1859119 FORMAT (' --> Galilei-Transform applied to horizontal advection:'/ &
1860            '     translation velocity = ',A/ &
1861            '     distance advected ',A,':  ',F8.3,' km(x)  ',F8.3,' km(y)')
1862120 FORMAT (' Accelerator boards: ',8X,I2)
1863122 FORMAT (' --> Time differencing scheme: ',A)
1864123 FORMAT (' --> Rayleigh-Damping active, starts ',A,' z = ',F8.2,' m'/ &
1865            '     maximum damping coefficient:',F6.3, ' 1/s')
1866129 FORMAT (' --> Additional prognostic equation for the specific humidity')
1867130 FORMAT (' --> Additional prognostic equation for the total water content')
1868131 FORMAT (' --> No pt-equation solved. Neutral stratification with pt = ', &
1869                  F6.2, ' K assumed')
1870132 FORMAT ('     Parameterization of long-wave radiation processes via'/ &
1871            '     effective emissivity scheme')
1872133 FORMAT ('     Precipitation parameterization via Kessler-Scheme')
1873134 FORMAT (' --> Additional prognostic equation for a passive scalar')
1874135 FORMAT (' --> Solve perturbation pressure via ',A,' method (', &
1875                  A,'-cycle)'/ &
1876            '     number of grid levels:                   ',I2/ &
1877            '     Gauss-Seidel red/black iterations:       ',I2)
1878136 FORMAT ('     gridpoints of coarsest subdomain (x,y,z): (',I3,',',I3,',', &
1879                  I3,')')
1880137 FORMAT ('     level data gathered on PE0 at level:     ',I2/ &
1881            '     gridpoints of coarsest subdomain (x,y,z): (',I3,',',I3,',', &
1882                  I3,')'/ &
1883            '     gridpoints of coarsest domain (x,y,z):    (',I3,',',I3,',', &
1884                  I3,')')
1885139 FORMAT (' --> Loop optimization method: ',A)
1886140 FORMAT ('     maximum residual allowed:                ',E10.3)
1887141 FORMAT ('     fixed number of multigrid cycles:        ',I4)
1888142 FORMAT ('     perturbation pressure is calculated at every Runge-Kutta ', &
1889                  'step')
1890143 FORMAT ('     Euler/upstream scheme is used for the SGS turbulent ', &
1891                  'kinetic energy')
1892144 FORMAT ('     masking method is used')
1893150 FORMAT (' --> Volume flow at the right and north boundary will be ', &
1894                  'conserved'/ &
1895            '     using the ',A,' mode')
1896151 FORMAT ('     with u_bulk = ',F7.3,' m/s and v_bulk = ',F7.3,' m/s')
1897152 FORMAT (' --> External pressure gradient directly prescribed by the user:',&
1898           /'     ',2(1X,E12.5),'Pa/m in x/y direction', &
1899           /'     starting from dp_level_b =', F8.3, 'm', A /)
1900160 FORMAT (//' Large scale forcing and nudging:'/ &
1901              ' -------------------------------'/)
1902161 FORMAT (' --> No large scale forcing from external is used (default) ')
1903162 FORMAT (' --> Large scale forcing from external file LSF_DATA is used: ')
1904163 FORMAT ('     - large scale advection tendencies ')
1905164 FORMAT ('     - large scale subsidence velocity w_subs ')
1906165 FORMAT ('     - large scale subsidence tendencies ')
1907167 FORMAT ('     - and geostrophic wind components ug and vg')
1908168 FORMAT (' --> Large-scale vertical motion is used in the ', &
1909                  'prognostic equation(s) for')
1910169 FORMAT ('     the scalar(s) only')
1911170 FORMAT (' --> Nudging is used')
1912171 FORMAT (' --> No nudging is used (default) ')
1913180 FORMAT ('     - prescribed surface values for temperature')
1914181 FORMAT ('     - prescribed surface fluxes for temperature')
1915182 FORMAT ('     - prescribed surface values for humidity')
1916183 FORMAT ('     - prescribed surface fluxes for humidity')
1917200 FORMAT (//' Run time and time step information:'/ &
1918             ' ----------------------------------'/)
1919201 FORMAT ( ' Timestep:             variable     maximum value: ',F6.3,' s', &
1920             '    CFL-factor:',F5.2)
1921202 FORMAT ( ' Timestep:          dt = ',F6.3,' s'/)
1922203 FORMAT ( ' Start time:          ',F9.3,' s'/ &
1923             ' End time:            ',F9.3,' s')
1924204 FORMAT ( A,F9.3,' s')
1925205 FORMAT ( A,F9.3,' s',5X,'restart every',17X,F9.3,' s')
1926206 FORMAT (/' Time reached:        ',F9.3,' s'/ &
1927             ' CPU-time used:       ',F9.3,' s     per timestep:               ', &
1928               '  ',F9.3,' s'/                                                    &
1929             '                                      per second of simulated tim', &
1930               'e: ',F9.3,' s')
1931207 FORMAT ( ' Coupling start time: ',F9.3,' s')
1932250 FORMAT (//' Computational grid and domain size:'/ &
1933              ' ----------------------------------'// &
1934              ' Grid length:      dx =    ',F7.3,' m    dy =    ',F7.3, &
1935              ' m    dz =    ',F7.3,' m'/ &
1936              ' Domain size:       x = ',F10.3,' m     y = ',F10.3, &
1937              ' m  z(u) = ',F10.3,' m'/)
1938252 FORMAT (' dz constant up to ',F10.3,' m (k=',I4,'), then stretched by', &
1939              ' factor:',F6.3/ &
1940            ' maximum dz not to be exceeded is dz_max = ',F10.3,' m'/)
1941254 FORMAT (' Number of gridpoints (x,y,z):  (0:',I4,', 0:',I4,', 0:',I4,')'/ &
1942            ' Subdomain size (x,y,z):        (  ',I4,',   ',I4,',   ',I4,')'/)
1943260 FORMAT (/' The model has a slope in x-direction. Inclination angle: ',F6.2,&
1944             ' degrees')
1945270 FORMAT (//' Topography information:'/ &
1946              ' ----------------------'// &
1947              1X,'Topography: ',A)
1948271 FORMAT (  ' Building size (x/y/z) in m: ',F5.1,' / ',F5.1,' / ',F5.1/ &
1949              ' Horizontal index bounds (l/r/s/n): ',I4,' / ',I4,' / ',I4, &
1950                ' / ',I4)
1951272 FORMAT (  ' Single quasi-2D street canyon of infinite length in ',A, &
1952              ' direction' / &
1953              ' Canyon height: ', F6.2, 'm, ch = ', I4, '.'      / &
1954              ' Canyon position (',A,'-walls): cxl = ', I4,', cxr = ', I4, '.')
1955278 FORMAT (' Topography grid definition convention:'/ &
1956            ' cell edge (staggered grid points'/  &
1957            ' (u in x-direction, v in y-direction))' /)
1958279 FORMAT (' Topography grid definition convention:'/ &
1959            ' cell center (scalar grid points)' /)
1960300 FORMAT (//' Boundary conditions:'/ &
1961             ' -------------------'// &
1962             '                     p                    uv             ', &
1963             '                     pt'// &
1964             ' B. bound.: ',A/ &
1965             ' T. bound.: ',A)
1966301 FORMAT (/'                     ',A// &
1967             ' B. bound.: ',A/ &
1968             ' T. bound.: ',A)
1969303 FORMAT (/' Bottom surface fluxes are used in diffusion terms at k=1')
1970304 FORMAT (/' Top surface fluxes are used in diffusion terms at k=nzt')
1971305 FORMAT (//'    Prandtl-Layer between bottom surface and first ', &
1972               'computational u,v-level:'// &
1973             '       zp = ',F6.2,' m   z0 =',F7.4,' m   z0h =',F8.5,&
1974             ' m   kappa =',F5.2/ &
1975             '       Rif value range:   ',F8.2,' <= rif <=',F6.2)
1976306 FORMAT ('       Predefined constant heatflux:   ',F9.6,' K m/s')
1977307 FORMAT ('       Heatflux has a random normal distribution')
1978308 FORMAT ('       Predefined surface temperature')
1979309 FORMAT ('       Predefined constant salinityflux:   ',F9.6,' psu m/s')
1980310 FORMAT (//'    1D-Model:'// &
1981             '       Rif value range:   ',F6.2,' <= rif <=',F6.2)
1982311 FORMAT ('       Predefined constant humidity flux: ',E10.3,' m/s')
1983312 FORMAT ('       Predefined surface humidity')
1984313 FORMAT ('       Predefined constant scalar flux: ',E10.3,' kg/(m**2 s)')
1985314 FORMAT ('       Predefined scalar value at the surface')
1986315 FORMAT ('       Humidity / scalar flux at top surface is 0.0')
1987316 FORMAT ('       Sensible heatflux and momentum flux from coupled ', &
1988                    'atmosphere model')
1989317 FORMAT (//' Lateral boundaries:'/ &
1990            '       left/right:  ',A/    &
1991            '       north/south: ',A)
1992318 FORMAT (/'       use_cmax: ',L1 / &
1993            '       pt damping layer width = ',F8.2,' m, pt ', &
1994                    'damping factor =',F7.4)
1995319 FORMAT ('       turbulence recycling at inflow switched on'/ &
1996            '       width of recycling domain: ',F7.1,' m   grid index: ',I4/ &
1997            '       inflow damping height: ',F6.1,' m   width: ',F6.1,' m')
1998320 FORMAT ('       Predefined constant momentumflux:  u: ',F9.6,' m**2/s**2'/ &
1999            '                                          v: ',F9.6,' m**2/s**2')
2000321 FORMAT (//' Initial profiles:'/ &
2001              ' ----------------')
2002322 FORMAT ('       turbulence recycling at inflow switched on'/ &
2003            '       y shift of the recycled inflow turbulence switched on'/ &
2004            '       width of recycling domain: ',F7.1,' m   grid index: ',I4/ &
2005            '       inflow damping height: ',F6.1,' m   width: ',F6.1,' m'/)
2006325 FORMAT (//' List output:'/ &
2007             ' -----------'//  &
2008            '    1D-Profiles:'/    &
2009            '       Output every             ',F8.2,' s')
2010326 FORMAT ('       Time averaged over       ',F8.2,' s'/ &
2011            '       Averaging input every    ',F8.2,' s')
2012330 FORMAT (//' Data output:'/ &
2013             ' -----------'/)
2014331 FORMAT (/'    1D-Profiles:')
2015332 FORMAT (/'       ',A)
2016333 FORMAT ('       Output every             ',F8.2,' s',/ &
2017            '       Time averaged over       ',F8.2,' s'/ &
2018            '       Averaging input every    ',F8.2,' s')
2019334 FORMAT (/'    2D-Arrays',A,':')
2020335 FORMAT (/'       ',A2,'-cross-section  Arrays: ',A/ &
2021            '       Output every             ',F8.2,' s  ',A/ &
2022            '       Cross sections at ',A1,' = ',A/ &
2023            '       scalar-coordinates:   ',A,' m'/)
2024336 FORMAT (/'    3D-Arrays',A,':')
2025337 FORMAT (/'       Arrays: ',A/ &
2026            '       Output every             ',F8.2,' s  ',A/ &
2027            '       Upper output limit at    ',F8.2,' m  (GP ',I4,')'/)
2028339 FORMAT ('       No output during initial ',F8.2,' s')
2029340 FORMAT (/'    Time series:')
2030341 FORMAT ('       Output every             ',F8.2,' s'/)
2031342 FORMAT (/'       ',A2,'-cross-section  Arrays: ',A/ &
2032            '       Output every             ',F8.2,' s  ',A/ &
2033            '       Time averaged over       ',F8.2,' s'/ &
2034            '       Averaging input every    ',F8.2,' s'/ &
2035            '       Cross sections at ',A1,' = ',A/ &
2036            '       scalar-coordinates:   ',A,' m'/)
2037343 FORMAT (/'       Arrays: ',A/ &
2038            '       Output every             ',F8.2,' s  ',A/ &
2039            '       Time averaged over       ',F8.2,' s'/ &
2040            '       Averaging input every    ',F8.2,' s'/ &
2041            '       Upper output limit at    ',F8.2,' m  (GP ',I4,')'/)
2042344 FORMAT ('       Output format: ',A/)
2043345 FORMAT (/'    Scaling lengths for output locations of all subsequent mask IDs:',/ &
2044            '       mask_scale_x (in x-direction): ',F9.3, ' m',/ &
2045            '       mask_scale_y (in y-direction): ',F9.3, ' m',/ &
2046            '       mask_scale_z (in z-direction): ',F9.3, ' m' )
2047346 FORMAT (/'    Masked data output',A,' for mask ID ',I2, ':')
2048347 FORMAT ('       Variables: ',A/ &
2049            '       Output every             ',F8.2,' s')
2050348 FORMAT ('       Variables: ',A/ &
2051            '       Output every             ',F8.2,' s'/ &
2052            '       Time averaged over       ',F8.2,' s'/ &
2053            '       Averaging input every    ',F8.2,' s')
2054349 FORMAT (/'       Output locations in ',A,'-direction in multiples of ', &
2055            'mask_scale_',A,' predefined by array mask_',I2.2,'_',A,':'/ &
2056            13('       ',8(F8.2,',')/) )
2057350 FORMAT (/'       Output locations in ',A,'-direction: ', &
2058            'all gridpoints along ',A,'-direction (default).' )
2059351 FORMAT (/'       Output locations in ',A,'-direction in multiples of ', &
2060            'mask_scale_',A,' constructed from array mask_',I2.2,'_',A,'_loop:'/ &
2061            '          loop begin:',F8.2,', end:',F8.2,', stride:',F8.2 )
2062352 FORMAT  (/'       Number of output time levels allowed: ',I3 /)
2063353 FORMAT  (/'       Number of output time levels allowed: unlimited' /)
2064354 FORMAT ('       Output format: ',A, '   compressed with level: ',I1/)
2065#if defined( __dvrp_graphics )
2066360 FORMAT ('    Plot-Sequence with dvrp-software:'/ &
2067            '       Output every      ',F7.1,' s'/ &
2068            '       Output mode:      ',A/ &
2069            '       Host / User:      ',A,' / ',A/ &
2070            '       Directory:        ',A// &
2071            '       The sequence contains:')
2072361 FORMAT (/'       Isosurface of "',A,'"    Threshold value: ', E12.3/ &
2073            '          Isosurface color: (',F4.2,',',F4.2,',',F4.2,') (R,G,B)')
2074362 FORMAT (/'       Slicer plane ',A/ &
2075            '       Slicer limits: [',F6.2,',',F6.2,']')
2076365 FORMAT (/'       Groundplate color: (',F4.2,',',F4.2,',',F4.2,') (R,G,B)'/ &
2077            '       Superelevation along (x,y,z): (',F4.1,',',F4.1,',',F4.1, &
2078                     ')'/ &
2079            '       Clipping limits: from x = ',F9.1,' m to x = ',F9.1,' m'/ &
2080            '                        from y = ',F9.1,' m to y = ',F9.1,' m')
2081366 FORMAT (/'       Topography color: (',F4.2,',',F4.2,',',F4.2,') (R,G,B)')
2082367 FORMAT ('       Polygon reduction for topography: cluster_size = ', I1)
2083#endif
2084400 FORMAT (//' Physical quantities:'/ &
2085              ' -------------------'/)
2086410 FORMAT ('    Geograph. latitude  :   phi    = ',F4.1,' degr'/   &
2087            '    Angular velocity    :   omega  =',E10.3,' rad/s'/  &
2088            '    Coriolis parameter  :   f      = ',F9.6,' 1/s'/    &
2089            '                            f*     = ',F9.6,' 1/s')
2090411 FORMAT (/'    Gravity             :   g      = ',F4.1,' m/s**2')
2091412 FORMAT (/'    Reference state used in buoyancy terms: ',A)
2092413 FORMAT ('       Reference density in buoyancy terms: ',F8.3,' kg/m**3')
2093414 FORMAT ('       Reference temperature in buoyancy terms: ',F8.4,' K')
2094415 FORMAT (/' Cloud physics parameters:'/ &
2095             ' ------------------------'/)
2096416 FORMAT ('    Surface pressure   :   p_0   = ',F7.2,' hPa'/      &
2097            '    Gas constant       :   R     = ',F5.1,' J/(kg K)'/ &
2098            '    Density of air     :   rho_0 =',F6.3,' kg/m**3'/  &
2099            '    Specific heat cap. :   c_p   = ',F6.1,' J/(kg K)'/ &
2100            '    Vapourization heat :   L_v   =',E9.2,' J/kg')
2101417 FORMAT ('    Geograph. longitude :   lambda = ',F4.1,' degr')
2102418 FORMAT (/'    Day of the year at model start :   day_init      =     ',I3 &
2103            /'    UTC time at model start        :   time_utc_init = ',F7.1' s')
2104420 FORMAT (/'    Characteristic levels of the initial temperature profile:'// &
2105            '       Height:        ',A,'  m'/ &
2106            '       Temperature:   ',A,'  K'/ &
2107            '       Gradient:      ',A,'  K/100m'/ &
2108            '       Gridpoint:     ',A)
2109421 FORMAT (/'    Characteristic levels of the initial humidity profile:'// &
2110            '       Height:      ',A,'  m'/ &
2111            '       Humidity:    ',A,'  kg/kg'/ &
2112            '       Gradient:    ',A,'  (kg/kg)/100m'/ &
2113            '       Gridpoint:   ',A)
2114422 FORMAT (/'    Characteristic levels of the initial scalar profile:'// &
2115            '       Height:                  ',A,'  m'/ &
2116            '       Scalar concentration:    ',A,'  kg/m**3'/ &
2117            '       Gradient:                ',A,'  (kg/m**3)/100m'/ &
2118            '       Gridpoint:               ',A)
2119423 FORMAT (/'    Characteristic levels of the geo. wind component ug:'// &
2120            '       Height:      ',A,'  m'/ &
2121            '       ug:          ',A,'  m/s'/ &
2122            '       Gradient:    ',A,'  1/100s'/ &
2123            '       Gridpoint:   ',A)
2124424 FORMAT (/'    Characteristic levels of the geo. wind component vg:'// &
2125            '       Height:      ',A,'  m'/ &
2126            '       vg:          ',A,'  m/s'/ &
2127            '       Gradient:    ',A,'  1/100s'/ &
2128            '       Gridpoint:   ',A)
2129425 FORMAT (/'    Characteristic levels of the initial salinity profile:'// &
2130            '       Height:     ',A,'  m'/ &
2131            '       Salinity:   ',A,'  psu'/ &
2132            '       Gradient:   ',A,'  psu/100m'/ &
2133            '       Gridpoint:  ',A)
2134426 FORMAT (/'    Characteristic levels of the subsidence/ascent profile:'// &
2135            '       Height:      ',A,'  m'/ &
2136            '       w_subs:      ',A,'  m/s'/ &
2137            '       Gradient:    ',A,'  (m/s)/100m'/ &
2138            '       Gridpoint:   ',A)
2139427 FORMAT (/'    Initial wind profiles (u,v) are interpolated from given'// &
2140                  ' profiles')
2141428 FORMAT (/'    Initial profiles (u, v, pt, q) are taken from file '/ &
2142             '    NUDGING_DATA')
2143430 FORMAT (//' Cloud physics quantities / methods:'/ &
2144              ' ----------------------------------'/)
2145431 FORMAT ('    Humidity is treated as purely passive scalar (no condensati', &
2146                 'on)')
2147432 FORMAT ('    Bulk scheme with liquid water potential temperature and'/ &
2148            '    total water content is used.'/ &
2149            '    Condensation is parameterized via 0% - or 100% scheme.')
2150433 FORMAT ('    Cloud droplets treated explicitly using the Lagrangian part', &
2151                 'icle model')
2152434 FORMAT ('    Curvature and solution effecs are considered for growth of', &
2153                 ' droplets < 1.0E-6 m')
2154435 FORMAT ('    Droplet collision is handled by ',A,'-kernel')
2155436 FORMAT ('       Fast kernel with fixed radius- and dissipation classes ', &
2156                    'are used'/ &
2157            '          number of radius classes:       ',I3,'    interval ', &
2158                       '[1.0E-6,2.0E-4] m'/ &
2159            '          number of dissipation classes:   ',I2,'    interval ', &
2160                       '[0,1000] cm**2/s**3')
2161437 FORMAT ('    Droplet collision is switched off')
2162450 FORMAT (//' LES / Turbulence quantities:'/ &
2163              ' ---------------------------'/)
2164451 FORMAT ('    Diffusion coefficients are constant:'/ &
2165            '    Km = ',F6.2,' m**2/s   Kh = ',F6.2,' m**2/s   Pr = ',F5.2)
2166453 FORMAT ('    Mixing length is limited to',F5.2,' * z')
2167454 FORMAT ('    TKE is not allowed to fall below ',E9.2,' (m/s)**2')
2168455 FORMAT ('    initial TKE is prescribed as ',E9.2,' (m/s)**2')
2169470 FORMAT (//' Actions during the simulation:'/ &
2170              ' -----------------------------'/)
2171471 FORMAT ('    Disturbance impulse (u,v) every :   ',F6.2,' s'/            &
2172            '    Disturbance amplitude           :    ',F5.2, ' m/s'/       &
2173            '    Lower disturbance level         : ',F8.2,' m (GP ',I4,')'/  &
2174            '    Upper disturbance level         : ',F8.2,' m (GP ',I4,')')
2175472 FORMAT ('    Disturbances continued during the run from i/j =',I4, &
2176                 ' to i/j =',I4)
2177473 FORMAT ('    Disturbances cease as soon as the disturbance energy exceeds',&
2178                 F6.3, ' m**2/s**2')
2179474 FORMAT ('    Random number generator used    : ',A/)
2180475 FORMAT ('    The surface temperature is increased (or decreased, ', &
2181                 'respectively, if'/ &
2182            '    the value is negative) by ',F5.2,' K at the beginning of the',&
2183                 ' 3D-simulation'/)
2184476 FORMAT ('    The surface humidity is increased (or decreased, ',&
2185                 'respectively, if the'/ &
2186            '    value is negative) by ',E8.1,' kg/kg at the beginning of', &
2187                 ' the 3D-simulation'/)
2188477 FORMAT ('    The scalar value is increased at the surface (or decreased, ',&
2189                 'respectively, if the'/ &
2190            '    value is negative) by ',E8.1,' kg/m**3 at the beginning of', &
2191                 ' the 3D-simulation'/)
2192480 FORMAT ('    Particles:'/ &
2193            '    ---------'// &
2194            '       Particle advection is active (switched on at t = ', F7.1, &
2195                    ' s)'/ &
2196            '       Start of new particle generations every  ',F6.1,' s'/ &
2197            '       Boundary conditions: left/right: ', A, ' north/south: ', A/&
2198            '                            bottom:     ', A, ' top:         ', A/&
2199            '       Maximum particle age:                 ',F9.1,' s'/ &
2200            '       Advection stopped at t = ',F9.1,' s'/)
2201481 FORMAT ('       Particles have random start positions'/)
2202482 FORMAT ('          Particles are advected only horizontally'/)
2203485 FORMAT ('       Particle data are written on file every ', F9.1, ' s')
2204486 FORMAT ('       Particle statistics are written on file'/)
2205487 FORMAT ('       Number of particle groups: ',I2/)
2206488 FORMAT ('       SGS velocity components are used for particle advection'/ &
2207            '          minimum timestep for advection:', F8.5/)
2208489 FORMAT ('       Number of particles simultaneously released at each ', &
2209                    'point: ', I5/)
2210490 FORMAT ('       Particle group ',I2,':'/ &
2211            '          Particle radius: ',E10.3, 'm')
2212491 FORMAT ('          Particle inertia is activated'/ &
2213            '             density_ratio (rho_fluid/rho_particle) =',F6.3/)
2214492 FORMAT ('          Particles are advected only passively (no inertia)'/)
2215493 FORMAT ('          Boundaries of particle source: x:',F8.1,' - ',F8.1,' m'/&
2216            '                                         y:',F8.1,' - ',F8.1,' m'/&
2217            '                                         z:',F8.1,' - ',F8.1,' m'/&
2218            '          Particle distances:  dx = ',F8.1,' m  dy = ',F8.1, &
2219                       ' m  dz = ',F8.1,' m'/)
2220494 FORMAT ('       Output of particle time series in NetCDF format every ', &
2221                    F8.2,' s'/)
2222495 FORMAT ('       Number of particles in total domain: ',I10/)
2223496 FORMAT ('       Initial vertical particle positions are interpreted ', &
2224                    'as relative to the given topography')
2225500 FORMAT (//' 1D-Model parameters:'/                           &
2226              ' -------------------'//                           &
2227            '    Simulation time:                   ',F8.1,' s'/ &
2228            '    Run-controll output every:         ',F8.1,' s'/ &
2229            '    Vertical profile output every:     ',F8.1,' s'/ &
2230            '    Mixing length calculation:         ',A/         &
2231            '    Dissipation calculation:           ',A/)
2232502 FORMAT ('    Damping layer starts from ',F7.1,' m (GP ',I4,')'/)
2233503 FORMAT (' --> Momentum advection via Wicker-Skamarock-Scheme 5th order')
2234504 FORMAT (' --> Scalar advection via Wicker-Skamarock-Scheme 5th order')
2235505 FORMAT ('    Precipitation parameterization via Seifert-Beheng-Scheme')
2236506 FORMAT ('    Cloud water sedimentation parameterization via Stokes law')
2237507 FORMAT ('    Turbulence effects on precipitation process')
2238508 FORMAT ('    Ventilation effects on evaporation of rain drops')
2239509 FORMAT ('    Slope limiter used for sedimentation process')
2240510 FORMAT ('    Droplet density    :   N_c   = ',F6.1,' 1/cm**3')
2241511 FORMAT ('    Sedimentation Courant number:                  '/&
2242            '                               C_s   =',F4.1,'        ')
2243512 FORMAT (/' Date:                 ',A8,6X,'Run:       ',A20/      &
2244            ' Time:                 ',A8,6X,'Run-No.:   ',I2.2/     &
2245            ' Run on host:        ',A10,6X,'En-No.:    ',I2.2)
2246513 FORMAT (' --> Scalar advection via Wicker-Skamarock-Scheme 5th order ' // & 
2247            '+ monotonic adjustment')
2248600 FORMAT (/' Nesting informations:'/ &
2249            ' --------------------'/ &
2250            ' Nesting mode:                     ',A/ &
2251            ' Nesting-datatransfer mode:        ',A// &
2252            ' Nest id  parent  number   lower left coordinates   name'/ &
2253            ' (*=me)     id    of PEs      x (m)     y (m)' )
2254601 FORMAT (2X,A1,1X,I2.2,6X,I2.2,5X,I5,5X,F8.2,2X,F8.2,5X,A)
2255
2256 END SUBROUTINE header
Note: See TracBrowser for help on using the repository browser.