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

Last change on this file since 1307 was 1300, checked in by heinze, 10 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 80.4 KB
Line 
1 SUBROUTINE header
2
3!--------------------------------------------------------------------------------!
4! This file is part of PALM.
5!
6! PALM is free software: you can redistribute it and/or modify it under the terms
7! of the GNU General Public License as published by the Free Software Foundation,
8! either version 3 of the License, or (at your option) any later version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 1997-2012  Leibniz University Hannover
18!--------------------------------------------------------------------------------!
19!
20! Current revisions:
21! -----------------
22!
23!
24! Former revisions:
25! -----------------
26! $Id: header.f90 1300 2014-03-06 13:18:34Z raasch $
27!
28! 1299 2014-03-06 13:15:21Z heinze
29! output for using large_scale subsidence in combination
30! with large_scale_forcing
31! reformatting, more detailed explanations
32!
33! 1241 2013-10-30 11:36:58Z heinze
34! output for nudging + large scale forcing from external file
35!
36! 1216 2013-08-26 09:31:42Z raasch
37! output for transpose_compute_overlap
38!
39! 1212 2013-08-15 08:46:27Z raasch
40! output for poisfft_hybrid removed
41!
42! 1179 2013-06-14 05:57:58Z raasch
43! output of reference_state, use_reference renamed use_single_reference_value
44!
45! 1159 2013-05-21 11:58:22Z fricke
46! +use_cmax
47!
48! 1115 2013-03-26 18:16:16Z hoffmann
49! descriptions for Seifert-Beheng-cloud-physics-scheme added
50!
51! 1111 2013-03-08 23:54:10Z raasch
52! output of accelerator board information
53! ibc_p_b = 2 removed
54!
55! 1108 2013-03-05 07:03:32Z raasch
56! bugfix for r1106
57!
58! 1106 2013-03-04 05:31:38Z raasch
59! some format changes for coupled runs
60!
61! 1092 2013-02-02 11:24:22Z raasch
62! unused variables removed
63!
64! 1036 2012-10-22 13:43:42Z raasch
65! code put under GPL (PALM 3.9)
66!
67! 1031 2012-10-19 14:35:30Z raasch
68! output of netCDF data format modified
69!
70! 1015 2012-09-27 09:23:24Z raasch
71! output of Aajustment of mixing length to the Prandtl mixing length at first
72! grid point above ground removed
73!
74! 1003 2012-09-14 14:35:53Z raasch
75! output of information about equal/unequal subdomain size removed
76!
77! 1001 2012-09-13 14:08:46Z raasch
78! all actions concerning leapfrog- and upstream-spline-scheme removed
79!
80! 978 2012-08-09 08:28:32Z fricke
81! -km_damp_max, outflow_damping_width
82! +pt_damping_factor, pt_damping_width
83! +z0h
84!
85! 964 2012-07-26 09:14:24Z raasch
86! output of profil-related quantities removed
87!
88! 940 2012-07-09 14:31:00Z raasch
89! Output in case of simulations for pure neutral stratification (no pt-equation
90! solved)
91!
92! 927 2012-06-06 19:15:04Z raasch
93! output of masking_method for mg-solver
94!
95! 868 2012-03-28 12:21:07Z raasch
96! translation velocity in Galilean transformation changed to 0.6 * ug
97!
98! 833 2012-02-22 08:55:55Z maronga
99! Adjusted format for leaf area density
100!
101! 828 2012-02-21 12:00:36Z raasch
102! output of dissipation_classes + radius_classes
103!
104! 825 2012-02-19 03:03:44Z raasch
105! Output of cloud physics parameters/quantities complemented and restructured
106!
107! 767 2011-10-14 06:39:12Z raasch
108! Output of given initial u,v-profiles
109!
110! 759 2011-09-15 13:58:31Z raasch
111! output of maximum number of parallel io streams
112!
113! 707 2011-03-29 11:39:40Z raasch
114! bc_lr/ns replaced by bc_lr/ns_cyc
115!
116! 667 2010-12-23 12:06:00Z suehring/gryschka
117! Output of advection scheme.
118! Modified output of Prandtl-layer height.
119!
120! 580 2010-10-05 13:59:11Z heinze
121! Renaming of ws_vertical_gradient to subs_vertical_gradient,
122! ws_vertical_gradient_level to subs_vertical_gradient_level and
123! ws_vertical_gradient_level_ind to subs_vertical_gradient_level_i
124!
125! 493 2010-03-01 08:30:24Z raasch
126! NetCDF data output format extendend for NetCDF4/HDF5
127!
128! 449 2010-02-02 11:23:59Z raasch
129! +large scale vertical motion (subsidence/ascent)
130! Bugfix: index problem concerning gradient_level indices removed
131!
132! 410 2009-12-04 17:05:40Z letzel
133! Masked data output: + dt_domask, mask_01~20_x|y|z, mask_01~20_x|y|z_loop,
134! mask_scale|_x|y|z, masks, skip_time_domask
135!
136! 346 2009-07-06 10:13:41Z raasch
137! initializing_actions='read_data_for_recycling' renamed to 'cyclic_fill'
138! Coupling with independent precursor runs.
139! Output of messages replaced by message handling routine.
140! Output of several additional dvr parameters
141! +canyon_height, canyon_width_x, canyon_width_y, canyon_wall_left,
142! canyon_wall_south, conserve_volume_flow_mode, dp_external, dp_level_b,
143! dp_smooth, dpdxy, u_bulk, v_bulk
144! topography_grid_convention moved from user_header
145! small bugfix concerning 3d 64bit netcdf output format
146!
147! 206 2008-10-13 14:59:11Z raasch
148! Bugfix: error in zu index in case of section_xy = -1
149!
150! 198 2008-09-17 08:55:28Z raasch
151! Format adjustments allowing output of larger revision numbers
152!
153! 197 2008-09-16 15:29:03Z raasch
154! allow 100 spectra levels instead of 10 for consistency with
155! define_netcdf_header,
156! bugfix in the output of the characteristic levels of potential temperature,
157! geostrophic wind, scalar concentration, humidity and leaf area density,
158! output of turbulence recycling informations
159!
160! 138 2007-11-28 10:03:58Z letzel
161! Allow new case bc_uv_t = 'dirichlet_0' for channel flow.
162! Allow two instead of one digit to specify isosurface and slicer variables.
163! Output of sorting frequency of particles
164!
165! 108 2007-08-24 15:10:38Z letzel
166! Output of informations for coupled model runs (boundary conditions etc.)
167! + output of momentumfluxes at the top boundary
168! Rayleigh damping for ocean, e_init
169!
170! 97 2007-06-21 08:23:15Z raasch
171! Adjustments for the ocean version.
172! use_pt_reference renamed use_reference
173!
174! 87 2007-05-22 15:46:47Z raasch
175! Bugfix: output of use_upstream_for_tke
176!
177! 82 2007-04-16 15:40:52Z raasch
178! Preprocessor strings for different linux clusters changed to "lc",
179! routine local_flush is used for buffer flushing
180!
181! 76 2007-03-29 00:58:32Z raasch
182! Output of netcdf_64bit_3d, particles-package is now part of the default code,
183! output of the loop optimization method, moisture renamed humidity,
184! output of subversion revision number
185!
186! 19 2007-02-23 04:53:48Z raasch
187! Output of scalar flux applied at top boundary
188!
189! RCS Log replace by Id keyword, revision history cleaned up
190!
191! Revision 1.63  2006/08/22 13:53:13  raasch
192! Output of dz_max
193!
194! Revision 1.1  1997/08/11 06:17:20  raasch
195! Initial revision
196!
197!
198! Description:
199! ------------
200! Writing a header with all important informations about the actual run.
201! This subroutine is called three times, two times at the beginning
202! (writing information on files RUN_CONTROL and HEADER) and one time at the
203! end of the run, then writing additional information about CPU-usage on file
204! header.
205!-----------------------------------------------------------------------------!
206
207    USE arrays_3d
208    USE control_parameters
209    USE cloud_parameters
210    USE cpulog
211    USE dvrp_variables
212    USE grid_variables
213    USE indices
214    USE model_1d
215    USE particle_attributes
216    USE pegrid
217    USE subsidence_mod
218    USE spectrum
219
220    IMPLICIT NONE
221
222    CHARACTER (LEN=1)  ::  prec
223    CHARACTER (LEN=2)  ::  do2d_mode
224    CHARACTER (LEN=5)  ::  section_chr
225    CHARACTER (LEN=10) ::  coor_chr, host_chr
226    CHARACTER (LEN=16) ::  begin_chr
227    CHARACTER (LEN=26) ::  ver_rev
228    CHARACTER (LEN=40) ::  output_format
229    CHARACTER (LEN=70) ::  char1, char2, dopr_chr, &
230                           do2d_xy, do2d_xz, do2d_yz, do3d_chr, &
231                           domask_chr, run_classification
232    CHARACTER (LEN=86) ::  coordinates, gradients, learde, slices,  &
233                           temperatures, ugcomponent, vgcomponent
234    CHARACTER (LEN=85) ::  roben, runten
235
236    CHARACTER (LEN=1), DIMENSION(1:3) ::  dir = (/ 'x', 'y', 'z' /)
237
238    INTEGER ::  av, bh, blx, bly, bxl, bxr, byn, bys, ch, count, cwx, cwy,  &
239                cxl, cxr, cyn, cys, dim, i, io, j, l, ll, mpi_type
240    REAL    ::  cpuseconds_per_simulated_second
241
242!
243!-- Open the output file. At the end of the simulation, output is directed
244!-- to unit 19.
245    IF ( ( runnr == 0 .OR. force_print_header )  .AND. &
246         .NOT. simulated_time_at_begin /= simulated_time )  THEN
247       io = 15   !  header output on file RUN_CONTROL
248    ELSE
249       io = 19   !  header output on file HEADER
250    ENDIF
251    CALL check_open( io )
252
253!
254!-- At the end of the run, output file (HEADER) will be rewritten with
255!-- new informations
256    IF ( io == 19 .AND. simulated_time_at_begin /= simulated_time ) REWIND( 19 )
257
258!
259!-- Determine kind of model run
260    IF ( TRIM( initializing_actions ) == 'read_restart_data' )  THEN
261       run_classification = '3D - restart run'
262    ELSEIF ( TRIM( initializing_actions ) == 'cyclic_fill' )  THEN
263       run_classification = '3D - run with cyclic fill of 3D - prerun data'
264    ELSEIF ( INDEX( initializing_actions, 'set_constant_profiles' ) /= 0 )  THEN
265       run_classification = '3D - run without 1D - prerun'
266    ELSEIF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
267       run_classification = '3D - run with 1D - prerun'
268    ELSEIF ( INDEX( initializing_actions, 'by_user' ) /=0 )  THEN
269       run_classification = '3D - run initialized by user'
270    ELSE
271       message_string = ' unknown action(s): ' // TRIM( initializing_actions )
272       CALL message( 'header', 'PA0191', 0, 0, 0, 6, 0 )
273    ENDIF
274    IF ( ocean )  THEN
275       run_classification = 'ocean - ' // run_classification
276    ELSE
277       run_classification = 'atmosphere - ' // run_classification
278    ENDIF
279
280!
281!-- Run-identification, date, time, host
282    host_chr = host(1:10)
283    ver_rev = TRIM( version ) // '  ' // TRIM( revision )
284    WRITE ( io, 100 )  ver_rev, TRIM( run_classification )
285    IF ( TRIM( coupling_mode ) /= 'uncoupled' )  THEN
286#if defined( __mpi2 )
287       mpi_type = 2
288#else
289       mpi_type = 1
290#endif
291       WRITE ( io, 101 )  mpi_type, coupling_mode
292    ENDIF
293#if defined( __parallel )
294    IF ( coupling_start_time /= 0.0 )  THEN
295       IF ( coupling_start_time > simulated_time_at_begin )  THEN
296          WRITE ( io, 109 )
297       ELSE
298          WRITE ( io, 114 )
299       ENDIF
300    ENDIF
301#endif
302    WRITE ( io, 102 )  run_date, run_identifier, run_time, runnr, &
303                       ADJUSTR( host_chr )
304#if defined( __parallel )
305    IF ( npex == -1  .AND.  pdims(2) /= 1 )  THEN
306       char1 = 'calculated'
307    ELSEIF ( ( host(1:3) == 'ibm'  .OR.  host(1:3) == 'nec'  .OR.  &
308               host(1:2) == 'lc' )  .AND.                          &
309             npex == -1  .AND.  pdims(2) == 1 )  THEN
310       char1 = 'forced'
311    ELSE
312       char1 = 'predefined'
313    ENDIF
314    IF ( threads_per_task == 1 )  THEN
315       WRITE ( io, 103 )  numprocs, pdims(1), pdims(2), TRIM( char1 )
316    ELSE
317       WRITE ( io, 104 )  numprocs*threads_per_task, numprocs, &
318                          threads_per_task, pdims(1), pdims(2), TRIM( char1 )
319    ENDIF
320    IF ( num_acc_per_node /= 0 )  WRITE ( io, 117 )  num_acc_per_node   
321    IF ( ( host(1:3) == 'ibm'  .OR.  host(1:3) == 'nec'  .OR.    &
322           host(1:2) == 'lc'   .OR.  host(1:3) == 'dec' )  .AND. &
323         npex == -1  .AND.  pdims(2) == 1 )                      &
324    THEN
325       WRITE ( io, 106 )
326    ELSEIF ( pdims(2) == 1 )  THEN
327       WRITE ( io, 107 )  'x'
328    ELSEIF ( pdims(1) == 1 )  THEN
329       WRITE ( io, 107 )  'y'
330    ENDIF
331    IF ( use_seperate_pe_for_dvrp_output )  WRITE ( io, 105 )
332    IF ( numprocs /= maximum_parallel_io_streams )  THEN
333       WRITE ( io, 108 )  maximum_parallel_io_streams
334    ENDIF
335#else
336    IF ( num_acc_per_node /= 0 )  WRITE ( io, 120 )  num_acc_per_node
337#endif
338    WRITE ( io, 99 )
339
340!
341!-- Numerical schemes
342    WRITE ( io, 110 )
343    IF ( psolver(1:7) == 'poisfft' )  THEN
344       WRITE ( io, 111 )  TRIM( fft_method )
345       IF ( transpose_compute_overlap )  WRITE( io, 115 )
346    ELSEIF ( psolver == 'sor' )  THEN
347       WRITE ( io, 112 )  nsor_ini, nsor, omega_sor
348    ELSEIF ( psolver == 'multigrid' )  THEN
349       WRITE ( io, 135 )  cycle_mg, maximum_grid_level, ngsrb
350       IF ( mg_cycles == -1 )  THEN
351          WRITE ( io, 140 )  residual_limit
352       ELSE
353          WRITE ( io, 141 )  mg_cycles
354       ENDIF
355       IF ( mg_switch_to_pe0_level == 0 )  THEN
356          WRITE ( io, 136 )  nxr_mg(1)-nxl_mg(1)+1, nyn_mg(1)-nys_mg(1)+1, &
357                             nzt_mg(1)
358       ELSEIF (  mg_switch_to_pe0_level /= -1 )  THEN
359          WRITE ( io, 137 )  mg_switch_to_pe0_level,            &
360                             mg_loc_ind(2,0)-mg_loc_ind(1,0)+1, &
361                             mg_loc_ind(4,0)-mg_loc_ind(3,0)+1, &
362                             nzt_mg(mg_switch_to_pe0_level),    &
363                             nxr_mg(1)-nxl_mg(1)+1, nyn_mg(1)-nys_mg(1)+1, &
364                             nzt_mg(1)
365       ENDIF
366       IF ( masking_method )  WRITE ( io, 144 )
367    ENDIF
368    IF ( call_psolver_at_all_substeps  .AND. timestep_scheme(1:5) == 'runge' ) &
369    THEN
370       WRITE ( io, 142 )
371    ENDIF
372
373    IF ( momentum_advec == 'pw-scheme' )  THEN
374       WRITE ( io, 113 )
375    ELSEIF (momentum_advec == 'ws-scheme' )  THEN
376       WRITE ( io, 503 )
377    ENDIF
378    IF ( scalar_advec == 'pw-scheme' )  THEN
379       WRITE ( io, 116 )
380    ELSEIF ( scalar_advec == 'ws-scheme' )  THEN
381       WRITE ( io, 504 )
382    ELSE
383       WRITE ( io, 118 )
384    ENDIF
385
386    WRITE ( io, 139 )  TRIM( loop_optimization )
387
388    IF ( galilei_transformation )  THEN
389       IF ( use_ug_for_galilei_tr )  THEN
390          char1 = '0.6 * geostrophic wind'
391       ELSE
392          char1 = 'mean wind in model domain'
393       ENDIF
394       IF ( simulated_time_at_begin == simulated_time )  THEN
395          char2 = 'at the start of the run'
396       ELSE
397          char2 = 'at the end of the run'
398       ENDIF
399       WRITE ( io, 119 )  TRIM( char1 ), TRIM( char2 ), &
400                          advected_distance_x/1000.0, advected_distance_y/1000.0
401    ENDIF
402    WRITE ( io, 122 )  timestep_scheme
403    IF ( use_upstream_for_tke )  WRITE ( io, 143 )
404    IF ( rayleigh_damping_factor /= 0.0 )  THEN
405       IF ( .NOT. ocean )  THEN
406          WRITE ( io, 123 )  'above', rayleigh_damping_height, &
407               rayleigh_damping_factor
408       ELSE
409          WRITE ( io, 123 )  'below', rayleigh_damping_height, &
410               rayleigh_damping_factor
411       ENDIF
412    ENDIF
413    IF ( neutral )  WRITE ( io, 131 )  pt_surface
414    IF ( humidity )  THEN
415       IF ( .NOT. cloud_physics )  THEN
416          WRITE ( io, 129 )
417       ELSE
418          WRITE ( io, 130 )
419       ENDIF
420    ENDIF
421    IF ( passive_scalar )  WRITE ( io, 134 )
422    IF ( conserve_volume_flow )  THEN
423       WRITE ( io, 150 )  conserve_volume_flow_mode
424       IF ( TRIM( conserve_volume_flow_mode ) == 'bulk_velocity' )  THEN
425          WRITE ( io, 151 )  u_bulk, v_bulk
426       ENDIF
427    ELSEIF ( dp_external )  THEN
428       IF ( dp_smooth )  THEN
429          WRITE ( io, 152 )  dpdxy, dp_level_b, ', vertically smoothed.'
430       ELSE
431          WRITE ( io, 152 )  dpdxy, dp_level_b, '.'
432       ENDIF
433    ENDIF
434    IF ( large_scale_subsidence )  THEN
435        WRITE ( io, 153 )
436        WRITE ( io, 154 )
437    ENDIF
438    IF ( nudging )  THEN
439        WRITE ( io, 155 )
440        WRITE ( io, 156 )
441    ENDIF
442    IF ( large_scale_forcing )  THEN
443        WRITE ( io, 157 )
444        WRITE ( io, 158 )
445    ENDIF
446    WRITE ( io, 99 )
447
448!
449!-- Runtime and timestep informations
450    WRITE ( io, 200 )
451    IF ( .NOT. dt_fixed )  THEN
452       WRITE ( io, 201 )  dt_max, cfl_factor
453    ELSE
454       WRITE ( io, 202 )  dt
455    ENDIF
456    WRITE ( io, 203 )  simulated_time_at_begin, end_time
457
458    IF ( time_restart /= 9999999.9  .AND. &
459         simulated_time_at_begin == simulated_time )  THEN
460       IF ( dt_restart == 9999999.9 )  THEN
461          WRITE ( io, 204 )  ' Restart at:       ',time_restart
462       ELSE
463          WRITE ( io, 205 )  ' Restart at:       ',time_restart, dt_restart
464       ENDIF
465    ENDIF
466
467    IF ( simulated_time_at_begin /= simulated_time )  THEN
468       i = MAX ( log_point_s(10)%counts, 1 )
469       IF ( ( simulated_time - simulated_time_at_begin ) == 0.0 )  THEN
470          cpuseconds_per_simulated_second = 0.0
471       ELSE
472          cpuseconds_per_simulated_second = log_point_s(10)%sum / &
473                                            ( simulated_time -    &
474                                              simulated_time_at_begin )
475       ENDIF
476       WRITE ( io, 206 )  simulated_time, log_point_s(10)%sum, &
477                          log_point_s(10)%sum / REAL( i ),     &
478                          cpuseconds_per_simulated_second
479       IF ( time_restart /= 9999999.9  .AND.  time_restart < end_time )  THEN
480          IF ( dt_restart == 9999999.9 )  THEN
481             WRITE ( io, 204 )  ' Next restart at:     ',time_restart
482          ELSE
483             WRITE ( io, 205 )  ' Next restart at:     ',time_restart, dt_restart
484          ENDIF
485       ENDIF
486    ENDIF
487
488!
489!-- Start time for coupled runs, if independent precursor runs for atmosphere
490!-- and ocean are used or have been used. In this case, coupling_start_time
491!-- defines the time when the coupling is switched on.
492    IF ( coupling_start_time /= 0.0 )  THEN
493       WRITE ( io, 207 )  coupling_start_time
494    ENDIF
495
496!
497!-- Computational grid
498    IF ( .NOT. ocean )  THEN
499       WRITE ( io, 250 )  dx, dy, dz, (nx+1)*dx, (ny+1)*dy, zu(nzt+1)
500       IF ( dz_stretch_level_index < nzt+1 )  THEN
501          WRITE ( io, 252 )  dz_stretch_level, dz_stretch_level_index, &
502                             dz_stretch_factor, dz_max
503       ENDIF
504    ELSE
505       WRITE ( io, 250 )  dx, dy, dz, (nx+1)*dx, (ny+1)*dy, zu(0)
506       IF ( dz_stretch_level_index > 0 )  THEN
507          WRITE ( io, 252 )  dz_stretch_level, dz_stretch_level_index, &
508                             dz_stretch_factor, dz_max
509       ENDIF
510    ENDIF
511    WRITE ( io, 254 )  nx, ny, nzt+1, MIN( nnx, nx+1 ), MIN( nny, ny+1 ), &
512                       MIN( nnz+2, nzt+2 )
513    IF ( sloping_surface )  WRITE ( io, 260 )  alpha_surface
514
515!
516!-- Topography
517    WRITE ( io, 270 )  topography
518    SELECT CASE ( TRIM( topography ) )
519
520       CASE ( 'flat' )
521          ! no actions necessary
522
523       CASE ( 'single_building' )
524          blx = INT( building_length_x / dx )
525          bly = INT( building_length_y / dy )
526          bh  = INT( building_height / dz )
527
528          IF ( building_wall_left == 9999999.9 )  THEN
529             building_wall_left = ( nx + 1 - blx ) / 2 * dx
530          ENDIF
531          bxl = INT ( building_wall_left / dx + 0.5 )
532          bxr = bxl + blx
533
534          IF ( building_wall_south == 9999999.9 )  THEN
535             building_wall_south = ( ny + 1 - bly ) / 2 * dy
536          ENDIF
537          bys = INT ( building_wall_south / dy + 0.5 )
538          byn = bys + bly
539
540          WRITE ( io, 271 )  building_length_x, building_length_y, &
541                             building_height, bxl, bxr, bys, byn
542
543       CASE ( 'single_street_canyon' )
544          ch  = NINT( canyon_height / dz )
545          IF ( canyon_width_x /= 9999999.9 )  THEN
546!
547!--          Street canyon in y direction
548             cwx = NINT( canyon_width_x / dx )
549             IF ( canyon_wall_left == 9999999.9 )  THEN
550                canyon_wall_left = ( nx + 1 - cwx ) / 2 * dx
551             ENDIF
552             cxl = NINT( canyon_wall_left / dx )
553             cxr = cxl + cwx
554             WRITE ( io, 272 )  'y', canyon_height, ch, 'u', cxl, cxr
555
556          ELSEIF ( canyon_width_y /= 9999999.9 )  THEN
557!
558!--          Street canyon in x direction
559             cwy = NINT( canyon_width_y / dy )
560             IF ( canyon_wall_south == 9999999.9 )  THEN
561                canyon_wall_south = ( ny + 1 - cwy ) / 2 * dy
562             ENDIF
563             cys = NINT( canyon_wall_south / dy )
564             cyn = cys + cwy
565             WRITE ( io, 272 )  'x', canyon_height, ch, 'v', cys, cyn
566          ENDIF
567
568    END SELECT
569
570    IF ( TRIM( topography ) /= 'flat' )  THEN
571       IF ( TRIM( topography_grid_convention ) == ' ' )  THEN
572          IF ( TRIM( topography ) == 'single_building' .OR.  &
573               TRIM( topography ) == 'single_street_canyon' )  THEN
574             WRITE ( io, 278 )
575          ELSEIF ( TRIM( topography ) == 'read_from_file' )  THEN
576             WRITE ( io, 279 )
577          ENDIF
578       ELSEIF ( TRIM( topography_grid_convention ) == 'cell_edge' )  THEN
579          WRITE ( io, 278 )
580       ELSEIF ( TRIM( topography_grid_convention ) == 'cell_center' )  THEN
581          WRITE ( io, 279 )
582       ENDIF
583    ENDIF
584
585    IF ( plant_canopy )  THEN
586
587       WRITE ( io, 280 ) canopy_mode, pch_index, drag_coefficient
588       IF ( passive_scalar )  THEN
589          WRITE ( io, 281 ) scalar_exchange_coefficient,   &
590                            leaf_surface_concentration
591       ENDIF
592
593!
594!--    Heat flux at the top of vegetation
595       WRITE ( io, 282 ) cthf
596
597!
598!--    Leaf area density profile
599!--    Building output strings, starting with surface value
600       WRITE ( learde, '(F6.4)' )  lad_surface
601       gradients = '------'
602       slices = '     0'
603       coordinates = '   0.0'
604       i = 1
605       DO  WHILE ( lad_vertical_gradient_level_ind(i) /= -9999 )
606
607          WRITE (coor_chr,'(F7.2)')  lad(lad_vertical_gradient_level_ind(i))
608          learde = TRIM( learde ) // ' ' // TRIM( coor_chr )
609
610          WRITE (coor_chr,'(F7.2)')  lad_vertical_gradient(i)
611          gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
612
613          WRITE (coor_chr,'(I7)')  lad_vertical_gradient_level_ind(i)
614          slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
615
616          WRITE (coor_chr,'(F7.1)')  lad_vertical_gradient_level(i)
617          coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
618
619          i = i + 1
620       ENDDO
621
622       WRITE ( io, 283 )  TRIM( coordinates ), TRIM( learde ), &
623                          TRIM( gradients ), TRIM( slices )
624
625    ENDIF
626
627!
628!-- Boundary conditions
629    IF ( ibc_p_b == 0 )  THEN
630       runten = 'p(0)     = 0      |'
631    ELSEIF ( ibc_p_b == 1 )  THEN
632       runten = 'p(0)     = p(1)   |'
633    ENDIF
634    IF ( ibc_p_t == 0 )  THEN
635       roben  = 'p(nzt+1) = 0      |'
636    ELSE
637       roben  = 'p(nzt+1) = p(nzt) |'
638    ENDIF
639
640    IF ( ibc_uv_b == 0 )  THEN
641       runten = TRIM( runten ) // ' uv(0)     = -uv(1)                |'
642    ELSE
643       runten = TRIM( runten ) // ' uv(0)     = uv(1)                 |'
644    ENDIF
645    IF ( TRIM( bc_uv_t ) == 'dirichlet_0' )  THEN
646       roben  = TRIM( roben  ) // ' uv(nzt+1) = 0                     |'
647    ELSEIF ( ibc_uv_t == 0 )  THEN
648       roben  = TRIM( roben  ) // ' uv(nzt+1) = ug(nzt+1), vg(nzt+1)  |'
649    ELSE
650       roben  = TRIM( roben  ) // ' uv(nzt+1) = uv(nzt)               |'
651    ENDIF
652
653    IF ( ibc_pt_b == 0 )  THEN
654       runten = TRIM( runten ) // ' pt(0)   = pt_surface'
655    ELSEIF ( ibc_pt_b == 1 )  THEN
656       runten = TRIM( runten ) // ' pt(0)   = pt(1)'
657    ELSEIF ( ibc_pt_b == 2 )  THEN
658       runten = TRIM( runten ) // ' pt(0) = from coupled model'
659    ENDIF
660    IF ( ibc_pt_t == 0 )  THEN
661       roben  = TRIM( roben  ) // ' pt(nzt+1) = pt_top'
662    ELSEIF( ibc_pt_t == 1 )  THEN
663       roben  = TRIM( roben  ) // ' pt(nzt+1) = pt(nzt)'
664    ELSEIF( ibc_pt_t == 2 )  THEN
665       roben  = TRIM( roben  ) // ' pt(nzt+1) = pt(nzt) + dpt/dz_ini'
666
667    ENDIF
668
669    WRITE ( io, 300 )  runten, roben
670
671    IF ( .NOT. constant_diffusion )  THEN
672       IF ( ibc_e_b == 1 )  THEN
673          runten = 'e(0)     = e(1)'
674       ELSE
675          runten = 'e(0)     = e(1) = (u*/0.1)**2'
676       ENDIF
677       roben = 'e(nzt+1) = e(nzt) = e(nzt-1)'
678
679       WRITE ( io, 301 )  'e', runten, roben       
680
681    ENDIF
682
683    IF ( ocean )  THEN
684       runten = 'sa(0)    = sa(1)'
685       IF ( ibc_sa_t == 0 )  THEN
686          roben =  'sa(nzt+1) = sa_surface'
687       ELSE
688          roben =  'sa(nzt+1) = sa(nzt)'
689       ENDIF
690       WRITE ( io, 301 ) 'sa', runten, roben
691    ENDIF
692
693    IF ( humidity )  THEN
694       IF ( ibc_q_b == 0 )  THEN
695          runten = 'q(0)     = q_surface'
696       ELSE
697          runten = 'q(0)     = q(1)'
698       ENDIF
699       IF ( ibc_q_t == 0 )  THEN
700          roben =  'q(nzt)   = q_top'
701       ELSE
702          roben =  'q(nzt)   = q(nzt-1) + dq/dz'
703       ENDIF
704       WRITE ( io, 301 ) 'q', runten, roben
705    ENDIF
706
707    IF ( passive_scalar )  THEN
708       IF ( ibc_q_b == 0 )  THEN
709          runten = 's(0)     = s_surface'
710       ELSE
711          runten = 's(0)     = s(1)'
712       ENDIF
713       IF ( ibc_q_t == 0 )  THEN
714          roben =  's(nzt)   = s_top'
715       ELSE
716          roben =  's(nzt)   = s(nzt-1) + ds/dz'
717       ENDIF
718       WRITE ( io, 301 ) 's', runten, roben
719    ENDIF
720
721    IF ( use_surface_fluxes )  THEN
722       WRITE ( io, 303 )
723       IF ( constant_heatflux )  THEN
724          IF ( large_scale_forcing .AND. lsf_surf )  THEN
725             WRITE ( io, 306 )  shf(0,0)
726          ELSE
727             WRITE ( io, 306 )  surface_heatflux
728          ENDIF
729          IF ( random_heatflux )  WRITE ( io, 307 )
730       ENDIF
731       IF ( humidity  .AND.  constant_waterflux )  THEN
732          IF ( large_scale_forcing .AND. lsf_surf )  THEN
733             WRITE ( io, 311 ) qsws(0,0)
734          ELSE
735             WRITE ( io, 311 ) surface_waterflux
736          ENDIF
737       ENDIF
738       IF ( passive_scalar  .AND.  constant_waterflux )  THEN
739          WRITE ( io, 313 ) surface_waterflux
740       ENDIF
741    ENDIF
742
743    IF ( use_top_fluxes )  THEN
744       WRITE ( io, 304 )
745       IF ( coupling_mode == 'uncoupled' )  THEN
746          WRITE ( io, 320 )  top_momentumflux_u, top_momentumflux_v
747          IF ( constant_top_heatflux )  THEN
748             WRITE ( io, 306 )  top_heatflux
749          ENDIF
750       ELSEIF ( coupling_mode == 'ocean_to_atmosphere' )  THEN
751          WRITE ( io, 316 )
752       ENDIF
753       IF ( ocean  .AND.  constant_top_salinityflux )  THEN
754          WRITE ( io, 309 )  top_salinityflux
755       ENDIF
756       IF ( humidity  .OR.  passive_scalar )  THEN
757          WRITE ( io, 315 )
758       ENDIF
759    ENDIF
760
761    IF ( prandtl_layer )  THEN
762       WRITE ( io, 305 )  (zu(1)-zu(0)), roughness_length, &
763                          z0h_factor*roughness_length, kappa, &
764                          rif_min, rif_max
765       IF ( .NOT. constant_heatflux )  WRITE ( io, 308 )
766       IF ( humidity  .AND.  .NOT. constant_waterflux )  THEN
767          WRITE ( io, 312 )
768       ENDIF
769       IF ( passive_scalar  .AND.  .NOT. constant_waterflux )  THEN
770          WRITE ( io, 314 )
771       ENDIF
772    ELSE
773       IF ( INDEX(initializing_actions, 'set_1d-model_profiles') /= 0 )  THEN
774          WRITE ( io, 310 )  rif_min, rif_max
775       ENDIF
776    ENDIF
777
778    WRITE ( io, 317 )  bc_lr, bc_ns
779    IF ( .NOT. bc_lr_cyc  .OR.  .NOT. bc_ns_cyc )  THEN
780       WRITE ( io, 318 )  use_cmax, pt_damping_width, pt_damping_factor       
781       IF ( turbulent_inflow )  THEN
782          WRITE ( io, 319 )  recycling_width, recycling_plane, &
783                             inflow_damping_height, inflow_damping_width
784       ENDIF
785    ENDIF
786
787!
788!-- Listing of 1D-profiles
789    WRITE ( io, 325 )  dt_dopr_listing
790    IF ( averaging_interval_pr /= 0.0 )  THEN
791       WRITE ( io, 326 )  averaging_interval_pr, dt_averaging_input_pr
792    ENDIF
793
794!
795!-- DATA output
796    WRITE ( io, 330 )
797    IF ( averaging_interval_pr /= 0.0 )  THEN
798       WRITE ( io, 326 )  averaging_interval_pr, dt_averaging_input_pr
799    ENDIF
800
801!
802!-- 1D-profiles
803    dopr_chr = 'Profile:'
804    IF ( dopr_n /= 0 )  THEN
805       WRITE ( io, 331 )
806
807       output_format = ''
808       IF ( netcdf_output )  THEN
809          output_format = output_format_netcdf
810       ENDIF
811       WRITE ( io, 344 )  output_format
812
813       DO  i = 1, dopr_n
814          dopr_chr = TRIM( dopr_chr ) // ' ' // TRIM( data_output_pr(i) ) // ','
815          IF ( LEN_TRIM( dopr_chr ) >= 60 )  THEN
816             WRITE ( io, 332 )  dopr_chr
817             dopr_chr = '       :'
818          ENDIF
819       ENDDO
820
821       IF ( dopr_chr /= '' )  THEN
822          WRITE ( io, 332 )  dopr_chr
823       ENDIF
824       WRITE ( io, 333 )  dt_dopr, averaging_interval_pr, dt_averaging_input_pr
825       IF ( skip_time_dopr /= 0.0 )  WRITE ( io, 339 )  skip_time_dopr
826    ENDIF
827
828!
829!-- 2D-arrays
830    DO  av = 0, 1
831
832       i = 1
833       do2d_xy = ''
834       do2d_xz = ''
835       do2d_yz = ''
836       DO  WHILE ( do2d(av,i) /= ' ' )
837
838          l = MAX( 2, LEN_TRIM( do2d(av,i) ) )
839          do2d_mode = do2d(av,i)(l-1:l)
840
841          SELECT CASE ( do2d_mode )
842             CASE ( 'xy' )
843                ll = LEN_TRIM( do2d_xy )
844                do2d_xy = do2d_xy(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
845             CASE ( 'xz' )
846                ll = LEN_TRIM( do2d_xz )
847                do2d_xz = do2d_xz(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
848             CASE ( 'yz' )
849                ll = LEN_TRIM( do2d_yz )
850                do2d_yz = do2d_yz(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
851          END SELECT
852
853          i = i + 1
854
855       ENDDO
856
857       IF ( ( ( do2d_xy /= ''  .AND.  section(1,1) /= -9999 )  .OR.    &
858              ( do2d_xz /= ''  .AND.  section(1,2) /= -9999 )  .OR.    &
859              ( do2d_yz /= ''  .AND.  section(1,3) /= -9999 ) )  .AND. &
860            ( netcdf_output  .OR.  iso2d_output ) )  THEN
861
862          IF (  av == 0 )  THEN
863             WRITE ( io, 334 )  ''
864          ELSE
865             WRITE ( io, 334 )  '(time-averaged)'
866          ENDIF
867
868          IF ( do2d_at_begin )  THEN
869             begin_chr = 'and at the start'
870          ELSE
871             begin_chr = ''
872          ENDIF
873
874          output_format = ''
875          IF ( netcdf_output )  THEN
876             output_format = output_format_netcdf
877          ENDIF
878          IF ( iso2d_output )  THEN
879             IF ( netcdf_output )  THEN
880                output_format = TRIM( output_format_netcdf ) // ' and iso2d'
881             ELSE
882                output_format = 'iso2d'
883             ENDIF
884          ENDIF
885          WRITE ( io, 344 )  output_format
886
887          IF ( do2d_xy /= ''  .AND.  section(1,1) /= -9999 )  THEN
888             i = 1
889             slices = '/'
890             coordinates = '/'
891!
892!--          Building strings with index and coordinate informations of the
893!--          slices
894             DO  WHILE ( section(i,1) /= -9999 )
895
896                WRITE (section_chr,'(I5)')  section(i,1)
897                section_chr = ADJUSTL( section_chr )
898                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
899
900                IF ( section(i,1) == -1 )  THEN
901                   WRITE (coor_chr,'(F10.1)')  -1.0
902                ELSE
903                   WRITE (coor_chr,'(F10.1)')  zu(section(i,1))
904                ENDIF
905                coor_chr = ADJUSTL( coor_chr )
906                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
907
908                i = i + 1
909             ENDDO
910             IF ( av == 0 )  THEN
911                WRITE ( io, 335 )  'XY', do2d_xy, dt_do2d_xy, &
912                                   TRIM( begin_chr ), 'k', TRIM( slices ), &
913                                   TRIM( coordinates )
914                IF ( skip_time_do2d_xy /= 0.0 )  THEN
915                   WRITE ( io, 339 )  skip_time_do2d_xy
916                ENDIF
917             ELSE
918                WRITE ( io, 342 )  'XY', do2d_xy, dt_data_output_av, &
919                                   TRIM( begin_chr ), averaging_interval, &
920                                   dt_averaging_input, 'k', TRIM( slices ), &
921                                   TRIM( coordinates )
922                IF ( skip_time_data_output_av /= 0.0 )  THEN
923                   WRITE ( io, 339 )  skip_time_data_output_av
924                ENDIF
925             ENDIF
926
927          ENDIF
928
929          IF ( do2d_xz /= ''  .AND.  section(1,2) /= -9999 )  THEN
930             i = 1
931             slices = '/'
932             coordinates = '/'
933!
934!--          Building strings with index and coordinate informations of the
935!--          slices
936             DO  WHILE ( section(i,2) /= -9999 )
937
938                WRITE (section_chr,'(I5)')  section(i,2)
939                section_chr = ADJUSTL( section_chr )
940                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
941
942                WRITE (coor_chr,'(F10.1)')  section(i,2) * dy
943                coor_chr = ADJUSTL( coor_chr )
944                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
945
946                i = i + 1
947             ENDDO
948             IF ( av == 0 )  THEN
949                WRITE ( io, 335 )  'XZ', do2d_xz, dt_do2d_xz, &
950                                   TRIM( begin_chr ), 'j', TRIM( slices ), &
951                                   TRIM( coordinates )
952                IF ( skip_time_do2d_xz /= 0.0 )  THEN
953                   WRITE ( io, 339 )  skip_time_do2d_xz
954                ENDIF
955             ELSE
956                WRITE ( io, 342 )  'XZ', do2d_xz, dt_data_output_av, &
957                                   TRIM( begin_chr ), averaging_interval, &
958                                   dt_averaging_input, 'j', TRIM( slices ), &
959                                   TRIM( coordinates )
960                IF ( skip_time_data_output_av /= 0.0 )  THEN
961                   WRITE ( io, 339 )  skip_time_data_output_av
962                ENDIF
963             ENDIF
964          ENDIF
965
966          IF ( do2d_yz /= ''  .AND.  section(1,3) /= -9999 )  THEN
967             i = 1
968             slices = '/'
969             coordinates = '/'
970!
971!--          Building strings with index and coordinate informations of the
972!--          slices
973             DO  WHILE ( section(i,3) /= -9999 )
974
975                WRITE (section_chr,'(I5)')  section(i,3)
976                section_chr = ADJUSTL( section_chr )
977                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
978
979                WRITE (coor_chr,'(F10.1)')  section(i,3) * dx
980                coor_chr = ADJUSTL( coor_chr )
981                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
982
983                i = i + 1
984             ENDDO
985             IF ( av == 0 )  THEN
986                WRITE ( io, 335 )  'YZ', do2d_yz, dt_do2d_yz, &
987                                   TRIM( begin_chr ), 'i', TRIM( slices ), &
988                                   TRIM( coordinates )
989                IF ( skip_time_do2d_yz /= 0.0 )  THEN
990                   WRITE ( io, 339 )  skip_time_do2d_yz
991                ENDIF
992             ELSE
993                WRITE ( io, 342 )  'YZ', do2d_yz, dt_data_output_av, &
994                                   TRIM( begin_chr ), averaging_interval, &
995                                   dt_averaging_input, 'i', TRIM( slices ), &
996                                   TRIM( coordinates )
997                IF ( skip_time_data_output_av /= 0.0 )  THEN
998                   WRITE ( io, 339 )  skip_time_data_output_av
999                ENDIF
1000             ENDIF
1001          ENDIF
1002
1003       ENDIF
1004
1005    ENDDO
1006
1007!
1008!-- 3d-arrays
1009    DO  av = 0, 1
1010
1011       i = 1
1012       do3d_chr = ''
1013       DO  WHILE ( do3d(av,i) /= ' ' )
1014
1015          do3d_chr = TRIM( do3d_chr ) // ' ' // TRIM( do3d(av,i) ) // ','
1016          i = i + 1
1017
1018       ENDDO
1019
1020       IF ( do3d_chr /= '' )  THEN
1021          IF ( av == 0 )  THEN
1022             WRITE ( io, 336 )  ''
1023          ELSE
1024             WRITE ( io, 336 )  '(time-averaged)'
1025          ENDIF
1026
1027          output_format = ''
1028          IF ( netcdf_output )  THEN
1029             output_format = output_format_netcdf
1030          ENDIF
1031          IF ( avs_output )  THEN
1032             IF ( netcdf_output )  THEN
1033                output_format = TRIM( output_format_netcdf ) // ' and avs'
1034             ELSE
1035                output_format = 'avs'
1036             ENDIF
1037          ENDIF
1038          WRITE ( io, 344 )  output_format
1039
1040          IF ( do3d_at_begin )  THEN
1041             begin_chr = 'and at the start'
1042          ELSE
1043             begin_chr = ''
1044          ENDIF
1045          IF ( av == 0 )  THEN
1046             WRITE ( io, 337 )  do3d_chr, dt_do3d, TRIM( begin_chr ), &
1047                                zu(nz_do3d), nz_do3d
1048          ELSE
1049             WRITE ( io, 343 )  do3d_chr, dt_data_output_av,           &
1050                                TRIM( begin_chr ), averaging_interval, &
1051                                dt_averaging_input, zu(nz_do3d), nz_do3d
1052          ENDIF
1053
1054          IF ( do3d_compress )  THEN
1055             do3d_chr = ''
1056             i = 1
1057             DO WHILE ( do3d(av,i) /= ' ' )
1058
1059                SELECT CASE ( do3d(av,i) )
1060                   CASE ( 'u' )
1061                      j = 1
1062                   CASE ( 'v' )
1063                      j = 2
1064                   CASE ( 'w' )
1065                      j = 3
1066                   CASE ( 'p' )
1067                      j = 4
1068                   CASE ( 'pt' )
1069                      j = 5
1070                END SELECT
1071                WRITE ( prec, '(I1)' )  plot_3d_precision(j)%precision
1072                do3d_chr = TRIM( do3d_chr ) // ' ' // TRIM( do3d(av,i) ) // &
1073                           ':' // prec // ','
1074                i = i + 1
1075
1076             ENDDO
1077             WRITE ( io, 338 )  do3d_chr
1078
1079          ENDIF
1080
1081          IF ( av == 0 )  THEN
1082             IF ( skip_time_do3d /= 0.0 )  THEN
1083                WRITE ( io, 339 )  skip_time_do3d
1084             ENDIF
1085          ELSE
1086             IF ( skip_time_data_output_av /= 0.0 )  THEN
1087                WRITE ( io, 339 )  skip_time_data_output_av
1088             ENDIF
1089          ENDIF
1090
1091       ENDIF
1092
1093    ENDDO
1094
1095!
1096!-- masked arrays
1097    IF ( masks > 0 )  WRITE ( io, 345 )  &
1098         mask_scale_x, mask_scale_y, mask_scale_z
1099    DO  mid = 1, masks
1100       DO  av = 0, 1
1101
1102          i = 1
1103          domask_chr = ''
1104          DO  WHILE ( domask(mid,av,i) /= ' ' )
1105             domask_chr = TRIM( domask_chr ) // ' ' //  &
1106                          TRIM( domask(mid,av,i) ) // ','
1107             i = i + 1
1108          ENDDO
1109
1110          IF ( domask_chr /= '' )  THEN
1111             IF ( av == 0 )  THEN
1112                WRITE ( io, 346 )  '', mid
1113             ELSE
1114                WRITE ( io, 346 )  ' (time-averaged)', mid
1115             ENDIF
1116
1117             output_format = ' '
1118             IF ( netcdf_output )  THEN
1119                output_format = output_format_netcdf
1120             ENDIF
1121             WRITE ( io, 344 )  output_format
1122
1123             IF ( av == 0 )  THEN
1124                WRITE ( io, 347 )  domask_chr, dt_domask(mid)
1125             ELSE
1126                WRITE ( io, 348 )  domask_chr, dt_data_output_av, &
1127                                   averaging_interval, dt_averaging_input
1128             ENDIF
1129
1130             IF ( av == 0 )  THEN
1131                IF ( skip_time_domask(mid) /= 0.0 )  THEN
1132                   WRITE ( io, 339 )  skip_time_domask(mid)
1133                ENDIF
1134             ELSE
1135                IF ( skip_time_data_output_av /= 0.0 )  THEN
1136                   WRITE ( io, 339 )  skip_time_data_output_av
1137                ENDIF
1138             ENDIF
1139!
1140!--          output locations
1141             DO  dim = 1, 3
1142                IF ( mask(mid,dim,1) >= 0.0 )  THEN
1143                   count = 0
1144                   DO  WHILE ( mask(mid,dim,count+1) >= 0.0 )
1145                      count = count + 1
1146                   ENDDO
1147                   WRITE ( io, 349 )  dir(dim), dir(dim), mid, dir(dim), &
1148                                      mask(mid,dim,:count)
1149                ELSEIF ( mask_loop(mid,dim,1) < 0.0 .AND.  &
1150                         mask_loop(mid,dim,2) < 0.0 .AND.  &
1151                         mask_loop(mid,dim,3) == 0.0 )  THEN
1152                   WRITE ( io, 350 )  dir(dim), dir(dim)
1153                ELSEIF ( mask_loop(mid,dim,3) == 0.0 )  THEN
1154                   WRITE ( io, 351 )  dir(dim), dir(dim), mid, dir(dim), &
1155                                      mask_loop(mid,dim,1:2)
1156                ELSE
1157                   WRITE ( io, 351 )  dir(dim), dir(dim), mid, dir(dim), &
1158                                      mask_loop(mid,dim,1:3)
1159                ENDIF
1160             ENDDO
1161          ENDIF
1162
1163       ENDDO
1164    ENDDO
1165
1166!
1167!-- Timeseries
1168    IF ( dt_dots /= 9999999.9 )  THEN
1169       WRITE ( io, 340 )
1170
1171       output_format = ''
1172       IF ( netcdf_output )  THEN
1173          output_format = output_format_netcdf
1174       ENDIF
1175       WRITE ( io, 344 )  output_format
1176       WRITE ( io, 341 )  dt_dots
1177    ENDIF
1178
1179#if defined( __dvrp_graphics )
1180!
1181!-- Dvrp-output
1182    IF ( dt_dvrp /= 9999999.9 )  THEN
1183       WRITE ( io, 360 )  dt_dvrp, TRIM( dvrp_output ), TRIM( dvrp_host ), &
1184                          TRIM( dvrp_username ), TRIM( dvrp_directory )
1185       i = 1
1186       l = 0
1187       m = 0
1188       DO WHILE ( mode_dvrp(i) /= ' ' )
1189          IF ( mode_dvrp(i)(1:10) == 'isosurface' )  THEN
1190             READ ( mode_dvrp(i), '(10X,I2)' )  j
1191             l = l + 1
1192             IF ( do3d(0,j) /= ' ' )  THEN
1193                WRITE ( io, 361 )  TRIM( do3d(0,j) ), threshold(l), &
1194                                   isosurface_color(:,l)
1195             ENDIF
1196          ELSEIF ( mode_dvrp(i)(1:6) == 'slicer' )  THEN
1197             READ ( mode_dvrp(i), '(6X,I2)' )  j
1198             m = m + 1
1199             IF ( do2d(0,j) /= ' ' )  THEN
1200                WRITE ( io, 362 )  TRIM( do2d(0,j) ), &
1201                                   slicer_range_limits_dvrp(:,m)
1202             ENDIF
1203          ELSEIF ( mode_dvrp(i)(1:9) == 'particles' )  THEN
1204             WRITE ( io, 363 )  dvrp_psize
1205             IF ( particle_dvrpsize /= 'none' )  THEN
1206                WRITE ( io, 364 )  'size', TRIM( particle_dvrpsize ), &
1207                                   dvrpsize_interval
1208             ENDIF
1209             IF ( particle_color /= 'none' )  THEN
1210                WRITE ( io, 364 )  'color', TRIM( particle_color ), &
1211                                   color_interval
1212             ENDIF
1213          ENDIF
1214          i = i + 1
1215       ENDDO
1216
1217       WRITE ( io, 365 )  groundplate_color, superelevation_x, &
1218                          superelevation_y, superelevation, clip_dvrp_l, &
1219                          clip_dvrp_r, clip_dvrp_s, clip_dvrp_n
1220
1221       IF ( TRIM( topography ) /= 'flat' )  THEN
1222          WRITE ( io, 366 )  topography_color
1223          IF ( cluster_size > 1 )  THEN
1224             WRITE ( io, 367 )  cluster_size
1225          ENDIF
1226       ENDIF
1227
1228    ENDIF
1229#endif
1230
1231#if defined( __spectra )
1232!
1233!-- Spectra output
1234    IF ( dt_dosp /= 9999999.9 )  THEN
1235       WRITE ( io, 370 )
1236
1237       output_format = ' '
1238       IF ( netcdf_output )  THEN
1239          output_format = output_format_netcdf
1240       ENDIF
1241       WRITE ( io, 344 )  output_format
1242       WRITE ( io, 371 )  dt_dosp
1243       IF ( skip_time_dosp /= 0.0 )  WRITE ( io, 339 )  skip_time_dosp
1244       WRITE ( io, 372 )  ( data_output_sp(i), i = 1,10 ),     &
1245                          ( spectra_direction(i), i = 1,10 ),  &
1246                          ( comp_spectra_level(i), i = 1,100 ), &
1247                          ( plot_spectra_level(i), i = 1,100 ), &
1248                          averaging_interval_sp, dt_averaging_input_pr
1249    ENDIF
1250#endif
1251
1252    WRITE ( io, 99 )
1253
1254!
1255!-- Physical quantities
1256    WRITE ( io, 400 )
1257
1258!
1259!-- Geostrophic parameters
1260    WRITE ( io, 410 )  omega, phi, f, fs
1261
1262!
1263!-- Other quantities
1264    WRITE ( io, 411 )  g
1265    WRITE ( io, 412 )  TRIM( reference_state )
1266    IF ( use_single_reference_value )  THEN
1267       IF ( ocean )  THEN
1268          WRITE ( io, 413 )  prho_reference
1269       ELSE
1270          WRITE ( io, 414 )  pt_reference
1271       ENDIF
1272    ENDIF
1273
1274!
1275!-- Cloud physics parameters
1276    IF ( cloud_physics )  THEN
1277       WRITE ( io, 415 )
1278       WRITE ( io, 416 ) surface_pressure, r_d, rho_surface, cp, l_v
1279       IF ( icloud_scheme == 0 )  THEN
1280          WRITE ( io, 510 ) 1.0E-6 * nc_const
1281          IF ( precipitation )  WRITE ( io, 511 ) c_sedimentation
1282       ENDIF
1283    ENDIF
1284
1285!-- Profile of the geostrophic wind (component ug)
1286!-- Building output strings
1287    WRITE ( ugcomponent, '(F6.2)' )  ug_surface
1288    gradients = '------'
1289    slices = '     0'
1290    coordinates = '   0.0'
1291    i = 1
1292    DO  WHILE ( ug_vertical_gradient_level_ind(i) /= -9999 )
1293     
1294       WRITE (coor_chr,'(F6.2,1X)')  ug(ug_vertical_gradient_level_ind(i))
1295       ugcomponent = TRIM( ugcomponent ) // '  ' // TRIM( coor_chr )
1296
1297       WRITE (coor_chr,'(F6.2,1X)')  ug_vertical_gradient(i)
1298       gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
1299
1300       WRITE (coor_chr,'(I6,1X)')  ug_vertical_gradient_level_ind(i)
1301       slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
1302
1303       WRITE (coor_chr,'(F6.1,1X)')  ug_vertical_gradient_level(i)
1304       coordinates = TRIM( coordinates ) // '  ' // TRIM( coor_chr )
1305
1306       IF ( i == 10 )  THEN
1307          EXIT
1308       ELSE
1309          i = i + 1
1310       ENDIF
1311
1312    ENDDO
1313
1314    IF ( .NOT. large_scale_forcing )  THEN
1315       WRITE ( io, 423 )  TRIM( coordinates ), TRIM( ugcomponent ), &
1316                          TRIM( gradients ), TRIM( slices )
1317    ELSE
1318       WRITE ( io, 429 ) 
1319    ENDIF
1320
1321!-- Profile of the geostrophic wind (component vg)
1322!-- Building output strings
1323    WRITE ( vgcomponent, '(F6.2)' )  vg_surface
1324    gradients = '------'
1325    slices = '     0'
1326    coordinates = '   0.0'
1327    i = 1
1328    DO  WHILE ( vg_vertical_gradient_level_ind(i) /= -9999 )
1329
1330       WRITE (coor_chr,'(F6.2,1X)')  vg(vg_vertical_gradient_level_ind(i))
1331       vgcomponent = TRIM( vgcomponent ) // '  ' // TRIM( coor_chr )
1332
1333       WRITE (coor_chr,'(F6.2,1X)')  vg_vertical_gradient(i)
1334       gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
1335
1336       WRITE (coor_chr,'(I6,1X)')  vg_vertical_gradient_level_ind(i)
1337       slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
1338
1339       WRITE (coor_chr,'(F6.1,1X)')  vg_vertical_gradient_level(i)
1340       coordinates = TRIM( coordinates ) // '  ' // TRIM( coor_chr )
1341
1342       IF ( i == 10 )  THEN
1343          EXIT
1344       ELSE
1345          i = i + 1
1346       ENDIF
1347 
1348    ENDDO
1349
1350    IF ( .NOT. large_scale_forcing )  THEN
1351       WRITE ( io, 424 )  TRIM( coordinates ), TRIM( vgcomponent ), &
1352                          TRIM( gradients ), TRIM( slices )
1353    ENDIF
1354
1355!
1356!-- Initial wind profiles
1357    IF ( u_profile(1) /= 9999999.9 )  WRITE ( io, 427 )
1358
1359!
1360!-- Initial temperature profile
1361!-- Building output strings, starting with surface temperature
1362    WRITE ( temperatures, '(F6.2)' )  pt_surface
1363    gradients = '------'
1364    slices = '     0'
1365    coordinates = '   0.0'
1366    i = 1
1367    DO  WHILE ( pt_vertical_gradient_level_ind(i) /= -9999 )
1368
1369       WRITE (coor_chr,'(F7.2)')  pt_init(pt_vertical_gradient_level_ind(i))
1370       temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr )
1371
1372       WRITE (coor_chr,'(F7.2)')  pt_vertical_gradient(i)
1373       gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
1374
1375       WRITE (coor_chr,'(I7)')  pt_vertical_gradient_level_ind(i)
1376       slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
1377
1378       WRITE (coor_chr,'(F7.1)')  pt_vertical_gradient_level(i)
1379       coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
1380
1381       IF ( i == 10 )  THEN
1382          EXIT
1383       ELSE
1384          i = i + 1
1385       ENDIF
1386
1387    ENDDO
1388
1389    IF ( .NOT. nudging )  THEN
1390       WRITE ( io, 420 )  TRIM( coordinates ), TRIM( temperatures ), &
1391                          TRIM( gradients ), TRIM( slices )
1392    ELSE
1393       WRITE ( io, 428 ) 
1394    ENDIF
1395
1396!
1397!-- Initial humidity profile
1398!-- Building output strings, starting with surface humidity
1399    IF ( humidity  .OR.  passive_scalar )  THEN
1400       WRITE ( temperatures, '(E8.1)' )  q_surface
1401       gradients = '--------'
1402       slices = '       0'
1403       coordinates = '     0.0'
1404       i = 1
1405       DO  WHILE ( q_vertical_gradient_level_ind(i) /= -9999 )
1406         
1407          WRITE (coor_chr,'(E8.1,4X)')  q_init(q_vertical_gradient_level_ind(i))
1408          temperatures = TRIM( temperatures ) // '  ' // TRIM( coor_chr )
1409
1410          WRITE (coor_chr,'(E8.1,4X)')  q_vertical_gradient(i)
1411          gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
1412         
1413          WRITE (coor_chr,'(I8,4X)')  q_vertical_gradient_level_ind(i)
1414          slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
1415         
1416          WRITE (coor_chr,'(F8.1,4X)')  q_vertical_gradient_level(i)
1417          coordinates = TRIM( coordinates ) // '  '  // TRIM( coor_chr )
1418
1419          IF ( i == 10 )  THEN
1420             EXIT
1421          ELSE
1422             i = i + 1
1423          ENDIF
1424
1425       ENDDO
1426
1427       IF ( humidity )  THEN
1428          IF ( .NOT. nudging )  THEN
1429             WRITE ( io, 421 )  TRIM( coordinates ), TRIM( temperatures ), &
1430                                TRIM( gradients ), TRIM( slices )
1431          ENDIF
1432       ELSE
1433          WRITE ( io, 422 )  TRIM( coordinates ), TRIM( temperatures ), &
1434                             TRIM( gradients ), TRIM( slices )
1435       ENDIF
1436    ENDIF
1437
1438!
1439!-- Initial salinity profile
1440!-- Building output strings, starting with surface salinity
1441    IF ( ocean )  THEN
1442       WRITE ( temperatures, '(F6.2)' )  sa_surface
1443       gradients = '------'
1444       slices = '     0'
1445       coordinates = '   0.0'
1446       i = 1
1447       DO  WHILE ( sa_vertical_gradient_level_ind(i) /= -9999 )
1448
1449          WRITE (coor_chr,'(F7.2)')  sa_init(sa_vertical_gradient_level_ind(i))
1450          temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr )
1451
1452          WRITE (coor_chr,'(F7.2)')  sa_vertical_gradient(i)
1453          gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
1454
1455          WRITE (coor_chr,'(I7)')  sa_vertical_gradient_level_ind(i)
1456          slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
1457
1458          WRITE (coor_chr,'(F7.1)')  sa_vertical_gradient_level(i)
1459          coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
1460
1461          IF ( i == 10 )  THEN
1462             EXIT
1463          ELSE
1464             i = i + 1
1465          ENDIF
1466
1467       ENDDO
1468
1469       WRITE ( io, 425 )  TRIM( coordinates ), TRIM( temperatures ), &
1470                          TRIM( gradients ), TRIM( slices )
1471    ENDIF
1472
1473!
1474!-- Profile for the large scale vertial velocity
1475!-- Building output strings, starting with surface value
1476    IF ( large_scale_subsidence )  THEN
1477       temperatures = '   0.0'
1478       gradients = '------'
1479       slices = '     0'
1480       coordinates = '   0.0'
1481       i = 1
1482       DO  WHILE ( subs_vertical_gradient_level_i(i) /= -9999 )
1483
1484          WRITE (coor_chr,'(E10.2,7X)')  &
1485                                w_subs(subs_vertical_gradient_level_i(i))
1486          temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr )
1487
1488          WRITE (coor_chr,'(E10.2,7X)')  subs_vertical_gradient(i)
1489          gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
1490
1491          WRITE (coor_chr,'(I10,7X)')  subs_vertical_gradient_level_i(i)
1492          slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
1493
1494          WRITE (coor_chr,'(F10.2,7X)')  subs_vertical_gradient_level(i)
1495          coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
1496
1497          IF ( i == 10 )  THEN
1498             EXIT
1499          ELSE
1500             i = i + 1
1501          ENDIF
1502
1503       ENDDO
1504
1505 
1506       IF ( .NOT. large_scale_forcing )  THEN
1507          WRITE ( io, 426 )  TRIM( coordinates ), TRIM( temperatures ), &
1508                             TRIM( gradients ), TRIM( slices )
1509        ELSE
1510           WRITE ( io, 460 ) 
1511       ENDIF
1512
1513
1514    ENDIF
1515
1516!
1517!-- Cloud physcis parameters / quantities / numerical methods
1518    WRITE ( io, 430 )
1519    IF ( humidity .AND. .NOT. cloud_physics .AND. .NOT. cloud_droplets)  THEN
1520       WRITE ( io, 431 )
1521    ELSEIF ( humidity  .AND.  cloud_physics )  THEN
1522       WRITE ( io, 432 )
1523       IF ( radiation )  WRITE ( io, 132 )
1524       IF ( icloud_scheme == 1 )  THEN
1525          IF ( precipitation )  WRITE ( io, 133 )
1526       ELSEIF ( icloud_scheme == 0 )  THEN
1527          IF ( drizzle )  WRITE ( io, 506 )
1528          IF ( precipitation )  THEN
1529             WRITE ( io, 505 )
1530             IF ( turbulence )  WRITE ( io, 507 )
1531             IF ( ventilation_effect )  WRITE ( io, 508 )
1532             IF ( limiter_sedimentation )  WRITE ( io, 509 )
1533          ENDIF
1534       ENDIF
1535    ELSEIF ( humidity  .AND.  cloud_droplets )  THEN
1536       WRITE ( io, 433 )
1537       IF ( curvature_solution_effects )  WRITE ( io, 434 )
1538       IF ( collision_kernel /= 'none' )  THEN
1539          WRITE ( io, 435 )  TRIM( collision_kernel )
1540          IF ( collision_kernel(6:9) == 'fast' )  THEN
1541             WRITE ( io, 436 )  radius_classes, dissipation_classes
1542          ENDIF
1543       ELSE
1544          WRITE ( io, 437 )
1545       ENDIF
1546    ENDIF
1547
1548!
1549!-- LES / turbulence parameters
1550    WRITE ( io, 450 )
1551
1552!--
1553! ... LES-constants used must still be added here
1554!--
1555    IF ( constant_diffusion )  THEN
1556       WRITE ( io, 451 )  km_constant, km_constant/prandtl_number, &
1557                          prandtl_number
1558    ENDIF
1559    IF ( .NOT. constant_diffusion)  THEN
1560       IF ( e_init > 0.0 )  WRITE ( io, 455 )  e_init
1561       IF ( e_min > 0.0 )  WRITE ( io, 454 )  e_min
1562       IF ( wall_adjustment )  WRITE ( io, 453 )  wall_adjustment_factor
1563    ENDIF
1564
1565!
1566!-- Special actions during the run
1567    WRITE ( io, 470 )
1568    IF ( create_disturbances )  THEN
1569       WRITE ( io, 471 )  dt_disturb, disturbance_amplitude,                   &
1570                          zu(disturbance_level_ind_b), disturbance_level_ind_b,&
1571                          zu(disturbance_level_ind_t), disturbance_level_ind_t
1572       IF ( .NOT. bc_lr_cyc  .OR.  .NOT. bc_ns_cyc )  THEN
1573          WRITE ( io, 472 )  inflow_disturbance_begin, inflow_disturbance_end
1574       ELSE
1575          WRITE ( io, 473 )  disturbance_energy_limit
1576       ENDIF
1577       WRITE ( io, 474 )  TRIM( random_generator )
1578    ENDIF
1579    IF ( pt_surface_initial_change /= 0.0 )  THEN
1580       WRITE ( io, 475 )  pt_surface_initial_change
1581    ENDIF
1582    IF ( humidity  .AND.  q_surface_initial_change /= 0.0 )  THEN
1583       WRITE ( io, 476 )  q_surface_initial_change       
1584    ENDIF
1585    IF ( passive_scalar  .AND.  q_surface_initial_change /= 0.0 )  THEN
1586       WRITE ( io, 477 )  q_surface_initial_change       
1587    ENDIF
1588
1589    IF ( particle_advection )  THEN
1590!
1591!--    Particle attributes
1592       WRITE ( io, 480 )  particle_advection_start, dt_prel, bc_par_lr, &
1593                          bc_par_ns, bc_par_b, bc_par_t, particle_maximum_age, &
1594                          end_time_prel, dt_sort_particles
1595       IF ( use_sgs_for_particles )  WRITE ( io, 488 )  dt_min_part
1596       IF ( random_start_position )  WRITE ( io, 481 )
1597       IF ( particles_per_point > 1 )  WRITE ( io, 489 )  particles_per_point
1598       WRITE ( io, 495 )  total_number_of_particles
1599       IF ( use_particle_tails  .AND.  maximum_number_of_tailpoints /= 0 )  THEN
1600          WRITE ( io, 483 )  maximum_number_of_tailpoints
1601          IF ( minimum_tailpoint_distance /= 0 )  THEN
1602             WRITE ( io, 484 )  total_number_of_tails,      &
1603                                minimum_tailpoint_distance, &
1604                                maximum_tailpoint_age
1605          ENDIF
1606       ENDIF
1607       IF ( dt_write_particle_data /= 9999999.9 )  THEN
1608          WRITE ( io, 485 )  dt_write_particle_data
1609          output_format = ' '
1610          IF ( netcdf_output )  THEN
1611             IF ( netcdf_data_format > 1 )  THEN
1612                output_format = 'netcdf (64 bit offset) and binary'
1613             ELSE
1614                output_format = 'netcdf and binary'
1615             ENDIF
1616          ELSE
1617             output_format = 'binary'
1618          ENDIF
1619          WRITE ( io, 344 )  output_format
1620       ENDIF
1621       IF ( dt_dopts /= 9999999.9 )  WRITE ( io, 494 )  dt_dopts
1622       IF ( write_particle_statistics )  WRITE ( io, 486 )
1623
1624       WRITE ( io, 487 )  number_of_particle_groups
1625
1626       DO  i = 1, number_of_particle_groups
1627          IF ( i == 1  .AND.  density_ratio(i) == 9999999.9 )  THEN
1628             WRITE ( io, 490 )  i, 0.0
1629             WRITE ( io, 492 )
1630          ELSE
1631             WRITE ( io, 490 )  i, radius(i)
1632             IF ( density_ratio(i) /= 0.0 )  THEN
1633                WRITE ( io, 491 )  density_ratio(i)
1634             ELSE
1635                WRITE ( io, 492 )
1636             ENDIF
1637          ENDIF
1638          WRITE ( io, 493 )  psl(i), psr(i), pss(i), psn(i), psb(i), pst(i), &
1639                             pdx(i), pdy(i), pdz(i)
1640          IF ( .NOT. vertical_particle_advection(i) )  WRITE ( io, 482 )
1641       ENDDO
1642
1643    ENDIF
1644
1645
1646!
1647!-- Parameters of 1D-model
1648    IF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
1649       WRITE ( io, 500 )  end_time_1d, dt_run_control_1d, dt_pr_1d, &
1650                          mixing_length_1d, dissipation_1d
1651       IF ( damp_level_ind_1d /= nzt+1 )  THEN
1652          WRITE ( io, 502 )  zu(damp_level_ind_1d), damp_level_ind_1d
1653       ENDIF
1654    ENDIF
1655
1656!
1657!-- User-defined informations
1658    CALL user_header( io )
1659
1660    WRITE ( io, 99 )
1661
1662!
1663!-- Write buffer contents to disc immediately
1664    CALL local_flush( io )
1665
1666!
1667!-- Here the FORMATs start
1668
1669 99 FORMAT (1X,78('-'))
1670100 FORMAT (/1X,'******************************',6X,42('-')/        &
1671            1X,'* ',A,' *',6X,A/                               &
1672            1X,'******************************',6X,42('-'))
1673101 FORMAT (37X,'coupled run using MPI-',I1,': ',A/ &
1674            37X,42('-'))
1675102 FORMAT (/' Date:                 ',A8,6X,'Run:       ',A20/      &
1676            ' Time:                 ',A8,6X,'Run-No.:   ',I2.2/     &
1677            ' Run on host:        ',A10)
1678#if defined( __parallel )
1679103 FORMAT (' Number of PEs:',10X,I6,6X,'Processor grid (x,y): (',I3,',',I3, &
1680              ')',1X,A)
1681104 FORMAT (' Number of PEs:',8X,I5,9X,'Tasks:',I4,'   threads per task:',I4/ &
1682              37X,'Processor grid (x,y): (',I3,',',I3,')',1X,A)
1683105 FORMAT (37X,'One additional PE is used to handle'/37X,'the dvrp output!')
1684106 FORMAT (37X,'A 1d-decomposition along x is forced'/ &
1685            37X,'because the job is running on an SMP-cluster')
1686107 FORMAT (37X,'A 1d-decomposition along ',A,' is used')
1687108 FORMAT (37X,'Max. # of parallel I/O streams is ',I5)
1688109 FORMAT (37X,'Precursor run for coupled atmos-ocean run'/ &
1689            37X,42('-'))
1690114 FORMAT (37X,'Coupled atmosphere-ocean run following'/ &
1691            37X,'independent precursor runs'/             &
1692            37X,42('-'))
1693117 FORMAT (' Accelerator boards / node:  ',I2)
1694#endif
1695110 FORMAT (/' Numerical Schemes:'/ &
1696             ' -----------------'/)
1697111 FORMAT (' --> Solve perturbation pressure via FFT using ',A,' routines')
1698112 FORMAT (' --> Solve perturbation pressure via SOR-Red/Black-Schema'/ &
1699            '     Iterations (initial/other): ',I3,'/',I3,'  omega = ',F5.3)
1700113 FORMAT (' --> Momentum advection via Piascek-Williams-Scheme (Form C3)', &
1701                  ' or Upstream')
1702115 FORMAT ('     FFT and transpositions are overlapping')
1703116 FORMAT (' --> Scalar advection via Piascek-Williams-Scheme (Form C3)', &
1704                  ' or Upstream')
1705118 FORMAT (' --> Scalar advection via Bott-Chlond-Scheme')
1706119 FORMAT (' --> Galilei-Transform applied to horizontal advection:'/ &
1707            '     translation velocity = ',A/ &
1708            '     distance advected ',A,':  ',F8.3,' km(x)  ',F8.3,' km(y)')
1709120 FORMAT (' Accelerator boards: ',8X,I2)
1710122 FORMAT (' --> Time differencing scheme: ',A)
1711123 FORMAT (' --> Rayleigh-Damping active, starts ',A,' z = ',F8.2,' m'/ &
1712            '     maximum damping coefficient: ',F5.3, ' 1/s')
1713129 FORMAT (' --> Additional prognostic equation for the specific humidity')
1714130 FORMAT (' --> Additional prognostic equation for the total water content')
1715131 FORMAT (' --> No pt-equation solved. Neutral stratification with pt = ', &
1716                  F6.2, ' K assumed')
1717132 FORMAT ('     Parameterization of long-wave radiation processes via'/ &
1718            '     effective emissivity scheme')
1719133 FORMAT ('     Precipitation parameterization via Kessler-Scheme')
1720134 FORMAT (' --> Additional prognostic equation for a passive scalar')
1721135 FORMAT (' --> Solve perturbation pressure via multigrid method (', &
1722                  A,'-cycle)'/ &
1723            '     number of grid levels:                   ',I2/ &
1724            '     Gauss-Seidel red/black iterations:       ',I2)
1725136 FORMAT ('     gridpoints of coarsest subdomain (x,y,z): (',I3,',',I3,',', &
1726                  I3,')')
1727137 FORMAT ('     level data gathered on PE0 at level:     ',I2/ &
1728            '     gridpoints of coarsest subdomain (x,y,z): (',I3,',',I3,',', &
1729                  I3,')'/ &
1730            '     gridpoints of coarsest domain (x,y,z):    (',I3,',',I3,',', &
1731                  I3,')')
1732139 FORMAT (' --> Loop optimization method: ',A)
1733140 FORMAT ('     maximum residual allowed:                ',E10.3)
1734141 FORMAT ('     fixed number of multigrid cycles:        ',I4)
1735142 FORMAT ('     perturbation pressure is calculated at every Runge-Kutta ', &
1736                  'step')
1737143 FORMAT ('     Euler/upstream scheme is used for the SGS turbulent ', &
1738                  'kinetic energy')
1739144 FORMAT ('     masking method is used')
1740150 FORMAT (' --> Volume flow at the right and north boundary will be ', &
1741                  'conserved'/ &
1742            '     using the ',A,' mode')
1743151 FORMAT ('     with u_bulk = ',F7.3,' m/s and v_bulk = ',F7.3,' m/s')
1744152 FORMAT (' --> External pressure gradient directly prescribed by the user:',&
1745           /'     ',2(1X,E12.5),'Pa/m in x/y direction', &
1746           /'     starting from dp_level_b =', F8.3, 'm', A /)
1747153 FORMAT (' --> Large-scale vertical motion is used in the ', &
1748                  'prognostic equation(s) for')
1749154 FORMAT ('     the scalar(s) only')
1750155 FORMAT (' --> Nudging is used - initial profiles for u, v, pt and q ',& 
1751                  'correspond to the')
1752156 FORMAT ('     first profiles in NUDGING_DATA')
1753157 FORMAT (' --> Large scale forcing from external file (LSF_DATA) is used: ')
1754158 FORMAT ('     prescribed surfaces fluxes and geostrophic wind components')
1755200 FORMAT (//' Run time and time step information:'/ &
1756             ' ----------------------------------'/)
1757201 FORMAT ( ' Timestep:             variable     maximum value: ',F6.3,' s', &
1758             '    CFL-factor: ',F4.2)
1759202 FORMAT ( ' Timestep:          dt = ',F6.3,' s'/)
1760203 FORMAT ( ' Start time:          ',F9.3,' s'/ &
1761             ' End time:            ',F9.3,' s')
1762204 FORMAT ( A,F9.3,' s')
1763205 FORMAT ( A,F9.3,' s',5X,'restart every',17X,F9.3,' s')
1764206 FORMAT (/' Time reached:        ',F9.3,' s'/ &
1765             ' CPU-time used:       ',F9.3,' s     per timestep:               ', &
1766               '  ',F9.3,' s'/                                                    &
1767             '                                      per second of simulated tim', &
1768               'e: ',F9.3,' s')
1769207 FORMAT ( ' Coupling start time: ',F9.3,' s')
1770250 FORMAT (//' Computational grid and domain size:'/ &
1771              ' ----------------------------------'// &
1772              ' Grid length:      dx =    ',F7.3,' m    dy =    ',F7.3, &
1773              ' m    dz =    ',F7.3,' m'/ &
1774              ' Domain size:       x = ',F10.3,' m     y = ',F10.3, &
1775              ' m  z(u) = ',F10.3,' m'/)
1776252 FORMAT (' dz constant up to ',F10.3,' m (k=',I4,'), then stretched by', &
1777              ' factor: ',F5.3/ &
1778            ' maximum dz not to be exceeded is dz_max = ',F10.3,' m'/)
1779254 FORMAT (' Number of gridpoints (x,y,z):  (0:',I4,', 0:',I4,', 0:',I4,')'/ &
1780            ' Subdomain size (x,y,z):        (  ',I4,',   ',I4,',   ',I4,')'/)
1781260 FORMAT (/' The model has a slope in x-direction. Inclination angle: ',F6.2,&
1782             ' degrees')
1783270 FORMAT (//' Topography informations:'/ &
1784              ' -----------------------'// &
1785              1X,'Topography: ',A)
1786271 FORMAT (  ' Building size (x/y/z) in m: ',F5.1,' / ',F5.1,' / ',F5.1/ &
1787              ' Horizontal index bounds (l/r/s/n): ',I4,' / ',I4,' / ',I4, &
1788                ' / ',I4)
1789272 FORMAT (  ' Single quasi-2D street canyon of infinite length in ',A, &
1790              ' direction' / &
1791              ' Canyon height: ', F6.2, 'm, ch = ', I4, '.'      / &
1792              ' Canyon position (',A,'-walls): cxl = ', I4,', cxr = ', I4, '.')
1793278 FORMAT (' Topography grid definition convention:'/ &
1794            ' cell edge (staggered grid points'/  &
1795            ' (u in x-direction, v in y-direction))' /)
1796279 FORMAT (' Topography grid definition convention:'/ &
1797            ' cell center (scalar grid points)' /)
1798280 FORMAT (//' Vegetation canopy (drag) model:'/ &
1799              ' ------------------------------'// &
1800              ' Canopy mode: ', A / &
1801              ' Canopy top: ',I4 / &
1802              ' Leaf drag coefficient: ',F6.2 /)
1803281 FORMAT (/ ' Scalar_exchange_coefficient: ',F6.2 / &
1804              ' Scalar concentration at leaf surfaces in kg/m**3: ',F6.2 /)
1805282 FORMAT (' Predefined constant heatflux at the top of the vegetation: ',F6.2,' K m/s')
1806283 FORMAT (/ ' Characteristic levels of the leaf area density:'// &
1807              ' Height:              ',A,'  m'/ &
1808              ' Leaf area density:   ',A,'  m**2/m**3'/ &
1809              ' Gradient:            ',A,'  m**2/m**4'/ &
1810              ' Gridpoint:           ',A)
1811               
1812300 FORMAT (//' Boundary conditions:'/ &
1813             ' -------------------'// &
1814             '                     p                    uv             ', &
1815             '                   pt'// &
1816             ' B. bound.: ',A/ &
1817             ' T. bound.: ',A)
1818301 FORMAT (/'                     ',A// &
1819             ' B. bound.: ',A/ &
1820             ' T. bound.: ',A)
1821303 FORMAT (/' Bottom surface fluxes are used in diffusion terms at k=1')
1822304 FORMAT (/' Top surface fluxes are used in diffusion terms at k=nzt')
1823305 FORMAT (//'    Prandtl-Layer between bottom surface and first ', &
1824               'computational u,v-level:'// &
1825             '       zp = ',F6.2,' m   z0 = ',F6.4,' m   z0h = ',F7.5,&
1826             ' m   kappa = ',F4.2/ &
1827             '       Rif value range:   ',F6.2,' <= rif <=',F6.2)
1828306 FORMAT ('       Predefined constant heatflux:   ',F9.6,' K m/s')
1829307 FORMAT ('       Heatflux has a random normal distribution')
1830308 FORMAT ('       Predefined surface temperature')
1831309 FORMAT ('       Predefined constant salinityflux:   ',F9.6,' psu m/s')
1832310 FORMAT (//'    1D-Model:'// &
1833             '       Rif value range:   ',F6.2,' <= rif <=',F6.2)
1834311 FORMAT ('       Predefined constant humidity flux: ',E10.3,' m/s')
1835312 FORMAT ('       Predefined surface humidity')
1836313 FORMAT ('       Predefined constant scalar flux: ',E10.3,' kg/(m**2 s)')
1837314 FORMAT ('       Predefined scalar value at the surface')
1838315 FORMAT ('       Humidity / scalar flux at top surface is 0.0')
1839316 FORMAT ('       Sensible heatflux and momentum flux from coupled ', &
1840                    'atmosphere model')
1841317 FORMAT (//' Lateral boundaries:'/ &
1842            '       left/right:  ',A/    &
1843            '       north/south: ',A)
1844318 FORMAT (/'       use_cmax: ',L1 / &
1845            '       pt damping layer width = ',F8.2,' m, pt ', &
1846                    'damping factor = ',F6.4)
1847319 FORMAT ('       turbulence recycling at inflow switched on'/ &
1848            '       width of recycling domain: ',F7.1,' m   grid index: ',I4/ &
1849            '       inflow damping height: ',F6.1,' m   width: ',F6.1,' m')
1850320 FORMAT ('       Predefined constant momentumflux:  u: ',F9.6,' m**2/s**2'/ &
1851            '                                          v: ',F9.6,' m**2/s**2')
1852325 FORMAT (//' List output:'/ &
1853             ' -----------'//  &
1854            '    1D-Profiles:'/    &
1855            '       Output every             ',F8.2,' s')
1856326 FORMAT ('       Time averaged over       ',F8.2,' s'/ &
1857            '       Averaging input every    ',F8.2,' s')
1858330 FORMAT (//' Data output:'/ &
1859             ' -----------'/)
1860331 FORMAT (/'    1D-Profiles:')
1861332 FORMAT (/'       ',A)
1862333 FORMAT ('       Output every             ',F8.2,' s',/ &
1863            '       Time averaged over       ',F8.2,' s'/ &
1864            '       Averaging input every    ',F8.2,' s')
1865334 FORMAT (/'    2D-Arrays',A,':')
1866335 FORMAT (/'       ',A2,'-cross-section  Arrays: ',A/ &
1867            '       Output every             ',F8.2,' s  ',A/ &
1868            '       Cross sections at ',A1,' = ',A/ &
1869            '       scalar-coordinates:   ',A,' m'/)
1870336 FORMAT (/'    3D-Arrays',A,':')
1871337 FORMAT (/'       Arrays: ',A/ &
1872            '       Output every             ',F8.2,' s  ',A/ &
1873            '       Upper output limit at    ',F8.2,' m  (GP ',I4,')'/)
1874338 FORMAT ('       Compressed data output'/ &
1875            '       Decimal precision: ',A/)
1876339 FORMAT ('       No output during initial ',F8.2,' s')
1877340 FORMAT (/'    Time series:')
1878341 FORMAT ('       Output every             ',F8.2,' s'/)
1879342 FORMAT (/'       ',A2,'-cross-section  Arrays: ',A/ &
1880            '       Output every             ',F8.2,' s  ',A/ &
1881            '       Time averaged over       ',F8.2,' s'/ &
1882            '       Averaging input every    ',F8.2,' s'/ &
1883            '       Cross sections at ',A1,' = ',A/ &
1884            '       scalar-coordinates:   ',A,' m'/)
1885343 FORMAT (/'       Arrays: ',A/ &
1886            '       Output every             ',F8.2,' s  ',A/ &
1887            '       Time averaged over       ',F8.2,' s'/ &
1888            '       Averaging input every    ',F8.2,' s'/ &
1889            '       Upper output limit at    ',F8.2,' m  (GP ',I4,')'/)
1890344 FORMAT ('       Output format: ',A/)
1891345 FORMAT (/'    Scaling lengths for output locations of all subsequent mask IDs:',/ &
1892            '       mask_scale_x (in x-direction): ',F9.3, ' m',/ &
1893            '       mask_scale_y (in y-direction): ',F9.3, ' m',/ &
1894            '       mask_scale_z (in z-direction): ',F9.3, ' m' )
1895346 FORMAT (/'    Masked data output',A,' for mask ID ',I2, ':')
1896347 FORMAT ('       Variables: ',A/ &
1897            '       Output every             ',F8.2,' s')
1898348 FORMAT ('       Variables: ',A/ &
1899            '       Output every             ',F8.2,' s'/ &
1900            '       Time averaged over       ',F8.2,' s'/ &
1901            '       Averaging input every    ',F8.2,' s')
1902349 FORMAT (/'       Output locations in ',A,'-direction in multiples of ', &
1903            'mask_scale_',A,' predefined by array mask_',I2.2,'_',A,':'/ &
1904            13('       ',8(F8.2,',')/) )
1905350 FORMAT (/'       Output locations in ',A,'-direction: ', &
1906            'all gridpoints along ',A,'-direction (default).' )
1907351 FORMAT (/'       Output locations in ',A,'-direction in multiples of ', &
1908            'mask_scale_',A,' constructed from array mask_',I2.2,'_',A,'_loop:'/ &
1909            '          loop begin:',F8.2,', end:',F8.2,', stride:',F8.2 )
1910#if defined( __dvrp_graphics )
1911360 FORMAT ('    Plot-Sequence with dvrp-software:'/ &
1912            '       Output every      ',F7.1,' s'/ &
1913            '       Output mode:      ',A/ &
1914            '       Host / User:      ',A,' / ',A/ &
1915            '       Directory:        ',A// &
1916            '       The sequence contains:')
1917361 FORMAT (/'       Isosurface of "',A,'"    Threshold value: ', E12.3/ &
1918            '          Isosurface color: (',F4.2,',',F4.2,',',F4.2,') (R,G,B)')
1919362 FORMAT (/'       Slicer plane ',A/ &
1920            '       Slicer limits: [',F6.2,',',F6.2,']')
1921363 FORMAT (/'       Particles'/ &
1922            '          particle size:  ',F7.2,' m')
1923364 FORMAT ('          particle ',A,' controlled by "',A,'" with interval [', &
1924                       F6.2,',',F6.2,']')
1925365 FORMAT (/'       Groundplate color: (',F4.2,',',F4.2,',',F4.2,') (R,G,B)'/ &
1926            '       Superelevation along (x,y,z): (',F4.1,',',F4.1,',',F4.1, &
1927                     ')'/ &
1928            '       Clipping limits: from x = ',F9.1,' m to x = ',F9.1,' m'/ &
1929            '                        from y = ',F9.1,' m to y = ',F9.1,' m')
1930366 FORMAT (/'       Topography color: (',F4.2,',',F4.2,',',F4.2,') (R,G,B)')
1931367 FORMAT ('       Polygon reduction for topography: cluster_size = ', I1)
1932#endif
1933#if defined( __spectra )
1934370 FORMAT ('    Spectra:')
1935371 FORMAT ('       Output every ',F7.1,' s'/)
1936372 FORMAT ('       Arrays:     ', 10(A5,',')/                         &
1937            '       Directions: ', 10(A5,',')/                         &
1938            '       height levels  k = ', 20(I3,',')/                  &
1939            '                          ', 20(I3,',')/                  &
1940            '                          ', 20(I3,',')/                  &
1941            '                          ', 20(I3,',')/                  &
1942            '                          ', 19(I3,','),I3,'.'/           &
1943            '       height levels selected for standard plot:'/        &
1944            '                      k = ', 20(I3,',')/                  &
1945            '                          ', 20(I3,',')/                  &
1946            '                          ', 20(I3,',')/                  &
1947            '                          ', 20(I3,',')/                  &
1948            '                          ', 19(I3,','),I3,'.'/           &
1949            '       Time averaged over ', F7.1, ' s,' /                &
1950            '       Profiles for the time averaging are taken every ', &
1951                    F6.1,' s')
1952#endif
1953400 FORMAT (//' Physical quantities:'/ &
1954              ' -------------------'/)
1955410 FORMAT ('    Angular velocity    :   omega = ',E9.3,' rad/s'/  &
1956            '    Geograph. latitude  :   phi   = ',F4.1,' degr'/   &
1957            '    Coriolis parameter  :   f     = ',F9.6,' 1/s'/    &
1958            '                            f*    = ',F9.6,' 1/s')
1959411 FORMAT (/'    Gravity             :   g     = ',F4.1,' m/s**2')
1960412 FORMAT (/'    Reference state used in buoyancy terms: ',A)
1961413 FORMAT ('       Reference density in buoyancy terms: ',F8.3,' kg/m**3')
1962414 FORMAT ('       Reference temperature in buoyancy terms: ',F8.4,' K')
1963415 FORMAT (/'    Cloud physics parameters:'/ &
1964             '    ------------------------'/)
1965416 FORMAT ('        Surface pressure   :   p_0   = ',F7.2,' hPa'/      &
1966            '        Gas constant       :   R     = ',F5.1,' J/(kg K)'/ &
1967            '        Density of air     :   rho_0 = ',F5.3,' kg/m**3'/  &
1968            '        Specific heat cap. :   c_p   = ',F6.1,' J/(kg K)'/ &
1969            '        Vapourization heat :   L_v   = ',E8.2,' J/kg')
1970420 FORMAT (/'    Characteristic levels of the initial temperature profile:'// &
1971            '       Height:        ',A,'  m'/ &
1972            '       Temperature:   ',A,'  K'/ &
1973            '       Gradient:      ',A,'  K/100m'/ &
1974            '       Gridpoint:     ',A)
1975421 FORMAT (/'    Characteristic levels of the initial humidity profile:'// &
1976            '       Height:      ',A,'  m'/ &
1977            '       Humidity:    ',A,'  kg/kg'/ &
1978            '       Gradient:    ',A,'  (kg/kg)/100m'/ &
1979            '       Gridpoint:   ',A)
1980422 FORMAT (/'    Characteristic levels of the initial scalar profile:'// &
1981            '       Height:                  ',A,'  m'/ &
1982            '       Scalar concentration:    ',A,'  kg/m**3'/ &
1983            '       Gradient:                ',A,'  (kg/m**3)/100m'/ &
1984            '       Gridpoint:               ',A)
1985423 FORMAT (/'    Characteristic levels of the geo. wind component ug:'// &
1986            '       Height:      ',A,'  m'/ &
1987            '       ug:          ',A,'  m/s'/ &
1988            '       Gradient:    ',A,'  1/100s'/ &
1989            '       Gridpoint:   ',A)
1990424 FORMAT (/'    Characteristic levels of the geo. wind component vg:'// &
1991            '       Height:      ',A,'  m'/ &
1992            '       vg:          ',A,'  m/s'/ &
1993            '       Gradient:    ',A,'  1/100s'/ &
1994            '       Gridpoint:   ',A)
1995425 FORMAT (/'    Characteristic levels of the initial salinity profile:'// &
1996            '       Height:     ',A,'  m'/ &
1997            '       Salinity:   ',A,'  psu'/ &
1998            '       Gradient:   ',A,'  psu/100m'/ &
1999            '       Gridpoint:  ',A)
2000426 FORMAT (/'    Characteristic levels of the subsidence/ascent profile:'// &
2001            '       Height:      ',A,'  m'/ &
2002            '       w_subs:      ',A,'  m/s'/ &
2003            '       Gradient:    ',A,'  (m/s)/100m'/ &
2004            '       Gridpoint:   ',A)
2005427 FORMAT (/'    Initial wind profiles (u,v) are interpolated from given'// &
2006                  ' profiles')
2007428 FORMAT (/'    Initial profiles (u, v, pt, q) are taken from file '/ &
2008             '    NUDGING_DATA')
2009429 FORMAT (/'    Geostrophic wind profiles (ug, vg) are are taken from file '/ &
2010             '    LSF_DATA')
2011430 FORMAT (//' Cloud physics quantities / methods:'/ &
2012              ' ----------------------------------'/)
2013431 FORMAT ('    Humidity is treated as purely passive scalar (no condensati', &
2014                 'on)')
2015432 FORMAT ('    Bulk scheme with liquid water potential temperature and'/ &
2016            '    total water content is used.'/ &
2017            '    Condensation is parameterized via 0% - or 100% scheme.')
2018433 FORMAT ('    Cloud droplets treated explicitly using the Lagrangian part', &
2019                 'icle model')
2020434 FORMAT ('    Curvature and solution effecs are considered for growth of', &
2021                 ' droplets < 1.0E-6 m')
2022435 FORMAT ('    Droplet collision is handled by ',A,'-kernel')
2023436 FORMAT ('       Fast kernel with fixed radius- and dissipation classes ', &
2024                    'are used'/ &
2025            '          number of radius classes:       ',I3,'    interval ', &
2026                       '[1.0E-6,2.0E-4] m'/ &
2027            '          number of dissipation classes:   ',I2,'    interval ', &
2028                       '[0,1000] cm**2/s**3')
2029437 FORMAT ('    Droplet collision is switched off')
2030450 FORMAT (//' LES / Turbulence quantities:'/ &
2031              ' ---------------------------'/)
2032451 FORMAT ('    Diffusion coefficients are constant:'/ &
2033            '    Km = ',F6.2,' m**2/s   Kh = ',F6.2,' m**2/s   Pr = ',F5.2)
2034453 FORMAT ('    Mixing length is limited to ',F4.2,' * z')
2035454 FORMAT ('    TKE is not allowed to fall below ',E9.2,' (m/s)**2')
2036455 FORMAT ('    initial TKE is prescribed as ',E9.2,' (m/s)**2')
2037460 FORMAT (/'    Profiles for large scale vertical velocity are '/ &
2038             '    taken from file LSF_DATA')
2039470 FORMAT (//' Actions during the simulation:'/ &
2040              ' -----------------------------'/)
2041471 FORMAT ('    Disturbance impulse (u,v) every :   ',F6.2,' s'/            &
2042            '    Disturbance amplitude           :     ',F4.2, ' m/s'/       &
2043            '    Lower disturbance level         : ',F8.2,' m (GP ',I4,')'/  &
2044            '    Upper disturbance level         : ',F8.2,' m (GP ',I4,')')
2045472 FORMAT ('    Disturbances continued during the run from i/j =',I4, &
2046                 ' to i/j =',I4)
2047473 FORMAT ('    Disturbances cease as soon as the disturbance energy exceeds',&
2048                 1X,F5.3, ' m**2/s**2')
2049474 FORMAT ('    Random number generator used    : ',A/)
2050475 FORMAT ('    The surface temperature is increased (or decreased, ', &
2051                 'respectively, if'/ &
2052            '    the value is negative) by ',F5.2,' K at the beginning of the',&
2053                 ' 3D-simulation'/)
2054476 FORMAT ('    The surface humidity is increased (or decreased, ',&
2055                 'respectively, if the'/ &
2056            '    value is negative) by ',E8.1,' kg/kg at the beginning of', &
2057                 ' the 3D-simulation'/)
2058477 FORMAT ('    The scalar value is increased at the surface (or decreased, ',&
2059                 'respectively, if the'/ &
2060            '    value is negative) by ',E8.1,' kg/m**3 at the beginning of', &
2061                 ' the 3D-simulation'/)
2062480 FORMAT ('    Particles:'/ &
2063            '    ---------'// &
2064            '       Particle advection is active (switched on at t = ', F7.1, &
2065                    ' s)'/ &
2066            '       Start of new particle generations every  ',F6.1,' s'/ &
2067            '       Boundary conditions: left/right: ', A, ' north/south: ', A/&
2068            '                            bottom:     ', A, ' top:         ', A/&
2069            '       Maximum particle age:                 ',F9.1,' s'/ &
2070            '       Advection stopped at t = ',F9.1,' s'/ &
2071            '       Particles are sorted every ',F9.1,' s'/)
2072481 FORMAT ('       Particles have random start positions'/)
2073482 FORMAT ('          Particles are advected only horizontally'/)
2074483 FORMAT ('       Particles have tails with a maximum of ',I3,' points')
2075484 FORMAT ('            Number of tails of the total domain: ',I10/ &
2076            '            Minimum distance between tailpoints: ',F8.2,' m'/ &
2077            '            Maximum age of the end of the tail:  ',F8.2,' s')
2078485 FORMAT ('       Particle data are written on file every ', F9.1, ' s')
2079486 FORMAT ('       Particle statistics are written on file'/)
2080487 FORMAT ('       Number of particle groups: ',I2/)
2081488 FORMAT ('       SGS velocity components are used for particle advection'/ &
2082            '          minimum timestep for advection: ', F7.5/)
2083489 FORMAT ('       Number of particles simultaneously released at each ', &
2084                    'point: ', I5/)
2085490 FORMAT ('       Particle group ',I2,':'/ &
2086            '          Particle radius: ',E10.3, 'm')
2087491 FORMAT ('          Particle inertia is activated'/ &
2088            '             density_ratio (rho_fluid/rho_particle) = ',F5.3/)
2089492 FORMAT ('          Particles are advected only passively (no inertia)'/)
2090493 FORMAT ('          Boundaries of particle source: x:',F8.1,' - ',F8.1,' m'/&
2091            '                                         y:',F8.1,' - ',F8.1,' m'/&
2092            '                                         z:',F8.1,' - ',F8.1,' m'/&
2093            '          Particle distances:  dx = ',F8.1,' m  dy = ',F8.1, &
2094                       ' m  dz = ',F8.1,' m'/)
2095494 FORMAT ('       Output of particle time series in NetCDF format every ', &
2096                    F8.2,' s'/)
2097495 FORMAT ('       Number of particles in total domain: ',I10/)
2098500 FORMAT (//' 1D-Model parameters:'/                           &
2099              ' -------------------'//                           &
2100            '    Simulation time:                   ',F8.1,' s'/ &
2101            '    Run-controll output every:         ',F8.1,' s'/ &
2102            '    Vertical profile output every:     ',F8.1,' s'/ &
2103            '    Mixing length calculation:         ',A/         &
2104            '    Dissipation calculation:           ',A/)
2105502 FORMAT ('    Damping layer starts from ',F7.1,' m (GP ',I4,')'/)
2106503 FORMAT (' --> Momentum advection via Wicker-Skamarock-Scheme 5th order')
2107504 FORMAT (' --> Scalar advection via Wicker-Skamarock-Scheme 5th order')
2108505 FORMAT ('    Precipitation parameterization via Seifert-Beheng-Scheme')
2109506 FORMAT ('    Drizzle parameterization via Stokes law')
2110507 FORMAT ('    Turbulence effects on precipitation process')
2111508 FORMAT ('    Ventilation effects on evaporation of rain drops')
2112509 FORMAT ('    Slope limiter used for sedimentation process')
2113510 FORMAT ('        Droplet density    :   N_c   = ',F6.1,' 1/cm**3')
2114511 FORMAT ('        Sedimentation Courant number:                  '/&
2115            '                               C_s   = ',F3.1,'        ')
2116
2117 END SUBROUTINE header
Note: See TracBrowser for help on using the repository browser.