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

Last change on this file since 211 was 206, checked in by raasch, 15 years ago

ocean-atmosphere coupling realized with MPI-1, adjustments in mrun, mbuild, subjob for lcxt4

  • Property svn:keywords set to Id
File size: 63.2 KB
Line 
1 SUBROUTINE header
2
3!------------------------------------------------------------------------------!
4! Actual revisions:
5! -----------------
6! Bugfix: error in zu index in case of section_xy = -1
7!
8! Former revisions:
9! -----------------
10! $Id: header.f90 206 2008-10-13 14:59:11Z raasch $
11!
12! 198 2008-09-17 08:55:28Z raasch
13! Format adjustments allowing output of larger revision numbers
14!
15! 197 2008-09-16 15:29:03Z raasch
16! allow 100 spectra levels instead of 10 for consistency with
17! define_netcdf_header,
18! bugfix in the output of the characteristic levels of potential temperature,
19! geostrophic wind, scalar concentration, humidity and leaf area density,
20! output of turbulence recycling informations
21!
22! 138 2007-11-28 10:03:58Z letzel
23! Allow new case bc_uv_t = 'dirichlet_0' for channel flow.
24! Allow two instead of one digit to specify isosurface and slicer variables.
25! Output of sorting frequency of particles
26!
27! 108 2007-08-24 15:10:38Z letzel
28! Output of informations for coupled model runs (boundary conditions etc.)
29! + output of momentumfluxes at the top boundary
30! Rayleigh damping for ocean, e_init
31!
32! 97 2007-06-21 08:23:15Z raasch
33! Adjustments for the ocean version.
34! use_pt_reference renamed use_reference
35!
36! 87 2007-05-22 15:46:47Z raasch
37! Bugfix: output of use_upstream_for_tke
38!
39! 82 2007-04-16 15:40:52Z raasch
40! Preprocessor strings for different linux clusters changed to "lc",
41! routine local_flush is used for buffer flushing
42!
43! 76 2007-03-29 00:58:32Z raasch
44! Output of netcdf_64bit_3d, particles-package is now part of the default code,
45! output of the loop optimization method, moisture renamed humidity,
46! output of subversion revision number
47!
48! 19 2007-02-23 04:53:48Z raasch
49! Output of scalar flux applied at top boundary
50!
51! RCS Log replace by Id keyword, revision history cleaned up
52!
53! Revision 1.63  2006/08/22 13:53:13  raasch
54! Output of dz_max
55!
56! Revision 1.1  1997/08/11 06:17:20  raasch
57! Initial revision
58!
59!
60! Description:
61! ------------
62! Writing a header with all important informations about the actual run.
63! This subroutine is called three times, two times at the beginning
64! (writing information on files RUN_CONTROL and HEADER) and one time at the
65! end of the run, then writing additional information about CPU-usage on file
66! header.
67!------------------------------------------------------------------------------!
68
69    USE arrays_3d
70    USE control_parameters
71    USE cloud_parameters
72    USE cpulog
73    USE dvrp_variables
74    USE grid_variables
75    USE indices
76    USE model_1d
77    USE particle_attributes
78    USE pegrid
79    USE spectrum
80
81    IMPLICIT NONE
82
83    CHARACTER (LEN=1)  ::  prec
84    CHARACTER (LEN=2)  ::  do2d_mode
85    CHARACTER (LEN=5)  ::  section_chr
86    CHARACTER (LEN=9)  ::  time_to_string
87    CHARACTER (LEN=10) ::  coor_chr, host_chr
88    CHARACTER (LEN=16) ::  begin_chr
89    CHARACTER (LEN=23) ::  ver_rev
90    CHARACTER (LEN=40) ::  output_format
91    CHARACTER (LEN=70) ::  char1, char2, dopr_chr, &
92                           do2d_xy, do2d_xz, do2d_yz, do3d_chr, &
93                           run_classification
94    CHARACTER (LEN=86) ::  coordinates, gradients, learde, slices,  &
95                           temperatures, ugcomponent, vgcomponent
96    CHARACTER (LEN=85) ::  roben, runten
97
98    INTEGER ::  av, bh, blx, bly, bxl, bxr, byn, bys, i, ihost, io, j, l, ll
99    REAL    ::  cpuseconds_per_simulated_second
100
101!
102!-- Open the output file. At the end of the simulation, output is directed
103!-- to unit 19.
104    IF ( ( runnr == 0 .OR. force_print_header )  .AND. &
105         .NOT. simulated_time_at_begin /= simulated_time )  THEN
106       io = 15   !  header output on file RUN_CONTROL
107    ELSE
108       io = 19   !  header output on file HEADER
109    ENDIF
110    CALL check_open( io )
111
112!
113!-- At the end of the run, output file (HEADER) will be rewritten with
114!-- new informations
115    IF ( io == 19 .AND. simulated_time_at_begin /= simulated_time ) REWIND( 19 )
116
117!
118!-- Determine kind of model run
119    IF ( TRIM( initializing_actions ) == 'read_restart_data' )  THEN
120       run_classification = '3D - restart run'
121    ELSEIF ( TRIM( initializing_actions ) == 'read_data_for_recycling' )  THEN
122       run_classification = '3D - run using 3D - prerun data'
123    ELSEIF ( INDEX( initializing_actions, 'set_constant_profiles' ) /= 0 )  THEN
124       run_classification = '3D - run without 1D - prerun'
125    ELSEIF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
126       run_classification = '3D - run with 1D - prerun'
127    ELSEIF ( INDEX( initializing_actions, 'by_user' ) /=0 )  THEN
128       run_classification = '3D - run initialized by user'
129    ELSE
130       PRINT*,'+++ header:  unknown action(s): ',initializing_actions
131    ENDIF
132    IF ( ocean )  THEN
133       run_classification = 'ocean - ' // run_classification
134    ELSE
135       run_classification = 'atmosphere - ' // run_classification
136    ENDIF
137
138!
139!-- Run-identification, date, time, host
140    host_chr = host(1:10)
141    ver_rev = TRIM( version ) // '  ' // TRIM( revision )
142    WRITE ( io, 100 )  ver_rev, TRIM( run_classification )
143    IF ( coupling_mode /= 'uncoupled' )  WRITE ( io, 101 )  coupling_mode
144    WRITE ( io, 102 )  run_date, run_identifier, run_time, runnr, &
145                       ADJUSTR( host_chr )
146#if defined( __parallel )
147    IF ( npex == -1  .AND.  pdims(2) /= 1 )  THEN
148       char1 = 'calculated'
149    ELSEIF ( ( host(1:3) == 'ibm'  .OR.  host(1:3) == 'nec'  .OR.  &
150               host(1:2) == 'lc' )  .AND.                          &
151             npex == -1  .AND.  pdims(2) == 1 )  THEN
152       char1 = 'forced'
153    ELSE
154       char1 = 'predefined'
155    ENDIF
156    IF ( threads_per_task == 1 )  THEN
157       WRITE ( io, 103 )  numprocs, pdims(1), pdims(2), TRIM( char1 )
158    ELSE
159       WRITE ( io, 104 )  numprocs*threads_per_task, numprocs, &
160                          threads_per_task, pdims(1), pdims(2), TRIM( char1 )
161    ENDIF
162    IF ( ( host(1:3) == 'ibm'  .OR.  host(1:3) == 'nec'  .OR.    &
163           host(1:2) == 'lc'   .OR.  host(1:3) == 'dec' )  .AND. &
164         npex == -1  .AND.  pdims(2) == 1 )                      &
165    THEN
166       WRITE ( io, 106 )
167    ELSEIF ( pdims(2) == 1 )  THEN
168       WRITE ( io, 107 )  'x'
169    ELSEIF ( pdims(1) == 1 )  THEN
170       WRITE ( io, 107 )  'y'
171    ENDIF
172    IF ( use_seperate_pe_for_dvrp_output )  WRITE ( io, 105 )
173#endif
174    WRITE ( io, 99 )
175
176!
177!-- Numerical schemes
178    WRITE ( io, 110 )
179    IF ( psolver(1:7) == 'poisfft' )  THEN
180       WRITE ( io, 111 )  TRIM( fft_method )
181       IF ( psolver == 'poisfft_hybrid' )  WRITE ( io, 138 )
182    ELSEIF ( psolver == 'sor' )  THEN
183       WRITE ( io, 112 )  nsor_ini, nsor, omega_sor
184    ELSEIF ( psolver == 'multigrid' )  THEN
185       WRITE ( io, 135 )  cycle_mg, maximum_grid_level, ngsrb
186       IF ( mg_cycles == -1 )  THEN
187          WRITE ( io, 140 )  residual_limit
188       ELSE
189          WRITE ( io, 141 )  mg_cycles
190       ENDIF
191       IF ( mg_switch_to_pe0_level == 0 )  THEN
192          WRITE ( io, 136 )  nxr_mg(1)-nxl_mg(1)+1, nyn_mg(1)-nys_mg(1)+1, &
193                             nzt_mg(1)
194       ELSEIF (  mg_switch_to_pe0_level /= -1 )  THEN
195          WRITE ( io, 137 )  mg_switch_to_pe0_level,            &
196                             mg_loc_ind(2,0)-mg_loc_ind(1,0)+1, &
197                             mg_loc_ind(4,0)-mg_loc_ind(3,0)+1, &
198                             nzt_mg(mg_switch_to_pe0_level),    &
199                             nxr_mg(1)-nxl_mg(1)+1, nyn_mg(1)-nys_mg(1)+1, &
200                             nzt_mg(1)
201       ENDIF
202    ENDIF
203    IF ( call_psolver_at_all_substeps  .AND. timestep_scheme(1:5) == 'runge' ) &
204    THEN
205       WRITE ( io, 142 )
206    ENDIF
207
208    IF ( momentum_advec == 'pw-scheme' )  THEN
209       WRITE ( io, 113 )
210    ELSE
211       WRITE ( io, 114 )
212       IF ( cut_spline_overshoot )  WRITE ( io, 124 )
213       IF ( overshoot_limit_u /= 0.0  .OR.  overshoot_limit_v /= 0.0  .OR. &
214            overshoot_limit_w /= 0.0 )  THEN
215          WRITE ( io, 127 )  overshoot_limit_u, overshoot_limit_v, &
216                             overshoot_limit_w
217       ENDIF
218       IF ( ups_limit_u /= 0.0  .OR.  ups_limit_v /= 0.0  .OR. &
219            ups_limit_w /= 0.0 )                               &
220       THEN
221          WRITE ( io, 125 )  ups_limit_u, ups_limit_v, ups_limit_w
222       ENDIF
223       IF ( long_filter_factor /= 0.0 )  WRITE ( io, 115 )  long_filter_factor
224    ENDIF
225    IF ( scalar_advec == 'pw-scheme' )  THEN
226       WRITE ( io, 116 )
227    ELSEIF ( scalar_advec == 'ups-scheme' )  THEN
228       WRITE ( io, 117 )
229       IF ( cut_spline_overshoot )  WRITE ( io, 124 )
230       IF ( overshoot_limit_e /= 0.0  .OR.  overshoot_limit_pt /= 0.0 )  THEN
231          WRITE ( io, 128 )  overshoot_limit_e, overshoot_limit_pt
232       ENDIF
233       IF ( ups_limit_e /= 0.0  .OR.  ups_limit_pt /= 0.0 )  THEN
234          WRITE ( io, 126 )  ups_limit_e, ups_limit_pt
235       ENDIF
236    ELSE
237       WRITE ( io, 118 )
238    ENDIF
239
240    WRITE ( io, 139 )  TRIM( loop_optimization )
241
242    IF ( galilei_transformation )  THEN
243       IF ( use_ug_for_galilei_tr )  THEN
244          char1 = 'geostrophic wind'
245       ELSE
246          char1 = 'mean wind in model domain'
247       ENDIF
248       IF ( simulated_time_at_begin == simulated_time )  THEN
249          char2 = 'at the start of the run'
250       ELSE
251          char2 = 'at the end of the run'
252       ENDIF
253       WRITE ( io, 119 )  TRIM( char1 ), TRIM( char2 ), &
254                          advected_distance_x/1000.0, advected_distance_y/1000.0
255    ENDIF
256    IF ( timestep_scheme == 'leapfrog' )  THEN
257       WRITE ( io, 120 )
258    ELSEIF ( timestep_scheme == 'leapfrog+euler' )  THEN
259       WRITE ( io, 121 )
260    ELSE
261       WRITE ( io, 122 )  timestep_scheme
262    ENDIF
263    IF ( use_upstream_for_tke )  WRITE ( io, 143 )
264    IF ( rayleigh_damping_factor /= 0.0 )  THEN
265       IF ( .NOT. ocean )  THEN
266          WRITE ( io, 123 )  'above', rayleigh_damping_height, &
267               rayleigh_damping_factor
268       ELSE
269          WRITE ( io, 123 )  'below', rayleigh_damping_height, &
270               rayleigh_damping_factor
271       ENDIF
272    ENDIF
273    IF ( humidity )  THEN
274       IF ( .NOT. cloud_physics )  THEN
275          WRITE ( io, 129 )
276       ELSE
277          WRITE ( io, 130 )
278          WRITE ( io, 131 )
279          IF ( radiation )      WRITE ( io, 132 )
280          IF ( precipitation )  WRITE ( io, 133 )
281       ENDIF
282    ENDIF
283    IF ( passive_scalar )  WRITE ( io, 134 )
284    IF ( conserve_volume_flow )  WRITE ( io, 150 )
285    WRITE ( io, 99 )
286
287!
288!-- Runtime and timestep informations
289    WRITE ( io, 200 )
290    IF ( .NOT. dt_fixed )  THEN
291       WRITE ( io, 201 )  dt_max, cfl_factor
292    ELSE
293       WRITE ( io, 202 )  dt
294    ENDIF
295    WRITE ( io, 203 )  simulated_time_at_begin, end_time
296
297    IF ( time_restart /= 9999999.9  .AND. &
298         simulated_time_at_begin == simulated_time )  THEN
299       IF ( dt_restart == 9999999.9 )  THEN
300          WRITE ( io, 204 )  ' Restart at:       ',time_restart
301       ELSE
302          WRITE ( io, 205 )  ' Restart at:       ',time_restart, dt_restart
303       ENDIF
304    ENDIF
305
306    IF ( simulated_time_at_begin /= simulated_time )  THEN
307       i = MAX ( log_point_s(10)%counts, 1 )
308       IF ( ( simulated_time - simulated_time_at_begin ) == 0.0 )  THEN
309          cpuseconds_per_simulated_second = 0.0
310       ELSE
311          cpuseconds_per_simulated_second = log_point_s(10)%sum / &
312                                            ( simulated_time -    &
313                                              simulated_time_at_begin )
314       ENDIF
315       WRITE ( io, 206 )  simulated_time, log_point_s(10)%sum, &
316                          log_point_s(10)%sum / REAL( i ),     &
317                          cpuseconds_per_simulated_second
318       IF ( time_restart /= 9999999.9  .AND.  time_restart < end_time )  THEN
319          IF ( dt_restart == 9999999.9 )  THEN
320             WRITE ( io, 204 )  ' Next restart at:  ',time_restart
321          ELSE
322             WRITE ( io, 205 )  ' Next restart at:  ',time_restart, dt_restart
323          ENDIF
324       ENDIF
325    ENDIF
326
327!
328!-- Computational grid
329    IF ( .NOT. ocean )  THEN
330       WRITE ( io, 250 )  dx, dy, dz, (nx+1)*dx, (ny+1)*dy, zu(nzt+1)
331       IF ( dz_stretch_level_index < nzt+1 )  THEN
332          WRITE ( io, 252 )  dz_stretch_level, dz_stretch_level_index, &
333                             dz_stretch_factor, dz_max
334       ENDIF
335    ELSE
336       WRITE ( io, 250 )  dx, dy, dz, (nx+1)*dx, (ny+1)*dy, zu(0)
337       IF ( dz_stretch_level_index > 0 )  THEN
338          WRITE ( io, 252 )  dz_stretch_level, dz_stretch_level_index, &
339                             dz_stretch_factor, dz_max
340       ENDIF
341    ENDIF
342    WRITE ( io, 254 )  nx, ny, nzt+1, MIN( nnx, nx+1 ), MIN( nny, ny+1 ), &
343                       MIN( nnz+2, nzt+2 )
344    IF ( numprocs > 1 )  THEN
345       IF ( nxa == nx  .AND.  nya == ny  .AND.  nza == nz )  THEN
346          WRITE ( io, 255 )
347       ELSE
348          WRITE ( io, 256 )  nnx-(nxa-nx), nny-(nya-ny), nzt+2
349       ENDIF
350    ENDIF
351    IF ( sloping_surface )  WRITE ( io, 260 )  alpha_surface
352
353!
354!-- Topography
355    WRITE ( io, 270 )  topography
356    SELECT CASE ( TRIM( topography ) )
357
358       CASE ( 'flat' )
359          ! no actions necessary
360
361       CASE ( 'single_building' )
362          blx = INT( building_length_x / dx )
363          bly = INT( building_length_y / dy )
364          bh  = INT( building_height / dz )
365
366          IF ( building_wall_left == 9999999.9 )  THEN
367             building_wall_left = ( nx + 1 - blx ) / 2 * dx
368          ENDIF
369          bxl = INT ( building_wall_left / dx + 0.5 )
370          bxr = bxl + blx
371
372          IF ( building_wall_south == 9999999.9 )  THEN
373             building_wall_south = ( ny + 1 - bly ) / 2 * dy
374          ENDIF
375          bys = INT ( building_wall_south / dy + 0.5 )
376          byn = bys + bly
377
378          WRITE ( io, 271 )  building_length_x, building_length_y, &
379                             building_height, bxl, bxr, bys, byn
380
381    END SELECT
382
383    IF ( plant_canopy ) THEN
384
385       WRITE ( io, 280 ) canopy_mode, pch_index, drag_coefficient
386       IF ( passive_scalar ) THEN
387          WRITE ( io, 281 ) scalar_exchange_coefficient,   &
388                            leaf_surface_concentration
389       ENDIF
390
391!
392!--    Heat flux at the top of vegetation
393       WRITE ( io, 282 ) cthf
394
395!
396!--    Leaf area density profile
397!--    Building output strings, starting with surface value
398       WRITE ( learde, '(F6.2)' )  lad_surface
399       gradients = '------'
400       slices = '     0'
401       coordinates = '   0.0'
402       i = 1
403       DO  WHILE ( lad_vertical_gradient_level_ind(i) /= -9999 )
404
405          WRITE (coor_chr,'(F7.2)')  lad(lad_vertical_gradient_level_ind(i))
406          learde = TRIM( learde ) // ' ' // TRIM( coor_chr )
407
408          WRITE (coor_chr,'(F7.2)')  lad_vertical_gradient(i)
409          gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
410
411          WRITE (coor_chr,'(I7)')  lad_vertical_gradient_level_ind(i)
412          slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
413
414          WRITE (coor_chr,'(F7.1)')  lad_vertical_gradient_level(i)
415          coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
416
417          i = i + 1
418       ENDDO
419
420       WRITE ( io, 283 )  TRIM( coordinates ), TRIM( learde ), &
421                          TRIM( gradients ), TRIM( slices )
422
423    ENDIF
424
425!
426!-- Boundary conditions
427    IF ( ibc_p_b == 0 )  THEN
428       runten = 'p(0)     = 0      |'
429    ELSEIF ( ibc_p_b == 1 )  THEN
430       runten = 'p(0)     = p(1)   |'
431    ELSE
432       runten = 'p(0)     = p(1) +R|'
433    ENDIF
434    IF ( ibc_p_t == 0 )  THEN
435       roben  = 'p(nzt+1) = 0      |'
436    ELSE
437       roben  = 'p(nzt+1) = p(nzt) |'
438    ENDIF
439
440    IF ( ibc_uv_b == 0 )  THEN
441       runten = TRIM( runten ) // ' uv(0)     = -uv(1)                |'
442    ELSE
443       runten = TRIM( runten ) // ' uv(0)     = uv(1)                 |'
444    ENDIF
445    IF ( TRIM( bc_uv_t ) == 'dirichlet_0' )  THEN
446       roben  = TRIM( roben  ) // ' uv(nzt+1) = 0                     |'
447    ELSEIF ( ibc_uv_t == 0 )  THEN
448       roben  = TRIM( roben  ) // ' uv(nzt+1) = ug(nzt+1), vg(nzt+1)  |'
449    ELSE
450       roben  = TRIM( roben  ) // ' uv(nzt+1) = uv(nzt)               |'
451    ENDIF
452
453    IF ( ibc_pt_b == 0 )  THEN
454       runten = TRIM( runten ) // ' pt(0)   = pt_surface'
455    ELSEIF ( ibc_pt_b == 1 )  THEN
456       runten = TRIM( runten ) // ' pt(0)   = pt(1)'
457    ELSEIF ( ibc_pt_b == 2 )  THEN
458       runten = TRIM( runten ) // ' pt(0) = from coupled model'
459    ENDIF
460    IF ( ibc_pt_t == 0 )  THEN
461       roben  = TRIM( roben  ) // ' pt(nzt+1) = pt_top'
462    ELSEIF( ibc_pt_t == 1 )  THEN
463       roben  = TRIM( roben  ) // ' pt(nzt+1) = pt(nzt)'
464    ELSEIF( ibc_pt_t == 2 )  THEN
465       roben  = TRIM( roben  ) // ' pt(nzt+1) = pt(nzt) + dpt/dz_ini'
466    ENDIF
467
468    WRITE ( io, 300 )  runten, roben
469
470    IF ( .NOT. constant_diffusion )  THEN
471       IF ( ibc_e_b == 1 )  THEN
472          runten = 'e(0)     = e(1)'
473       ELSE
474          runten = 'e(0)     = e(1) = (u*/0.1)**2'
475       ENDIF
476       roben = 'e(nzt+1) = e(nzt) = e(nzt-1)'
477
478       WRITE ( io, 301 )  'e', runten, roben       
479
480    ENDIF
481
482    IF ( ocean )  THEN
483       runten = 'sa(0)    = sa(1)'
484       IF ( ibc_sa_t == 0 )  THEN
485          roben =  'sa(nzt+1) = sa_surface'
486       ELSE
487          roben =  'sa(nzt+1) = sa(nzt)'
488       ENDIF
489       WRITE ( io, 301 ) 'sa', runten, roben
490    ENDIF
491
492    IF ( humidity )  THEN
493       IF ( ibc_q_b == 0 )  THEN
494          runten = 'q(0)     = q_surface'
495       ELSE
496          runten = 'q(0)     = q(1)'
497       ENDIF
498       IF ( ibc_q_t == 0 )  THEN
499          roben =  'q(nzt)   = q_top'
500       ELSE
501          roben =  'q(nzt)   = q(nzt-1) + dq/dz'
502       ENDIF
503       WRITE ( io, 301 ) 'q', runten, roben
504    ENDIF
505
506    IF ( passive_scalar )  THEN
507       IF ( ibc_q_b == 0 )  THEN
508          runten = 's(0)     = s_surface'
509       ELSE
510          runten = 's(0)     = s(1)'
511       ENDIF
512       IF ( ibc_q_t == 0 )  THEN
513          roben =  's(nzt)   = s_top'
514       ELSE
515          roben =  's(nzt)   = s(nzt-1) + ds/dz'
516       ENDIF
517       WRITE ( io, 301 ) 's', runten, roben
518    ENDIF
519
520    IF ( use_surface_fluxes )  THEN
521       WRITE ( io, 303 )
522       IF ( constant_heatflux )  THEN
523          WRITE ( io, 306 )  surface_heatflux
524          IF ( random_heatflux )  WRITE ( io, 307 )
525       ENDIF
526       IF ( humidity  .AND.  constant_waterflux )  THEN
527          WRITE ( io, 311 ) surface_waterflux
528       ENDIF
529       IF ( passive_scalar  .AND.  constant_waterflux )  THEN
530          WRITE ( io, 313 ) surface_waterflux
531       ENDIF
532    ENDIF
533
534    IF ( use_top_fluxes )  THEN
535       WRITE ( io, 304 )
536       IF ( coupling_mode == 'uncoupled' )  THEN
537          WRITE ( io, 320 )  top_momentumflux_u, top_momentumflux_v
538          IF ( constant_top_heatflux )  THEN
539             WRITE ( io, 306 )  top_heatflux
540          ENDIF
541       ELSEIF ( coupling_mode == 'ocean_to_atmosphere' )  THEN
542          WRITE ( io, 316 )
543       ENDIF
544       IF ( ocean  .AND.  constant_top_salinityflux )  THEN
545          WRITE ( io, 309 )  top_salinityflux
546       ENDIF
547       IF ( humidity  .OR.  passive_scalar )  THEN
548          WRITE ( io, 315 )
549       ENDIF
550    ENDIF
551
552    IF ( prandtl_layer )  THEN
553       WRITE ( io, 305 )  0.5 * (zu(1)-zu(0)), roughness_length, kappa, &
554                          rif_min, rif_max
555       IF ( .NOT. constant_heatflux )  WRITE ( io, 308 )
556       IF ( humidity  .AND.  .NOT. constant_waterflux )  THEN
557          WRITE ( io, 312 )
558       ENDIF
559       IF ( passive_scalar  .AND.  .NOT. constant_waterflux )  THEN
560          WRITE ( io, 314 )
561       ENDIF
562    ELSE
563       IF ( INDEX(initializing_actions, 'set_1d-model_profiles') /= 0 )  THEN
564          WRITE ( io, 310 )  rif_min, rif_max
565       ENDIF
566    ENDIF
567
568    WRITE ( io, 317 )  bc_lr, bc_ns
569    IF ( bc_lr /= 'cyclic'  .OR.  bc_ns /= 'cyclic' )  THEN
570       WRITE ( io, 318 )  outflow_damping_width, km_damp_max
571       IF ( turbulent_inflow )  THEN
572          WRITE ( io, 319 )  recycling_width, recycling_plane, &
573                             inflow_damping_height, inflow_damping_width
574       ENDIF
575    ENDIF
576
577!
578!-- Listing of 1D-profiles
579    WRITE ( io, 325 )  dt_dopr_listing
580    IF ( averaging_interval_pr /= 0.0 )  THEN
581       WRITE ( io, 326 )  averaging_interval_pr, dt_averaging_input_pr
582    ENDIF
583
584!
585!-- DATA output
586    WRITE ( io, 330 )
587    IF ( averaging_interval_pr /= 0.0 )  THEN
588       WRITE ( io, 326 )  averaging_interval_pr, dt_averaging_input_pr
589    ENDIF
590
591!
592!-- 1D-profiles
593    dopr_chr = 'Profile:'
594    IF ( dopr_n /= 0 )  THEN
595       WRITE ( io, 331 )
596
597       output_format = ''
598       IF ( netcdf_output )  THEN
599          IF ( netcdf_64bit )  THEN
600             output_format = 'netcdf (64 bit offset)'
601          ELSE
602             output_format = 'netcdf'
603          ENDIF
604       ENDIF
605       IF ( profil_output )  THEN
606          IF ( netcdf_output )  THEN
607             output_format = TRIM( output_format ) // ' and profil'
608          ELSE
609             output_format = 'profil'
610          ENDIF
611       ENDIF
612       WRITE ( io, 345 )  output_format
613
614       DO  i = 1, dopr_n
615          dopr_chr = TRIM( dopr_chr ) // ' ' // TRIM( data_output_pr(i) ) // ','
616          IF ( LEN_TRIM( dopr_chr ) >= 60 )  THEN
617             WRITE ( io, 332 )  dopr_chr
618             dopr_chr = '       :'
619          ENDIF
620       ENDDO
621
622       IF ( dopr_chr /= '' )  THEN
623          WRITE ( io, 332 )  dopr_chr
624       ENDIF
625       WRITE ( io, 333 )  dt_dopr, averaging_interval_pr, dt_averaging_input_pr
626       IF ( skip_time_dopr /= 0.0 )  WRITE ( io, 339 )  skip_time_dopr
627    ENDIF
628
629!
630!-- 2D-arrays
631    DO  av = 0, 1
632
633       i = 1
634       do2d_xy = ''
635       do2d_xz = ''
636       do2d_yz = ''
637       DO  WHILE ( do2d(av,i) /= ' ' )
638
639          l = MAX( 2, LEN_TRIM( do2d(av,i) ) )
640          do2d_mode = do2d(av,i)(l-1:l)
641
642          SELECT CASE ( do2d_mode )
643             CASE ( 'xy' )
644                ll = LEN_TRIM( do2d_xy )
645                do2d_xy = do2d_xy(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
646             CASE ( 'xz' )
647                ll = LEN_TRIM( do2d_xz )
648                do2d_xz = do2d_xz(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
649             CASE ( 'yz' )
650                ll = LEN_TRIM( do2d_yz )
651                do2d_yz = do2d_yz(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
652          END SELECT
653
654          i = i + 1
655
656       ENDDO
657
658       IF ( ( ( do2d_xy /= ''  .AND.  section(1,1) /= -9999 )  .OR.    &
659              ( do2d_xz /= ''  .AND.  section(1,2) /= -9999 )  .OR.    &
660              ( do2d_yz /= ''  .AND.  section(1,3) /= -9999 ) )  .AND. &
661            ( netcdf_output  .OR.  iso2d_output ) )  THEN
662
663          IF (  av == 0 )  THEN
664             WRITE ( io, 334 )  ''
665          ELSE
666             WRITE ( io, 334 )  '(time-averaged)'
667          ENDIF
668
669          IF ( do2d_at_begin )  THEN
670             begin_chr = 'and at the start'
671          ELSE
672             begin_chr = ''
673          ENDIF
674
675          output_format = ''
676          IF ( netcdf_output )  THEN
677             IF ( netcdf_64bit )  THEN
678                output_format = 'netcdf (64 bit offset)'
679             ELSE
680                output_format = 'netcdf'
681             ENDIF
682          ENDIF
683          IF ( iso2d_output )  THEN
684             IF ( netcdf_output )  THEN
685                output_format = TRIM( output_format ) // ' and iso2d'
686             ELSE
687                output_format = 'iso2d'
688             ENDIF
689          ENDIF
690          WRITE ( io, 345 )  output_format
691
692          IF ( do2d_xy /= ''  .AND.  section(1,1) /= -9999 )  THEN
693             i = 1
694             slices = '/'
695             coordinates = '/'
696!
697!--          Building strings with index and coordinate informations of the
698!--          slices
699             DO  WHILE ( section(i,1) /= -9999 )
700
701                WRITE (section_chr,'(I5)')  section(i,1)
702                section_chr = ADJUSTL( section_chr )
703                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
704
705                IF ( section(i,1) == -1 )  THEN
706                   WRITE (coor_chr,'(F10.1)')  -1.0
707                ELSE
708                   WRITE (coor_chr,'(F10.1)')  zu(section(i,1))
709                ENDIF
710                coor_chr = ADJUSTL( coor_chr )
711                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
712
713                i = i + 1
714             ENDDO
715             IF ( av == 0 )  THEN
716                WRITE ( io, 335 )  'XY', do2d_xy, dt_do2d_xy, &
717                                   TRIM( begin_chr ), 'k', TRIM( slices ), &
718                                   TRIM( coordinates )
719                IF ( skip_time_do2d_xy /= 0.0 )  THEN
720                   WRITE ( io, 339 )  skip_time_do2d_xy
721                ENDIF
722             ELSE
723                WRITE ( io, 342 )  'XY', do2d_xy, dt_data_output_av, &
724                                   TRIM( begin_chr ), averaging_interval, &
725                                   dt_averaging_input, 'k', TRIM( slices ), &
726                                   TRIM( coordinates )
727                IF ( skip_time_data_output_av /= 0.0 )  THEN
728                   WRITE ( io, 339 )  skip_time_data_output_av
729                ENDIF
730             ENDIF
731
732          ENDIF
733
734          IF ( do2d_xz /= ''  .AND.  section(1,2) /= -9999 )  THEN
735             i = 1
736             slices = '/'
737             coordinates = '/'
738!
739!--          Building strings with index and coordinate informations of the
740!--          slices
741             DO  WHILE ( section(i,2) /= -9999 )
742
743                WRITE (section_chr,'(I5)')  section(i,2)
744                section_chr = ADJUSTL( section_chr )
745                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
746
747                WRITE (coor_chr,'(F10.1)')  section(i,2) * dy
748                coor_chr = ADJUSTL( coor_chr )
749                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
750
751                i = i + 1
752             ENDDO
753             IF ( av == 0 )  THEN
754                WRITE ( io, 335 )  'XZ', do2d_xz, dt_do2d_xz, &
755                                   TRIM( begin_chr ), 'j', TRIM( slices ), &
756                                   TRIM( coordinates )
757                IF ( skip_time_do2d_xz /= 0.0 )  THEN
758                   WRITE ( io, 339 )  skip_time_do2d_xz
759                ENDIF
760             ELSE
761                WRITE ( io, 342 )  'XZ', do2d_xz, dt_data_output_av, &
762                                   TRIM( begin_chr ), averaging_interval, &
763                                   dt_averaging_input, 'j', TRIM( slices ), &
764                                   TRIM( coordinates )
765                IF ( skip_time_data_output_av /= 0.0 )  THEN
766                   WRITE ( io, 339 )  skip_time_data_output_av
767                ENDIF
768             ENDIF
769          ENDIF
770
771          IF ( do2d_yz /= ''  .AND.  section(1,3) /= -9999 )  THEN
772             i = 1
773             slices = '/'
774             coordinates = '/'
775!
776!--          Building strings with index and coordinate informations of the
777!--          slices
778             DO  WHILE ( section(i,3) /= -9999 )
779
780                WRITE (section_chr,'(I5)')  section(i,3)
781                section_chr = ADJUSTL( section_chr )
782                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
783
784                WRITE (coor_chr,'(F10.1)')  section(i,3) * dx
785                coor_chr = ADJUSTL( coor_chr )
786                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
787
788                i = i + 1
789             ENDDO
790             IF ( av == 0 )  THEN
791                WRITE ( io, 335 )  'YZ', do2d_yz, dt_do2d_yz, &
792                                   TRIM( begin_chr ), 'i', TRIM( slices ), &
793                                   TRIM( coordinates )
794                IF ( skip_time_do2d_yz /= 0.0 )  THEN
795                   WRITE ( io, 339 )  skip_time_do2d_yz
796                ENDIF
797             ELSE
798                WRITE ( io, 342 )  'YZ', do2d_yz, dt_data_output_av, &
799                                   TRIM( begin_chr ), averaging_interval, &
800                                   dt_averaging_input, 'i', TRIM( slices ), &
801                                   TRIM( coordinates )
802                IF ( skip_time_data_output_av /= 0.0 )  THEN
803                   WRITE ( io, 339 )  skip_time_data_output_av
804                ENDIF
805             ENDIF
806          ENDIF
807
808       ENDIF
809
810    ENDDO
811
812!
813!-- 3d-arrays
814    DO  av = 0, 1
815
816       i = 1
817       do3d_chr = ''
818       DO  WHILE ( do3d(av,i) /= ' ' )
819
820          do3d_chr = TRIM( do3d_chr ) // ' ' // TRIM( do3d(av,i) ) // ','
821          i = i + 1
822
823       ENDDO
824
825       IF ( do3d_chr /= '' )  THEN
826          IF ( av == 0 )  THEN
827             WRITE ( io, 336 )  ''
828          ELSE
829             WRITE ( io, 336 )  '(time-averaged)'
830          ENDIF
831
832          output_format = ''
833          IF ( netcdf_output )  THEN
834             IF ( netcdf_64bit .AND. netcdf_64bit_3d )  THEN
835                output_format = 'netcdf (64 bit offset)'
836             ELSE
837                output_format = 'netcdf'
838             ENDIF
839          ENDIF
840          IF ( avs_output )  THEN
841             IF ( netcdf_output )  THEN
842                output_format = TRIM( output_format ) // ' and avs'
843             ELSE
844                output_format = 'avs'
845             ENDIF
846          ENDIF
847          WRITE ( io, 345 )  output_format
848
849          IF ( do3d_at_begin )  THEN
850             begin_chr = 'and at the start'
851          ELSE
852             begin_chr = ''
853          ENDIF
854          IF ( av == 0 )  THEN
855             WRITE ( io, 337 )  do3d_chr, dt_do3d, TRIM( begin_chr ), &
856                                zu(nz_do3d), nz_do3d
857          ELSE
858             WRITE ( io, 343 )  do3d_chr, dt_data_output_av,           &
859                                TRIM( begin_chr ), averaging_interval, &
860                                dt_averaging_input, zu(nz_do3d), nz_do3d
861          ENDIF
862
863          IF ( do3d_compress )  THEN
864             do3d_chr = ''
865             i = 1
866             DO WHILE ( do3d(av,i) /= ' ' )
867
868                SELECT CASE ( do3d(av,i) )
869                   CASE ( 'u' )
870                      j = 1
871                   CASE ( 'v' )
872                      j = 2
873                   CASE ( 'w' )
874                      j = 3
875                   CASE ( 'p' )
876                      j = 4
877                   CASE ( 'pt' )
878                      j = 5
879                END SELECT
880                WRITE ( prec, '(I1)' )  plot_3d_precision(j)%precision
881                do3d_chr = TRIM( do3d_chr ) // ' ' // TRIM( do3d(av,i) ) // &
882                           ':' // prec // ','
883                i = i + 1
884
885             ENDDO
886             WRITE ( io, 338 )  do3d_chr
887
888          ENDIF
889
890          IF ( av == 0 )  THEN
891             IF ( skip_time_do3d /= 0.0 )  THEN
892                WRITE ( io, 339 )  skip_time_do3d
893             ENDIF
894          ELSE
895             IF ( skip_time_data_output_av /= 0.0 )  THEN
896                WRITE ( io, 339 )  skip_time_data_output_av
897             ENDIF
898          ENDIF
899
900       ENDIF
901
902    ENDDO
903
904!
905!-- Timeseries
906    IF ( dt_dots /= 9999999.9 )  THEN
907       WRITE ( io, 340 )
908
909       output_format = ''
910       IF ( netcdf_output )  THEN
911          IF ( netcdf_64bit )  THEN
912             output_format = 'netcdf (64 bit offset)'
913          ELSE
914             output_format = 'netcdf'
915          ENDIF
916       ENDIF
917       IF ( profil_output )  THEN
918          IF ( netcdf_output )  THEN
919             output_format = TRIM( output_format ) // ' and profil'
920          ELSE
921             output_format = 'profil'
922          ENDIF
923       ENDIF
924       WRITE ( io, 345 )  output_format
925       WRITE ( io, 341 )  dt_dots
926    ENDIF
927
928#if defined( __dvrp_graphics )
929!
930!-- Dvrp-output
931    IF ( dt_dvrp /= 9999999.9 )  THEN
932       WRITE ( io, 360 )  dt_dvrp, TRIM( dvrp_output ), TRIM( dvrp_host ), &
933                          TRIM( dvrp_username ), TRIM( dvrp_directory )
934       i = 1
935       l = 0
936       DO WHILE ( mode_dvrp(i) /= ' ' )
937          IF ( mode_dvrp(i)(1:10) == 'isosurface' )  THEN
938             READ ( mode_dvrp(i), '(10X,I2)' )  j
939             l = l + 1
940             IF ( do3d(0,j) /= ' ' )  THEN
941                WRITE ( io, 361 )  TRIM( do3d(0,j) ), threshold(l)
942             ENDIF
943          ELSEIF ( mode_dvrp(i)(1:6) == 'slicer' )  THEN
944             READ ( mode_dvrp(i), '(6X,I2)' )  j
945             IF ( do2d(0,j) /= ' ' )  WRITE ( io, 362 )  TRIM( do2d(0,j) )
946          ELSEIF ( mode_dvrp(i)(1:9) == 'particles' )  THEN
947             WRITE ( io, 363 )
948          ENDIF
949          i = i + 1
950       ENDDO
951    ENDIF
952#endif
953
954#if defined( __spectra )
955!
956!-- Spectra output
957    IF ( dt_dosp /= 9999999.9 ) THEN
958       WRITE ( io, 370 )
959
960       output_format = ''
961       IF ( netcdf_output )  THEN
962          IF ( netcdf_64bit )  THEN
963             output_format = 'netcdf (64 bit offset)'
964          ELSE
965             output_format = 'netcdf'
966          ENDIF
967       ENDIF
968       IF ( profil_output )  THEN
969          IF ( netcdf_output )  THEN
970             output_format = TRIM( output_format ) // ' and profil'
971          ELSE
972             output_format = 'profil'
973          ENDIF
974       ENDIF
975       WRITE ( io, 345 )  output_format
976       WRITE ( io, 371 )  dt_dosp
977       IF ( skip_time_dosp /= 0.0 )  WRITE ( io, 339 )  skip_time_dosp
978       WRITE ( io, 372 )  ( data_output_sp(i), i = 1,10 ),     &
979                          ( spectra_direction(i), i = 1,10 ),  &
980                          ( comp_spectra_level(i), i = 1,100 ), &
981                          ( plot_spectra_level(i), i = 1,100 ), &
982                          averaging_interval_sp, dt_averaging_input_pr
983    ENDIF
984#endif
985
986    WRITE ( io, 99 )
987
988!
989!-- Physical quantities
990    WRITE ( io, 400 )
991
992!
993!-- Geostrophic parameters
994    WRITE ( io, 410 )  omega, phi, f, fs
995
996!
997!-- Other quantities
998    WRITE ( io, 411 )  g
999    IF ( use_reference )  THEN
1000       IF ( ocean )  THEN
1001          WRITE ( io, 412 )  prho_reference
1002       ELSE
1003          WRITE ( io, 413 )  pt_reference
1004       ENDIF
1005    ENDIF
1006
1007!
1008!-- Cloud physics parameters
1009    IF ( cloud_physics ) THEN
1010       WRITE ( io, 415 )
1011       WRITE ( io, 416 ) surface_pressure, r_d, rho_surface, cp, l_v
1012    ENDIF
1013
1014!-- Profile of the geostrophic wind (component ug)
1015!-- Building output strings
1016    WRITE ( ugcomponent, '(F6.2)' )  ug_surface
1017    gradients = '------'
1018    slices = '     0'
1019    coordinates = '   0.0'
1020    i = 1
1021    DO  WHILE ( ug_vertical_gradient_level_ind(i) /= -9999 )
1022     
1023       WRITE (coor_chr,'(F6.2,1X)')  ug(ug_vertical_gradient_level_ind(i))
1024       ugcomponent = TRIM( ugcomponent ) // '  ' // TRIM( coor_chr )
1025
1026       WRITE (coor_chr,'(F6.2,1X)')  ug_vertical_gradient(i)
1027       gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
1028
1029       WRITE (coor_chr,'(I6,1X)')  ug_vertical_gradient_level_ind(i)
1030       slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
1031
1032       WRITE (coor_chr,'(F6.1,1X)')  ug_vertical_gradient_level(i)
1033       coordinates = TRIM( coordinates ) // '  ' // TRIM( coor_chr )
1034
1035       i = i + 1
1036    ENDDO
1037
1038    WRITE ( io, 423 )  TRIM( coordinates ), TRIM( ugcomponent ), &
1039                       TRIM( gradients ), TRIM( slices )
1040
1041!-- Profile of the geostrophic wind (component vg)
1042!-- Building output strings
1043    WRITE ( vgcomponent, '(F6.2)' )  vg_surface
1044    gradients = '------'
1045    slices = '     0'
1046    coordinates = '   0.0'
1047    i = 1
1048    DO  WHILE ( vg_vertical_gradient_level_ind(i) /= -9999 )
1049
1050       WRITE (coor_chr,'(F6.2,1X)')  vg(vg_vertical_gradient_level_ind(i))
1051       vgcomponent = TRIM( vgcomponent ) // '  ' // TRIM( coor_chr )
1052
1053       WRITE (coor_chr,'(F6.2,1X)')  vg_vertical_gradient(i)
1054       gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
1055
1056       WRITE (coor_chr,'(I6,1X)')  vg_vertical_gradient_level_ind(i)
1057       slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
1058
1059       WRITE (coor_chr,'(F6.1,1X)')  vg_vertical_gradient_level(i)
1060       coordinates = TRIM( coordinates ) // '  ' // TRIM( coor_chr )
1061
1062       i = i + 1 
1063    ENDDO
1064
1065    WRITE ( io, 424 )  TRIM( coordinates ), TRIM( vgcomponent ), &
1066                       TRIM( gradients ), TRIM( slices )
1067
1068!
1069!-- Initial temperature profile
1070!-- Building output strings, starting with surface temperature
1071    WRITE ( temperatures, '(F6.2)' )  pt_surface
1072    gradients = '------'
1073    slices = '     0'
1074    coordinates = '   0.0'
1075    i = 1
1076    DO  WHILE ( pt_vertical_gradient_level_ind(i) /= -9999 )
1077
1078       WRITE (coor_chr,'(F7.2)')  pt_init(pt_vertical_gradient_level_ind(i))
1079       temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr )
1080
1081       WRITE (coor_chr,'(F7.2)')  pt_vertical_gradient(i)
1082       gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
1083
1084       WRITE (coor_chr,'(I7)')  pt_vertical_gradient_level_ind(i)
1085       slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
1086
1087       WRITE (coor_chr,'(F7.1)')  pt_vertical_gradient_level(i)
1088       coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
1089
1090       i = i + 1
1091    ENDDO
1092
1093    WRITE ( io, 420 )  TRIM( coordinates ), TRIM( temperatures ), &
1094                       TRIM( gradients ), TRIM( slices )
1095
1096!
1097!-- Initial humidity profile
1098!-- Building output strings, starting with surface humidity
1099    IF ( humidity  .OR.  passive_scalar )  THEN
1100       WRITE ( temperatures, '(E8.1)' )  q_surface
1101       gradients = '--------'
1102       slices = '       0'
1103       coordinates = '     0.0'
1104       i = 1
1105       DO  WHILE ( q_vertical_gradient_level_ind(i) /= -9999 )
1106         
1107          WRITE (coor_chr,'(E8.1,4X)')  q_init(q_vertical_gradient_level_ind(i))
1108          temperatures = TRIM( temperatures ) // '  ' // TRIM( coor_chr )
1109
1110          WRITE (coor_chr,'(E8.1,4X)')  q_vertical_gradient(i)
1111          gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
1112         
1113          WRITE (coor_chr,'(I8,4X)')  q_vertical_gradient_level_ind(i)
1114          slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
1115         
1116          WRITE (coor_chr,'(F8.1,4X)')  q_vertical_gradient_level(i)
1117          coordinates = TRIM( coordinates ) // '  '  // TRIM( coor_chr )
1118
1119          i = i + 1
1120       ENDDO
1121
1122       IF ( humidity )  THEN
1123          WRITE ( io, 421 )  TRIM( coordinates ), TRIM( temperatures ), &
1124                             TRIM( gradients ), TRIM( slices )
1125       ELSE
1126          WRITE ( io, 422 )  TRIM( coordinates ), TRIM( temperatures ), &
1127                             TRIM( gradients ), TRIM( slices )
1128       ENDIF
1129    ENDIF
1130
1131!
1132!-- Initial salinity profile
1133!-- Building output strings, starting with surface salinity
1134    IF ( ocean )  THEN
1135       WRITE ( temperatures, '(F6.2)' )  sa_surface
1136       gradients = '------'
1137       slices = '     0'
1138       coordinates = '   0.0'
1139       i = 1
1140       DO  WHILE ( sa_vertical_gradient_level_ind(i) /= -9999 )
1141
1142          WRITE (coor_chr,'(F7.2)')  sa_init(sa_vertical_gradient_level_ind(i))
1143          temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr )
1144
1145          WRITE (coor_chr,'(F7.2)')  sa_vertical_gradient(i)
1146          gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
1147
1148          WRITE (coor_chr,'(I7)')  sa_vertical_gradient_level_ind(i)
1149          slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
1150
1151          WRITE (coor_chr,'(F7.1)')  sa_vertical_gradient_level(i)
1152          coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
1153
1154          i = i + 1
1155       ENDDO
1156
1157       WRITE ( io, 425 )  TRIM( coordinates ), TRIM( temperatures ), &
1158                          TRIM( gradients ), TRIM( slices )
1159    ENDIF
1160
1161!
1162!-- LES / turbulence parameters
1163    WRITE ( io, 450 )
1164
1165!--
1166! ... LES-constants used must still be added here
1167!--
1168    IF ( constant_diffusion )  THEN
1169       WRITE ( io, 451 )  km_constant, km_constant/prandtl_number, &
1170                          prandtl_number
1171    ENDIF
1172    IF ( .NOT. constant_diffusion)  THEN
1173       IF ( e_init > 0.0 )  WRITE ( io, 455 )  e_init
1174       IF ( e_min > 0.0 )  WRITE ( io, 454 )  e_min
1175       IF ( wall_adjustment )  WRITE ( io, 453 )  wall_adjustment_factor
1176       IF ( adjust_mixing_length  .AND.  prandtl_layer )  WRITE ( io, 452 )
1177    ENDIF
1178
1179!
1180!-- Special actions during the run
1181    WRITE ( io, 470 )
1182    IF ( create_disturbances )  THEN
1183       WRITE ( io, 471 )  dt_disturb, disturbance_amplitude,                   &
1184                          zu(disturbance_level_ind_b), disturbance_level_ind_b,&
1185                          zu(disturbance_level_ind_t), disturbance_level_ind_t
1186       IF ( bc_lr /= 'cyclic'  .OR.  bc_ns /= 'cyclic' )  THEN
1187          WRITE ( io, 472 )  inflow_disturbance_begin, inflow_disturbance_end
1188       ELSE
1189          WRITE ( io, 473 )  disturbance_energy_limit
1190       ENDIF
1191       WRITE ( io, 474 )  TRIM( random_generator )
1192    ENDIF
1193    IF ( pt_surface_initial_change /= 0.0 )  THEN
1194       WRITE ( io, 475 )  pt_surface_initial_change
1195    ENDIF
1196    IF ( humidity  .AND.  q_surface_initial_change /= 0.0 )  THEN
1197       WRITE ( io, 476 )  q_surface_initial_change       
1198    ENDIF
1199    IF ( passive_scalar  .AND.  q_surface_initial_change /= 0.0 )  THEN
1200       WRITE ( io, 477 )  q_surface_initial_change       
1201    ENDIF
1202
1203    IF ( particle_advection )  THEN
1204!
1205!--    Particle attributes
1206       WRITE ( io, 480 )  particle_advection_start, dt_prel, bc_par_lr, &
1207                          bc_par_ns, bc_par_b, bc_par_t, particle_maximum_age, &
1208                          end_time_prel, dt_sort_particles
1209       IF ( use_sgs_for_particles )  WRITE ( io, 488 )  dt_min_part
1210       IF ( random_start_position )  WRITE ( io, 481 )
1211       IF ( particles_per_point > 1 )  WRITE ( io, 489 )  particles_per_point
1212       WRITE ( io, 495 )  total_number_of_particles
1213       IF ( .NOT. vertical_particle_advection )  WRITE ( io, 482 )
1214       IF ( maximum_number_of_tailpoints /= 0 )  THEN
1215          WRITE ( io, 483 )  maximum_number_of_tailpoints
1216          IF ( minimum_tailpoint_distance /= 0 )  THEN
1217             WRITE ( io, 484 )  total_number_of_tails,      &
1218                                minimum_tailpoint_distance, &
1219                                maximum_tailpoint_age
1220          ENDIF
1221       ENDIF
1222       IF ( dt_write_particle_data /= 9999999.9 )  THEN
1223          WRITE ( io, 485 )  dt_write_particle_data
1224          output_format = ''
1225          IF ( netcdf_output )  THEN
1226             IF ( netcdf_64bit )  THEN
1227                output_format = 'netcdf (64 bit offset) and binary'
1228             ELSE
1229                output_format = 'netcdf and binary'
1230             ENDIF
1231          ELSE
1232             output_format = 'binary'
1233          ENDIF
1234          WRITE ( io, 345 )  output_format
1235       ENDIF
1236       IF ( dt_dopts /= 9999999.9 )  WRITE ( io, 494 )  dt_dopts
1237       IF ( write_particle_statistics )  WRITE ( io, 486 )
1238
1239       WRITE ( io, 487 )  number_of_particle_groups
1240
1241       DO  i = 1, number_of_particle_groups
1242          IF ( i == 1  .AND.  density_ratio(i) == 9999999.9 )  THEN
1243             WRITE ( io, 490 )  i, 0.0
1244             WRITE ( io, 492 )
1245          ELSE
1246             WRITE ( io, 490 )  i, radius(i)
1247             IF ( density_ratio(i) /= 0.0 )  THEN
1248                WRITE ( io, 491 )  density_ratio(i)
1249             ELSE
1250                WRITE ( io, 492 )
1251             ENDIF
1252          ENDIF
1253          WRITE ( io, 493 )  psl(i), psr(i), pss(i), psn(i), psb(i), pst(i), &
1254                             pdx(i), pdy(i), pdz(i)
1255       ENDDO
1256
1257    ENDIF
1258
1259
1260!
1261!-- Parameters of 1D-model
1262    IF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
1263       WRITE ( io, 500 )  end_time_1d, dt_run_control_1d, dt_pr_1d, &
1264                          mixing_length_1d, dissipation_1d
1265       IF ( damp_level_ind_1d /= nzt+1 )  THEN
1266          WRITE ( io, 502 )  zu(damp_level_ind_1d), damp_level_ind_1d
1267       ENDIF
1268    ENDIF
1269
1270!
1271!-- User-defined informations
1272    CALL user_header( io )
1273
1274    WRITE ( io, 99 )
1275
1276!
1277!-- Write buffer contents to disc immediately
1278    CALL local_flush( io )
1279
1280!
1281!-- Here the FORMATs start
1282
1283 99 FORMAT (1X,78('-'))
1284100 FORMAT (/1X,'***************************',9X,42('-')/        &
1285            1X,'* ',A,' *',9X,A/                               &
1286            1X,'***************************',9X,42('-'))
1287101 FORMAT (37X,'coupled run: ',A/ &
1288            37X,42('-'))
1289102 FORMAT (/' Date:              ',A8,9X,'Run:       ',A20/      &
1290            ' Time:              ',A8,9X,'Run-No.:   ',I2.2/     &
1291            ' Run on host:     ',A10)
1292#if defined( __parallel )
1293103 FORMAT (' Number of PEs:',8X,I5,9X,'Processor grid (x,y): (',I3,',',I3, &
1294              ')',1X,A)
1295104 FORMAT (' Number of PEs:',8X,I5,9X,'Tasks:',I4,'   threads per task:',I4/ &
1296              37X,'Processor grid (x,y): (',I3,',',I3,')',1X,A)
1297105 FORMAT (37X,'One additional PE is used to handle'/37X,'the dvrp output!')
1298106 FORMAT (37X,'A 1d-decomposition along x is forced'/ &
1299            37X,'because the job is running on an SMP-cluster')
1300107 FORMAT (37X,'A 1d-decomposition along ',A,' is used')
1301#endif
1302110 FORMAT (/' Numerical Schemes:'/ &
1303             ' -----------------'/)
1304111 FORMAT (' --> Solve perturbation pressure via FFT using ',A,' routines')
1305112 FORMAT (' --> Solve perturbation pressure via SOR-Red/Black-Schema'/ &
1306            '     Iterations (initial/other): ',I3,'/',I3,'  omega = ',F5.3)
1307113 FORMAT (' --> Momentum advection via Piascek-Williams-Scheme (Form C3)', &
1308                  ' or Upstream')
1309114 FORMAT (' --> Momentum advection via Upstream-Spline-Scheme')
1310115 FORMAT ('     Tendencies are smoothed via Long-Filter with factor ',F5.3) 
1311116 FORMAT (' --> Scalar advection via Piascek-Williams-Scheme (Form C3)', &
1312                  ' or Upstream')
1313117 FORMAT (' --> Scalar advection via Upstream-Spline-Scheme')
1314118 FORMAT (' --> Scalar advection via Bott-Chlond-Scheme')
1315119 FORMAT (' --> Galilei-Transform applied to horizontal advection', &
1316            '     Translation velocity = ',A/ &
1317            '     distance advected ',A,':  ',F8.3,' km(x)  ',F8.3,' km(y)')
1318120 FORMAT (' --> Time differencing scheme: leapfrog only (no euler in case', &
1319                  ' of timestep changes)')
1320121 FORMAT (' --> Time differencing scheme: leapfrog + euler in case of', &
1321                  ' timestep changes')
1322122 FORMAT (' --> Time differencing scheme: ',A)
1323123 FORMAT (' --> Rayleigh-Damping active, starts ',A,' z = ',F8.2,' m'/ &
1324            '     maximum damping coefficient: ',F5.3, ' 1/s')
1325124 FORMAT ('     Spline-overshoots are being suppressed')
1326125 FORMAT ('     Upstream-Scheme is used if Upstream-differences fall short', &
1327                  ' of'/                                                       &
1328            '     delta_u = ',F6.4,4X,'delta_v = ',F6.4,4X,'delta_w = ',F6.4)
1329126 FORMAT ('     Upstream-Scheme is used if Upstream-differences fall short', &
1330                  ' of'/                                                       &
1331            '     delta_e = ',F6.4,4X,'delta_pt = ',F6.4)
1332127 FORMAT ('     The following absolute overshoot differences are tolerated:'/&
1333            '     delta_u = ',F6.4,4X,'delta_v = ',F6.4,4X,'delta_w = ',F6.4)
1334128 FORMAT ('     The following absolute overshoot differences are tolerated:'/&
1335            '     delta_e = ',F6.4,4X,'delta_pt = ',F6.4)
1336129 FORMAT (' --> Additional prognostic equation for the specific humidity')
1337130 FORMAT (' --> Additional prognostic equation for the total water content')
1338131 FORMAT (' --> Parameterization of condensation processes via (0%-or100%)')
1339132 FORMAT (' --> Parameterization of long-wave radiation processes via'/ &
1340            '     effective emissivity scheme')
1341133 FORMAT (' --> Precipitation parameterization via Kessler-Scheme')
1342134 FORMAT (' --> Additional prognostic equation for a passive scalar')
1343135 FORMAT (' --> Solve perturbation pressure via multigrid method (', &
1344                  A,'-cycle)'/ &
1345            '     number of grid levels:                   ',I2/ &
1346            '     Gauss-Seidel red/black iterations:       ',I2)
1347136 FORMAT ('     gridpoints of coarsest subdomain (x,y,z): (',I3,',',I3,',', &
1348                  I3,')')
1349137 FORMAT ('     level data gathered on PE0 at level:     ',I2/ &
1350            '     gridpoints of coarsest subdomain (x,y,z): (',I3,',',I3,',', &
1351                  I3,')'/ &
1352            '     gridpoints of coarsest domain (x,y,z):    (',I3,',',I3,',', &
1353                  I3,')')
1354138 FORMAT ('     Using hybrid version for 1d-domain-decomposition')
1355139 FORMAT (' --> Loop optimization method: ',A)
1356140 FORMAT ('     maximum residual allowed:                ',E10.3)
1357141 FORMAT ('     fixed number of multigrid cycles:        ',I4)
1358142 FORMAT ('     perturbation pressure is calculated at every Runge-Kutta ', &
1359                  'step')
1360143 FORMAT ('     Euler/upstream scheme is used for the SGS turbulent ', &
1361                  'kinetic energy')
1362150 FORMAT (' --> Volume flow at the right and north boundary will be ', &
1363                  'conserved')
1364200 FORMAT (//' Run time and time step information:'/ &
1365             ' ----------------------------------'/)
1366201 FORMAT ( ' Timestep:          variable     maximum value: ',F6.3,' s', &
1367             '    CFL-factor: ',F4.2)
1368202 FORMAT ( ' Timestep:       dt = ',F6.3,' s'/)
1369203 FORMAT ( ' Start time:       ',F9.3,' s'/ &
1370             ' End time:         ',F9.3,' s')
1371204 FORMAT ( A,F9.3,' s')
1372205 FORMAT ( A,F9.3,' s',5X,'restart every',17X,F9.3,' s')
1373206 FORMAT (/' Time reached:     ',F9.3,' s'/ &
1374             ' CPU-time used:    ',F9.3,' s     per timestep:               ', &
1375               '  ',F9.3,' s'/                                                 &
1376             '                                   per second of simulated tim', &
1377               'e: ',F9.3,' s')
1378250 FORMAT (//' Computational grid and domain size:'/ &
1379              ' ----------------------------------'// &
1380              ' Grid length:      dx =    ',F7.3,' m    dy =    ',F7.3, &
1381              ' m    dz =    ',F7.3,' m'/ &
1382              ' Domain size:       x = ',F10.3,' m     y = ',F10.3, &
1383              ' m  z(u) = ',F10.3,' m'/)
1384252 FORMAT (' dz constant up to ',F10.3,' m (k=',I4,'), then stretched by', &
1385              ' factor: ',F5.3/ &
1386            ' maximum dz not to be exceeded is dz_max = ',F10.3,' m'/)
1387254 FORMAT (' Number of gridpoints (x,y,z):  (0:',I4,', 0:',I4,', 0:',I4,')'/ &
1388            ' Subdomain size (x,y,z):        (  ',I4,',   ',I4,',   ',I4,')'/)
1389255 FORMAT (' Subdomains have equal size')
1390256 FORMAT (' Subdomains at the upper edges of the virtual processor grid ', &
1391              'have smaller sizes'/                                          &
1392            ' Size of smallest subdomain:    (  ',I4,',   ',I4,',   ',I4,')')
1393260 FORMAT (/' The model has a slope in x-direction. Inclination angle: ',F6.2,&
1394             ' degrees')
1395270 FORMAT (//' Topography informations:'/ &
1396              ' -----------------------'// &
1397              1X,'Topography: ',A)
1398271 FORMAT (  ' Building size (x/y/z) in m: ',F5.1,' / ',F5.1,' / ',F5.1/ &
1399              ' Horizontal index bounds (l/r/s/n): ',I4,' / ',I4,' / ',I4, &
1400                ' / ',I4)
1401280 FORMAT (//' Vegetation canopy (drag) model:'/ &
1402              ' ------------------------------'// &
1403              ' Canopy mode: ', A / &
1404              ' Canopy top: ',I4 / &
1405              ' Leaf drag coefficient: ',F6.2 /)
1406281 FORMAT (/ ' Scalar_exchange_coefficient: ',F6.2 / &
1407              ' Scalar concentration at leaf surfaces in kg/m**3: ',F6.2 /)
1408282 FORMAT (' Predefined constant heatflux at the top of the vegetation: ',F6.2,' K m/s')
1409283 FORMAT (/ ' Characteristic levels of the leaf area density:'// &
1410              ' Height:              ',A,'  m'/ &
1411              ' Leaf area density:   ',A,'  m**2/m**3'/ &
1412              ' Gradient:            ',A,'  m**2/m**4'/ &
1413              ' Gridpoint:           ',A)
1414               
1415300 FORMAT (//' Boundary conditions:'/ &
1416             ' -------------------'// &
1417             '                     p                    uv             ', &
1418             '                   pt'// &
1419             ' B. bound.: ',A/ &
1420             ' T. bound.: ',A)
1421301 FORMAT (/'                     ',A// &
1422             ' B. bound.: ',A/ &
1423             ' T. bound.: ',A)
1424303 FORMAT (/' Bottom surface fluxes are used in diffusion terms at k=1')
1425304 FORMAT (/' Top surface fluxes are used in diffusion terms at k=nzt')
1426305 FORMAT (//'    Prandtl-Layer between bottom surface and first ', &
1427               'computational u,v-level:'// &
1428             '       zp = ',F6.2,' m   z0 = ',F6.4,' m   kappa = ',F4.2/ &
1429             '       Rif value range:   ',F6.2,' <= rif <=',F6.2)
1430306 FORMAT ('       Predefined constant heatflux:   ',F9.6,' K m/s')
1431307 FORMAT ('       Heatflux has a random normal distribution')
1432308 FORMAT ('       Predefined surface temperature')
1433309 FORMAT ('       Predefined constant salinityflux:   ',F9.6,' psu m/s')
1434310 FORMAT (//'    1D-Model:'// &
1435             '       Rif value range:   ',F6.2,' <= rif <=',F6.2)
1436311 FORMAT ('       Predefined constant humidity flux: ',E10.3,' m/s')
1437312 FORMAT ('       Predefined surface humidity')
1438313 FORMAT ('       Predefined constant scalar flux: ',E10.3,' kg/(m**2 s)')
1439314 FORMAT ('       Predefined scalar value at the surface')
1440315 FORMAT ('       Humidity / scalar flux at top surface is 0.0')
1441316 FORMAT ('       Sensible heatflux and momentum flux from coupled ', &
1442                    'atmosphere model')
1443317 FORMAT (//' Lateral boundaries:'/ &
1444            '       left/right:  ',A/    &
1445            '       north/south: ',A)
1446318 FORMAT (/'       outflow damping layer width: ',I3,' gridpoints with km_', &
1447                    'max =',F5.1,' m**2/s')
1448319 FORMAT ('       turbulence recycling at inflow switched on'/ &
1449            '       width of recycling domain: ',F7.1,' m   grid index: ',I4/ &
1450            '       inflow damping height: ',F6.1,' m   width: ',F6.1,' m')
1451320 FORMAT ('       Predefined constant momentumflux:  u: ',F9.6,' m**2/s**2'/ &
1452            '                                          v: ',F9.6,' m**2/s**2')
1453325 FORMAT (//' List output:'/ &
1454             ' -----------'//  &
1455            '    1D-Profiles:'/    &
1456            '       Output every             ',F8.2,' s')
1457326 FORMAT ('       Time averaged over       ',F8.2,' s'/ &
1458            '       Averaging input every    ',F8.2,' s')
1459330 FORMAT (//' Data output:'/ &
1460             ' -----------'/)
1461331 FORMAT (/'    1D-Profiles:')
1462332 FORMAT (/'       ',A)
1463333 FORMAT ('       Output every             ',F8.2,' s',/ &
1464            '       Time averaged over       ',F8.2,' s'/ &
1465            '       Averaging input every    ',F8.2,' s')
1466334 FORMAT (/'    2D-Arrays',A,':')
1467335 FORMAT (/'       ',A2,'-cross-section  Arrays: ',A/ &
1468            '       Output every             ',F8.2,' s  ',A/ &
1469            '       Cross sections at ',A1,' = ',A/ &
1470            '       scalar-coordinates:   ',A,' m'/)
1471336 FORMAT (/'    3D-Arrays',A,':')
1472337 FORMAT (/'       Arrays: ',A/ &
1473            '       Output every             ',F8.2,' s  ',A/ &
1474            '       Upper output limit at    ',F8.2,' m  (GP ',I4,')'/)
1475338 FORMAT ('       Compressed data output'/ &
1476            '       Decimal precision: ',A/)
1477339 FORMAT ('       No output during initial ',F8.2,' s')
1478340 FORMAT (/'    Time series:')
1479341 FORMAT ('       Output every             ',F8.2,' s'/)
1480342 FORMAT (/'       ',A2,'-cross-section  Arrays: ',A/ &
1481            '       Output every             ',F8.2,' s  ',A/ &
1482            '       Time averaged over       ',F8.2,' s'/ &
1483            '       Averaging input every    ',F8.2,' s'/ &
1484            '       Cross sections at ',A1,' = ',A/ &
1485            '       scalar-coordinates:   ',A,' m'/)
1486343 FORMAT (/'       Arrays: ',A/ &
1487            '       Output every             ',F8.2,' s  ',A/ &
1488            '       Time averaged over       ',F8.2,' s'/ &
1489            '       Averaging input every    ',F8.2,' s'/ &
1490            '       Upper output limit at    ',F8.2,' m  (GP ',I4,')'/)
1491345 FORMAT ('       Output format: ',A/)
1492#if defined( __dvrp_graphics )
1493360 FORMAT ('    Plot-Sequence with dvrp-software:'/ &
1494            '       Output every      ',F7.1,' s'/ &
1495            '       Output mode:      ',A/ &
1496            '       Host / User:      ',A,' / ',A/ &
1497            '       Directory:        ',A// &
1498            '       The sequence contains:')
1499361 FORMAT ('       Isosurface of ',A,'  Threshold value: ', E12.3)
1500362 FORMAT ('       Sectional plane ',A)
1501363 FORMAT ('       Particles')
1502#endif
1503#if defined( __spectra )
1504370 FORMAT ('    Spectra:')
1505371 FORMAT ('       Output every ',F7.1,' s'/)
1506372 FORMAT ('       Arrays:     ', 10(A5,',')/                         &
1507            '       Directions: ', 10(A5,',')/                         &
1508            '       height levels  k = ', 20(I3,',')/                  &
1509            '                          ', 20(I3,',')/                  &
1510            '                          ', 20(I3,',')/                  &
1511            '                          ', 20(I3,',')/                  &
1512            '                          ', 19(I3,','),I3,'.'/           &
1513            '       height levels selected for standard plot:'/        &
1514            '                      k = ', 20(I3,',')/                  &
1515            '                          ', 20(I3,',')/                  &
1516            '                          ', 20(I3,',')/                  &
1517            '                          ', 20(I3,',')/                  &
1518            '                          ', 19(I3,','),I3,'.'/           &
1519            '       Time averaged over ', F7.1, ' s,' /                &
1520            '       Profiles for the time averaging are taken every ', &
1521                    F6.1,' s')
1522#endif
1523400 FORMAT (//' Physical quantities:'/ &
1524              ' -------------------'/)
1525410 FORMAT ('    Angular velocity    :   omega = ',E9.3,' rad/s'/  &
1526            '    Geograph. latitude  :   phi   = ',F4.1,' degr'/   &
1527            '    Coriolis parameter  :   f     = ',F9.6,' 1/s'/    &
1528            '                            f*    = ',F9.6,' 1/s')
1529411 FORMAT (/'    Gravity             :   g     = ',F4.1,' m/s**2')
1530412 FORMAT (/'    Reference density in buoyancy terms: ',F8.3,' kg/m**3')
1531413 FORMAT (/'    Reference temperature in buoyancy terms: ',F8.4,' K')
1532415 FORMAT (/'    Cloud physics parameters:'/ &
1533             '    ------------------------'/)
1534416 FORMAT ('        Surface pressure   :   p_0   = ',F7.2,' hPa'/      &
1535            '        Gas constant       :   R     = ',F5.1,' J/(kg K)'/ &
1536            '        Density of air     :   rho_0 = ',F5.3,' kg/m**3'/  &
1537            '        Specific heat cap. :   c_p   = ',F6.1,' J/(kg K)'/ &
1538            '        Vapourization heat :   L_v   = ',E8.2,' J/kg')
1539420 FORMAT (/'    Characteristic levels of the initial temperature profile:'// &
1540            '       Height:        ',A,'  m'/ &
1541            '       Temperature:   ',A,'  K'/ &
1542            '       Gradient:      ',A,'  K/100m'/ &
1543            '       Gridpoint:     ',A)
1544421 FORMAT (/'    Characteristic levels of the initial humidity profile:'// &
1545            '       Height:      ',A,'  m'/ &
1546            '       Humidity:    ',A,'  kg/kg'/ &
1547            '       Gradient:    ',A,'  (kg/kg)/100m'/ &
1548            '       Gridpoint:   ',A)
1549422 FORMAT (/'    Characteristic levels of the initial scalar profile:'// &
1550            '       Height:                  ',A,'  m'/ &
1551            '       Scalar concentration:    ',A,'  kg/m**3'/ &
1552            '       Gradient:                ',A,'  (kg/m**3)/100m'/ &
1553            '       Gridpoint:               ',A)
1554423 FORMAT (/'    Characteristic levels of the geo. wind component ug:'// &
1555            '       Height:      ',A,'  m'/ &
1556            '       ug:          ',A,'  m/s'/ &
1557            '       Gradient:    ',A,'  1/100s'/ &
1558            '       Gridpoint:   ',A)
1559424 FORMAT (/'    Characteristic levels of the geo. wind component vg:'// &
1560            '       Height:      ',A,'  m'/ &
1561            '       vg:          ',A,'  m/s'/ &
1562            '       Gradient:    ',A,'  1/100s'/ &
1563            '       Gridpoint:   ',A)
1564425 FORMAT (/'    Characteristic levels of the initial salinity profile:'// &
1565            '       Height:     ',A,'  m'/ &
1566            '       Salinity:   ',A,'  psu'/ &
1567            '       Gradient:   ',A,'  psu/100m'/ &
1568            '       Gridpoint:  ',A)
1569450 FORMAT (//' LES / Turbulence quantities:'/ &
1570              ' ---------------------------'/)
1571451 FORMAT ('   Diffusion coefficients are constant:'/ &
1572            '   Km = ',F6.2,' m**2/s   Kh = ',F6.2,' m**2/s   Pr = ',F5.2)
1573452 FORMAT ('   Mixing length is limited to the Prandtl mixing lenth.')
1574453 FORMAT ('   Mixing length is limited to ',F4.2,' * z')
1575454 FORMAT ('   TKE is not allowed to fall below ',E9.2,' (m/s)**2')
1576455 FORMAT ('   initial TKE is prescribed as ',E9.2,' (m/s)**2')
1577470 FORMAT (//' Actions during the simulation:'/ &
1578              ' -----------------------------'/)
1579471 FORMAT ('    Disturbance impulse (u,v) every :   ',F6.2,' s'/            &
1580            '    Disturbance amplitude           :     ',F4.2, ' m/s'/       &
1581            '    Lower disturbance level         : ',F8.2,' m (GP ',I4,')'/  &
1582            '    Upper disturbance level         : ',F8.2,' m (GP ',I4,')')
1583472 FORMAT ('    Disturbances continued during the run from i/j =',I4, &
1584                 ' to i/j =',I4)
1585473 FORMAT ('    Disturbances cease as soon as the disturbance energy exceeds',&
1586                 1X,F5.3, ' m**2/s**2')
1587474 FORMAT ('    Random number generator used    : ',A/)
1588475 FORMAT ('    The surface temperature is increased (or decreased, ', &
1589                 'respectively, if'/ &
1590            '    the value is negative) by ',F5.2,' K at the beginning of the',&
1591                 ' 3D-simulation'/)
1592476 FORMAT ('    The surface humidity is increased (or decreased, ',&
1593                 'respectively, if the'/ &
1594            '    value is negative) by ',E8.1,' kg/kg at the beginning of', &
1595                 ' the 3D-simulation'/)
1596477 FORMAT ('    The scalar value is increased at the surface (or decreased, ',&
1597                 'respectively, if the'/ &
1598            '    value is negative) by ',E8.1,' kg/m**3 at the beginning of', &
1599                 ' the 3D-simulation'/)
1600480 FORMAT ('    Particles:'/ &
1601            '    ---------'// &
1602            '       Particle advection is active (switched on at t = ', F7.1, &
1603                    ' s)'/ &
1604            '       Start of new particle generations every  ',F6.1,' s'/ &
1605            '       Boundary conditions: left/right: ', A, ' north/south: ', A/&
1606            '                            bottom:     ', A, ' top:         ', A/&
1607            '       Maximum particle age:                 ',F9.1,' s'/ &
1608            '       Advection stopped at t = ',F9.1,' s'/ &
1609            '       Particles are sorted every ',F9.1,' s'/)
1610481 FORMAT ('       Particles have random start positions'/)
1611482 FORMAT ('       Particles are advected only horizontally'/)
1612483 FORMAT ('       Particles have tails with a maximum of ',I3,' points')
1613484 FORMAT ('            Number of tails of the total domain: ',I10/ &
1614            '            Minimum distance between tailpoints: ',F8.2,' m'/ &
1615            '            Maximum age of the end of the tail:  ',F8.2,' s')
1616485 FORMAT ('       Particle data are written on file every ', F9.1, ' s')
1617486 FORMAT ('       Particle statistics are written on file'/)
1618487 FORMAT ('       Number of particle groups: ',I2/)
1619488 FORMAT ('       SGS velocity components are used for particle advection'/ &
1620            '          minimum timestep for advection: ', F7.5/)
1621489 FORMAT ('       Number of particles simultaneously released at each ', &
1622                    'point: ', I5/)
1623490 FORMAT ('       Particle group ',I2,':'/ &
1624            '          Particle radius: ',E10.3, 'm')
1625491 FORMAT ('          Particle inertia is activated'/ &
1626            '             density_ratio (rho_fluid/rho_particle) = ',F5.3/)
1627492 FORMAT ('          Particles are advected only passively (no inertia)'/)
1628493 FORMAT ('          Boundaries of particle source: x:',F8.1,' - ',F8.1,' m'/&
1629            '                                         y:',F8.1,' - ',F8.1,' m'/&
1630            '                                         z:',F8.1,' - ',F8.1,' m'/&
1631            '          Particle distances:  dx = ',F8.1,' m  dy = ',F8.1, &
1632                       ' m  dz = ',F8.1,' m'/)
1633494 FORMAT ('       Output of particle time series in NetCDF format every ', &
1634                    F8.2,' s'/)
1635495 FORMAT ('       Number of particles in total domain: ',I10/)
1636500 FORMAT (//' 1D-Model parameters:'/                           &
1637              ' -------------------'//                           &
1638            '    Simulation time:                   ',F8.1,' s'/ &
1639            '    Run-controll output every:         ',F8.1,' s'/ &
1640            '    Vertical profile output every:     ',F8.1,' s'/ &
1641            '    Mixing length calculation:         ',A/         &
1642            '    Dissipation calculation:           ',A/)
1643502 FORMAT ('    Damping layer starts from ',F7.1,' m (GP ',I4,')'/)
1644
1645
1646 END SUBROUTINE header
Note: See TracBrowser for help on using the repository browser.