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

Last change on this file since 2550 was 2550, checked in by boeske, 7 years ago

enable simulations with complex terrain

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