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

Last change on this file since 121 was 117, checked in by raasch, 17 years ago

further small changes added to last revision

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